第四章 数值计算优秀PPT.ppt

《第四章 数值计算优秀PPT.ppt》由会员分享，可在线阅读，更多相关《第四章 数值计算优秀PPT.ppt（50页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、第四章 数值计算1第一页，本课件共有50页数值计算内容数值计算内容:数值计算包括数值计算包括:多项式运算多项式运算,线性方程组求解线性方程组求解,矩阵特征值问题的矩阵特征值问题的解解,卷积卷积,数据分析数据分析,泛函指令的使用泛函指令的使用,信号处理和系统分析信号处理和系统分析.4.1 多项式运算多项式运算多项式运算是数学中最基本的运算之一多项式运算是数学中最基本的运算之一.多项式一般可多项式一般可以表示为以表示为:f(x)=a0 xn+a1xn-1+a2xn-2+an-1x+an对于这种表示形式，很容易用一个行向量表示，即：对于这种表示形式，很容易用一个行向量表示，即：T=a0,a1,a2,

2、an-1,an在在MATLAB中中,多项式正是用这样一个行向量表示的多项式正是用这样一个行向量表示的,系数系数是按降序排列的是按降序排列的幂幂2第二页，本课件共有50页4.1.1 多项式构造多项式构造多项式可以直接用向量表示多项式可以直接用向量表示,因此因此,构造多项式最简单的方法是构造多项式最简单的方法是直接输入向量直接输入向量.例例4.1.1-1 直接输入向量构造直接输入向量构造 f(x)=2x5+5x4+4x2+x+4T=2,5,0,4,1,4;fx=poly2sym(T)函数函数poly2sym是符号工具箱中的函数是符号工具箱中的函数,在用此种方式在用此种方式构造多项式时构造多项式时,

3、必须把多项式各项的系数写完整必须把多项式各项的系数写完整,而不而不管此项的系数是否为管此项的系数是否为0fx=2*x5+5*x4+4*x2+x+43第三页，本课件共有50页r=1,2,3,4;T1=poly(r);y=poly2sym(T1)y_class=class(y)例例4.1.1-2 用多项式的根构造多项式用多项式的根构造多项式,根为根为r=1,2,3,4T1=1 -10 35 -50 24y=x4-10*x3+35*x2-50*x+24y_class=sym4第四页，本课件共有50页4.1.2 多项式的运算多项式的运算多项式的运算主要包括多项式的四则运算多项式的运算主要包括多项式的四

4、则运算,导数运算导数运算,估值估值运算运算,求根运算以及多项式的拟合等求根运算以及多项式的拟合等1 多项式的四则运算多项式的四则运算多项式四则运算主要是多项式的加多项式四则运算主要是多项式的加,减减,乘乘,除除.其中其中,多项多项式的加减运算要求两个相加式的加减运算要求两个相加,减的多项式的大小必须相等减的多项式的大小必须相等,此时加此时加,减才有效减才有效.当两个相加当两个相加,减的多项式阶次不同时必减的多项式阶次不同时必须通过补须通过补0使其相同使其相同.加加/减减-+/-乘/除除-conv,deconvT=deconv(T1,T3)T,r=deconv(T1,T3)商多项式商多项式余式余

5、式T3为分母为分母5第五页，本课件共有50页T1=2,5,0,4,1,4;T2=0,0,5,1,3,2;T3=5,1,3,2;%除法运算中分母多项式第一个系数不能为除法运算中分母多项式第一个系数不能为0T=T1+T2;%必须是同维的才能相加必须是同维的才能相加T_add=poly2sym(T)T=T1-T2;T_sub=poly2sym(T)T=conv(T1,T2);%乘法不要求同维乘法不要求同维T_mul=poly2sym(T)A_coe,A_r=deconv(T1,T3);T_coe=poly2sym(A_coe)T_rem=poly2sym(A_r)例例4.1.1-3 多项式的加减乘除

6、运算多项式的加减乘除运算f1(x)=2x5+5x4+4x2+x+4,f2(x)=5x3+x2+3x+26第六页，本课件共有50页例例4.1.1-4 多项式求值多项式求值,求上式求上式f1(x)在在x=0.5处的函数值处的函数值T1=2,5,0,4,1,4;x=0.5;y=polyval(T1,x)y=5.87507第七页，本课件共有50页3 多项式的导数运算多项式的导数运算-polyder例例4.1.1-6 求多项式求多项式f1(x)=2x5+5x4+4x2+x+4的导数的导数 T1=2,5,0,4,1,4;h=polyder(T1);poly2sym(h)2 多项式求根多项式求根-roots

