線性方程組數(shù)值求解

1、 數(shù)值分析課程設(shè)計(jì)線性方程組數(shù)值求解算法研究院（系）名稱 專 業(yè) 班 級(jí) 學(xué) 號(hào) 學(xué) 生 姓 名 指 導(dǎo) 教 師 摘 要本課程設(shè)計(jì)用matlab就線性方程組數(shù)值方法，Jacobi迭代法，Gauss-Seidel迭代法，超松弛法對(duì)所設(shè)計(jì)的問題進(jìn)行求解，并編寫程序在Matlab中實(shí)現(xiàn)，在文章中對(duì)各種迭代法進(jìn)行了收斂性分析。接著用幾種不同方法對(duì)線性方程組進(jìn)行求解及結(jié)果分析，最后對(duì)此次課程設(shè)計(jì)進(jìn)行了總結(jié)。 關(guān)鍵詞：線性方程組，迭代，Matlab，結(jié)果分析一、 提出問題直接法得到的解是理論上準(zhǔn)確的，但是我們可以看得出，它們的計(jì)算量都是n3數(shù)量級(jí)，存儲(chǔ)量為n2量級(jí)，這在n比較小的時(shí)候還比較合適（n<

2、1000 ），但是對(duì)于現(xiàn)在的很多實(shí)際問題，往往要我們求解很大的n的矩陣，而且這些矩陣往往是系數(shù)矩陣就是這些矩陣含有大量的0元素。對(duì)于這類的矩陣，在用直接法時(shí)就會(huì)耗費(fèi)大量的時(shí)間和存儲(chǔ)單元。因此我們有必要引入一類新的方法：迭代法。迭代法具有的特點(diǎn)是速度快。與非線性方程的迭代方法一樣，需要我們構(gòu)造一個(gè)等價(jià)的方程，從而構(gòu)造一個(gè)收斂序列，序列的極限值就是方程組的根,如何建立迭代格式，向量序列Xk是否收斂條件等。二、 算法的基本思想線性方程組的數(shù)值法a) 直接法直接法是指在無舍入誤差存在的情況下，經(jīng)過有限步運(yùn)算即可求得方程組解的算法，因此，又稱為精確法。這種方法是通過矩陣約化將原方程組化為與之等價(jià)的三角形

3、方程組或其他形式的可直接求解的方程組而實(shí)現(xiàn)的。由于實(shí)際中存在舍入誤差的影響，這種方法只能求得線性方程組的近似解。通常用來解低階稠密矩陣方程組及某些大型稀疏矩陣方程組比較著名的方法有：Gauss消去法，主元消去法，LU分解法，QR分解法，Cholesky分解法等。b) 迭代法迭代法是用某種極限過程去逐步逼近線性方程組精確解的方法。迭代法具有需要計(jì)算機(jī)的存儲(chǔ)單元較少，程序設(shè)計(jì)簡單和原始系數(shù)矩陣在計(jì)算過程中始終不變等優(yōu)點(diǎn)，但存在收斂性與收斂速度問題。迭代法是解大型稀疏矩陣方程組，尤其是由微分方程離散后得到的大型方程組的重用方法。比較著名的方法有：jacobi迭代法，Gauss-seidel迭代法，超

4、松弛迭代法，梯度法，理查森迭代法。與解f(x)=0 的不動(dòng)點(diǎn)迭代相類似,將AX=b改寫為X=BX+f 的形 式,建立雅可比方法的迭代格式： k+1=BX(k)+f ,其中，B稱為迭代矩陣.其計(jì)算精度可控,特別適用于求解系數(shù)為大型稀疏矩陣(sparse matrices)的方程組.三、 算法的推導(dǎo)及步驟對(duì)方程組做等價(jià)變換則，我們可以構(gòu)造序列若同時(shí)：所以，序列收斂與初值的選取無關(guān)定義：（收斂矩陣）定理： 矩陣G為收斂矩陣，當(dāng)且僅當(dāng)G的譜半徑<1由知，若有某種范數(shù) 則，迭代收斂1. Jacobi迭代格式很簡單：2. GaussSeidel迭代在Jacobi迭代中，使用最新計(jì)算出的分量值3. 松

5、弛迭代(SOR)記則可以看作在前一步上加一個(gè)修正量。若在修正量前乘以一個(gè)因子 有對(duì)Gauss-Seidel迭代格式寫成分量形式，有四、 算法分析（穩(wěn)定.收斂.誤差）1. 算法及收斂性：1.1. Jacobi迭代矩陣記易知，Jacobi迭代有收斂條件：迭代格式收斂的充要條件是G的譜半徑<1。對(duì)于Jacobi迭代，我們有一些保證收斂的充分條件若A滿足下列條件之一，則Jacobi迭代收斂。 A為行對(duì)角占優(yōu)陣 A為列對(duì)角占優(yōu)陣 A滿足1.2. GaussSeidel迭代矩陣GaussSeidel收斂條件：迭代格式收斂的充要條件是G的譜半徑<1。我們看一些充分條件若A滿足下列條件之一，則Ga

