matlab求解代數(shù)方程組

1、第三講 Matlab求解代數(shù)方程組理論介紹：直接法+迭代法，簡單介紹相關(guān)知識和應(yīng)用條件及注意事項軟件求解：各種求解程序討論如下表示含有個未知數(shù)、由個方程構(gòu)成的線性方程組： （1）一、直接法1.高斯消元法：高斯消元法的基本原理：在（1）中設(shè)將第一行乘以加到第得： （2）其中再設(shè)將（2）式的第二行乘以加到第行，如此進行下去最終得到： （3）從（3）式最后一個方程解出，代入它上面的一個方程解出，并如此進行下去，即可依次將全部解出，這樣在的假設(shè)下，由上而下的消元由下而上的回代，構(gòu)成了方程組的高斯消元法.高斯消元法的矩陣表示：若記，則（1）式可表為于是高斯消元法的過程可用矩陣表示為：其中：高斯消元法的M

2、atlab程序：%順序gauss消去法,gauss函數(shù)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;練習(xí)和分析與思考： 用高斯消元法解方程組：2.列主元素消元法 在高斯消元法中進行到第步時，不論是否為0，都按列選擇中最大的一個，稱為列主元，將列主元所在行與第行交換再按高斯消元法進行下去稱為列主元素消元法。列主元素消元法的matlab程序%列主元guass消去函數(shù)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;練習(xí)和分析與思考： 用列主元消去法重新求解gauss消元法求解的上述問題。二、迭代法1.迭代法的總體思想（1）迭代公式的構(gòu)造：對線性方程組，可以構(gòu)造迭代公式給出由迭代公式的如果收斂于則就是原方程組的解。（2）矩陣的分解：設(shè)線性方程組，其中非奇異，則可以把矩陣分解：，于是化為與之對應(yīng)的迭代公式為：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;%迭代的次數(shù)n=1;%當(dāng)誤差沒有滿足要求時繼續(xù)迭代while norm(y-x0)>=epsx0=y;y=B*x0+f;n=n+1; end練習(xí)和分析與思考：利用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;%迭代的次數(shù)n=1;%當(dāng)誤差沒有滿足要求時繼續(xù)迭代while norm(y-

8、x0)>=epsx0=y;y=G*x0+f;n=n+1; end4.超松弛（Successive Over Relaxation method）迭代法（SOR）三、求解非線性方程的方法1.根的隔離二分法利用函數(shù)的性質(zhì)確定根的大致范圍取的中點若，則即是根，否則如果令如果令則在內(nèi)至少有一根，且再取的中點，如此進行下去，包含根的區(qū)間的長度每次縮小一半，足夠大時即可獲得滿意精度。2.切線法對于方程將在作Taylor展式保留線性項，即設(shè)然后用代替右端的就得到迭代公式：，這種方法稱為牛頓（Newton）切線法。3.割線法在牛頓切線法中用差商代替，則（割線法）練習(xí)和分析與思考：求方程的一個正根。四、求

9、解非線性方程組的牛頓法方程組記作其中設(shè)且方程組的第步近似解，與牛頓切線法類似，在作Taylor展式，線性化后用代替可得：記則有，若可逆，我們得到解非線性方程組的牛頓迭代公式：.注：實際上在計算過程的第步，往往先計算和，再解方程組得到后，令即可.迭代過程需要分析其收斂性、分支與混沌.練習(xí)和分析與思考：用牛頓迭代法解非線性方程組：五、Matlab求解（非）線性方程組的內(nèi)置工具箱1.線性方程組求解（1）直接求解：線性方程組，解法（注意：此處不是/符號）若為方陣，若不是方陣，練習(xí)和分析與思考：求方程組的解.（2）利用矩陣分解求解線性方程組矩陣分解是指根據(jù)一定的原理用某種算法將一個矩陣分解成若干個矩陣的

10、乘積.常見的矩陣分解有LU分解、QR分解、Cholesky分解以及Schur分解、Hessenberg分解、奇異分解等. LU分解:矩陣的LU分解就是將一個矩陣表示為一個交換下三角矩陣和一個上三角矩陣的乘積形式.只要方陣非奇異，矩陣的LU分解總是可以進行的.Matlab提供的lu函數(shù)用于對矩陣進行LU分解,其調(diào)用格式為：L,U=lu(A):產(chǎn)生一個上三角陣U和一個變換形式的下三角矩陣L（行交換）使之滿足A=LU，注意：這里的矩陣必須是方陣.L,U,P=lu(A): 產(chǎn)生一個上三角陣U和一個變換形式的下三角矩陣L以及一個置換矩陣P使之滿足PA=LU，注意：這里的矩陣必須是方陣.實現(xiàn)LU分解后，線

11、性方程組的解為：QR分解：矩陣的QR分解就是將一個矩陣A分解成一個正交矩陣Q和一個上三角矩陣R的乘積形式. QR分解只能對方陣進行. Matlab的函數(shù)qr用于對矩陣進行QR分解,其調(diào)用格式為：Q,R=qr(A):產(chǎn)生一個正交矩陣Q和一個上三角矩陣R使之滿足A=QR. Q,R,E=qr(A): 產(chǎn)生一個正交矩陣Q和一個上三角矩陣R以及一個置換矩陣E使之滿足AE=QR方陣.實現(xiàn)QR分解后，線性方程組的解為：Cholesky分解：如果矩陣是對稱正定的,則Cholesky分解將矩陣分解成一個下三角矩陣R和上三角矩陣的乘積R(R的轉(zhuǎn)置),即A= RR. Matlab的函數(shù)chol(A)用于對矩陣進行C

12、holesky分解,其調(diào)用格式為：R=chol(A):產(chǎn)生一個上三角矩陣R使之滿足A=RR.若矩陣不是對稱正定的，則輸出一個出錯信息. R,p=chol(A): 這個命令格式不輸出出錯信息. 當(dāng)是對稱正定的,則p=0,R與上述格式得到的結(jié)果相同，否則p為一個正整數(shù)，如果A為滿秩矩陣，則R為一個階數(shù)為q=p-1的上三角陣，且滿足RR=A(1:q,1:q).實現(xiàn)Cholesky分解后， 變?yōu)椋?.單變量非線性方程求解Matlab的函數(shù)fzero可用于對矩陣求單變量非線性方程的根,其調(diào)用格式為：z=fzero(fname,x0,tol,trace)其中fname是待求根的函數(shù)文件名，x0為搜索起

13、點，一個函數(shù)可能有多個根，但fzero函數(shù)只給出離x0最近的那個根.tol控制結(jié)果的相對精度,缺省時取tol=eps,trace指定迭代信息是否在運算中顯示，為1時顯示，為0時不顯示，缺省時trace=0.練習(xí)和分析與思考：求在附近的根.3.非線性方程組的求解對于非線性方程組用fsolve函數(shù)求其數(shù)值解.fsolve函數(shù)的調(diào)用格式為：x=fsolve(fun,x0,option)其中x為返回的解，fun用于定義需求解的非線性方程組的函數(shù)文件名，x0是求根過程的初值，option為最優(yōu)化工具箱的選項設(shè)定.最優(yōu)化工具箱提供了20多個選項，用戶可以使用optimset命令將它們顯示出來.如果想改變其中某個選項，可以調(diào)用optimset()函數(shù)來完成.例如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