7、例例4.1.1-5 求多项式求多项式f1(x)=2x5+5x4+4x2+x+4的根的根 T1=2,5,0,4,1,4;root=roots(T1);root=-2.7709 0.5611+0.7840i 0.5611-0.7840i -0.4257+0.7716i -0.4257-0.7716i10*x4+20*x3+8*x+18第八页，本课件共有50页4 拟合和插值拟合和插值-polyfit,interp1例例4.1.1-7 对下列数据对对下列数据对(x0,y0)求三次拟合多项式并绘图求三次拟合多项式并绘图.x0=0:0.1:1;y0=-0.447,1.978,3.11,5.25,5.02,

8、4.66,4.01,4.58,3.45,5.35,9.22;p=polyfit(x0,y0,n);p,s=polyfit(x0,y0,n);x0,y0-给定数据对给定数据对n-拟合出的多项式次数拟合出的多项式次数p-多项式向量多项式向量s-偏差信息偏差信息yi=interp1(x0,y0,xi,cubic);xi,yi-得到的新的插值点得到的新的插值点cubic-插值方式插值方式9第九页，本课件共有50页x0=0:0.1:1;y0=-0.447,1.978,3.11,5.25,5.02,4.66,4.01,4.58,3.45,5.35,9.22;n=3;%设定拟合次数为设定拟合次数为3p,s=

9、polyfit(x0,y0,n);%得到拟合多项式向量和相关偏差信得到拟合多项式向量和相关偏差信息息 xx=0:0.01:1;yy=polyval(p,xx);%按拟合曲线计算采样值按拟合曲线计算采样值n1=6;%设定拟合次数为设定拟合次数为6p1,s1=polyfit(x0,y0,n1);yy1=polyval(p1,xx);plot(xx,yy,-b,xx,yy1,-m,x0,y0,.r,MarkerSize,20);y3=polyval(p,0.5);y6=polyval(p1,0.5);text(0.5,y3-0.3,n=3);text(0.5,y6+0.2,n=6);10第十页，本课

10、件共有50页y=51.48x3-77.74x2+35.06x-0.20y=348.21x6-1060.23x5+1297.68x4-758.90 x3+181.00 x2+1.00 x+0.4811第十一页，本课件共有50页拟合只能在给定数据所限定的区间内使用拟合只能在给定数据所限定的区间内使用,不要任意向外拓展不要任意向外拓展例例4.1.1-8 按上例所给数据研究插值按上例所给数据研究插值,并观察插值和拟合的区别并观察插值和拟合的区别.x0=0:0.1:1;y0=-0.447,1.978,3.11,5.25,5.02,4.66,4.01,4.58,3.45,5.35,9.22;xi=0:0.

11、02:1;yi=interp1(x0,y0,xi,cubic);subplot(1,2,1);plot(xi,yi,-b,x0,y0,.r,MarkerSize,20);xlabel(x),ylabel(y);p,s=polyfit(x0,y0,3);xx=0:0.01:1;yy=polyval(p,xx);subplot(1,2,2);plot(xx,yy,-r,x0,y0,.r,MarkerSize,20);xlabel(x);ylabel(y);12第十二页，本课件共有50页插值插值拟合拟合13第十三页，本课件共有50页插值方式插值方式:linear-线性插值线性插值nearest-最近

12、插值最近插值cubic-三次多项式插值三次多项式插值spline-样条插值样条插值pchip-分段三次分段三次Hermite插值插值v5cubic-MATLAB 5.0版的三次多项式插值版的三次多项式插值.拟合和插值的区别拟合和插值的区别:曲线拟合是研究如何寻找曲线拟合是研究如何寻找“平滑平滑”曲线曲线最好地表现带噪声的最好地表现带噪声的“测量数据测量数据”,但并不要求拟合曲线通但并不要求拟合曲线通过这些过这些“测试数据测试数据”点点.插值是在认定所给数据完全正确的插值是在认定所给数据完全正确的情况下情况下,研究如何平滑地估算出研究如何平滑地估算出“基准数据基准数据”之间其他点的之间其他点的函

