matlab求解代数方程组

1、第三讲 Matlab求解代数方程组理论介绍：直接法+迭代法，简单介绍相关知识和应用条件及注意事项软件求解：各种求解程序讨论如下表示含有个未知数、由个方程构成的线性方程组： （1）一、直接法1.高斯消元法：高斯消元法的基本原理：在（1）中设将第一行乘以加到第得： （2）其中再设将（2）式的第二行乘以加到第行，如此进行下去最终得到： （3）从（3）式最后一个方程解出，代入它上面的一个方程解出，并如此进行下去，即可依次将全部解出，这样在的假设下，由上而下的消元由下而上的回代，构成了方程组的高斯消元法.高斯消元法的矩阵表示：若记，则（1）式可表为于是高斯消元法的过程可用矩阵表示为：其中：高斯消元法的M

2、atlab程序：%顺序gauss消去法,gauss函数functionA,u=gauss(a,n)for k=1:n-1%消去过程for i=k+1:n for j=k+1:n+1 %如果a(k,k)=0,则不能削去 if abs(a(k,k)>1e-6 %计算第k步的增广矩阵 a(i,j)=a(i,j)-a(i,k)/a(k,k)*a(k,j); else %a(k,k)=0,顺序gauss消去失败 disp(顺序gauss消去失败)； pause; exit; end end endend%回代过程x(n)=a(n,n+1)/a(n,n);for i=n-1:-1:1 s=0;for

3、 j=i+1:n s=s+a(i,j)*x(j);endx(i)=(a(i,n+1)-s)/a(i,i);end%返回gauss消去后的增广矩阵A=triu(a);%返回方程组的解u=x;练习和分析与思考： 用高斯消元法解方程组：2.列主元素消元法 在高斯消元法中进行到第步时，不论是否为0，都按列选择中最大的一个，称为列主元，将列主元所在行与第行交换再按高斯消元法进行下去称为列主元素消元法。列主元素消元法的matlab程序%列主元guass消去函数functionA,u=gauss(a,n)%消去过程for k=1:n-1 %选主元 c=0; for q=k:n if abs(a(q,k)&g

4、t;c c=a(q,k); l=q; end end%如果主元为0，则矩阵A不可逆 if abs(c)<1e-10 disp(error); pause; exit end%如果l不等于k，则交换第l行和第k行 if l=k for q=k:n+1 temp=a(k,q); a(k,q)=a(l,q); a(l,q)=temp; end end%计算第k步的元素值 for i=k+1:n for j=k+1:n a(i,j)=a(i,j)-a(i,k)/a(k,k)*a(k,j); end endend%回代过程 x(n)=a(n,n+1)/a(n,n); for i=n-1:-1:1

5、s=0; for j=i+1:n s=a+a(i,j)*x(j); end x(i)=(a(i,n+1)-s)/a(i,i); end%返回列主元gauss消去后的增广矩阵 A=triu(a);%返回方程组的解 u=x;练习和分析与思考： 用列主元消去法重新求解gauss消元法求解的上述问题。二、迭代法1.迭代法的总体思想（1）迭代公式的构造：对线性方程组，可以构造迭代公式给出由迭代公式的如果收敛于则就是原方程组的解。（2）矩阵的分解：设线性方程组，其中非奇异，则可以把矩阵分解：，于是化为与之对应的迭代公式为：2.雅可比（Jacobian）迭代法公式：用矩阵的元素表示为：Jacobian迭代法

6、的Matlab程序function y,n=jacobi(A,b,x0,eps)%误差if nargin=3eps=1.0e-6;elseif nargin<3 error returnend%求A的对角矩阵,下三角阵，上三角阵D=diag(A);diag(diag(A)?L=-tril(A,-1);U=-triu(A,1);B=D(L+U);f=Db;y=B*x0+f;%迭代的次数n=1;%当误差没有满足要求时继续迭代while norm(y-x0)>=epsx0=y;y=B*x0+f;n=n+1; end练习和分析与思考：利用Jacobian迭代法求解方程组：3.高斯-塞德尔（

7、Gauss-Seidel）迭代法公式： 将Jacobi迭代公式改进为，于是得到.用矩阵的元素表示为：Gauss-Seidel迭代法的Matlab程序function y,n=gauseidel(A,b,x0,eps)%误差if nargin=3eps=1.0e-6;elseif nargin<3 error returnend%求A的对角矩阵,下三角阵，上三角阵D=diag(A);diag(diag(A)?L=-tril(A,-1);U=-triu(A,1);G=(D-L)U;f=(D-L)b;y=G*x0+f;%迭代的次数n=1;%当误差没有满足要求时继续迭代while norm(y-

8、x0)>=epsx0=y;y=G*x0+f;n=n+1; end4.超松弛（Successive Over Relaxation method）迭代法（SOR）三、求解非线性方程的方法1.根的隔离二分法利用函数的性质确定根的大致范围取的中点若，则即是根，否则如果令如果令则在内至少有一根，且再取的中点，如此进行下去，包含根的区间的长度每次缩小一半，足够大时即可获得满意精度。2.切线法对于方程将在作Taylor展式保留线性项，即设然后用代替右端的就得到迭代公式：，这种方法称为牛顿（Newton）切线法。3.割线法在牛顿切线法中用差商代替，则（割线法）练习和分析与思考：求方程的一个正根。四、求

9、解非线性方程组的牛顿法方程组记作其中设且方程组的第步近似解，与牛顿切线法类似，在作Taylor展式，线性化后用代替可得：记则有，若可逆，我们得到解非线性方程组的牛顿迭代公式：.注：实际上在计算过程的第步，往往先计算和，再解方程组得到后，令即可.迭代过程需要分析其收敛性、分支与混沌.练习和分析与思考：用牛顿迭代法解非线性方程组：五、Matlab求解（非）线性方程组的内置工具箱1.线性方程组求解（1）直接求解：线性方程组，解法（注意：此处不是/符号）若为方阵，若不是方阵，练习和分析与思考：求方程组的解.（2）利用矩阵分解求解线性方程组矩阵分解是指根据一定的原理用某种算法将一个矩阵分解成若干个矩阵的

10、乘积.常见的矩阵分解有LU分解、QR分解、Cholesky分解以及Schur分解、Hessenberg分解、奇异分解等. LU分解:矩阵的LU分解就是将一个矩阵表示为一个交换下三角矩阵和一个上三角矩阵的乘积形式.只要方阵非奇异，矩阵的LU分解总是可以进行的.Matlab提供的lu函数用于对矩阵进行LU分解,其调用格式为：L,U=lu(A):产生一个上三角阵U和一个变换形式的下三角矩阵L（行交换）使之满足A=LU，注意：这里的矩阵必须是方阵.L,U,P=lu(A): 产生一个上三角阵U和一个变换形式的下三角矩阵L以及一个置换矩阵P使之满足PA=LU，注意：这里的矩阵必须是方阵.实现LU分解后，线

11、性方程组的解为：QR分解：矩阵的QR分解就是将一个矩阵A分解成一个正交矩阵Q和一个上三角矩阵R的乘积形式. QR分解只能对方阵进行. Matlab的函数qr用于对矩阵进行QR分解,其调用格式为：Q,R=qr(A):产生一个正交矩阵Q和一个上三角矩阵R使之满足A=QR. Q,R,E=qr(A): 产生一个正交矩阵Q和一个上三角矩阵R以及一个置换矩阵E使之满足AE=QR方阵.实现QR分解后，线性方程组的解为：Cholesky分解：如果矩阵是对称正定的,则Cholesky分解将矩阵分解成一个下三角矩阵R和上三角矩阵的乘积R(R的转置),即A= RR. Matlab的函数chol(A)用于对矩阵进行C

12、holesky分解,其调用格式为：R=chol(A):产生一个上三角矩阵R使之满足A=RR.若矩阵不是对称正定的，则输出一个出错信息. R,p=chol(A): 这个命令格式不输出出错信息. 当是对称正定的,则p=0,R与上述格式得到的结果相同，否则p为一个正整数，如果A为满秩矩阵，则R为一个阶数为q=p-1的上三角阵，且满足RR=A(1:q,1:q).实现Cholesky分解后， 变为，所以2.单变量非线性方程求解Matlab的函数fzero可用于对矩阵求单变量非线性方程的根,其调用格式为：z=fzero(fname,x0,tol,trace)其中fname是待求根的函数文件名，x0为搜索起

13、点，一个函数可能有多个根，但fzero函数只给出离x0最近的那个根.tol控制结果的相对精度,缺省时取tol=eps,trace指定迭代信息是否在运算中显示，为1时显示，为0时不显示，缺省时trace=0.练习和分析与思考：求在附近的根.3.非线性方程组的求解对于非线性方程组用fsolve函数求其数值解.fsolve函数的调用格式为：x=fsolve(fun,x0,option)其中x为返回的解，fun用于定义需求解的非线性方程组的函数文件名，x0是求根过程的初值，option为最优化工具箱的选项设定.最优化工具箱提供了20多个选项，用户可以使用optimset命令将它们显示出来.如果想改变其中某个选项，可以调用optimset()函数来完成.例如optimset(display,off).六、实例赏析a=1 -1 4 -2; 1 -1 -1 2;3 1 7 -2;1 -3 -12 6; rref(a)b=2 3 1;1 -2 4;3 8 -2;4 -1 9b1=4 -5