數(shù)值分析第六章課件

1、數(shù)值計(jì)算方法任課教師: 徐昱 工程學(xué)院海洋工程系 1201第五章解線性方程組的數(shù)值解法引言 在自然科學(xué)和工程技術(shù)中,很多問題歸結(jié)為解線性方程組.例如電學(xué)中的網(wǎng)絡(luò)問題，船體數(shù)學(xué)放樣中建立三次樣條函數(shù)問題，用最小二乘法求實(shí)驗(yàn)數(shù)據(jù)的曲線擬合問題，解非線性方程組問題，用差分法或者有限元方法解常微分方程、偏微分方程邊值問題等都導(dǎo)致求解線性方程組，而這些方程組的系數(shù)矩陣大致分為兩種，一種是低階稠密矩陣(例如，階數(shù)不超過150). 另一種是大型稀疏矩陣(即矩陣階數(shù)高且零元素較多). 有的問題的數(shù)學(xué)模型中雖不直接表現(xiàn)為含線性方程組,但它的數(shù)值解法中將問題“離散化”或“線性化”為線性方程組.因此線性方程組的求解

2、是數(shù)值分析課程中最基本的內(nèi)容之一. 關(guān)于線性方程組的解法一般有兩大類： 但實(shí)際計(jì)算中由于舍入誤差的存在和影響，這種方法也只能求得線性方程組的近似解. 本章將闡述這類算法中最基本的和具有代表性的算法就是高斯消元法,以及它的某些變形和應(yīng)用.這類方法是解低階稠密矩陣方程組及某些大型稀疏矩陣方程組(例如，大型帶狀方程組)的有效方法. 1. 直接法 經(jīng)過有限次的算術(shù)運(yùn)算,可以求得方程組的精確解(假定計(jì)算過程沒有舍入誤差).如線性代數(shù)課程中提到的克萊姆算法就是一種直接法.但該法對(duì)高階方程組計(jì)算量太大,不是一種實(shí)用的算法. 就是用某種極限過程去逐步逼近方程組精確解的方法. 迭代法具有計(jì)算機(jī)的存儲(chǔ)單元較少、程

3、序設(shè)計(jì)簡(jiǎn)單、原始系數(shù)矩陣在計(jì)算過程中始終不變等優(yōu)點(diǎn)，但存在收斂條件和收斂速度問題.迭代法是解大型稀疏矩陣方程組(尤其是由微分方程離散后得到的大型方程組)的重要方法. 2. 迭代法5.2 高斯消去法 本節(jié)介紹高斯消去法(逐次消去法)及消去法和矩陣三角分解之間的關(guān)系. 雖然高斯消去法是一種古老的求解線性方程組的方法(早在公元前250年我國(guó)就掌握了解方程組的消去法)，但由它改進(jìn)、變形得到的選主元素消去法、三角分解法仍然是目前計(jì)算機(jī)上常用的有效方法.我們?cè)谥袑W(xué)學(xué)過消去法,高斯消去法就是它的標(biāo)準(zhǔn)化的、適合在計(jì)算機(jī)上自動(dòng)計(jì)算的一種方法.5.2.1 高斯消去法 設(shè)有線性方程組或?qū)懗删仃囆问胶?jiǎn)記為Ax=b.如

4、: 上三角矩陣所對(duì)應(yīng)的線性方程組解為 當(dāng)m=n時(shí)，對(duì)三角形方程組的解非常容易，如下三角矩陣所對(duì)應(yīng)的線性方程組計(jì)算量（乘除法的主要部分）都為 n2/2.解為 因此，我們將一般的線性方程組化成等價(jià)的三角形方程組來求解. 首先舉一個(gè)簡(jiǎn)單的例子來說明消去法的基本思想. 例1 用消去法解方程組 解 第1步，將方程(2.2)乘上-2加到方程(2.4)上去，消去(2.4)中的未知數(shù)x1，得到 第2步，將方程(2.3) 加到方程(2.5)上去，消去(2.5)中的未知數(shù)x2，得到與原方程組等價(jià)的三角形方程組顯然，方程組是(2.6)是容易求解的，解為 上述過程相當(dāng)于對(duì)方程的增廣陣做初等行變換其中ri用表示矩陣的第