13、数值函数值,因此插值一定通过因此插值一定通过“基准数据基准数据”14第十四页，本课件共有50页通俗地讲通俗地讲,拟合就是由已知点得到一条曲线拟合就是由已知点得到一条曲线,使该曲线使该曲线最能反映点所代表的规律最能反映点所代表的规律.比如做欧姆定理的实验的时候比如做欧姆定理的实验的时候,由于实验中存在误差由于实验中存在误差,最后拟合得到的曲线是一条直线最后拟合得到的曲线是一条直线,而且肯定只有部分点落在拟合的直线上而且肯定只有部分点落在拟合的直线上,但此时该直线和但此时该直线和测试点的方差最小测试点的方差最小.由拟合直线的斜率就可以知道电阻的由拟合直线的斜率就可以知道电阻的阻值阻值.拟合是探测事

14、物变化规律的办法拟合是探测事物变化规律的办法.插值就是根据函数上某些已知点插值就是根据函数上某些已知点(或实验数据或实验数据),按一定规律按一定规律(插值方法插值方法)寻求未知的点寻求未知的点,比如已知一个常用对数比如已知一个常用对数y=log(x)表表,是按照是按照x=0.1:0.1:10制表的制表的,如果按已知数据求如果按已知数据求y=log(2.897)就可以用插值得到就可以用插值得到.表制得越密表制得越密,插值越准确插值越准确15第十五页，本课件共有50页例例4.1.1-9 已知通过实验得到电压和电流的数值如下：已知通过实验得到电压和电流的数值如下：i0=0:0.5:10;v0=0,9

15、8,202,306,395,506,592,708,789,890,1009,1108,1186,1310,1391,1510,1601,1716,1782,1920,2014;按一次按一次,二次二次,三次多项式对数值关系进行拟合三次多项式对数值关系进行拟合说明：说明：为保证较好的拟合效果，多项式阶数要取得适当，过低，残为保证较好的拟合效果，多项式阶数要取得适当，过低，残差较大，过高，拟合模型将包含噪声影响，通常要保证拟合差较大，过高，拟合模型将包含噪声影响，通常要保证拟合阶数小于数据对的数目。阶数小于数据对的数目。16第十六页，本课件共有50页i0=0:0.5:10;v0=0,98,202,

16、306,395,506,592,708,789,890,1009,1108,.1186,1310,1391,1510,1612,1716,1782,1920,2014;n1=1;p1,s=polyfit(i0,v0,n1);n2=2;p2,s=polyfit(i0,v0,n2);n3=3;p3,s=polyfit(i0,v0,n3);ii=0:0.1:10;v1=polyval(p1,ii);v2=polyval(p2,ii);v3=polyval(p3,ii);subplot(1,3,1),plot(ii,v1,-b,i0,v0,.r,MarkerSize,8),xlabel(一次拟合一次拟

17、合);subplot(1,3,2),plot(ii,v2,-b,i0,v0,.r,MarkerSize,8),xlabel(二次拟合二次拟合);subplot(1,3,3),plot(ii,v3,-b,i0,v0,.r,MarkerSize,8),xlabel(三三次拟合次拟合);p1=201.00,-2.89P2=0.31,197.90,2.02P3=0.03,-0.10,199.49,0.8517第十七页，本课件共有50页p1=201.00,-2.89P2=0.31,197.90,2.02P3=0.03,-0.10,199.49,0.8518第十八页，本课件共有50页对于含对于含n个未知数

18、的个未知数的n个方程构成的方程组个方程构成的方程组Ax=b,在线性代数教科在线性代数教科书中书中,最常介绍的解法有最常介绍的解法有:Cramer法法;逆阵法逆阵法,即即x=A-1b;高斯消元法高斯消元法;LU法法 对于维数不高对于维数不高,条件数不大的矩阵条件数不大的矩阵,以上四种所得结果不会有明显以上四种所得结果不会有明显的差别的差别.但前三种解法的更多的意义在理论上但前三种解法的更多的意义在理论上,而不在实际的数值而不在实际的数值计算上计算上.在在MATLAB中中,方程采用方程采用LU法求解法求解,并且出于算法稳定性并且出于算法稳定性的考虑的考虑,行列式和逆的计算也都是在行列式和逆的计算也

19、都是在LU分解的基础上进行的分解的基础上进行的MATLAB中的四条指令中的四条指令:L,U,P=lu(A)-矩阵的矩阵的LU分解分解,使使LU=PA.多种格式多种格式.det(A)-求矩阵求矩阵A的行列式的行列式,它由它由detA=detUuii 算得算得inv(A)-求矩阵求矩阵A的逆的逆,由由A-1=U-1L-1P算得算得x=Ab-除法指令求方程的解除法指令求方程的解(对恰定方程对恰定方程,采用采用LU法执行法执行)4.2.1 方程组求解方程组求解4.2 线性方程组的解线性方程组的解19第十九页，本课件共有50页对于方程组对于方程组Ax=b,采用采用x=Ab计算，如果方程组为计算，如果方程

