數(shù)值線性代數(shù)課程設計報告

燕山大學課程設計題目：Gauss消去法及其變形方法求解線性方程組學院〔系〕：年級專業(yè)：學號：學生姓名：指導教師：教師職稱：燕山大學課程設計〔論文〕任務書院〔系〕：理學院基層教學單位：信息與計算科學系學號學生姓名專業(yè)〔班級〕設計題目Gauss消去法及其變形方法求解線性方程組設計技術參數(shù)Gauss消去法是一個古老的求解線性方程組的方法，但由于它的改良、變形得到的平方根法和改良的平方根法仍是在計算機上求解系數(shù)矩陣為中、低階稠密矩陣的線性方程組常用的有效方法，本文對這幾種方法進行了比擬。

用高斯消去法求解n階線性方程組Ax=b的根本思想是在逐步消元的過程中，把方程組的系數(shù)矩陣化為上三角矩陣，從而將原方程組約化為容易求解的等價三角方程組，然后進行回代求解。

設計要求通過查閱資料，了解課題在實際中的應用；運用所學知識，建立問題的數(shù)學模型；利用MATLAB軟件進行計算，并給出計算機程序；運用所學知識，對計算結果進行分析；要閱讀不少于四篇的相關文獻。工作量本課程設計要在一周內(nèi)完成；報告的正文在4000字左右。工作計劃第1-2天:閱讀參考文獻,收集資料.第3-5天:構建數(shù)學模型,編制計算機程序.第5-6天:上機運算,調試程序,計算結果分析.第7天:撰寫報告.參考資料[1]張德豐.MATLAB數(shù)值分析與應用[M].北京：國防工業(yè)出版社.2007.1[2]李慶揚，王能超，易大義.數(shù)值分析[M].武漢：華中科技大學出版社.2006.7[3]蕭樹鐵.數(shù)學實驗[M].北京：高等教育出版社.1990[4]萊〔LayD.C.〕，劉深泉.數(shù)值線性代數(shù)及其應用[M].機械工業(yè)出版社.2005.8指導教師簽字基層教學單位主任簽字說明：此表一式四份，學生、指導教師、基層教學單位、系部各一份。年月日燕山大學課程設計評審意見表指導教師評語：成績：指導教師：年月日辯論小組評語：成績：評閱人：年月日課程設計總成績：辯論小組成員簽字：年月日目錄摘要………………2一、問題的提出………………2二、數(shù)學模型的建立………3三、計算機程序………………3四、計算結果及分析………5五、參考文獻…………………20六、心得體會…………………20

