方程及方程组解法.ppt
第八讲第八讲 方程及方程组解法方程及方程组解法2n本讲教学目标本讲教学目标 掌握线性方程解法掌握线性方程解法 掌握线性方程组的直接解法掌握线性方程组的直接解法 掌握线性方程组的迭代解法掌握线性方程组的迭代解法 了解非线性方程和方程组的函数解法了解非线性方程和方程组的函数解法 了解非线性方程和方程组的迭代解法了解非线性方程和方程组的迭代解法3n线性方程组的求解不仅在工程技术领域涉及，线性方程组的求解不仅在工程技术领域涉及，而且在其他许多领域也经常碰到，因此它的应而且在其他许多领域也经常碰到，因此它的应用相当广泛。用相当广泛。n线性方程组的解法一般分为两类：线性方程组的解法一般分为两类：一类是直接一类是直接法，另一类是迭代法。法，另一类是迭代法。n非线性方程组的解法主要采用迭代法求解非线性方程组的解法主要采用迭代法求解，常，常用的有不动点迭代法、牛顿迭代法和用的有不动点迭代法、牛顿迭代法和Broyden迭代法等几种。迭代法等几种。4n8.1 线性方程及方程组的解法线性方程及方程组的解法 8.1.1 线性方程的解法线性方程的解法n通过调用函数通过调用函数 roots 求解。求解。例例1 求解方程求解方程 的根。的根。n p=1,-6,-72,-27n r=roots(p)p=1 -6 -72 -27 r=12.1229 -5.7345 -0.3884 5n8.1.2 线性方程组的解法线性方程组的解法 1.直接解法直接解法n(1)利用左除和右除运算符。利用左除和右除运算符。n通过这两个运算符，程序会自动根据输入的系数矩通过这两个运算符，程序会自动根据输入的系数矩阵判断选用哪种方法进行求解。阵判断选用哪种方法进行求解。n对线性方程组对线性方程组 Ax=b，利用左除运算符，利用左除运算符“”求解：求解：x=Abn右除同样。右除同样。6n例例2：求解以下方程组。求解以下方程组。a=1 2;2 3;b=8;13;x=ab x=2.0000 3.00007n(2)利用矩阵分解方法利用矩阵分解方法 LU分解分解n将非奇异的方阵分解为一个下三角矩阵将非奇异的方阵分解为一个下三角矩阵L和一个上三和一个上三角矩阵角矩阵U的乘积，即的乘积，即A=L*U。Ax=b 可变换为可变换为 x=U(Lb)或或 x=U(LPb)n分解格式：分解格式：L,U=lu(A)产生三角阵产生三角阵U和下三角阵和下三角阵L。L,U,P=lu(A)产生三角阵产生三角阵U和下三角阵和下三角阵L以以 及一个置换矩阵及一个置换矩阵P。8n QR分解：分解：把矩阵把矩阵A分解为一个正交矩阵分解为一个正交矩阵Q和一个上三角和一个上三角矩阵矩阵R的乘积，即的乘积，即A=Q*R。Ax=b 可变换为可变换为 x=R(Qb)或或 x=E(R(Qb)分解格式：分解格式：nQ,R=qr(A)产生正交阵产生正交阵Q和上三角阵和上三角阵R。nQ,R,E=qr(A)产生正交阵产生正交阵Q、上三角阵、上三角阵R及及置换阵置换阵E。9n Cholesky 分解：分解：把对称正定的矩阵把对称正定的矩阵A分解成一个上三角矩阵分解成一个上三角矩阵R和和其转置矩阵其转置矩阵R的乘积，即的乘积，即A=RR。Ax=b 可变成可变成 x=R(Rb)分解格式：分解格式：nR=chol(A)产生一个上三角阵产生一个上三角阵R。nR,p=chol(A)A对称正定时，对称正定时，p=0；否则；否则p为为一个正整数。如果一个正整数。如果A满秩，则满秩，则R为一个阶为为一个阶为p-1的上的上三角阵，且满足三角阵，且满足RR=A(1:p-1,1:p-1)。10n例例3：求解方程组。求解方程组。a=12-3 3;-16 3-1;1 1 1;b=15;-13;6;L,U=lu(a);x=U(Lb)x=1.0000 2.0000 3.0000 11n2.迭代解法迭代解法迭代解法适合求解大型系数矩阵的方程组。迭代解法适合求解大型系数矩阵的方程组。迭代解法主要包括：迭代解法主要包括：n Jacobi迭代法迭代法n Gauss-Serdel迭代法迭代法n 超松弛迭代法超松弛迭代法n 两步迭代法两步迭代法12n(1)Jacobi迭代法迭代法对线性方程组对线性方程组Ax=b，若，若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如果序列如果序列 x(k+1)收敛于收敛于 x，则，则 x 必是方程必是方程Ax=b 的解。的解。13nfunction y,n=jacobi(A,b,x0,eps)nif nargin=3n eps=1.0e-6;nelseif nargin=epsn x0=y;n y=B*x0+f;n n=n+1;nendnJacobi 迭代法的迭代法的Matlab 函数文件函数文件Jacobi.m如下：如下：14n例例4：用用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)15n(2)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迭代相比，迭代相比，Gauss-Serdel迭代用新迭代用新分量代替旧分量，精度会高些。分量代替旧分量，精度会高些。16nfunction y,n=gauseidel(A,b,x0,eps)nif nargin=3n eps=1.0e-6;nelseif nargin=epsn x0=y;n y=G*x0+f;n n=n+1;nendnGauss-Serdel迭代法迭代法的的Matlab函数文件函数文件gauseidel.m如下：如下：17n例例5：用用Gauss-Serdel迭代法求解线性方程组。迭代法求解线性方程组。设迭代初值为设迭代初值为0，迭代精度为，迭代精度为10-6。在命令中调用函数文件在命令中调用函数文件gauseidel.m，命令如下：，命令如下：A=10,-1,0;-1,10,-2;0,-2,10;b=9,7,6;x,n=gauseidel(A,b,0,0,0,1.0e-6)18n8.2 非线性方程及方程组的解法非线性方程及方程组的解法 8.2.1 调用函数求解调用函数求解n1.使用使用fzero函数求解，其格式为：函数求解，其格式为：z=fzero(fname,x0,tol,trace)n其中其中fname是函数文件名，是函数文件名，x0为搜索起点。一个函为搜索起点。一个函数可能有多个根，但只给出离数可能有多个根，但只给出离 x0 最近的根。最近的根。tol 控控制结果的相对精度，缺省取制结果的相对精度，缺省取eps，trace指定迭代信指定迭代信息是否在运算中显示，为息是否在运算中显示，为1时显示，为时显示，为0时不显示，时不显示，缺省时取缺省时取0。19n例例6：求求 f(x)=x-10 x+2=0 在在 0.5 附近的根。附近的根。步骤如下：步骤如下：x0=0.5;%初值初值 z=fzero(x-10 x+2,x0)z=0.3758 x1=0 1;%定义包含根的区间定义包含根的区间 z=fzero(x-10 x+2,x1)z=0.3758n函数调用式函数调用式如下：如下：n 建立函数文件建立函数文件funx.m。function fx=funx(x)fx=x-10 x+2;n 调用调用fzero函数求根。函数求根。z=fzero(funx,0.5)z=0.375820n2.使用使用fsolve函数求解函数求解 调用格式：调用格式：x=fsolve(fun,x0,option)n其中其中fun是函数文件名，是函数文件名，x0是初值，是初值，option为最优为最优化工具箱的选项设定，默认缺省。化工具箱的选项设定，默认缺省。n最优化工具箱提供了最优化工具箱提供了 20 多个选项，用户可以使用多个选项，用户可以使用optimset 命令将它们显示出来。命令将它们显示出来。n如果想改变其中某个选项，则可以调用如果想改变其中某个选项，则可以调用optimset()函数来完成。函数来完成。21n例例7：求解方程求解方程 在在0,10内内的解。的解。首先作图如下：首先作图如下：fplot(5*x2*sin(x)-exp(-x),0,0,10)发现在发现在0,10区间中有区间中有4个解，分别在个解，分别在0,3,6,9附近，附近，所以所以 用命令：用命令：fun=inline(5*x.2.*sin(x)-exp(-x);fsolve(fun,0 3 6 9,1e-6)ans=0.5018 3.1407 6.2832 9.424822n例例8：求非线性方程组在求非线性方程组在(0.5,0.5)附近的解。附近的解。首先编制函数文件首先编制函数文件fc.m如下：如下：function y=fc(x)y=x(1)-0.7*sin(x(1)-0.2*cos(x(2);x(2)-0.7*cos(x(1)+0.2*sin(x(2);然后在命令窗口输入：然后在命令窗口输入：x0=0.5 0.5;fsolve(fc,x0)或或 fsolve(fc,x0)n结果结果如下：如下：ans=0.5265 0.507923n8.2.2 其他解法其他解法 非线性方程非线性方程n (1)对分法对分法n (2)迭代法迭代法 非线性方程组非线性方程组n (1)不动点迭代法不动点迭代法n (2)牛顿迭代法牛顿迭代法n (3)Broyden迭代法迭代法n具体使用格式具体使用格式和方法参考教和方法参考教材材P86-94。24n例例 9 编制程序，对任何函数求根，在运行该编制程序，对任何函数求根，在运行该程序时将显示函数与零点的图形，用户在可能程序时将显示函数与零点的图形，用户在可能的零点附近用鼠标点击采点，程序通过循环计的零点附近用鼠标点击采点，程序通过循环计算每个采集点的根，并输出结果。算每个采集点的根，并输出结果。以函数以函数 求零点为求零点为例，综合叙述相关命令的用法。例，综合叙述相关命令的用法。25ny=inline(sin(t)2*exp(-a*t)-b*abs(t),t,a,b);na=0.1;b=0.5;nt=-10:0.01:10;%对自变量采样对自变量采样ny_char=vectorize(y);nY=feval(y_char,t,a,b);nclf%作人机对话界面作人机对话界面nplot(t,Y,r);hold on,plot(t,zeros(size(t),k);%画坐标横轴画坐标横轴nxlabel(t);ylabel(y(t),hold off ntitle(用鼠标在函数零点附近连续点五点用鼠标在函数零点附近连续点五点)26%循环使用命令循环使用命令fzero，求鼠标采集点附近的函数根。，求鼠标采集点附近的函数根。nfor i=1:5 t1,y1=ginput(1);ttt(i),yyy(i),exitflag=fzero(y,t1,a,b)%在在ttt(i)附近搜索附近搜索0点点nendnA(:,1)=ttt;nA(:,2)=yyy;nA n计算结果为：计算结果为：A=x f(x)-2.0074 0.0000 -2.0074 0.0000 -0.5198 0.0000 1.6738 0.0000 1.6738 -0.000027 本节完，谢谢！本节完，谢谢！
