《MATLAB解方程和函数极值.ppt》由会员分享，可在线阅读，更多相关《MATLAB解方程和函数极值.ppt（52页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、第三讲解方程与函数极值第三讲解方程与函数极值 方程问题和极值问题始终是数学问题中的核心问题！学习内容学习内容1多项式运算多项式运算2线性方程组求解线性方程组求解3非线性方程数值求解非线性方程数值求解4常微分方程初值问题的数值解法常微分方程初值问题的数值解法5函数极值函数极值6线性插值线性插值（一）多项式的表示方法（一）多项式的表示方法n对于多项式的表达式约定如下对于多项式的表达式约定如下n对于多项式对于多项式n对于上述多项式一般用以下行向量表示对于上述多项式一般用以下行向量表示（二）多项式的创建（二）多项式的创建1.系数向量直接输入法系数向量直接输入法由于MATLAB中多项式是以向量形式存储的
2、,因此最简单的多项式输入即向量输入。例:输入多项式p=1-5 6-33；poly2sym(p)%此函数将多项式转换换为符号多项式（二）多项式的创建（二）多项式的创建2.特征多项式输入法特征多项式输入法 多项式创建的另一途径是从矩阵求其特征多项式获得，由函数poly实现例如：例如：a=1 2 3;2 3 4;3 4 5p1=poly(a)%矩阵a对应的特征多项式poly2sym(p1)%将多项式p1转换为符号多 项式（二）多项式的创建（二）多项式的创建3.由根创立多项式由根创立多项式 给定的根也可产生其相应的多项式,此功能还是由函数poly实现例：root=-5-3+4i-3-4i;p=poly
3、(root)poly2sym(p)注：若要生成实系数多项式,则根中的复数比为对称共轭复数（三）多项式运算（三）多项式运算1.求多项式的值求多项式的值求多项式的值可以由两种形式(1)输入变量值代入多项式计算时是以数组为单元的,此时函数为polyval(2)以矩阵为计算单元,进行矩阵式运算以求的多项式的值,此时的函数为polyvalm注：这两种计算在数值上由很大差别。当进行矩阵运算时,变量矩阵须为方阵。实例演示实例演示p=1 11 55 125;b=1 1;1 1;c=5poly2sym(p)polyval(p,b)poly(p,c)polyvalm(p,b)（三）多项式运算（三）多项式运算2.求
4、多项式的根求多项式的根求多项式的根可以由两种方法(1)直接调用MATLAB的函数roots求解多项式的所有根(2)通过建立多项式的伴随矩阵再求其特征值的方法得到多项式的所有根实例演示实例演示n用两种方法解方程 的所有根p=2-5 6-1 9;roots(p)a=compan(p)%求多项式的特征矩阵eig(a)%求特征矩阵a的特征根（三）多项式运算（三）多项式运算3.多项式的乘除法多项式的乘除法(1)多项式的乘法由函数多项式的乘法由函数conv来实现来实现,此函数同用此函数同用于向量的卷积于向量的卷积(2)多项式的除法由函数多项式的除法由函数deconv来实现来实现,此函数与此函数与向量的解卷
5、函数相同向量的解卷函数相同实例演示实例演示n计算两多项式的乘除法np=2-5 6-1 9;npoly2sym(p);nd=3-90-18;nploy2sym(d)npd=conv(p,d)%多项式p与d相乘npoly2sym(pd)npl=deconv(pd,d)%多项式pd除以d（三）多项式运算（三）多项式运算4.多项式微分多项式微分多项式的微分可以用函数polyer进行例：p=2-5 6-1 9；poly2sym(p);Dp=polyer(p)（三）多项式运算（三）多项式运算5.多项式的拟合 多项式拟合是多项式运算的一个重要组成部分，在工程及科研工作中得到了广泛的应用，其一方面可以由矩阵的
6、除法解超定方程来进行;另一方面在MATLAB中还提供了专门的拟合函数polyfit,其调用格式如下：(1)polyfit(X,Y,n)%其中X,Y为拟合数据,n为你和多项式的阶数(2)p,s=polyfit(X,Y,n)%其中p为拟合多项式的系数向量，s为拟合多项式系数向量的结构信息实例演示实例演示例:用五阶多项式对0,pi/2上的正弦函数进行最小二乘拟合x=0:pi/20:pi/2;y=sin(x);a=ployfit(x,y,5);%用五阶多项式拟合y=sin(x),a为你和多项式的系数x1=0:pi/30:pi/2;y1=sin(x1);y2=a(1)*x1.5+a(2)*x1.4+a(
7、3)*x1.3+a(4)*x1.2+a(5)*x1.+a(6)plot(x1,y1,b-,x2,y2,r*)legend(原曲线，拟合曲线)axis(0,7,-1.2,4)%坐标轴的限制学习内容学习内容1多项式运算多项式运算2线性方程组求解线性方程组求解3非线性方程数值求解非线性方程数值求解4常微分方程初值问题的数值解法常微分方程初值问题的数值解法5函数极值函数极值6线性插值线性插值线性方程组解的结构n一、齐次线性方程组解的结构n齐次线性方程组的矩阵形式为AX=0n其中A是mn阶矩阵；X为未知向量二、线性方程组求解二、线性方程组求解1.直接解法直接解法（1）利用左除运算符的直接解法）利用左除运
8、算符的直接解法 对于线性方程组对于线性方程组Ax=b，可以利用左除运算符，可以利用左除运算符“”求解求解:x=Ab例例:用直接解法求解下列线性方程组。用直接解法求解下列线性方程组。命令如下：命令如下：A=2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4;b=13,-9,6,0;x=Ab二、线性方程组求解二、线性方程组求解 (2)利用矩阵的分解求解线性方程组利用矩阵的分解求解线性方程组 矩阵分解是指根据一定的原理用某种算法将一矩阵分解是指根据一定的原理用某种算法将一个矩阵分解成若干个矩阵的乘积。个矩阵分解成若干个矩阵的乘积。常见的矩阵分解有常见的矩阵分解有:特征值分解特征
9、值分解,奇异值分解奇异值分解,LU分解分解,QR分解分解,Cholesky分解分解,Schur分解分解 Hessenberg分解分解 等等 下面将对他们一一介绍下面将对他们一一介绍a)特征值分解特征值分解n矩阵的特征值分解也调用函数eig,其调用格式为n(1)V,D=eig(X)%得到矩阵X的特征值对角矩阵D和其列为对应特征值的特征向量矩阵V,矩阵的特征值分解为XV=VDn(2)V,D=eig(X,nobalance)%此形式为关闭平衡算法的求解方法.平衡算法对于某些问题可以得到更高的准确度。n(3)V,D=eig(A,B)%对矩阵A和B做广义特征值分解,即AV=BVD实例演示实例演示-矩阵的
10、特征值分解矩阵的特征值分解n单位阵的特征值分解na=-149-50-154;537 180 546;-27-9-25nv,d=eig(a)n双矩阵的特征值分解nb=2 10 2;10 5-8;2-8 11;nv,d=eig(a,b)附附:复数特征值对角阵与实数特征值对复数特征值对角阵与实数特征值对角阵的转换角阵的转换n即使对于实阵,其特征值也可能出现复数.而在实际使用中,常需要把这些共轭复数特征值转换为实数块.nMATLAB提供两个函数进行转换n(1)V,D=cdf2rdf(V,D)%将复数对角型转换为实数块对角型n(2)U,T=rdf2cdf(U,T)%将实数对角块转换为复数对角型实例演示实
11、例演示-复数对角型矩阵与复数对角型矩阵与实数对角型矩阵的转换实数对角型矩阵的转换na=1-3;2 2/3;nv,d=eig(a)%将矩阵a进行特征值分解,其结果复对角型nvs,ds=cdf2rdf(v,d)%将复对角型转换为实对角型nb=vs*ds/vs%其结果与a相等b)奇异值分解奇异值分解n设矩阵A为一个mn阶的(实)矩阵(mn),则存在正交矩阵V和U,使得:此式称为A的奇异值分解n矩阵奇异值分解可由函数svd实现,其调用格式为n(1)U,S,V=svd(X)n(2)U,S,V=svd(X,0)n作用:生成U,S,和V使得X=USV实例演示实例演示例：对矩阵a进行奇异值分解a=1;1;U,
12、S,V=svd(a)注：进行奇异值分解时,矩阵的维数没有限制c)LU分解分解LU分解是高斯消去法的基础分解是高斯消去法的基础,它将一矩阵表示为一个交它将一矩阵表示为一个交换下三角矩阵和一个上三角矩阵的乘积形式换下三角矩阵和一个上三角矩阵的乘积形式.线性代数线性代数中已经证明中已经证明,只要方阵只要方阵A是非奇异的是非奇异的,LU分解总是可以进分解总是可以进行的。行的。lu函数用于对矩阵进行函数用于对矩阵进行LU分解，其调用格式为分解，其调用格式为(1)L,U=lu(X)-产生一个上三角阵产生一个上三角阵U和一个变换形和一个变换形式的下三角阵式的下三角阵L(行交换行交换)，使之满足，使之满足X=
13、LU。注意。注意,这里这里的矩阵的矩阵X必须是方阵必须是方阵。(2)L,U,P=lu(X)-产生一个上三角阵产生一个上三角阵U和一个下三和一个下三角阵角阵L以及一个置换矩阵以及一个置换矩阵P，使之满足，使之满足PX=LU。当然矩。当然矩阵阵X同样必须是方阵同样必须是方阵。实现实现LU分解后，线性方程组分解后，线性方程组Ax=b的解的解x=U(Lb)或或x=U(LPb)，这样可以大大提高运算速度。，这样可以大大提高运算速度。实例演示实例演示例例:用用LU分解求解例分解求解例1中的线性方程组。中的线性方程组。命令如下：命令如下：A=2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1
14、,-4;b=13,-9,6,0;L,U=lu(A);X=U(Lb)或采用或采用LU分解的第分解的第2种格式，命令如下：种格式，命令如下：L,U,P=lu(A);X=U(LP*b)d)QR分解分解n在求解矩阵的特征值时在求解矩阵的特征值时,引入一种分解方法引入一种分解方法,即实阵即实阵A可以写成可以写成A=QR,其中其中Q为一个正交矩阵为一个正交矩阵和和R为一为一个上三角矩阵个上三角矩阵.QR分解只能对方阵进行分解只能对方阵进行.nQR分解由函数分解由函数qr实现，其调用格式为实现，其调用格式为:(1)Q,R=qr(X):产生一个一个正交矩阵产生一个一个正交矩阵Q和一个和一个上三角矩阵上三角矩阵
15、R,使之满足使之满足X=QR。(2)Q,R,E=qr(X):产生一个一个正交矩阵产生一个一个正交矩阵Q,一个一个上三角矩阵上三角矩阵R以及一个置换矩阵以及一个置换矩阵E,使之满足使之满足XE=QR。n实现实现QR分解后，线性方程组分解后，线性方程组Ax=b的解的解x=R(Qb)或或x=E(R(Qb)。实例演示实例演示例例:用用QR分解求解前例中的线性方程组。分解求解前例中的线性方程组。命令如下：命令如下：A=2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4;b=13,-9,6,0;Q,R=qr(A);x=R(Qb)或采用或采用QR分解的第分解的第2种格式，命令如下：种格
16、式，命令如下：Q,R,E=qr(A);x=E*(R(Qb)e)Cholesky分解分解n如果矩阵如果矩阵X是对称正定的是对称正定的,则则Cholesky分解将矩阵分解将矩阵X分解成一个分解成一个对角元素为正的下三角矩阵对角元素为正的下三角矩阵R和上三角矩阵和上三角矩阵R的乘积的乘积.使得使得X=RR.n函数函数chol(X)用于对矩阵用于对矩阵X进行进行Cholesky分解分解,其调用格式为其调用格式为:(1)R=chol(X):产生一个上三角阵产生一个上三角阵R,使使RR=X.若若X为非对称正定为非对称正定,则输出一个出错信息。则输出一个出错信息。(2)R,p=chol(X):这个命令格式将
17、不输出出错信息这个命令格式将不输出出错信息.当当X为对称为对称正定的正定的,则则p=0,R与上述格式得到的结果相同与上述格式得到的结果相同;否则否则p为一个正整为一个正整数数.如果如果X为满秩矩阵为满秩矩阵,则则R为一个阶数为为一个阶数为q=p-1阶的上三角阵阶的上三角阵,且满足且满足RR=X(1:q,1:q).n实现实现Cholesky分解后分解后,线性方程组线性方程组Ax=b变成变成RRx=b,所以所以x=R(Rb)。实例演示实例演示例例:用用Cholesky分解求解前例中的线性方程组。分解求解前例中的线性方程组。命令如下：命令如下：A=2,1,-5,1;1,-5,0,7;0,2,1,-1
18、;1,6,-1,-4;b=13,-9,6,0;R=chol(A)结果为结果为?Error using=cholMatrix must be positive definite%说明说明A为非正定为非正定矩阵矩阵?可以执行命令可以执行命令R,p=chol(X)x=R(Rb)2.迭代解法迭代解法n迭代解法非常适合求解大型系数矩阵的方程组迭代解法非常适合求解大型系数矩阵的方程组.在数值在数值分析中分析中,迭代解法主要包括迭代解法主要包括:Jacobi迭代法、迭代法、Gauss-Serdel迭代法、超松弛迭代法和两步迭代法。迭代法、超松弛迭代法和两步迭代法。1Jacobi迭代法迭代法:对于线性方程组对
19、于线性方程组Ax=b,如果如果A为非奇异为非奇异方阵方阵,即即aii0(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-1bn这就是这就是Jacobi迭代公式迭代公式.如果序列如果序列x(k+1)收敛于收敛于x，则，则x必是方程必是方程Ax=b的解的解附：附：Jacobi迭代法的迭代法的MATLAB函数文件函数文件
20、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例例:用用Jacobi迭代法求解下列线性方程组。设迭迭代法求解下列线性方程组。设迭代初值为代初值为0，迭代精度为，迭代精度为10-6。在命令中调用函数文件在命令中调用函数文件Jacobi.m，命令如下：，命令如下：A=10,-1,0;-1,10,-2;0,-2,10;b=9,7,6;x,n=jacobi(A,b,0,0,0,1.0e-6)！其他迭代法在此一律不在详细说明！其他迭代法
21、在此一律不在详细说明！学习内容1线性方程组求解线性方程组求解2非线性方程数值求解非线性方程数值求解3常微分方程初值问题的数值解法常微分方程初值问题的数值解法4函数极值函数极值5线性插值线性插值二、非线性方程组的求解二、非线性方程组的求解n对于非线性方程组对于非线性方程组F(X)=0,用用fsolve函数求其解函数求其解.fsolve函数的调用格式为函数的调用格式为:X=fsolve(fun,X0,option)n其中其中X为返回的解为返回的解,fun是用于定义需求解的非线性是用于定义需求解的非线性方程组的函数文件名方程组的函数文件名,X0是求根过程的初值是求根过程的初值,option为最优化工
22、具箱的选项设定为最优化工具箱的选项设定.最优化工具箱提供了最优化工具箱提供了20多个选项多个选项.用户可以使用用户可以使用optimset命令将它们显示出命令将它们显示出来来.如果想改变其中某个选项如果想改变其中某个选项.则可以调用则可以调用optimset()函数来完成函数来完成.例如例如,Display选项决定函数调用时中间选项决定函数调用时中间结果的显示方式结果的显示方式,其中其中off为不显示为不显示,iter表表示每步都显示示每步都显示,final只显示最终结果只显示最终结果.optimset(Display,off)将设定将设定Display选选项为项为off.例例:求下列非线性方
23、程组在求下列非线性方程组在(0.5,0.5)附近的数值解。附近的数值解。(1)建立函数文件建立函数文件myfun.m。function q=myfun(p)x=p(1);y=p(2);q(1)=x-0.6*sin(x)-0.3*cos(y);q(2)=y-0.6*cos(x)+0.3*sin(y);(2)在给定的初值在给定的初值x0=0.5,y0=0.5下下,调用调用fsolve函数求方函数求方程的根。程的根。x=fsolve(myfun,0.5,0.5,optimset(Display,off)x=0.6354 0.3734学习内容1线性方程组求解线性方程组求解2非线性方程数值求解非线性方程
24、数值求解3常微分方程初值问题的数值解法常微分方程初值问题的数值解法4函数极值函数极值5线性插值线性插值三、常微分方程初值问题的解法三、常微分方程初值问题的解法n一、解析解一、解析解 求解微风方程（组）的解析解的MATLAB命令为dsolve(eqn1,eqn2,x)其中eqni表示第i个方程与初始条件等式x表示微风方程（组）中的自变量，默认是自变量为t实例演示实例演示n例1：求解一阶微分方程 的通解及x=0,y=1时的特解 n求通解：dsolve(Dy=1+y2,x)n求特解：dsolve(Dy=1+y2,y(0)=1,x)n例2：求解二阶微分方程实例演示实例演示n输入命令：输入命令：ndso
25、lve(x2*D2y+x*Dy+(x2-(1/2)2)*y=0,y(pi/2)=2,Dy(pi/2)=-2/pi,x)n结果为：npi(1/2)*2(1/2)/x(1/2)*sin(x)n上机练习：求方程组二、数值解二、数值解n微分方程数值解的求法有多种,常用的主要是Runge Kutta(龙格-库塔法)和欧拉方法n基于龙格库塔法，MATLAB提供了使用2阶（3阶）或4阶（5阶）龙格库塔公式，一般调用格式为：n t,y=ode23(fname,tspan,y0，options)n t,y=ode45(fname,tspan,y0,options)n其中fname是定义f(t,y)的函数文件名，
26、ntspan的取法有几种，当tspan=t0,tf时,t0,tf分别表示自变量的初值和终值；当tspan=t0,t1,tf时则输出在指定时刻t0,t1,tf处给出；当tspan=t0:k:tf时则输出在t0,tf的等分点给出ny0是函数的初值。t,y为输出矩阵，分别表示自变量t和因变量y的取值实例演示实例演示n例求方程 的数值解n首先建立M函数文件：jie3.mnfunction f=jie3(x,y)nf=x2*D2y+x*Dy+(x2-(1/2)2)*y=0nx,y=ode23(jie3,pi/2,pi),2,-2/pi实例演示实例演示n用经典的R-K方法求解n解：编制函数文件：fun.m
27、nfunction f=fun（x，y）nf=-2*y+2*x.2+2*xn在MATLAB中输入nx,y=ode23(fun,0,0.5,1);学习内容1线性方程组求解线性方程组求解2非线性方程数值求解非线性方程数值求解3常微分方程初值问题的数值解法常微分方程初值问题的数值解法4函数极值函数极值5线性插值线性插值四、函数极值四、函数极值 MATLAB提供了求解无约束条件函数极值的命令提供了求解无约束条件函数极值的命令fminunc和和fminuncs，它们分别用于单变量函数和它们分别用于单变量函数和多变量函数的最小值，其调用格式为：多变量函数的最小值，其调用格式为：x=fminunc(fnam
28、e,x1,x2)x=fminuncs(fname,x0)注：这两个函数的调用格式相似。其中注：这两个函数的调用格式相似。其中fminunc函数函数用于求单变量函数的最小值点。用于求单变量函数的最小值点。fname是被最小化是被最小化的目标函数名，的目标函数名，x1和和x2限定自变量的取值范围。限定自变量的取值范围。fminuncs函数用于求多变量函数的最小值点，函数用于求多变量函数的最小值点，x0是是求解的初始值向量。求解的初始值向量。实例演示实例演示例：例：求求f(x)=x3-2x-5在在0,5内的最小值点。内的最小值点。(1)建立函数文件建立函数文件mymin.m。function fx=
29、mymin(x)fx=x.3-2*x-5;(2)调用调用fmin函数求最小值点。函数求最小值点。x=fminunc(mymin,0,5)x=0.8165学习内容1线性方程组求解线性方程组求解2非线性方程数值求解非线性方程数值求解3常微分方程初值问题的数值解法常微分方程初值问题的数值解法4函数极值函数极值5线性插值线性插值插值插值n插值的定义是对某些集合给定的数据点之间函数的估值方法。n当不能很快地求出所需中间点的函数时，插值是一个非常有价值的工具。nMatlab提供了一维、二维、三次样条等许多插值选择插值函数插值函数ntable1 ntable2 nintep1 ninterp2 nspline n利用已知点确定未知点n粗糙 精确n集合大的 简化的插值函数插值函数