20、组为yC=d,要要使用右除，即指令为使用右除，即指令为y=d/CAx=bxA=byC=dx=Abx=b/Ay=d/C例例4.2.1-1 解下列方程组解下列方程组2x1+2x2+3x3=34x1+7x2+7x3=1-2x1+4x2+5x3=-7A=2,2,3;4,7,7;-2,4,5,b=3;1;-7x=Ab,y=b/A求解方程时求解方程时,尽量不要使用指令尽量不要使用指令inv(A)*b进行进行,它不仅计算速度它不仅计算速度没有除法快没有除法快,而且计算精度也没除法高而且计算精度也没除法高.此外此外,除法还适用于除法还适用于“超定超定”和和“欠定欠定”方程方程x=9.00 -36.00 31.

21、00y=9.00 -36.00 31.0020第二十页，本课件共有50页方程组的解的三种情况方程组的解的三种情况:对于方程对于方程Ax=b,A为为Amn矩阵矩阵,有三种情况有三种情况:当当m=n时时,此方程成为此方程成为恰定恰定方程方程,求解精确解求解精确解 当当mn时时,此方程成为此方程成为“超定超定”方程方程,寻求寻求最小二乘解最小二乘解(直线拟合直线拟合)当当mn时此时无解时此时无解,但实用中可以求最小二乘解即但实用中可以求最小二乘解即:方程解方程解(AA)x=A b x=(AA)-1Ab-求逆法求逆法 x=Ab-matlab用最小二乘法找一个近似解用最小二乘法找一个近似解.x1+2x2

22、=1 2x1+3x2=23x1+4x2=3例例4.2.1-3 求下列超定方程的解求下列超定方程的解 求逆法：求逆法：x=inv(A*A)*A*b 使用使用Ab 的方式的方式A=1,2;2,3;3,4;b=1;2;3;x1=inv(A*A)*A*b;x2=Abx1=1.00 0.00 x2=1.00 022第二十二页，本课件共有50页 3)3)欠定方程欠定方程当方程数少于未知量个数时当方程数少于未知量个数时,即不定情况即不定情况,有无穷多个解存在有无穷多个解存在.matlab可求出两个解：可求出两个解：a.用除法求的解用除法求的解x是具有最多零元素的解是具有最多零元素的解b.具有最小长度或范数的

23、解具有最小长度或范数的解,这个解是基于伪逆这个解是基于伪逆pinv求得的求得的.例例4.2.1-4求下列欠定方程组的解求下列欠定方程组的解x1+2x2+3x3=12x1+3x2+4x3=2A=1,2,3;2,3,4;b=1;2;x1=Ab;x2=pinv(A)*bx1=1 0 0 x2=0.8333 0.3333 -0.166723第二十三页，本课件共有50页4.2.2 奇异值分解和矩阵结构奇异值分解和矩阵结构相应的相应的MATLAB指令如下指令如下:U,S,V=svd(A)-矩阵的奇异值分解三对组阵，使矩阵的奇异值分解三对组阵，使A=USVTnorm(A,flag)-计算矩阵计算矩阵A的范数

24、，范数类型由的范数，范数类型由flag指定指定.cond(A)-计算矩阵计算矩阵A的条件数的条件数.pinv(A)-求矩阵的广义逆求矩阵的广义逆flag-1,2,inf,fro24第二十四页，本课件共有50页4.2.3 线性二乘问题的解线性二乘问题的解对于超定方程对于超定方程,求其最小二乘解有三种方法：求其最小二乘解有三种方法：正则方程法得解正则方程法得解x=(A A)-1A b,数值精度差数值精度差广义逆法得解广义逆法得解x=A+b,所得解可靠所得解可靠.用矩阵除法得解用矩阵除法得解x=Ab.可靠性稍差，但速度快可靠性稍差，但速度快.例例4.2.3-1 对于超定方程对于超定方程y=Ax,进行

25、三种解法比较进行三种解法比较.(1)生成矩阵生成矩阵A及及y,并用三种方法求解并用三种方法求解A=gallery(5),A(:,1)=;y=1.7,7.5,6.3,0.83,-0.082;x=inv(A*A)*A*y,xx=pinv(A)*y,xxx=Ayx=3.4811 5.1595 0.9534 -0.0466xx=3.4759 5.1948 0.7121 -0.1101xxx=3.4605 5.2987 0 -0.2974由于由于阵缺一个列秩阵缺一个列秩,除法给出的最小二乘除法给出的最小二乘解就只有三个非零元解就只有三个非零元素素.它是只有三个非零它是只有三个非零元素的所有最小二乘元素的