6、uss-Seidel迭代收斂。 A為行或列對(duì)角占優(yōu)陣 A對(duì)稱正定陣1.3. 松弛迭代(SOR)松弛迭代收斂A對(duì)稱正定，則松弛迭代收斂SOR方法收斂的快慢與松弛因子w的選擇有密切關(guān)系.但是如何選取最佳松弛因子,即選取w=w*,使r(Gw )達(dá)到最小,是一個(gè)尚未很好解決的問題.實(shí)際上可采用試算的方法來確定較好的松弛因子.經(jīng)驗(yàn)上可取1.4<w<1.6.2. 穩(wěn)定性分析?。篴*x=d a=5 2 1;-1 4 2;2 -3 10 d=-12;20;3(1)考察用J、G-S及SOR解此方程組的收斂性.(2)分別用雅可比迭代法、高斯-賽德爾迭代法及逐次超松弛迭代法解此方程組。要求時(shí)迭代中止，觀

7、察收斂速度。(1) J迭代矩陣的譜半徑：max(abs(eig(D*(L+U)= 0.506G-S迭代矩陣的譜半徑：max(abs(eig(inv(D-L)*U)= 0.2000SOR迭代矩陣的譜半徑：max(abs(eig(inv(D-w*L)*(1-w)*D+w*U)=0.9756所以用J、G-S及SOR均收斂（且有 ）。(2)Jacbio： X =-4.0000 3.0000 2.0000 n =18GaussSeidel： X =-4.0000 3.0000 2.0000 n =8SOR( )： X =-4.0000 3.0000 2.0000 n =482可見高斯-賽德爾迭代法要比雅

8、可比迭代公式的收斂速度快。同時(shí)，如果超松弛迭代法的 選取不當(dāng)不但不會(huì)加速收斂還會(huì)對(duì)收斂速度造成很惡劣的結(jié)果.3. 誤差分析設(shè)方程組 AX=b,其中， A=(aij)n 為非奇異陣， X=(x1,x2,xn)T, b=(b1,b2,bn)T.由于原始數(shù)據(jù)aij, bi往往是觀測(cè)數(shù)據(jù)，難免帶有誤差，因此， 我們下面討論原始數(shù)據(jù)的微小變化對(duì)方程組的影響.設(shè)A精確b有誤差，得到解即 而于是 說明右端相對(duì)誤差 在解放大 倍.設(shè)b精確A有誤差，得到解即 或者 只要|A|充分小使得上式表明：當(dāng)|A|充分小，矩陣 A 的相對(duì)誤差在解中可能放大了倍.五、 算法的實(shí)現(xiàn)（代碼+算例）1. Jacobi迭代代碼fun

9、ction x,n=jacbio(a,d,x)stop=1.0e-4; %迭代的精度L=-tril(a,-1);U=-triu(a,1);D=inv(diag(diag(a);X=D*(L+U)*x+D*d; % J迭代公式n=1;while norm(X-x,inf)>=stop % 時(shí)迭代中止否則繼續(xù)x=X;X=D*(L+U)*x+D*d;n=n+1;endXn2. GaussSeidel迭代代碼function x,n=Gauss_Seidel_iterative(a,d,x)stop=1.0e-4;%迭代的精度L=-tril(a,-1);U=-triu(a,1);D=(diag(

10、diag(a);X=inv(D-L)*U*x+inv(D-L)*d; % G-S迭代公式n=1;while norm(X-x,inf)>= stop % 時(shí)迭代中止否則繼續(xù)x=X;X=inv(D-L)*U*x+inv(D-L)*d;n=n+1;endXn3. 松弛迭代(SOR)迭代代碼function x,n=successive(a,d,x)stop=1.0e-4;%迭代的精度w=1.44; %松弛因子L=-tril(a,-1);U=-triu(a,1);D=(diag(diag(a);X=inv(D-w*L)*(1-w)*D+w*U)*x+w*inv(D-w*L)*d;% SOR迭代

11、公n=1;while norm(X-x,inf)>=stop % 時(shí)迭代中止否則繼續(xù)x=X;X=inv(D-w*L)*(1-w)*D+w*U)*x+w*inv(D-w*L)*d;n=n+1;endXn4. 算法實(shí)例4.1 jacbio迭代a=5 2 1;-1 4 2;2 -3 10;>> x=0;0;0; >> d=-12;20;3;>> x,n=jacbio(a,d,x)X =-4.0000;3.0000;2.0000n =184.2 Gauss迭代a=5 2 1;-1 4 2;2 -3 10;>> x=0;0;0; >> d

12、=-12;20;3;>> x,n=Gauss_Seidel_iterative(a,d,x)X =-4.0000;3.0000;2.0000n = 84.3 SOR迭代a=5 2 1;-1 4 2;2 -3 10;>> x=0;0;0; >> d=-12;20;3;>> x,n=jacbio(a,d,x)X =-4.0000;3.0000;2.0000n =482六、 應(yīng)用實(shí)例1. 道路交通流量2. 電網(wǎng)中的電流流量3. 營養(yǎng)減肥食譜七、 知識(shí)拓展線性方程組是近似求解微分方程的基礎(chǔ)，因而是人類揭示自然規(guī)律的基礎(chǔ)；邊界問題的內(nèi)點(diǎn)個(gè)數(shù)決定著線性方程組的階數(shù)，在實(shí)際問題中我們遇到的求解區(qū)域往往很大，如計(jì)算天氣預(yù)報(bào)的大氣層區(qū)域、計(jì)算海洋動(dòng)力的海域、計(jì)算油氣儲(chǔ)量的含油區(qū)域等等，為提高精度就要將網(wǎng)格劃