matlab实验指导书.docx
- 文档编号:18188684
- 上传时间:2023-04-24
- 格式:DOCX
- 页数:123
- 大小:2.09MB
matlab实验指导书.docx
《matlab实验指导书.docx》由会员分享,可在线阅读,更多相关《matlab实验指导书.docx(123页珍藏版)》请在冰豆网上搜索。
matlab实验指导书
实验一MATLAB基础训练
一、实验目的
本次上机实验主要练习使用Matlab的基本操作和基础知识,包括数组(复数、向量、矩阵、结构体数组等)的创建和数组元素的操作和运算、矩阵的运算、Matlab的运算符(尤其是点运算‘.’)、脚本M文件和函数M文件的编写、Matlab文件的编程(基本的流程控制结构)、基本的二维和三维绘图方法以及图形的标注等。
希望通过本次实验使大家尽量在短时间内(4学时)掌握Matlab的基本操作和基础知识,为后面的实验项目奠定基础。
二、实验原理
参见PPT中有关内容。
三、实验内容
1.上机练习课件中的例子。
2.设两个复数a=1+2i,b=3-4i,计算a+b,a-b,a*b,a/b,a和b的模。
3.计算下式的结果,其中x=-3.5°,y=6.7°
(提示:
①应将角度单位由度转换为数学函数所能处理的弧度值;②求根函数sqrt,取绝对值函数abs,具体用法用help查询)
4.对矩阵
实现下列操作:
(1)左右翻转(fliplr命令)
(2)上下翻转(flipud命令)
(3)利用cat命令分别将A扩展成3×6和6×3的矩阵
(3)分别提取A的第2行,第2列,对角线元素
(4)删除A的第2行2列的元素
(提示:
将矩阵元素赋空阵[]可以删除元素,注意此时元素的访问只能使用单下标的方式。
观察删除元素后,A中元素的排列方式的变化)
5.创建[0,2]区间上拥有100个等间隔元素的列向量x(Matlab默认是行向量),并绘制y=sin(x1/3)的函数图像。
6.创建如下图所示的单结构体数组。
7.编程训练:
下图所示电路中,R1=2,R2=4,R3=12,R4=4,R5=12,R6=4,R7=2,us=10V,求i3。
(要求:
以脚本M文件方式建立程序。
此项练习M脚本文件的建立以及使用Matlab解线性方程组。
)
(提示:
先建立电路方程如下
该线性方程组等价于矩阵形式:
其中,A是包含R1~R7的系数矩阵,B=(us,0,0)T是常数矩阵,X=(ia,ib,ic)T=A-1B=A\B为解向量。
因此使用矩阵除法(注意区别左除和右除)可以迅速解出ia,ib,ic。
最后i3=ia-ib=0.3704。
)
8.编程训练:
已知某电路的电流i(t)的表达式为
画出[-,]区间上的电流波形。
(要求:
①以脚本M文件方式建立程序;②用12号字体给图形加上标注‘t’和‘i(t)’,红色曲线+圆形标记‘o’,如下图所示。
此项练习二维绘图方法、图形的标注和循环控制结构)
9.编程训练:
已知
的取值范围是
,若a=﹣8,b=8,画出
所表示的三维曲面。
(要求:
①以脚本M文件方式建立程序;②使用surf和mesh函数分别绘出曲面图和网格图,参考图形如下。
此项练习三维绘图方法)
10.以函数M文件的方式重做第9个项目,其中函数的输入参数为a和b。
(此项练习M函数文件的建立)
11.建立如下图所示的用户界面菜单。
要求:
(A)把用户菜单'Option'设置为顶层的第3菜单项;(B)下拉菜单被两条分隔线分为三个菜单区;(C)最下菜单项又有两个子菜单组成。
(此项练习用户菜单的建立,写出M文件)
四、实验要求
以上训练项目要求在4个实验学时内完成,并提交项目7-11的程序。
五、思考题
1.举例说明如何创建三维数组?
2.第9个训练项目中,函数z的表达式中,分子可能被0除,导致绘图数据点中出现非数NaN,绘出的图形不正确。
如何解决这一问题?
3.如何将第9个项目中的两幅图放在同一个图形窗口中显示?
(提示:
使用subplot命令将图形窗口分割成1×2的子窗口。
subplot(m,n,i)将图形窗口分割成m×n的子窗口,并指定第i个子窗口为图形的绘制区域)
4.surf函数绘出的三维图一般有黑色格线,如何消掉黑色格线?
以第9个项目为例说明。
实验二 MATLAB数值计算:
二阶电路的时域分析
一、实验目的
在物理学和工程技术上,很多问题都可以用一个或一组常微分方程来描述,因此要解决相应的实际问题往往需要首先求解对应的微分方程(组)。
在大多数情况下这些微分方程(组)通常是非线性的或者是超越方程(比如范德堡方程,波导本征值方程等),很难解析地求解(精确解),因此往往需要使用计算机数值求解(近似解)。
MATLAB作为一种强大的科学计算语言,其在数值计算和数据的可视化方面具有无以伦比的优势。
在解决常微分方程(组)问题上,MATLAB就提供了多种可适用于不同场合(如刚性和非刚性问题)下的求解器(Solver),例如ode45,ode15s,ode23,ode23s等等。
本次实验将以二阶线性电路-RLC电路和二阶非线性电路-范德堡电路的时域计算为例,了解和学习使用MATLAB作为计算工具来解算复杂的微分方程,以期达到如下几个目的:
1.熟练使用dsolve函数解析求解常微分方程;
2.熟练运用ode45求解器数值求解常微分方程;
3.了解状态方程的概念,能使用MATLAB对二阶电路进行计算和分析;
二、实验预备知识
1.微分方程的概念
未知的函数以及它的某些阶的导数连同自变量都由一已知方程联系在一起的方程称为微分方程。
如果未知函数是一元函数,称为常微分方程(Ordinarydifferentialequations,简称odes)。
n阶常微分方程的一般形式(隐式)为:
(1)
其中t为自变量。
若方程中未知函数及其各阶导数都是一次的,称为线性常微分方程,否则就是非线性微分方程,例如方程
就是非线性的。
2.常微分方程的解及MATLAB指令
一阶常微分方程与高阶微分方程可以互化,已知一个n阶常微分方程(显式):
(2)
若令
,可将上式化为n个一阶常微分方程组:
(3)
(3)式称为状态方程,y1,y2,…,yn(即y,y,y,…,y(n-1))称为状态变量,其中y1(即y)就是常微分方程
(2)式的解。
(3)式中右边的函数f1、f2、…、fn代表各个状态变量的一阶导数的函数表达式,对于具体的方程它们有具体的形式,例如下列二阶非线性微分方程:
若令
,可将其改写成2个一阶微分方程组(状态方程)的形式:
因此
。
●解析解
只有少部分的线性常微分方程可以解析地求解(即可以算出精确的解表达式),例如一阶常系数常微分方程
可以通过直接积分解出,而多数微分方程尤其是非线性方程则很难得到解析解。
有解析解的方程虽然可以手算解出,但是MATLAB也提供了dslove指令来求方程的解析解,其使用格式:
S=dsolve(‘方程1’,‘方程2’,…,’初始条件1’,’初始条件2’…,’自变量’)
方程用字符串表示,自变量缺省值为t。
1阶导数用D表示,2阶导数用D2表示,以此类推。
S用于返回方程解析解的表达式。
如果是求解方程组,则S为一个结构体数组,它的每个域存放方程组每一个解的表达式。
例1:
求下列微分方程的解析解
>>s=dsolve(’D2y=sin(2*x)-y’,’y(0)=0’,’Dy(0)=1’,’x’);
>>simplify(s) %以最简形式显示s
ans=
-1/3*sin(x)*(-5+2*cos(x))%方程的解(符号表达式)
●数值解
对于没有解析解的方程主要依靠计算机进行数值求解(得到的是近似解),例如方程
就须通过计算机数值求解(结果是一系列解的数值而非表达式)。
考虑n阶微分方程
(2)式的数值求解,它等价于一阶常微分方程组(3)式。
现将(3)式写成矩阵形式:
其中
为n个分量的列向量(Columnvector),也称状态向量
n个分量的列向量,其每个元素分别为(3)式右边的函数表达式
我们知道,微分方程要有唯一的确定解,必须给定初值条件。
因此方程(4)式要有确定的解必须给定初值条件(t0为初始时刻):
(7)
Matlab提供了ode45指令(ode是常微分方程的英文缩写)来求解方程(4)的数值解(近似解)。
基本使用格式:
[tout,Yout]=ode45(odefun,tspan,Y0,options)
其参数说明如下:
①odefun一般是用M文件编写的函数,odefun代表函数名,由用户自己定义。
函数返回值为(4)式右边的F(t,Y)=(f1,f2,…,fn)T。
故odefun函数的返回值应是列向量,其最简单的编写格式为:
functionF=odefun(t,Y)【作用是计算并返回4式中的F(t,Y)】
其中t时间变量,为标量,代表计算进程中的某时刻点
Y代表状态变量的列向量(即5式)
F返回值F(t,Y),为列向量(见6式)
ode45求解指令在计算时将会不断地在各个时间点调用odefun函数,并自动给输入参数t和Y赋值。
②tspan指定方程的求解区间[t0,tf],t0是初始时刻。
③Y0用户给定的初值条件,为n个分量的列向量,见(7)式。
④options可选项。
一般情况下可缺省即可,若用户有特殊要求则须使用odeset指令设置options选项,具体用法可使用helpodeset命令查询,此处不做要求。
⑤tout列向量,输出求解过程中区间[t0,tf]上各个计算点的时刻,即
[t0,tf]上计算点的数目是由Matlab自动生成。
⑥Yout输出矩阵,其排列格式如下:
Yout的第1列代表的是状态变量y1在
(依次为tout的每个元素)各个时刻的值,由于y1=y,所以Yout的第一列就是待求方程的数值解,它显然是一系列离散的y(t)值:
;将Yout第一列提取出来并用plot(tout,Yout(:
1))指令即可绘制解y(t)的函数图像。
Yout的第2列是y的一阶导数
的数值解,第3列是y的二阶导数
的数值解…。
三、实验内容
1.二阶线性电路—RLC回路的零输入响应
当电路中含有二个动态元件(如电感、电容)时,建立的电路方程为二阶微分方程,这样的电路称为二阶电路,如果微分方程是线性的,则为线性电路,若为非线性方程,则是非线性电路,如范德堡电路。
所谓零输入响应指的是电路中无外加的激励源,仅由动态元件初始储能所产生的响应。
考虑如图所示的RLC电路,假设电容原已充电,根据电路理论,此二阶电路的零输入响应可用如下二阶线性微分方程描述:
其中uc代表电容电压。
给定初值条件:
元件参数L=0.5H,R=12.5,C=0.02F。
要得到t≥0时的零输入响应,就必须求解(9)式。
A.解析解
方程(9)可以直接求解,因此有解析解。
下面改用MATLAB中的dsolve指令来求方程(9)的解析解,程序如下:
S=dsolve('D2u=-R/L*Du-1/L/C*u','u(0)=1','Du(0)=0','t');
%S为字符型数组(字串),其值为方程9的解表达式
L=0.5;C=0.02;R=12.5;%元件参数
t=0:
0.01:
1;%定义区间[0,1]上的时间序列
y=eval(S);%eval串演算函数,计算字串S(即9式的解)在t时刻的值
plot(t,y)%绘制电压波形,即uc(t)的零输入响应
B.数值解
方程(9)虽然可以解析求解,但也可以使用ode45指令来近似计算。
令y1=uc,y2=duc/dt,则方程9改写成如下状态方程:
上式也可以写成(4)式的形式,因此状态向量
Y=(y1,y2)T,F=(f1,f2)T=(y2,﹣1/L/C*y1﹣R/L*y2)T(12)
第一步:
通过M文件创建ode函数,函数名circuit_2order_odefun。
functionF=circuit_2order_odefun(t,Y)
globalLCR%定义全局变量L、C、R,以实现参数在MATLAB的基本工作空%间和函数的专用空间之间数据的传递
F=[Y
(2);-1/L/C*Y
(1)-R/L*Y
(2)];%函数返回值F,列向量,见(12)式
注意,函数创建完后,必须存盘,存储文件名和函数名应一致!
第二步:
创建M脚本文件circuit_2order.m
globalLCR
L=0.5;C=0.02;R=12.5;%元件参数
Y0=[1;0];%列向量,初值条件
tspan=[0,1];%定义求解区间[0,1]
[tout,Yout]=ode45('circuit_2order_odefun',tspan,Y0);
plot(tout,Yout(:
1));%绘电压波形,即uc(t)的零输入响应
2.二阶非线性电路—范德堡(VandePol)电路
范德堡电路由一个线性电感、一个线性电容和一个非线性电阻构成,如图(a)所示。
非线性电阻的伏安特性如图(b)所示,可很明显看出电阻是非线性的。
范德堡电路的特性可由一个二阶非线性微分方程描述(参见邱关源的《电路》第5版P463-P465):
其中
。
上式就是著名的范德堡方程。
该方程很难解析求解,必须借助计算机数值求解。
令
,将范德堡方程改写成一阶微分方程组(状态方程)的形式:
上式也可以写成(4)式的形式。
因此状态向量Y和F分别为
,
(15)
假定参数=0.1,求解区间[0,100],以及初值条件
请仿照实验内容1,利用ode45指令求解范德堡方程,并绘制电流iL的波形(即电流的零输入响应图)。
四、实验任务
1.输入实验内容1中提供的程序上机练习,熟悉dslove指令以及ode45指令的使用方法,并尝试解答思考题1。
2.在任务1的基础上,独立编制实验内容2的计算程序,并尝试解答思考题3。
五、思考题
1.在实验内容1中,若使用ode45指令解方程9式,如何得到电流的零输入响应波形?
2.在实验内容1中,试研究R=1,2,3,…,10时电压Uc(t)的响应波形,并把它们绘制在一张图上。
3.从任务2的电流波形图(如下图)中,你能得到什么结论?
如何绘制范德堡电路的uc~iL函数图像(称为相图)?
(提示:
)
附:
实验任务2的参考图
实验三 牛顿环实验的MATLAB仿真
一、实验目的
MATLAB在光学实验的计算机仿真方面有着重要的应用。
使用MATLAB可以仿真大多数光学实验,例如杨氏双缝干涉实验、牛顿环实验、夫琅和费衍射等,使得原本抽象的必须借助实验仪器才能感知的光学现象可以直观而且动态的显示在计算机上,从而获得对某一特定光学现象充分的感性认知,加深我们对抽象光学现象的理解和认识。
这里将以牛顿环实验的MATLAB仿真为例,向大家介绍光学实验的MATLAB动态仿真的基本方法并通过上机练习以达到如下几个目的:
1.掌握用imshow实现光强度二维分布的可视化显示方法;
2.掌握MATLAB动态仿真技术—影片动画技术;
二、实验原理
1.牛顿环干涉原理
右图所示为牛顿环装置的示意图。
R为牛顿环透镜的曲率半径,d为空气膜的厚度(
)。
垂直入射光经空气膜的上下两表面反射后产生干涉,干涉后的光强
其中I1和I2是两束相干光的光强,可近似认为
I1=I2=I0。
为两束光相遇时的位相差
由图中的几何关系以及
条件可得
(2)、(3)式代入
(1)式后有
上式中为了方便取系数2I0=1。
在直角坐标系中
,(x,y)代表光强的二维分布点的坐标。
(4)式是实验仿真的基础,对于任意给定点(x,y),如果该点的光强I取最大值1,则该点为明条纹所在;若光强I取最小值0,则是暗纹所在;其他值则介于两者之间。
(4)式给出了牛顿环干涉光强的二维平面分布,那么如何将光强的平面分布(数值)可视化显示出来?
下面介绍一种简单的实现方法。
2.光强分布的可视化实现
对于数据的可视化,MATLAB提供了很多实现方法,比如前面介绍的plot、plot3、surf函数等。
对于牛顿环实验来说,虽可以使用surf函数将光强I在xy平面上的分布表现出来,但是得到的是3维曲面图,和实验观察到的2维环状的干涉条纹图形不一致,达不到仿真的目的。
所谓仿真,也就是把实验观察屏上二维的干涉条纹图像通过计算机再现出来,这就是图像显示技术。
图像显示是一种特殊的图形绘制,MATLAB提供了一系列创建和显示图像对象的命令,例如image,pcolor,imshow等。
考虑到干涉条纹的明与暗可用白与黑来显示,而强度介于两者之间的条纹则可用不同层次的灰度来显示,因此选择灰度图像来显示牛顿环的干涉条纹比用彩色图像显示效果更逼真;此外,相对来说用imshow创建灰度图像比用image和pcolor在语句上要简单些,所以这里我们选择imshow指令来进行光强的可视化操作(即干涉图像的再现),至于image和pcolor指令亦可实现不过语句稍复杂些,感兴趣的同学可以参考有关的资料。
使用imshow建立灰度图像的使用格式:
imshow(I,N)
其中参数N为正整数,指定灰度的层次,当缺省该参数时,系统默认为256级的灰度级;参数I为数值矩阵,imshow的作用就是将数值矩阵I的元素值用N个灰度级的黑白图像可视化显示出来。
实际上是在数值矩阵I和N个灰色调之间建立了一种颜色映射关系:
I当中元素值最大者映射为白色(将该元素值作为白色显示),元素值最小者映射为黑色(将该元素值作为黑色显示),元素值介于最大和最小之间的则按照某种约定的规则映射到其它的灰度级(显示为不同灰度的灰色)。
为了方便叙述,假定再现的图像尺寸2mm×2mm,使用上述指令可以很方便的将牛顿环干涉条纹在该区域内再现:
●首先,利用(4)式获取干涉光强I在该区域(
,假定观察屏是xy平面)的数值分布
x=linspace(-0.001,0.001,200);
y=linspace(-0.001,0.001,200);
[X,Y]=meshgrid(x,y);%将xy平面2mm×2mm的区域分割为200×200的网
%格(像素),矩阵X、Y分别输出格点的x和y坐标
r2=X.^2+Y.^2 ;
I=abs(sin(pi*r2/R/)).^2;
%200×200的数值矩阵;计算格点上光强I,得到光强的数值分布
●接下来使用imshow指令将数值矩阵I可视化
imshow(I)
3.动态仿真技术
牛顿环实验(也包括其他光学实验)的仿真有2个环节,其一是将观察屏(xy平面)上干涉光强的分布可视化显示;其二是动态仿真,比如当空气薄膜的厚度连续改变时(通过向上移动牛顿环中的透镜),干涉条纹也会随之移动,采用动态仿真可以再现这一过程。
前者我们在第2小节中做了介绍,下面介绍一种MATLAB动态仿真技术—影片动画技术。
顾名思义,影片动画技术类似于电影的制作,其原理是首先对仿真的过程按时间次序进行“拍照”,获得一帧一帧的画面(称为帧),并将之存档,然后再按时间顺序以高于视觉暂留的帧频率播放帧,即可获得类似于电影的动画效果。
这种动画技术适用于难以实时快速绘制的复杂画面,但计算量大,占用内存较多。
在MATLAB中实现影片动画依次要用到下列几个函数:
①moviein函数该函数将产生一个结构体数组(structure,以下称帧结构体)来存放动画的帧(即所拍摄的一幅幅画面),每帧画面作为结构体的一个元素保存。
调用格式
fmat=moviein(N)
作用是创建一个能存放N个帧的(1×N)结构体数组fmat。
该结构体包含两个域cdata和colormap,前者用于存放帧的图像数据,后者存放帧使用到的颜色表。
②getframe函数该函数作用是对当前的图像进行快照(“抓拍”),通常有两种使用格式:
A.getframe该格式不带参数,“抓拍”当前坐标轴里的内容;
B.getframe(h)“抓拍”某个图形窗口或坐标轴里的内容,该图形窗口或坐标轴以句柄h标识(图形窗口和坐标轴都是一种图形对象,每一种图形对象都有自己特有的句柄handle,即标识,类似于“身份证”)。
例如
>>fmat
(1)=getframe(gcf)
抓拍当前图形窗口下的内容,并将该帧画面存放于帧结构体fmat的第一个元素中;gcf为Gethandletocurrentfigure的缩写,意思是获取当前图形窗口的句柄。
在命令窗口中输入gcf可显示当前图形窗口的句柄值,是个整数。
③movie函数
作用是按顺序回放帧结构体fmat中存放的各帧画面以产生动画感。
一般格式:
movie(h,famt,n,fps)
h是播放动画的图形窗口或坐标轴的句柄,缺省时表示在当前的坐标轴中播放动画;
famt是帧结构体名,不可省;
n是重复播放的次数,缺省时,只播放一次;
fps代表每秒播放的帧数(即帧频),一般应快于视觉暂留,缺省时系统默认fps=12。
在动画播放前,movie函数首先将每帧图像的数据载入内存(此时图像会一帧一帧的显示出来),然后再按照用户设定的参数(重复次数n、帧频fps等)播放动画。
除了movie函数,MATLAB还提供了一个函数movie2avi,该函数能够在当前的工作目录下创建一个avi视频格式的动画文件。
一般的调用格式:
movie2avi(fmat,filename)
famt为前述的帧结构体名,filename是字符串,指定avi格式的文件名。
例1:
Z=peaks;%MATLAB提供的三维函数
surf(Z);
TheAxis=axis;%保存坐标值,使得下面所有帧都在同一坐标系。
变量TheAxis
%为6个元素的向量,分别代表x、y、z轴的最小、最大值
F=moviein(20);%创建可以存放20帧的帧结构体
forj=1:
20%该循环“抓拍”20帧画面并存放到F中
surf(sin(2*pi*j/20)*Z,Z);%画出每一步的曲面(帧)
axis(TheAxis);%使用相同的坐标系。
F(j)=getframe;%“抓拍”帧,并存在到帧结构体
end
movie(F,10)%回放保存在帧结构体F中的画面,重复10次,帧频12(缺省值)
movie2avi(F,’example’);%制作avi格式的视频动画
注意:
“抓拍”到的画面总帧数对动画的播放速度、连续感会产生影响,请将帧数修改为10,重新运行程序观察结果。
另,帧数越大帧结构体的内存开销也跟着增加。
四、实验内容及要求
1.绘制牛顿环干涉条纹图
在牛顿环实验中,假定透镜的曲率半径R=0.855m,入射光波长为589.3nm的钠黄光,实验观察到一幅6mm×6mm大小的牛顿环干涉条纹图。
请结合(4)式使用imshow指令对该图像进行静态仿真。
要求:
以M文件的形式编写程序,并调试通过,得到的干涉图像存为.jpg文件格式。
2.空气膜厚度连续变化时的动态仿真
若将图1中的透镜缓慢的向上移动,则每个点处空气薄膜厚度d将连续增加,假设空气膜厚度的变化用d表示,则薄膜上下表面的反射光之间的位相差由
(2)式变化为(5)式:
其中,
代表空气膜厚度变化导致的附加位相差。
相应地光强表达式也要修正为:
假设d的初终值分别为
。
空气膜厚度每变化0.04m时
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- matlab 实验 指导书