5、i行. 由此看出，用消去法解方程組的基本思想是用逐次消去未知數(shù)的方法把原方程組Ax=b化為與其等價(jià)的三角形方程組，而求解三角形方程組可用回代的方法求解. 換句話說，上述過程就是用初等行變換將原方程組系數(shù)矩陣化為簡(jiǎn)單形式(上三角矩陣)，從而求解原方程組(2.1)的問題轉(zhuǎn)化為求解簡(jiǎn)單方程組的問題. 或者說，對(duì)系數(shù)矩陣A施行一些行變換(用一些簡(jiǎn)單矩陣左乘A)將其約化為上三角矩陣. 這就是高斯消去法. 下面討論求解一般線性方程組的高斯消去法.由 將(2.1)記為A(1)x=b(1)，其中 (1) 第1步(k=1). 設(shè)a11(1)0，首先計(jì)算乘數(shù) mi1= ai1(1) /a11(1) (i=2,3,

6、m).用-mi1乘(2.1)的第一個(gè)方程，加到第i個(gè)(i=2,3, ,m)方程上，消去(2.1)的從第二個(gè)方程到第m個(gè)方程中的未知數(shù)x1，得到與(2.1)等價(jià)的方程組 (2.7)簡(jiǎn)記為 A(2)x=b(2)，其中A(2), b(2)的元素計(jì)算公式為 (2) 第k次消元(k=1,2,s, s=min(m-1,n). 設(shè)上述第1步， ，第k-1步消元過程計(jì)算已經(jīng)完成，即已計(jì)算好與(2.1)等價(jià)的方程組，簡(jiǎn)記為A(k)x=b(k)，其中 設(shè)akk(k)0，計(jì)算乘數(shù) mik= aik(k) /akk(k) (i=k+1, ,m).用-mik乘(2.8)的第k個(gè)方程加到第 i個(gè)(i= k+1, , m)

7、方程上，消去從第k+1個(gè)方程到第m個(gè)方程中的未知數(shù)xk，得到與(2.1)等價(jià)的方程組A(k+1)x=b(k+1)，其中A(k+1), b(k+1)的元素計(jì)算公式為，顯然A(k+1)中從第1行到第k行與A(k)相同. (3) 繼續(xù)上述過程，且設(shè)aii(i)0 (i=1,2, ,s)，直到完成第s步消元計(jì)算. 最后得到與原方程組等價(jià)的簡(jiǎn)單方程組A(s+1)x=b(s+1) ，其中A(s+1)為上階梯. 特別當(dāng)m=n時(shí)，與原方程組等價(jià)的簡(jiǎn)單方程組為A(n)x=b(n)，即由(2.1)約化為(2.10)的過程稱為消元過程. 如果A是非奇異矩陣，且akk(k)0 (k=1,2,n-1)，求解三角形方程組

8、(2.10)，得到求解公式(2.10)的求解過程(2.11) 稱為回代過程. 注意：設(shè)Ax=b，其中ARnn為非奇異矩陣，如果a11(1)=0，由于A為非奇異矩陣，所以A的第1列一定有元素不等于零，例如al10，于是可交換兩行元素(即r1rl)，將al1 調(diào)到(1,1)位置，然后進(jìn)行消元計(jì)算，這時(shí)A(2)右下角矩陣為n-1階非奇異矩陣，繼續(xù)這過程，高斯消去法照樣可進(jìn)行計(jì)算. 總結(jié)上述討論即有 定理5 設(shè)Ax=b，其中ARnn. (1) 如果akk(k)0 (k=1,2,n)，則可通過高斯消去法將Ax=b約化為等價(jià)的三角形方程組(2.10)，且計(jì)算公式為 (a) 消元計(jì)算(k=1,2, ,n-1

9、) (b) 回代計(jì)算 (2) 如果A為非奇異矩陣，則可通過高斯消去法(及交換兩行的初等變換)將方程組Ax=b約化為等價(jià)的三角形方程組(2.10).例子 解方程組解：消元回代得5.3 高斯主元素消元法 由高斯消去法知道, 在消元過程中可能有akk(k)0的情況，這時(shí)消去法將無法進(jìn)行；即使主元素akk(k)0但很小時(shí)，用其作除數(shù)，會(huì)導(dǎo)致其它元素?cái)?shù)量級(jí)的嚴(yán)重增長(zhǎng)和舍入誤差的擴(kuò)散，最后也使得計(jì)算解不可靠. 高斯消去法也稱主元素消去法 (條件det A0) 即當(dāng)akk(k)=0 時(shí)，高斯消元法無法進(jìn)行；或 | akk(k) |1時(shí)，帶來舍入誤差的擴(kuò)散。（小主元）5.3.1 列主元素消去法 設(shè)方程組(2.

10、1)的增廣矩陣為 首先在A的第1列中選取絕對(duì)值最大的元素作為主元素，例如交換消元法的應(yīng)用1. 求行列式 高斯消元法2. 求逆矩陣(用高斯-若當(dāng)消去法) 例子 分別用高斯消元法和列主元消元法求矩陣的行列式.解：高斯消元法大數(shù)“吃掉”了小數(shù)！列主元消元法5.3.2 高斯-若當(dāng)消去法 高斯-若當(dāng)消去法是高斯消去法的一種變形. 高斯消去法是消去對(duì)角元素下方的元素. 如果同時(shí)消去對(duì)角元素上方和下方的元素，而且將對(duì)角元素化為1，這種方法就稱為高斯-若當(dāng)(Gauss-Jordan)消去法. 設(shè)高斯-若當(dāng)消去法已完成k-1步，于是Ax=b化為等價(jià)方程組A(k)x=b(k)，其中增廣陣 在第k步計(jì)算時(shí)(k=1,

11、2,n),考慮將第k行上、下的第k列元素都進(jìn)行消元計(jì)算化為零，且akk化為1. 1. 按列選主元素，即確定ik使 2. 換行(當(dāng)ikk時(shí))，交換第k行與第ik行元素(j=k,k+1,n), 3. 計(jì)算乘數(shù) mik=-aik/akk(i=1,2,n, ik) mkk=1/akk . (mik可保存在存放aik的單元中). 4. 消元計(jì)算 aijaij+mikakj(i=1,2,n,且ik,j=k+1,n) bibi+mikbk (i=1,2,n,且ik) 5.計(jì)算主行 akjakjmkk (j=k,k+1,n)， bkbkmkk .上述過程結(jié)束后，即當(dāng)k=n時(shí)有顯然xi=bi,i=1,2,n 就

