《第四章 数值计算精选PPT.ppt》由会员分享,可在线阅读,更多相关《第四章 数值计算精选PPT.ppt(50页珍藏版)》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、第四章 数值计算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第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第3页,此课件共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第4页,此课件共50页哦4.1.2 多项式的运算多项式的运算多项式的运算主要包括多项式的四则运算多项式的运算主要包括多项式的四
4、则运算,导数运算导数运算,估估值运算值运算,求根运算以及多项式的拟合等求根运算以及多项式的拟合等1 多项式的四则运算多项式的四则运算多项式四则运算主要是多项式的加多项式四则运算主要是多项式的加,减减,乘乘,除除.其中其中,多项多项式的加减运算要求两个相加式的加减运算要求两个相加,减的多项式的大小必须相等减的多项式的大小必须相等,此时加此时加,减才有效减才有效.当两个相加当两个相加,减的多项式阶次不同时必减的多项式阶次不同时必须通过补须通过补0使其相同使其相同.加加/减减-+/-乘/除除-conv,deconvT=deconv(T1,T3)T,r=deconv(T1,T3)商多项式商多项式余式余
5、式T3为分母为分母5第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第6页,此课件共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第7页,此课件共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第8页,此课件共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第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页,此
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第11页,此课件共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第12页,此课件共50页哦插值插值拟合拟合13第13页,此课件共50页哦插值方式插值方式:linear-线性插值线性插值nearest-最
12、近插值最近插值cubic-三次多项式插值三次多项式插值spline-样条插值样条插值pchip-分段三次分段三次Hermite插值插值v5cubic-MATLAB 5.0版的三次多项式插值版的三次多项式插值.拟合和插值的区别拟合和插值的区别:曲线拟合是研究如何寻找曲线拟合是研究如何寻找“平滑平滑”曲线最曲线最好地表现带噪声的好地表现带噪声的“测量数据测量数据”,但并不要求拟合曲线通过这但并不要求拟合曲线通过这些些“测试数据测试数据”点点.插值是在认定所给数据完全正确的情况下插值是在认定所给数据完全正确的情况下,研究如何平滑地估算出研究如何平滑地估算出“基准数据基准数据”之间其他点的函数值之间其
13、他点的函数值,因因此插值一定通过此插值一定通过“基准数据基准数据”14第14页,此课件共50页哦通俗地讲通俗地讲,拟合就是由已知点得到一条曲线拟合就是由已知点得到一条曲线,使该曲线最能反使该曲线最能反映点所代表的规律映点所代表的规律.比如做欧姆定理的实验的时候比如做欧姆定理的实验的时候,由于实验中由于实验中存在误差存在误差,最后拟合得到的曲线是一条直线最后拟合得到的曲线是一条直线,而且肯定只有部分而且肯定只有部分点落在拟合的直线上点落在拟合的直线上,但此时该直线和测试点的方差最小但此时该直线和测试点的方差最小.由拟由拟合直线的斜率就可以知道电阻的阻值合直线的斜率就可以知道电阻的阻值.拟合是探测
14、事物变化规律拟合是探测事物变化规律的办法的办法.插值就是根据函数上某些已知点插值就是根据函数上某些已知点(或实验数据或实验数据),按一定规律按一定规律(插插值方法值方法)寻求未知的点寻求未知的点,比如已知一个常用对数比如已知一个常用对数y=log(x)表表,是是按照按照x=0.1:0.1:10制表的制表的,如果按已知数据求如果按已知数据求y=log(2.897)就可以用插值得到就可以用插值得到.表制得越密表制得越密,插值越准确插值越准确15第15页,此课件共50页哦例例4.1.1-9 已知通过实验得到电压和电流的数值如下:已知通过实验得到电压和电流的数值如下:i0=0:0.5:10;v0=0,
15、98,202,306,395,506,592,708,789,890,1009,1108,1186,1310,1391,1510,1601,1716,1782,1920,2014;按一次按一次,二次二次,三次多项式对数值关系进行拟合三次多项式对数值关系进行拟合说明:说明:为保证较好的拟合效果,多项式阶数要取得适当,过低,残为保证较好的拟合效果,多项式阶数要取得适当,过低,残差较大,过高,拟合模型将包含噪声影响,通常要保证拟合差较大,过高,拟合模型将包含噪声影响,通常要保证拟合阶数小于数据对的数目。阶数小于数据对的数目。16第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第17页,此课件共50页哦p1=201.00,-2.89P2=0.31,197.90,2.02P3=0.03,-0.10,199.49,0.8518第18页,此课件共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第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第20页,此课件共50页哦方程组的解的三种情况方程组的解的三种情况:对于方程对于方程Ax=b,A为为Amn矩阵矩阵,有三种情况有三种情况:当当m=n时时,此方程成为此方程成为恰定恰定方程方程,求解精确解求解精确解 当当mn时时,此方程成为此方程成为“超定超定”方程方程,寻求寻求最小二乘解最小二乘解(直线拟合直线拟合)当当mn时此时无解时此时无解,但实用中可以求最小二乘解即但实用中可以求最小二乘解即:方程解方程解(AA)x=A b x=(AA)-1Ab-求逆法求逆法 x=Ab-matlab用最小二乘法找一个近似解用最小二乘法找一个近似解.x1+2x
22、2=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第22页,此课件共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第23页,此课件共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第24页,此课件共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第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第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分解分解,使使A*
28、V=V*J,J是对角阵是对角阵c=condeig(A)-计算各特征值的条件数计算各特征值的条件数V,D,C=condeig(A)-相当于相当于V,D=eig(A)和和c=condeig(A)两两条指令条指令27第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第28页,此课件
29、共50页哦4.3 数据分析函数数据分析函数(1)若输入宗量若输入宗量x是向量是向量,那么不管是行向量还是列向量那么不管是行向量还是列向量,运算是对整个向量运算是对整个向量进行的进行的(2)若输入宗量若输入宗量x是二维数组,那么指令运算是按列进行的是二维数组,那么指令运算是按列进行的MATLAB在进行数据分析时的约定在进行数据分析时的约定:29第29页,此课件共50页哦4.3.1 随机数发生器和统计分析指令随机数发生器和统计分析指令函数函数含义含义rand(m,n)产生产生mxn维的维的0,1区间均匀分布随机数组区间均匀分布随机数组randn(m,n)产生产生mxn维的均值为维的均值为0标准差为
30、标准差为1的正态分布随的正态分布随机数组机数组min(x)对矩阵各列分别求最小值对矩阵各列分别求最小值max(x)对矩阵各列分别求最大值对矩阵各列分别求最大值median(x)对矩阵各列分别求中位数对矩阵各列分别求中位数mean(x)对矩阵各列分别求均值对矩阵各列分别求均值std(x)对矩阵各列分别求标准差对矩阵各列分别求标准差var(x)对矩阵各列分别求方差对矩阵各列分别求方差C=conv(x)给出矩阵各列的协方差阵给出矩阵各列的协方差阵P=corrcoef(x)给出矩阵各列间的相关系数给出矩阵各列间的相关系数30第30页,此课件共50页哦例例4.3-1 基本统计示例基本统计示例randn(
31、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)各列方差各列方差,标准差:标准差:var(A),st
32、d(A)A=1,2,3;4,5,6;7,8,9 7 8 931第31页,此课件共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)沿第沿第n维求累计
33、乘积维求累计乘积cumsum(X,n)沿第沿第n维求累计和维求累计和cumtrapz(x,Y,n)用梯形法沿第用梯形法沿第n维求函数维求函数Y关于关于x的累计积分的累计积分XS,KK=sort(X,n)沿第沿第n维对维对X元素按模增大排列元素按模增大排列32第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(b)0 0 0 0
34、1.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第33页,此课件共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第34页,此课件共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.0060s4=8.0000
36、0135第35页,此课件共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第36页,此课件共50页哦q=quad(fun,a,b)-采用递推自适应采用递推自适应simpson法计算积分法计算积分q=quadl(fun,a,b)-采用递推自适应采用递推自适应Lobatto法求数值积分法求数值积分SS=dblquad(fun,inmin,inmax,outmin,outmax)-二重闭环数值积分二重闭环数值积分t,YY=o
38、de45(fun,tspan,Y0)-采用采用4,5阶阶Runge-Kutta方程解算方程解算ODE初初值问题值问题指令中被处理的函数指令中被处理的函数fun,可以取可以取:字符串表达式字符串表达式,内联函数内联函数,M函数文件的函函数文件的函数句柄数句柄37第37页,此课件共50页哦4.4.1 求函数零点求函数零点对于任意函数对于任意函数f(x)=0,零点情况复杂零点情况复杂,没有一个通用解法没有一个通用解法.一般来说一般来说,零点的数值计算过程是零点的数值计算过程是:先猜测一个初始零点或该零点所在的区间先猜测一个初始零点或该零点所在的区间,然后通过计算然后通过计算,使猜测值不断精确化使猜测
39、值不断精确化,或使猜测区间不断收缩或使猜测区间不断收缩,直到直到达到预定的精度达到预定的精度,终止计算终止计算.求函数的零点的步骤为求函数的零点的步骤为:(1)利用利用matlab作图指令获取初步近似值作图指令获取初步近似值具体做法具体做法:先确定一个零点可能存在的自变量区间先确定一个零点可能存在的自变量区间,然后利用然后利用plot指令画指令画出出f(x)在该区间中的图形在该区间中的图形,观察观察f(x)与横轴的交点坐标与横轴的交点坐标,也可以用也可以用zoom对对交点处局部放大再读数交点处局部放大再读数.借助借助ginput指令获得更精确的交点坐标值指令获得更精确的交点坐标值38第38页,
40、此课件共50页哦(2)利用利用fzero指令求精确解指令求精确解z=fzero(fun,x0)-求一元函数零点指令的最简格式求一元函数零点指令的最简格式z,f_z,exitflag=fzero(fun,x0,options,p1,p2,)-最完整格式最完整格式 说明:说明:由于由于fzero是根据函数是否越横轴来决定零点的是根据函数是否越横轴来决定零点的,因此本程序无法确因此本程序无法确定函数曲线仅触及横轴和不穿越的零点定函数曲线仅触及横轴和不穿越的零点 第二个输入量第二个输入量x0表示零点的初始猜测表示零点的初始猜测.他可以是标量或二元向量他可以是标量或二元向量.当当x0取标量时取标量时,该
41、指令将在它两侧寻找一个与之最靠近的零点该指令将在它两侧寻找一个与之最靠近的零点;当取二元当取二元向量向量a,b时时,该指令将在区间该指令将在区间a,b内寻找一个零点内寻找一个零点.输入宗量输入宗量options是优化迭代采用的参数是优化迭代采用的参数.P1,P2是向函数是向函数fun传递传递的参数的参数.z是所求零点的自变量值是所求零点的自变量值.f_z是函数值是函数值,exitflag表明程序终止计表明程序终止计算的条件算的条件.若若exitflag0,表明找到零点后退出表明找到零点后退出39第39页,此课件共50页哦例例4.4.1-1 求函数求函数f(x)=2x-10 x的零点的零点(1)
42、绘制函数图形绘制函数图形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第40页,此课件共50页哦func=2x-10*x;x01,err1,exitflag1=fzero(func,tt(1),);x02,err2,exitflag2=fzero(fun
43、c,tt(2),);(3)利用利用fzero指令求精确解指令求精确解str=方程的根为方程的根为:x1=,num2str(x01),和和 x2=,num2str(x02);disp(str);41第41页,此课件共50页哦4.4.2 求函数极值点求函数极值点许多科学研究和工程计算问题都可以归结为一个极值问题许多科学研究和工程计算问题都可以归结为一个极值问题,如能量最小,时间最如能量最小,时间最短短,最佳拟合最佳拟合,最短路径等最短路径等.又因为又因为f(x)的极小值问题等价于的极小值问题等价于-f(x)的极大值问题的极大值问题,所所以以matlab只有处理极小值的指令只有处理极小值的指令确切地
44、说确切地说,这里讨论的只是这里讨论的只是“局域极值局域极值”问题问题.“全域最小全域最小”问题要复杂得问题要复杂得多多.至今没有一个系统性的方法可求解一般的至今没有一个系统性的方法可求解一般的“全域最小全域最小”问题问题.对于一元对于一元二元函数二元函数,作图观察对确定全域最小值有很好的应用价值作图观察对确定全域最小值有很好的应用价值.但更多元的函数但更多元的函数,就很难使用作图法就很难使用作图法.matlab求函数极值的三条指令求函数极值的三条指令x,fval,exitflag,output=fminbnd(fun,x1,x2,opt,p1,p2,)-求一元求一元函数在区间函数在区间(x1,
45、x2)中极小值的指令的最完整格式中极小值的指令的最完整格式x,fval,exitflag,output=fminsearch(fun,x0,opt,p1,p2,)-单纯形法单纯形法求多元函数极值点指令的最完整格式求多元函数极值点指令的最完整格式x,fval,exitflag,output,grad,hessian=fminunc(fun,x0,opt,p1,p2,)-拟牛顿法求多元函数极值点指令的最完整格式拟牛顿法求多元函数极值点指令的最完整格式42第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第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-被积函数被积函数,a,b-积分上下限积分上下限tol
48、-绝对误差控制量,缺省时积分的绝对精度为绝对误差控制量,缺省时积分的绝对精度为1x10-6trace-取非时,将随积分的进程逐点画出被积函数取非时,将随积分的进程逐点画出被积函数44第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第45页,此课件共50页哦4.5 解常微分方程解常微分方程matlab为解常微分方程初值问题提供了一组配套齐全为解常微分方程初值问题提供了一组配套齐全,结构严整的指令结构严整的指令,包包括括:ode45
49、,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第46页,此课件共50页哦例例4.5-1 解微分方程:解微分方程:y=3x2,y(0)=1,微分区间微分
50、区间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第47页,此课件共50页哦4.6 数字滤波数字滤波y=filter(B,A,x,dim)-运用数字滤波对输入信号运用数字滤波对输入信号x滤波滤波滤波器的作用滤波器的作用:通过对采样数据信号进行数学运算处理来达到频域滤波的目的通过对采样数据信号进行数学运算处理来达到频域滤波的目的.例例 设计一个低通滤波器设计一个低通滤波器,从受噪声干扰的多频