《Matlab基础与应用.ppt》由会员分享,可在线阅读,更多相关《Matlab基础与应用.ppt(85页珍藏版)》请在三一办公上搜索。
1、Matlab基础与应用,报告人:曹富军单 位:内蒙古科技大学Office:秋实楼A507Email:,内容提纲,MATLAB 语言的简洁高效性MATLAB 语言的科学运算功能MATLAB 语言的绘图功能MATLAB 庞大的工具箱与模块集MATLAB 强大的动态系统仿真功能,Matlab简介,MATLAB语言是当前国际上自动控制领域的首选计算机语言,也是很多理工科专业最适合的计算机数学语言.MATLAB作为线性系统的一种分析和仿真工具,是理工科大学生应该掌握的技术工具,它作为一种编程语言和可视化工具,可解决工程、科学计算和数学学科中许多问题。MATLAB是一种交互式的以矩阵为基础的系统计算平台,
2、它用于科学和工程的计算与可视化。它的优点在于快速开发计算方法,而不在于计算速度。,Matlab的特点,高度适应性、开放性:MATLAB的工具箱可以任意增减,任何人可以自己生成MATLAB工具箱可扩充性:MATLAB的函数大多为ASCII文件,可以直接编辑、修改基于矩阵运算的工作平台。多版本:windows/unix/dos/Macintosh极多的工具箱。,能做什么?,基本的数据处理优化和解方程动态过程仿真:实时的和非实时的数据来源:Excel、数据库、A/D等等嵌入式的控制:Pc/104和DSP神经元网络、小波分析、GA等等虚拟现实仿真,如何学习?,help 命令:已知命令不知道用法look
3、for命令:知道命令的关键词Matlab Help:命令查找,索引,说明书Demo有效的利用帮助,Matlab变量,基本运算 变量逗号、分号的意义行内编辑、注释、运算规则变量要求变量区分字母的大小写名字不能超过19个字符,第19个字符后的字符被忽略变量必须以字母开头,之后可以是任意字母、数字或者下滑线变量中不能含有标点符号。,特殊变量,显示格式,format命令改变显示格式,常用的的格式有:long(16位),bank(2个十进制位),hex(十六进制),short(缺省),short e(5位加指数)+(符号)long e(16位加指数)rat(有理数近似),数组,矩阵是MATLAB的核心创
4、建简单的数组,数组的访问,访问一块元素:x(a:b:c)表示访问数组x的从第a个元素开始,以步长为b到第c个元素(但不超过c),b可以为负数,b缺损时为1.矩阵建立 逗号或空格用于分隔某一行的元素,分号用于区分不同的行.除了分号,在输入矩阵时,按Enter键也表示开始一新行.输入矩阵时,严格要求所有行有相同的列.,访问一个元素:x(i)表示访问数组x的第i个元素,Matlab科学运算功能,直接赋值语句,矩阵定义,Matlab科学运算功能,特殊矩阵,Matlab科学运算功能,矩阵表示,矩阵加减法,矩阵乘法,Matlab科学运算功能,矩阵除法,AX=B,求 X,MATLAB 求解:X=AB,非奇异
5、:,矩阵右除:XA=B,求 X,非奇异:,MATLAB求解:X=B/A,Matlab科学运算功能,行列式 Determinent,矩阵的秩 Rank,矩阵的范数Norm,Matlab科学运算功能,矩阵特征值 Eigenvalue,E=eig(A):求矩阵A的全部特征值,构成向量E。,矩阵特征向量 Eigenvector,V,D=eig(A):求矩阵A的全部特征值,构成对角阵D,并求A的特征向量构成V的列向量。,Matlab科学运算功能,矩阵的逆 Inverse,线性方程组的求解,Matlab符号运算,建立符号对象 1建立符号变量和符号常量 MATLAB提供了两个建立符号对象的函数:sym和sy
6、ms,两个函数的用法不同。函数sym一次只能定义一个符号变量 sym a函数syms一次可以定义多个符号变量 syms a b c d A=a b;c d;det(A)=a*c-bd,符号运算,极限(limits)limit函数的调用格式为:(1)limit(f,x,a):求符号函数f(x)的极限值。即计算当变量x趋近于常数a时,f(x)函数的极限值。(2)limit(f,x,a,right):求符号函数f的极限值。right表示变量x从右边趋近于a。(3)limit(f,x,a,left):求符号函数f的极限值。left表示变量x从左边趋近于a。,符号极限,【例5】求极限 syms x;%定
7、义符号变量 f=(x*(exp(sin(x)+1)-2*(exp(tan(x)-)/sin(x)3;%确定符号表达式 w=limit(f)%求函数的极限,默认趋于0 w=-1/2,例:求下列极限。Problem 1:syms a m x;f=(x*(exp(sin(x)+1)-2*(exp(tan(x)-1)/(x+a);limit(f,x,a)ans=(1/2*a*exp(sin(a)+1/2*a-exp(tan(a)+1)/aProblem 2:syms x t;limit(1+2*t/x)(3*x),x,inf)ans=exp(6*t),符号运算,符号导数(differential)(1
8、)diff函数用于对符号表达式求导数。该函数的一般调用格式为:(2)diff(s,v):以v为自变量,对符号表达式s求一阶导数。(3)diff(s,v,n):以v为自变量,对符号表达式s求n阶导数。,符号导数,【例6】求导数:x=sym(x);%定义符号变量t=sym(t);diff(sin(x2)%求导运算ans=2*cos(x2)*x,符号运算,符号积分(integral)(1)符号积分由函数int来实现。该函数的一般调用格式为:(2)int(s):没有指定积分变量和积分阶数时,系统按findsym函数指示的默认变量对被积函数或符号表达式s求不定积分。(3)int(s,v):以v为自变量,
9、对被积函数或符号表达式s求不定积分。(4)int(s,v,a,b):求定积分运算。a,b分别表示定积分的下限和上限。该函数求被积函数在区间a,b上的定积分,符号积分,【例7】求下述积分。求积分:syms x int(1/(1+x2)ans=atan(x),符号运算,因式分解(factorization)factor函数的功能为:把多项式 S 分解为多个因式,各多项式的系数均为有理数。格式为:factor(s),【例9】将表达式(x9-1)分解为多个因式。syms x factor(x9-1)ans=(x-1)*(x2+x+1)*(x6+x3+1),级数求和,4.级数(级数求和)级数求和运算是数
10、学中常见的一种运算。例 如:f(x)=a0+a1x+a2x2+a3x3+anxn 函数symsum可以用于此类对符号函数f的求和运算。该函数的引用时,应确定级数的通项式S,变量的变化范围a和b。该函数的引用格式为:symsum(s,a,b),级数求和,【例8】求级数的和:1+1/22+1/32+1/42+syms k symsum(1/k2,1,Inf)%k值为1到无穷大 ans=1/6*pi2 其结果为:1/12+1/22+1/32+1/42+=2/6,二维图形功能Plot,基本绘图功能,Matlab提供了强大的图形绘制功能。在大多数情况下,用户只需要指定绘图的方式,提供绘图数据,利用Mat
11、lab提供的丰富的二维,三维图形函数,就可以绘制出所需的图形。1绘制二维连续函数Matlab中最常用的绘图函数是plot,plot的命令格式有以下几种:,Plot绘图函数,(1)plot(y)当y为一向量时,以y的序号作为x轴,按 向量y的值绘制曲线。(2)plot(x,y)x,y均为向量时,以x向量作为X轴,向量y 作为Y轴绘制曲线。注意:x和y种元素的个数必须相同!,Plot绘图,例一:画出衰减震荡曲线 及其包 络线 t的取值范围是 程序如下:,Plot绘图,Plot绘图,例二.用图形表示连续调制波形 及其包络。,程序如下:,Plot绘图,函数 m 文件,我们可以根据自己的需要建立自己的函
12、数文件,它与库文件一样方便调用,从而极大地扩展了Matlab的功能函数m文件的第一行有特殊的要求,如下文件名必须是.m,M函数举例,inline 函数和匿名函数,inline 函数,可以免去文件MATLAB 7.0,循环结构,for 结构while 结构,循环结构,利用循环求和,多项式插值,用MATLAB作插值计算,yi=interp1(x,y,xi,method),nearest:最邻近插值linear:线性插值;spline:三次样条插值;cubic:立方插值。缺省时:分段线性插值。,插值,例:在1-12的11小时内,每隔1小时测量一次温度,测得的温度依次为:5,8,9,15,25,29,
13、31,30,22,25,27,24。试估计每隔1/10小时的温度值。,hours=1:12;temps=5 8 9 15 25 29 31 30 22 25 27 24;h=1:0.1:12;t=interp1(hours,temps,h,spline);plot(hours,temps,+,h,t,hours,temps,r:)%作图xlabel(Hour),ylabel(Degrees Celsius),二维数据插值,在MATLAB中,二维插值问题的函数interp2,其调用格式为:Z1=interp2(X,Y,Z,X1,Y1,method)其中X,Y是两个向量,分别描述两个参数的采样点,
14、Z是与参数采样点对应的函数值,X1,Y1是两个向量或标量,描述欲插值的点。Z1是根据相应的插值方法得到的插值结果。method的取值与一维插值函数相同。X,Y,Z也可以是矩阵形式。注:X1,Y1的取值范围不能超出X,Y的给定范围,否则,会给出“NaN”错误,拟合,数据拟合是求一个简单的函数,例如是一个低次多项式,不要求通过已知的这些点,而是要求在整体上“尽量好”的逼近原函数。用插值的方法对一函数进行近似,要求所得到的插值多项式经过已知插值节点;在n比较大的情况下,插值多项式往往是高次多项式,这也就容易出现振荡现象(龙格现象),即虽然在插值节点上没有误差,但在插值节点之外插值误差变得很大,从“整
15、体”上看,插值逼近效果将变得“很差”。,拟 合 问 题 1,设 R=at+ba,b为待定系数,多项式拟合polyfit,用MATLAB作线性最小二乘拟合,1.作多项式f(x)=a1xm+amx+am+1拟合,可利用已有程序:,3.多项式在x处的值y可用以下命令计算:y=polyval(a,x),x0=-1+2*0:10/10;y0=1./(1+25*x0.2);x=-1:.01:1;ya=1./(1+25*x.2);p3=polyfit(x0,y0,3);y1=polyval(p3,x);p5=polyfit(x0,y0,5);y2=polyval(p5,x);p8=polyfit(x0,y0
16、,8);y3=polyval(p8,x);p10=polyfit(x0,y0,10);y4=polyval(p10,x);plot(x,ya,x,y1,x,y2,-.,x,y3,-,x,y4,:),多项式拟合的效果并不一定总是很精确的。,最小二乘法拟合,例,x=0.1:0.1:1;y=2.3201,2.6470,2.9707,3.2885,3.6008,3.9090,4.2147,4.5191,4.8232,5.1275;,function y=c8f3(a,x)y=a(1)*x+a(2)*x.2.*exp(-a(3)*x)+a(4);,a=lsqcurvefit(c8f3,1;2;2;3,x
17、,y);aMaximum number of function evaluations exceeded;increase options.MaxFunEvalsans=2.4575 2.4557 1.4437 2.0720,绘制曲线:y1=c8f3(a,x);plot(x,y,x,y1),Optimization terminated successfully:Relative function value changing by less than OPTIONS.TolFunans=0.1197 0.2125 0.5404 0.1702 1.2300res=7.1637e-007,绘制曲
18、线:x1=0:0.01:10;y1=f(xx,x1);plot(x1,y1,x,y,o),lsqnonlin用以求含参量x(向量)的向量值函数 f(x)=(f1(x),f2(x),fn(x)T 中的参量x,使得 最小。其中 fi(x)=f(x,xdatai,ydatai)=F(x,xdatai)-ydatai,非线性最小二乘拟合 lsqnonlin,已知数据点:xdata=(xdata1,xdata2,xdatan)ydata=(ydata1,ydata2,ydatan),用下面一组数据拟合参数a,b,k,该问题即解最优化问题:,数值积分与微分,数值积分数值积分基本原理 求解定积分的数值方法多
19、种多样,如简单的梯形法、辛普生(Simpson)法、牛顿柯特斯(Newton-Cotes)法等都是经常采用的方法。基本思想都是将整个积分区间a,b分成n个子区间xi,xi+1,i=1,2,n,其中,x1=a,xn+1=b。这样求定积分问题就分解为求和问题。,数值积分,变步长辛普生法 基于变步长辛普生法,MATLAB给出了quad函数来求定积分。该函数的调用格式为:I,n=quad(fname,a,b,tol,trace)fname是被积函数名 a和b分别是定积分的下限和上限。tol用来控制积分精度,缺省时取tol=0.001。trace控制是否展现积分过程,若取非0则展现积分过程,取0则不展现
20、,缺省时取trace=0。返回参数I即定积分值,n为被积函数的调用次数。,数值积分,例:求定积分(1)建立被积函数文件fesin.m。function f=fesin(x)f=exp(-0.5*x).*sin(x+pi/6);(2)调用数值积分函数quad求定积分 S,n=quad(fesin,0,3*pi)S=0.9008 n=77,牛顿柯特斯法,基于牛顿柯特斯法,MATLAB给出了quad8函数来求定积分。该函数的调用格式为:I,n=quad8(fname,a,b,tol,trace)该函数可以更精确地求出定积分的值,且一般情况下函数调用的步数明显小于quad函数,从而保证能以更高的效率求
21、出所需的定积分值.,二重积分,二重定积分的数值求解使用MATLAB提供的dblquad函数就可以直接求出上述二重定积分的数值解。该函数的调用格式为:I=dblquad(f,a,b,c,d,tol,trace)该函数求f(x,y)在a,bc,d区域上的二重定积分。参数tol,trace的用法与函数quad完全相同。,二重积分,计算二重定积分(1)建立一个函数文件fxy.m:function f=fxy(x,y)global ki;ki=ki+1;%ki用于统计被积函数的调用次数 f=exp(-x.2/2).*sin(x.2+y);(2)调用dblquad函数求解 global ki;ki=0;I
22、=dblquad(fxy,-2,2,-1,1)ki(3)I=1.57449318974494 ki=1038,数值微分,数值差分与差商在MATLAB中,没有直接提供求数值导数的函数,只有计算向前差分的函数diff,其调用格式为:DX=diff(X):计算向量X的向前差分,DX(i)=X(i+1)-X(i),i=1,2,n-1。DX=diff(X,n):计算X的n阶向前差分。,线性方程组求解,直接解法左除运算符 对于线性方程组Ax=b,利用左除运算符 求解:x=AbLU分解 L,U=lu(X):产生一个上三角阵U和一个变换形式的下三角阵L(行交换),使之满足X=LU。注意,这里的矩阵X必须是方阵
23、。实现LU分解后,线性方程组Ax=b的解 x=U(Lb)这样可以大大提高运算速度,线性方程组求解,迭代法求解雅可比(Jacobi)高斯赛德尔(Gauss-Sedial)超松弛迭代法(SOR),迭代法求解,雅可比(Jacobi)对于线性方程组Ax=b,如果A为非奇异方阵,即det(A)0(i=1,2,n),则可将A分解为A=D-L-U,其中D为对角阵,其元素为A的对角元素,L与U为A的下三角阵和上三角阵,于是Ax=b化为:x=D-1(L+U)x+D-1b与之对应的迭代公式为:x(k+1)=D-1(L+U)x(k)+D-1b 这就是Jacobi迭代公式。如果序列x(k+1)收敛于x,则x必是方程A
24、x=b的解。,Jacobi迭代法的MATLAB函数文件Jacobi.m如下:function y,n=jacobi(A,b,x0,eps)if nargin=3 eps=1.0e-6;elseif nargin=eps x0=y;y=B*x0+f;n=n+1;end,Gauss-Sedial,Gauss-Serdel迭代法在Jacobi迭代过程中,迭代公式为 Dx(k+1)=(L+U)x(k)+b可以改进为 Dx(k+1)=Lx(k+1)+Ux(k)+b,于是得到:x(k+1)=(D-L)-1Ux(k)+(D-L)-1b 该式即为Gauss-Serdel迭代公式。和Jacobi迭代相比,Gau
25、ss-Serdel迭代用新分量代替旧分量,精度会高些,Gauss-Serdel迭代法的MATLAB函数文件gauseidel.m如下:function y,n=gauseidel(A,b,x0,eps)if nargin=3 eps=1.0e-6;elseif nargin=eps x0=y;y=G*x0+f;n=n+1;end,非线性方程数值求解,单变量非线性方程求解在MATLAB中提供了一个fzero函数,可以用来求单变量非线性方程的根。该函数的调用格式为:z=fzero(fname,x0,tol,trace)其中fname是待求根的函数文件名,x0为搜索的起点。一个函数可能有多个根,但f
26、zero函数只给出离x0最近的那个根。tol控制结果的相对精度,缺省时取tol=eps,trace指定迭代信息是否在运算中显示,为1时显示,为0时不显示,缺省时取trace=0。,单变量非线性方程求解,例6.33 求f(x)在x0=-5和x0=1作为迭代初值的零点。先建立函数文件fz.m:function f=fz(x)f=x-1/x+5;然后调用fzero函数求根。:fzero(fz,-5)%以-5作为迭代初值 ans=-5.1926fzero(fz,1)%以1作为迭代初值 ans=0.1926,非线性方程组的求解,非线性方程组的求解对于非线性方程组F(X)=0,用fsolve函数求其数值解
27、。fsolve函数的调用格式为:X=fsolve(fun,X0,option)其中X为返回的解,fun是用于定义需求解的非线性方程组的函数文件名,X0是求根过程的初值,option为最优化工具箱的选项设定。最优化工具箱提供了20多个选项,用户可以使用optimset命令将它们显示出来。如果想改变其中某个选项,则可以调用optimset()函数来完成。例如,Display选项决定函数调用时中间结果的显示方式,其中off为不显示,iter表示每步都显示,final只显示最终结果。optimset(Display,off)将设定Display选项为off。,例6.34 求下列方程组在(1,1,1)附
28、近的解并对结果进行验证。首先建立函数文件myfun.m。function F=myfun(X)x=X(1);y=X(2);z=X(3);F(1)=sin(x)+y+z2*exp(x);F(2)=x+y+z;F(3)=x*y*z;在给定的初值x0=1,y0=1,z0=1下,调用fsolve函数求方程的根。X=fsolve(myfun,1,1,1,optimset(Display,off)X=0.0224-0.0224-0.0000,最优化问题求解,无约束最优化问题求解 在实际应用中,许多科学研究和工程计算问题都可以归结为一个最小化问题,如能量最小、时间最短等。MATLAB提供了3个求最小值的函数
29、,它们的调用格式为:求一元函数在(xl,x2)区间中的极小值点x和最小值fval x,fval=fminbnd(filename,x1,x2,option):单纯形算法求多元函数的极小值点x和最小值fval x,fval=fminsearch(filename,x0,option)基于拟牛顿法求多元函数的极小值点x和最小值 fval x,fval=fminunc(filename,x0,option):,例6.36 求函数在区间(-10,1)和(1,10)上的最小值点。首先建立函数文件fx.m:function f=f(x)f=x-1/x+5;上述函数文件也可用一个语句函数代替:f=inlin
30、e(x-1/x+5)再在MATLAB命令窗口,输入命令:fminbnd(fx,-10,-1)%求函数在(-10,-1)内的最小值点和最小值 fminbnd(f,1,10)%求函数在(1,10)内的最小值点。注意函数名f不用加例6.37 求函数f在(0.5,0.5,0.5)附近的最小值。建立函数文件fxyz.m:function f=fxyz(u)x=u(1);y=u(2);z=u(3);f=x+y.2./x/4+z.2./y+2./z;在MALAB命令窗口,输入命令:U,fmin=fminsearch(fxyz,0.5,0.5,0.5)%求函数的最小值点和最小值,有约束最优化问题求解,MATL
31、AB最优化工具箱提供了一个fmincon函数,专门用于求解各种约束下的最优化问题。该函数的调用格式为:x,fval=fmincon(filename,x0,A,b,Aeq,beq,Lbnd,Ubnd,NonF,option)其中 x、fval、filename、x0和option的含义与求最小值函数相同。其余参数为约束条件,参数NonF为非线性约束函数的M文件名。如果某个约束不存在,则用空矩阵来表示。,有约束最优化问题求解,求解有约束最优化问题。首先编写目标函数M文件fop.m。function f=fop(x)f=0.4*x(2)+x(1)2+x(2)2-x(1)*x(2)+1/30*x(1
32、)3;设定约束条件,调用fmincon函数求解此约束最优化问题 x0=0.5;0.5;A=-1,-0.5;-0.5,-1;b=-0.4;-0.5;lb=0;0;option=optimset;option.LargeScale=off;option.Display=off;x,f=fmincon(fop,x0,A,b,lb,option),常微分方程的数值求解,龙格库塔法基于龙格库塔法,MATLAB提供了求常微分方程数值解的函数,一般调用格式为:t,y=ode23(fname,tspan,y0)t,y=ode45(fname,tspan,y0)其中fname是定义f(t,y)的函数文件名,该函
33、数文件必须返回一个列向量。tspan形式为t0,tf,表示求解区间。y0是初始状态列向量。t和y分别给出时间向量和相应的状态向量,设有初值问题,试求其数值解,并与精确解相比较。(1)建立函数文件funt.m。function yp=funt(t,y)yp=(y2-t-2)/4/(t+1);(2)求解微分方程。t0=0;tf=10;y0=2;t,y=ode23(funt,t0,tf,y0);%求数值解 y1=sqrt(t+1)+1;%求精确解y为数值解,y1为精确值,显然两者近似。,已知一个二阶线性系统的微分方程,绘制系统的时间响应曲线和相平面图。函数ode23和ode45是对一阶常微分方程组设计的,因此对高阶常微分方程,需先将它转化为一阶常微分方程组,即状态方程。建立一个函数文件sys.m:function xdot=sys(t,x)xdot=-2*x(2);x(1);取t0=0,tf=20,求微分方程的解:t0=0;tf=20;t,x=ode45(sys,t0,tf,1,0);t,x subplot(1,2,1);plot(t,x(:,2);%解的曲线,即t-x subplot(1,2,2);plot(x(:,2),x(:,1)%相平面曲线,即x-x axis equal,
链接地址:https://www.31ppt.com/p-5439066.html