matlab求解代数方程组 matlab求解代数方程组.doc

matlab求解代数方程组.doc

1第三讲Matlab求解代数方程组理论介绍：直接法+迭代法，简单介绍相关知识和应用条件及注意事项软件求解：各种求解程序讨论如下表示含有个未知数、由个方程构成的线性方程组：nn121212nnnaxaxb(1)一、直接法1.高斯消元法：高斯消元法的基本原理：在(1)中设将第一行乘以加到第得：10,a1,ka(2,3),kn()(1)(1)(1)122(2)(2)()nnnnxaxba(2)其中再设将(2)式的第二行乘以加(1)(1),.kab()0,a(2),3,)ka到第行，如此进行下去最终得到：(1)()(1)(1)221()(1)(1),,nnnnaxaxbaxb(3)从(3)式最后一个方程解出，代入它上面的一个方程解出，并如此进行nx1nx下去，即可依次将全部解出，这样在的假设121,,n()0,2)ka2下，由上而下的消元由下而上的回代，构成了方程组的高斯消元法.高斯消元法的矩阵表示：若记，则(1)式可表为于是高11(),(,),(,)TTijnnnAaxb.Axb斯消元法的过程可用矩阵表示为：121121.nnMAxMb其中：(1)21(1)naMa(2)32(2)11naMa高斯消元法的Matlab程序：%顺序gauss消去法,gauss函数function[A,u]=gauss(a,n)fork=1:n-1%消去过程fori=k+1:nforj=k+1:n+1%如果a(k,k)=0,则不能削去ifabs(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;3endendendend%回代过程x(n)=a(n,n+1)/a(n,n);fori=n-1:-1:1s=0;forj=i+1:ns=s+a(i,j)*x(j);endx(i)=(a(i,n+1)-s)/a(i,i);end%返回gauss消去后的增广矩阵A=triu(a);%返回方程组的解u=x;练习和分析与思考：用高斯消元法解方程组：1234512457184693xxxx2.列主元素消元法在高斯消元法中进行到第步时，不论是否为0，都按列选择k()kia中最大的一个，称为列主元，将列主元所在行与第行交换再()|,)kiank按高斯消元法进行下去称为列主元素消元法。列主元素消元法的matlab程序4%列主元guass消去函数function[A,u]=gauss(a,n)%消去过程fork=1:n-1%选主元c=0;forq=k:nifabs(a(q,k))>cc=a(q,k);l=q;endend%如果主元为0，则矩阵A不可逆ifabs(c)=epsx0=y;y=B*x0+f;n=n+1;end练习和分析与思考：利用Jacobian迭代法求解方程组：1324658x3.高斯-塞德尔(Gauss-Seidel)迭代法公式：将Jacobi迭代公式改进为，(1)()kkDxLUxb(1)(1)()kkkDxLUxb于是得到.(1)()1kx1,()Bf用矩阵的元素表示为：A1()()1(1),2,)inkkjijjkiibaxxnGauss-Seidel迭代法的Matlab程序function[y,n]=gauseidel(A,b,x0,eps)%误差ifnargin==3eps=1.0e-6;elseifnargin=epsx0=y;y=G*x0+f;n=n+1;end4.超松弛(SuccessiveOverRelaxation)迭代法(SOR)三、求解非线性方程的方法1.根的隔离——二分法利用函数的性质确定根的大致范围取的中点若(,)ab(,)0,2abx，则即是根，否则如果令如果0()fx0 x0fx1,令则在内至少有一根，且再,b10,a1,b1(,)1(),)ab取的中点，如此进行下去，包含根的区间的1(,)a2x,2n长度每次缩小一半，足够大时即可获得满意精度。n2.切线法对于方程将在作Taylor展式保留线性项，即()0fx()fxk )(),kkffx设然后用代替右端的就得到迭代公式：，这种 ()0,kfx1kx1 ()kkfxx9方法称为牛顿(Newton)切线法。3.割线法在牛顿切线法中用差商代替，则1())kkfxf ()kfx(割线法)11())kkkfxf练习和分析与思考：求方程的一个正根。2()4fx四、求解非线性方程组的牛顿法方程组记作其中设()0,F112(,),(),(,(),TTnnxFxfxf且方程组的第步近似解，与牛顿切线法类似，在()()()1,kkkTnxx)0k作Taylor展式，线性化后用代替可得：()k(1)kx()()(1)()()()(1)()1,2)kkkkkkiiiinffxfxfxin记1122 112()nnnffxxFffxx则有，若可逆，我们得到解非线(1)() ()()()kkkkFxF ()kF性方程组的牛顿迭代公式：.(1)() ()1()[]xx注：实际上在计算过程的第步，往往先计算和，再解方程k()k ()kx组得到后，令即可.迭代过程需要 ()()()kkFxFxA()A(1)()()kxA分析其收敛性、分支与混沌.练习和分析与思考：10用牛顿迭代法解非线性方程组：22132070 xx五、Matlab求解(非)线性方程组的内置工具箱1.线性方程组求解(1)直接求解：线性方程组，解法(注意：此处\不是/符号)Axb,xAb若为方阵，，若不是方阵，A()*xinv()*.pinv练习和分析与思考：求方程组的解.12324354601xx(2)利用矩阵分解求解线性方程组矩阵分解是指根据一定的原理用某种算法将一个矩阵分解成若干个矩阵的乘积.常见的矩阵分解有LU分解、QR分解、Cholesky分解以及Schur分解、Hessenberg分解、奇异分解等.LU分解:矩阵的LU分解就是将一个矩阵表示为一个交换下三角矩阵和一个上三角矩阵的乘积形式.只要方阵非奇异，矩阵的LU分解总是可以进行的.AMatlab提供的lu函数用于对矩阵进行LU分解,其调用格式为：[L,U]=lu(A):产生一个上三角阵U和一个变换形式的下三角矩阵L(行交换)使之满足A=LU，注意：这里的矩阵必须是方阵.[L,U,P]=lu(A):产生一个上三角阵U和一个变换形式的下三角矩阵L以及一个置换矩阵P使之满足PA=LU，注意：这里的矩阵必须是方阵.A实现LU分解后，线性方程组的解为：Axb\()\().xLbxUPb或QR分解：矩阵的QR分解就是将一