26、所有最小二乘解中范数最小的解中范数最小的Gallery-MATLAB 设置设置的特殊矩阵的特殊矩阵.25第二十五页，本课件共有50页(2)计算三个解的范数计算三个解的范数nx=norm(x),nxx=norm(xx),nxxx=norm(xxx)nx=6.2968nxx=6.2918nxxx=6.3356用广义逆所得的最小二乘范数最小(3)比较三种解法的方程误差比较三种解法的方程误差e=norm(y-A*x),ee=norm(y-A*xx),eee=norm(y-A*xxx)e=0.6986ee=0.0474eee=0.0474误差的平方根最大，精度最差误差的平方根最大，精度最差26第二十六页

27、，本课件共有50页4.2.4 特征值分解和矩阵函数特征值分解和矩阵函数涉及矩阵特征值分解的常用的涉及矩阵特征值分解的常用的MATLAB指令如下指令如下:d=eig(A)-计算矩阵计算矩阵A的特征值的特征值,以向量形式存放以向量形式存放V,D=eig(A)-计算矩阵的特征向量计算矩阵的特征向量V和特征值对角阵和特征值对角阵D,使使A*V=V*DVR,DR=cdf2rdf(VC,DC)-把复数对角形转化成实数块对角形把复数对角形转化成实数块对角形.VC,DC=rsf2csf(VR,DR)-把实数块对角形转换成复数对角形把实数块对角形转换成复数对角形V,J=jordan(A)-jordan分解分解,

28、使使A*V=V*J,J是对角阵是对角阵c=condeig(A)-计算各特征值的条件数计算各特征值的条件数V,D,C=condeig(A)-相当于相当于V,D=eig(A)和和c=condeig(A)两条指令两条指令27第二十七页，本课件共有50页例例4.2.4-1 验证验证det(A)=iA=2,-1,-1;3,4,-2;3,-2,4det(A)d=eig(A)ans=60d=2.0000+2.4495i 2.0000-2.4495i 6.0000 d(1)*d(2)*d(3)ans=60.0000行列式行列式特征值特征值%矩阵行列式等于其所有特征值的积矩阵行列式等于其所有特征值的积28第二十

29、八页，本课件共有50页4.3 数据分析函数数据分析函数(1)若输入宗量若输入宗量x是向量是向量,那么不管是行向量还是列向量那么不管是行向量还是列向量,运算是对运算是对整个向量进行的整个向量进行的(2)若输入宗量若输入宗量x是二维数组，那么指令运算是按列进行的是二维数组，那么指令运算是按列进行的MATLAB在进行数据分析时的约定在进行数据分析时的约定:29第二十九页，本课件共有50页4.3.1 随机数发生器和统计分析指令随机数发生器和统计分析指令函数函数含义含义rand(m,n)产生产生mxn维的维的0,1区间均匀分布随机数组区间均匀分布随机数组randn(m,n)产生产生mxn维的均值为维的均

30、值为0标准差为标准差为1的正态分布的正态分布随机数组随机数组min(x)对矩阵各列分别求最小值对矩阵各列分别求最小值max(x)对矩阵各列分别求最大值对矩阵各列分别求最大值median(x)对矩阵各列分别求中位数对矩阵各列分别求中位数mean(x)对矩阵各列分别求均值对矩阵各列分别求均值std(x)对矩阵各列分别求标准差对矩阵各列分别求标准差var(x)对矩阵各列分别求方差对矩阵各列分别求方差C=conv(x)给出矩阵各列的协方差阵给出矩阵各列的协方差阵P=corrcoef(x)给出矩阵各列间的相关系数给出矩阵各列间的相关系数30第三十页，本课件共有50页例例4.3-1 基本统计示例基本统计示

31、例randn(state,0);A=randn(1000,4);AMAX=max(A),AMIN=min(A)AMED=median(A);AMEAN=mean(A);ASTD=std(A)已知已知A=1,2,3;4,5,6;7,8,9,演示演示max,min,median,mean,std,var,var,conv,corrcoefconv,corrcoef各列最大各列最大,最小值：最小值：max(A),min(A)各行最大各行最大,最小值：最小值：max(A),min(A)各列中间值：各列中间值：median(A)各列平均值：各列平均值：mean(A)各列方差各列方差,标准差：标准差：va