12、是 Ax=b的解. 說明用高斯-若當(dāng)消去法雖比高斯消去法略復(fù)雜些，但將A約化為單位矩陣，計(jì)算解就在常數(shù)項(xiàng)位置得到，因此用不著回代求解，故也稱為無回代的高斯消元法. 它的計(jì)算量大約需要n3/2次乘除法，要比高斯消去法大，但用高斯-若當(dāng)消去法求一個(gè)矩陣的逆矩陣還是比較合適的. 定理9(高斯-若當(dāng)法求逆矩陣) 設(shè)A為非奇異矩陣，方程組AX=In的增廣矩陣為C=(A|In)，如果對(duì)C應(yīng)用高斯-若當(dāng)方法化為(In|T) ，則A-1=T. 事實(shí)上, 求A的逆矩陣A-1, 即求n階矩陣X使AX=In, 其中In為單位矩陣，將X按列分塊X=(x1,x2,xn)， I=(e1,e2,en)，于是求解AX=In等

13、價(jià)于求解n個(gè)方程組Axj=ej (j=1,2,n).我們可用高斯-若當(dāng)方法求解AX=In. 例4 用高斯-若當(dāng)方法求下列矩陣的逆矩陣. 解 由分塊矩陣所以得小節(jié) 比較而言,Gauss順序消去法條件苛刻,且數(shù)值不穩(wěn)定; Gauss全主元消去法工作量偏大，需要比較 個(gè)元素及行列交換工作，算法復(fù)雜；對(duì)于Gauss-Jordan消去法形式上比其他消元法簡(jiǎn)單，且無回代求解，但計(jì)算量大,比Gauss順序消去法多 計(jì)算量。因此從算法優(yōu)化的角度考慮， Gauss列主元消去法比較好。5.4 矩陣的三角分解法 我們知道對(duì)矩陣進(jìn)行一次初等變換，就相當(dāng)于用相應(yīng)的初等矩陣去左乘原來的矩陣。因此我們這個(gè)觀點(diǎn)來考察Gaus