摘要： Gauss消去法是求解線性方程組的最根本的方法之一。實際計算中，線性方程組的系數(shù)矩陣常常具有對稱正定性，即其各階順序主子式及全部特征值均大于零。矩陣的這一特性使它的三角分解也有更簡單的形式，從而導出一些特殊的解法，如平方根法和改良的平方根法。關鍵詞：線性方程組，Gauss消去法，平方根法，改良的平方根法，Hilbert矩陣。一、問題的提出1、先用你所熟悉的計算機語言將不選主元和列主元Gauss消去法編寫成通用程序，然后用你編寫的程序求解下面的84階方程組最后，將你的計算結果與方程組的精確解進行比擬，并就此談談你對Gauss消去法的看法。2、先用你所熟悉的計算機語言將平方根法和改良的平方根法編寫成通用程序，然后用你編寫的程序求解對稱正定方程組，其中〔1〕隨機地選取，系數(shù)矩陣為100階矩陣〔2〕系數(shù)矩陣為40階Hilbert矩陣，即系數(shù)矩陣的第行第列元素為，向量的第個分量為。3、用1中的程序求解2中的兩個方程組并比擬所有的計算結果，然后評論各個方法的優(yōu)劣。。二、數(shù)學模型的建立１、Gauss消去法：用消去法解方程組的根本思想是用逐次消去未知數(shù)的方法把原來方程組Ax=b化為與其等價的三角方程組，而求解三角方程組就容易了。換句話說，上述過程就是用行的初等變換將原方程組系數(shù)矩陣化為簡單形式，從而將求解原方程組的問題轉化為求解簡單方程組的問題。２、平方根法：平方根法就是利用對稱正定矩陣的三角分解而得到的求解對稱正定方程組的一種有效方法。平方根法遞推公式可以證明對于對稱正定矩陣A,可以唯一地分解成A=LLT,其中L是非奇異下三角形矩陣。３、改良的平方根法：平方根算法的計算量約為普通三角分解算法的一半，但由于要用到開平方，效率不是很高，所以，便有了為效率而存在的改良版平方根算法。對Ａ進行ＬＵ分解，即A=LU，Ａ對稱正定，可得ｕｉｉ＞０，那么可得A=LDLT分解計算公式，Ｌ是對角元素為１的單位下三角矩陣，Ｄ為對角矩陣。三、計算機程序模型一：利用Gauss消去法對線性方程組Ax=b進行求解。用MATLAB建立m文件DelGauss.m,程序如下：functionx=DelGauss(a,b)[n,m]=size(a);nb=length(b);det=1;x=zeros(n,1);fork=1:n-1fori=k+1:nifa(k,k)==0returnendm=a(i,k)/a(k,k);forj=k+1:na(i,j)=a(i,j)-m*a(k,j);endb(i)=b(i)-m*b(k);enddet=det*a(k,k);enddet=det*a(n,n);fork=n:-1:1forj=k+1:nb(k)=b(k)-a(k,j)*x(j);endx(k)=b(k)/a(k,k);end模型二：利用平方根法對線性方程組Ax=b進行求解。用MATLAB建立m文件pingfg.m,程序如下：function[x]=pingfg(A,b)%喬累斯基分解[n,n]=size(A);L=zeros(n,n);%實際上不用為L申請空間，使用A即可L(1,1)=sqrt(A(1,1));fork=2:nL(k,1)=A(k,1)/L(1,1);endfork=2:n-1L(k,k)=sqrt(A(k,k)-sum(L(k,1:k-1).^2));fori=k+1:nL(i,k)=(A(i,k)-sum(L(i,1:k-1).*L(k,1:k-1)))/L(k,k);endendL(n,n)=sqrt(A(n,n)-sum(L(n,1:n-1).^2));%解下三角方程組Ly=by=zeros(n,1);fork=1:nj=1:k-1;y(k)=(b(k)-L(k,j)*y(j))/L(k,k);end%解上三角方程組L'x=yx=zeros(n,1);U=L';fork=n:-1:1j=k+1:n;x(k)=(y(k)-U(k,j)*x(j))/U(k,k);end模型三：利用改良的平方根法對線性方程組Ax=b進行求解。用MATLAB建立m文件ave.m,程序如下：function[x]=ave(A,b,n)%用改良平方根法求解Ax=bL=zeros(n,n);%L為n*n矩陣D=diag(n,0);%D為n*n的主對角矩陣S=L*D;fori=1:n%L的主對角元素均為1L(i,i)=1;endfori=1:nforj=1:n%驗證A是否為對稱正定矩陣if(eig(A)<=0)|(A(i,j)~=A(j,i))%A的特征值小于0或A非對稱時，輸出wrongdisp('wrong');break;endendendHilbert矩陣用MATLAB建立m文件Hil.m,程序如下：functionb=Hil()fork=1:40form=1:40s=0;t=s+1/(k+m-1);s=t;endb(k,1)=s;end四、計算結果及分析問題1：>>A=ones(1,83);>>B=6*ones(1,84);>>C=8*ones(1,83);>>D=diag(A,1)+diag(B,0)+diag(C,-1);>>b=[7;15*ones(82,1);14];>>DelGauss(D,b)ans=1.0e+008*0.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00010.0002-0.00030.0007-0.00130.0026-0.00520.0105-0.02090.0419-0.08360.1665-0.33030.6501-1.25822.3487-4.02635.3684方程組的精確解為x1=x2=…=x84=1.0000，與Gauss消去法求得的解差距很大，所得結果不夠準確，計算簡單但其消元過程有時不能進行到底而使求解出現(xiàn)解失真的情況。問題2：>>A=ones(1,99);B=10*ones(1,100);C=diag(A,1)+diag(B,0)+diag(A,-1);b=ones(100,1);pingfg(C,b)ans=0.09180.08250.08340.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08340.08250.0918>>A=ones(1,99);B=10*ones(1,100);C=diag(A,1)+diag(B,0)+diag(A,-1);b=ones(100,1);>>ave(C,b,100)ans=0.09180.08250.08340.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08340.08250.0918A=hilb(40);b=Hil();pingfg(A,b)ans=1.0e+007*0.0000-0.00000.0001-0.00190.0171-0.09390.3305-0.76851.1892-1.20770.7761-0.30260.06840.0135-0.0220-0.04500.0954-0.06750.0217-0.0032-0.00060.0013-0.00090.00020.0009-0.00150.0009-0.00060.0005-0.00020.0001-0.00010.00000.0000-0.0000-0.00000.0000-0.0000-0.00000.0000>>A=hilb(40);b=Hil();ave(A,b,40)ans=-0.00000.0000-0.00010.0010-0.00880.0459-0.14160.2546-0.25990.2041-0.28120.2680-0.10590.5017-0.7878-0.13990.36920.5230-0.23950.0245-0.4819-0.20150.6467-0.2634-0.02930.4175-0.43350.4124-0.4021-0.14970.15800.5615-0.2600-0.72070.34570.09950.27410.2579-0.79171.3333問題3：>>A=ones(1,99);B=10*ones(1,100);C=diag(A,1)+diag(B,0)+diag(A,-1);b=ones(100,1);>>DelGauss(C,b)ans=0.09180.08250.08340.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.08330.0833