32、r(A),std(A)A=1,2,3;4,5,6;7,8,9 7 8 931第三十一页，本课件共有50页4.3.2 差分和累计指令差分和累计指令 函数函数含义含义del2(U,hx,hy)五点五点Laplaciandiff(X,m,n)沿第沿第n维求维求m次差分次差分(二维：二维：n=1;列列,n=2,行行)DZx,DZy=gradient(z,dx,dy)对对Z求求x,y方向梯度方向梯度prod(X,n)沿第沿第n维求乘积维求乘积sum(X,n)沿第沿第n维求和维求和trapz(x,Y,n)用梯形法沿第用梯形法沿第n维求函数关于自变量维求函数关于自变量x的积分的积分cumprod(X,n)沿

33、第沿第n维求累计乘积维求累计乘积cumsum(X,n)沿第沿第n维求累计和维求累计和cumtrapz(x,Y,n)用梯形法沿第用梯形法沿第n维求函数维求函数Y关于关于x的累计积分的累计积分XS,KK=sort(X,n)沿第沿第n维对维对X元素按模增大排列元素按模增大排列32第三十二页，本课件共有50页diff,trapz,cumtrapz指令的演示指令的演示a=1,2,3,4,5;b=a,a+1,a+2,a+3diff(a)=diff(a,1)=1,1,1,1diff(b)=1,1,1,1;1,1,1,1;1,1,1,1;1,1,1,1trapz(b)=12 16 20 24cumtrapz(

34、b)0 0 0 01.50 2.50 3.50 4.504.00 6.00 8.00 10.007.50 10.50 13.50 16.5012.00 16.00 20.00 24.001 2 3 42 3 4 53 4 5 64 5 6 75 6 7 833第三十三页，本课件共有50页例例4.3.2-1求求1+2+3+100以及以及50!x=1:1:100;sum_x=sum(x);sum_x=5050a=1:1:50;prod_a=prod(a)prod_a=3.0414e+064cumsum,cumprod的用法的用法上两个函数分别是求向量的累计和和累计乘积上两个函数分别是求向量的累计和

35、和累计乘积.如果如果a=1:1:n,则则cumsum(a)=1,1+2,1+2+3,.,1+2+3+ncumprod(a)=1,1*2,1*2*3,1*2*3*n34第三十四页，本课件共有50页例例4.3.2-1 求求f(x)=3x2在区间在区间0,2的积分的积分dt=0.001;t=(0:dt:2);y=3*t.2;s1=dt*sum(y)s2=dt*trapz(y)s=dt*cumsum(y);s3=s(end)s=dt*cumtrapz(y);s4=s(end)matlab还有更精良的积分指令还有更精良的积分指令:quad,quda8s1=8.0060s2=8.000001s3=8.00

36、60s4=8.00000135第三十五页，本课件共有50页4.4 MATLAB泛函指令泛函指令在在MATLAB中中,凡以函数为输入宗量的指令凡以函数为输入宗量的指令,都被统称为泛函指令都被统称为泛函指令.最常见的有最常见的有:z=fzero(fun,x0)-求一元函数零点指令的最简单格式求一元函数零点指令的最简单格式x=fsolve(fun,x0)-解非线性方程组的最简单格式解非线性方程组的最简单格式x=fminbnd(fun,x1,x2)-求函数在区间求函数在区间(x1,x2)中极小值的指令中极小值的指令 最简格最简格式式x=fminsearch(fun,x0)-单纯形法求多元函数极值点指令

37、的最单纯形法求多元函数极值点指令的最简格式简格式.x=fminunc(fun,x0)-拟牛顿法求多元函数极值点指令的最简拟牛顿法求多元函数极值点指令的最简格式格式a=lsqnonlin(fun,a0)-解非线性最小二乘问题指令的最简格式解非线性最小二乘问题指令的最简格式.36第三十六页，本课件共有50页q=quad(fun,a,b)-采用递推自适应采用递推自适应simpson法计算积分法计算积分q=quadl(fun,a,b)-采用递推自适应采用递推自适应Lobatto法求数值积分法求数值积分SS=dblquad(fun,inmin,inmax,outmin,outmax)-二重闭环数二重闭环