14、s消元法用矩陣乘法來表示，即可得到求解線性方程組的另一種直接法：矩陣的三角分解。 5.4.1 Gauss消元法的矩陣形式5.4.2 Doolittle分解Doolittle分解若矩陣A有分解：A=LU，其中L為單位下三角陣，U為上三角陣，則稱該分解為Doolittle分解，可以證明，當(dāng)A的各階順序主子式均不為零時(shí)，Doolittle分解可以實(shí)現(xiàn)并且唯一。A的各階順序主子式均不為零，即Doolittle分解Doolittle分解Doolittle分解(看看就可)Doolittle分解Doolittle分解Doolittle分解例題例題例題例題例題上機(jī)P27 例題第六章 解線性方程組的迭代方法6.

15、1 引言6.2 基本迭代法6.3 迭代法的收斂性 例1 求解方程組記為Ax=b，其中方程組的精確解是x*=(3,2,1)T. 現(xiàn)將改寫為或?qū)憺閤=B0 x+f，其中 我們?nèi)稳〕跏贾?，例如取x(0)=(0, 0, 0)T. 將這些值代入(1.3)式右邊(若(1.3)式為等式即求得方程組的解，但一般不滿足)，得到新的值x(1)=(x1(1), x2(1), x3(1)T=(3.5, 3, 3)T ，再將x(1)分量代入(1.3)式右邊得到 x(2)，反復(fù)利用這個(gè)計(jì)算程序，得到一向量序列和一般的計(jì)算公式(迭代公式)簡(jiǎn)寫為 x(k+1)=B0 x(k) +f，其中k表示迭代次數(shù)(k=0,1,2,).

16、迭代到第10次有從此例看出，由迭代法產(chǎn)生的向量序列x(k)逐步逼近方程組的精確解是x*=(3,2,1)T. 即有 對(duì)于任何一個(gè)方程組x=Bx+f(由Ax=b變形得到的等價(jià)方程組)，由迭代法產(chǎn)生的向量序列x(k)是否一定逐步逼近方程組的解x*呢？回答是不一定. 請(qǐng)同學(xué)們考慮用迭代法解下述方程組但 x(k)并不是所有的都收斂到解x*！6.2 基本迭代法其中，A=(aij)Rnn為非奇異矩陣，下面研究任何建立Ax=b的各種迭代法. 設(shè)線性方程組 Ax=b， (2.1) 其中，M為可選擇的非奇異矩陣，且使Mx=d容易求解，一般選擇A的某種近似，稱M為分裂矩陣. 將A分裂為 A=M-N. (2.2) 于

17、是，求解Ax=b轉(zhuǎn)化為求解Mx=Nx+b ，即求解可構(gòu)造一階定常迭代法其中 B=M-1N=M-1(M-A)=I-M-1A ， f=M-1b. 稱 B=I-M-1A為迭代法的迭代矩陣，選取M矩陣，就得到解Ax=b的各種迭代法. 設(shè)aii0 (i=1,2,n)，并將A寫成三部分即A=D-L-U6.2.1 雅可比(Jacobi)迭代法 設(shè)aii0 (i=1,2,n)，選取M為A的對(duì)角元素部分,即選取M=D(對(duì)角陣)，A=D-N，由(2.3)式得到解方程組Ax=b的雅可比(Jacobi)迭代法. 又稱簡(jiǎn)單迭代法.其中B=I-D-1A=D-1(L+U)=J, f=D-1b. 稱J為解Ax=b的雅可比迭代

18、法的迭代矩陣.于是雅可比迭代法可寫為矩陣形式其Jacobi迭代矩陣為等價(jià)方程組其中 aii(i)0 (i=1,2,n)即由方程組Ax=b得到的建立的雅可比迭代格式為于是，解Ax=b的雅可比迭代法的計(jì)算公式為 由(2.6)式可知，雅可比迭代法計(jì)算公式簡(jiǎn)單，每迭代一次只需計(jì)算一次矩陣和向量的乘法且計(jì)算過程中原始矩陣A始終不變. 6.2.2 高斯賽德爾迭代法在 Jacobi 迭代中，計(jì)算 xi(k+1)(2 i n)時(shí)，使用xj(k+1)代替xj(k) (1 j i-1)，即有建立迭代格式或縮寫為稱為高斯塞德爾(Gauss Seidel)迭代法.解Ax=b的高斯塞德爾迭代法的計(jì)算公式為或 雅可比迭代法不使用變量的最新信息計(jì)算xi(k+1)，而由高斯塞德爾迭代公式(2