38、数值积分值积分t,YY=ode45(fun,tspan,Y0)-采用采用4,5阶阶Runge-Kutta方程解算方程解算ODE初值问题初值问题指令中被处理的函数指令中被处理的函数fun,可以取可以取:字符串表达式字符串表达式,内联函数内联函数,M函数函数文件的函数句柄文件的函数句柄37第三十七页，本课件共有50页4.4.1 求函数零点求函数零点对于任意函数对于任意函数f(x)=0,零点情况复杂零点情况复杂,没有一个通用解法没有一个通用解法.一一般来说般来说,零点的数值计算过程是零点的数值计算过程是:先猜测一个初始零点或该先猜测一个初始零点或该零点所在的区间零点所在的区间,然后通过计算然后通过计

39、算,使猜测值不断精确化使猜测值不断精确化,或使或使猜测区间不断收缩猜测区间不断收缩,直到达到预定的精度直到达到预定的精度,终止计算终止计算.求函数的零点的步骤为求函数的零点的步骤为:(1)利用利用matlab作图指令获取初步近似值作图指令获取初步近似值具体做法具体做法:先确定一个零点可能存在的自变量区间先确定一个零点可能存在的自变量区间,然后利用然后利用plot指指令画出令画出f(x)在该区间中的图形在该区间中的图形,观察观察f(x)与横轴的交点坐标与横轴的交点坐标,也可以用也可以用zoom对交点处局部放大再读数对交点处局部放大再读数.借助借助ginput指令获得更精确的交点指令获得更精确的交

40、点坐标值坐标值38第三十八页，本课件共有50页(2)利用利用fzero指令求精确解指令求精确解z=fzero(fun,x0)-求一元函数零点指令的最简格式求一元函数零点指令的最简格式z,f_z,exitflag=fzero(fun,x0,options,p1,p2,)-最完整格式最完整格式 说明：说明：由于由于fzero是根据函数是否越横轴来决定零点的是根据函数是否越横轴来决定零点的,因此本程序无法确因此本程序无法确定函数曲线仅触及横轴和不穿越的零点定函数曲线仅触及横轴和不穿越的零点 第二个输入量第二个输入量x0表示零点的初始猜测表示零点的初始猜测.他可以是标量或二元向量他可以是标量或二元向量

41、.当当x0取标量时取标量时,该指令将在它两侧寻找一个与之最靠近的零点该指令将在它两侧寻找一个与之最靠近的零点;当取二元当取二元向量向量a,b时时,该指令将在区间该指令将在区间a,b内寻找一个零点内寻找一个零点.输入宗量输入宗量options是优化迭代采用的参数是优化迭代采用的参数.P1,P2是向函数是向函数fun传递传递的参数的参数.z是所求零点的自变量值是所求零点的自变量值.f_z是函数值是函数值,exitflag表明程序终止表明程序终止计算的条件计算的条件.若若exitflag0,表明找到零点后退出表明找到零点后退出39第三十九页，本课件共有50页例例4.4.1-1 求函数求函数f(x)=

42、2x-10 x的零点的零点(1)绘制函数图形绘制函数图形x=-10:0.1:10;y=2.x-10*x;hold on;plot(x,y,r,LineWidth,5);plot(x,zeros(size(x),k);%便于查看零点便于查看零点xlabel(x),ylabel(y);title(曲线方程式为：曲线方程式为：y=2x-10*x);(2)获取初步近似解获取初步近似解zoom on;tt,yy=ginput(2);zoom off;40第四十页，本课件共有50页func=2x-10*x;x01,err1,exitflag1=fzero(func,tt(1),);x02,err2,exi

43、tflag2=fzero(func,tt(2),);(3)利用利用fzero指令求精确解指令求精确解str=方程的根为方程的根为:x1=,num2str(x01),和和 x2=,num2str(x02);disp(str);41第四十一页，本课件共有50页4.4.2 求函数极值点求函数极值点许多科学研究和工程计算问题都可以归结为一个极值问题许多科学研究和工程计算问题都可以归结为一个极值问题,如能量最小，时间最短如能量最小，时间最短,最佳最佳拟合拟合,最短路径等最短路径等.又因为又因为f(x)的极小值问题等价于的极小值问题等价于-f(x)的极大值问题的极大值问题,所以所以matlab只有处只有处

44、理极小值的指令理极小值的指令确切地说确切地说,这里讨论的只是这里讨论的只是“局域极值局域极值”问题问题.“全域最小全域最小”问题要复杂问题要复杂得多得多.至今没有一个系统性的方法可求解一般的至今没有一个系统性的方法可求解一般的“全域最小全域最小”问题问题.对于对于一元二元函数一元二元函数,作图观察对确定全域最小值有很好的应用价值作图观察对确定全域最小值有很好的应用价值.但更多元但更多元的函数的函数,就很难使用作图法就很难使用作图法.matlab求函数极值的三条指令求函数极值的三条指令x,fval,exitflag,output=fminbnd(fun,x1,x2,opt,p1,p2,)-求一元

45、函数在区间求一元函数在区间(x1,x2)中极小值的指令的最完整格式中极小值的指令的最完整格式x,fval,exitflag,output=fminsearch(fun,x0,opt,p1,p2,)-单纯形法求多元函数极值点指令的最完整格式单纯形法求多元函数极值点指令的最完整格式x,fval,exitflag,output,grad,hessian=fminunc(fun,x0,opt,p1,p2,)-拟牛顿法求多元函数极值点指令的最完整格式拟牛顿法求多元函数极值点指令的最完整格式42第四十二页，本课件共有50页例例4.4.2-1 求例求例4.4.1-1 中的函数的最小值中的函数的最小值(1)函

46、数采用字符串表达式函数采用字符串表达式func=2x-10*x;x1=0,x2=8;%通过观察通过观察,发现最小值在此区间发现最小值在此区间Xmin,Ymin=fminbnd(func,x1,x2)Xmin=3.8507Ymin=-24.0800(2)函数采用函数采用M文件函数的句柄文件函数的句柄%M文件内容，保存为文件内容，保存为func.mfunction y=f(x)y=2x-10*x;%M文件内容文件内容x1=0;x2=8;Hf=func;Xmin,Ymin=fminbnd(Hf,x1,x2)43第四十三页，本课件共有50页4.4.3 数值积分数值积分数值积分有闭型和开型算法之分两者的

47、区别在于：是否需要计算积数值积分有闭型和开型算法之分两者的区别在于：是否需要计算积分区间端点处的函数值分区间端点处的函数值matlab中仅提供闭型数值积分指令中仅提供闭型数值积分指令:q=quad(fun,a,b,tol,trace,p1,p2,)-采用递推自采用递推自 适应适应Simpson法计算积分法计算积分q=qudal(fun,a,b,tol,trace,p1,p2,)-采用递推自适采用递推自适应应Lobatlo法求数值积分法求数值积分SS=dblquad(fun,inmin,inmax,outmin,outmax,tol,method)-二重积分指令二重积分指令fun-被积函数被积函

48、数,a,b-积分上下限积分上下限tol-绝对误差控制量，缺省时积分的绝对精度为绝对误差控制量，缺省时积分的绝对精度为1x10-6trace-取非时，将随积分的进程逐点画出被积函数取非时，将随积分的进程逐点画出被积函数44第四十四页，本课件共有50页例例4.4.3-1求求f(x)=3x2在区间在区间0,2的积分的积分y=3*x.2;a=0;b=2;q1=quad(y,a,b)q2=quadl(y,a,b)q1=8.00q2=8.0045第四十五页，本课件共有50页4.5 解常微分方程解常微分方程matlab为解常微分方程初值问题提供了一组配套齐全为解常微分方程初值问题提供了一组配套齐全,结构严整

49、的指令结构严整的指令,包括包括:ode45,ode23,ode113,ode23t,ode15s,ode23s,ode23tb.在此只介绍最常用的在此只介绍最常用的ode45的基本使用方法的基本使用方法.ode45使用方法使用方法:x,y=ode45(file_name,xspan,y0)file_name=file_namex-微分区间的点系列微分区间的点系列y-原函数在微分区间点系列上的函数值原函数在微分区间点系列上的函数值file_name-档名档名,保存的保存的M文件名文件名xspan-微分区间微分区间y0-初始值初始值46第四十六页，本课件共有50页例例4.5-1 解微分方程：解微分

50、方程：y=3x2,y(0)=1,微分区间微分区间0,5%M文件文件,保存为保存为dy_file.mfunction dy=f(x,y)dy=3*x2;y0=1;xspan=0,5;x,y=ode45(dy_file,xspan,y0);plot(x,y,r),xlabel(x),ylabel(y)47第四十七页，本课件共有50页4.6 数字滤波数字滤波y=filter(B,A,x,dim)-运用数字滤波对输入信号运用数字滤波对输入信号x滤波滤波滤波器的作用滤波器的作用:通过对采样数据信号进行数学运算处理来达到频域滤波的目的通过对采样数据信号进行数学运算处理来达到频域滤波的目的.例例 设计一个低