數(shù)值分析解線性方程組的直接方法演示文稿

數(shù)值分析解線性方程組的直接方法演示文稿目前一頁\總數(shù)九十一頁\編于二十一點優(yōu)選數(shù)值分析解線性方程組的直接方法Ppt目前二頁\總數(shù)九十一頁\編于二十一點

在自然科學(xué)和工程技術(shù)中有很多問題的解決常常歸結(jié)為解線性代數(shù)方程組．如三次樣條函數(shù)問題，用最小二乘法求實驗數(shù)據(jù)的曲線擬合問題，解非線性方程組問題，用差分法或者有限元方法解常微分方程、偏微分方程的邊值問題等都導(dǎo)致求解線性代數(shù)方程組，而這些方程組的系數(shù)矩陣大致分為兩種，一種是低階稠密矩陣，另一種是大型稀疏矩陣。

關(guān)于線性方程組的數(shù)值解法一般有兩類：1.直接法2.迭代法§5.1引言與預(yù)備知識5.1.1引言目前三頁\總數(shù)九十一頁\編于二十一點

本章討論n元線性方程組

（5.1）的直接解法。方程組(5.1)的矩陣形式為

Ax=b其中

若矩陣A非奇異,即det(A)≠0,則方程組(2.1)有唯一解。目前四頁\總數(shù)九十一頁\編于二十一點

所謂直接解法是指，若不考慮計算過程中的舍入誤差，經(jīng)過有限次算術(shù)運算就能求出線性方程組的精確解的方法。但由于實際計算中舍入誤差的存在，用直接解法一般也只能求出方程組的近似解。Cramer法則是一種不實用的直接法，本章將介紹幾種實用的直接法。目前五頁\總數(shù)九十一頁\編于二十一點5.1.2預(yù)備知識M行n列矩陣.n維列向量.目前六頁\總數(shù)九十一頁\編于二十一點矩陣的基本運算:(1)矩陣的加法(7)矩陣的行列式行列式性質(zhì):(a)det(AB)=det(A)det(B)(6)非奇異矩陣(5)單位矩陣(4)轉(zhuǎn)置矩陣(3)矩陣與矩陣的乘法(2)矩陣與標(biāo)量的乘法目前七頁\總數(shù)九十一頁\編于二十一點矩陣特征值與譜半徑定義1設(shè)若存在一個數(shù)λ(實數(shù)或復(fù)數(shù))和非零向量使(1.1)則稱λ為A的特征值,x為A對應(yīng)λ的特征向量,A的全體特征值稱為A的譜，記作稱為A的譜半徑.(1.2)目前八頁\總數(shù)九十一頁\編于二十一點由式(1.1)知,λ

可使齊次方程組有非零解,故系數(shù)行列式記稱為矩陣Ａ的特征多項式，方程(1.3)稱為Ａ的特征方程．(1.3)目前九頁\總數(shù)九十一頁\編于二十一點在復(fù)數(shù)域中有n個根故由行列式(1.3)展開可知：的n個特征值故是它的特征方程(1.3)的幾個根,并有(1.4)(1.5)A的跡.目前十頁\總數(shù)九十一頁\編于二十一點A的特征值λ和特征向量x還有以下性質(zhì):(1)AT與A有相同的特征值λ及相同的特征向量x.(2)若A非奇異，則A-1的特征值為λ-1,特征向量為x.

(3)相似矩陣B=S-1AS有相同的特征多項式.目前十一頁\總數(shù)九十一頁\編于二十一點例1求的特征值及譜半徑.解:A的特征方程為故A的特征值為A的譜半徑為目前十二頁\總數(shù)九十一頁\編于二十一點5.1.4特殊矩陣目前十三頁\總數(shù)九十一頁\編于二十一點目前十四頁\總數(shù)九十一頁\編于二十一點定理1.目前十五頁\總數(shù)九十一頁\編于二十一點定理2.目前十六頁\總數(shù)九十一頁\編于二十一點定理3.定理4(Jordan標(biāo)準(zhǔn)型)設(shè)A為n階矩陣,則存在一個非奇異矩陣P使得目前十七頁\總數(shù)九十一頁\編于二十一點其中:(1)當(dāng)A的若當(dāng)標(biāo)準(zhǔn)型中所有若當(dāng)塊Ji均為一階時,此標(biāo)準(zhǔn)型變成對角矩陣.返回主頁目前十八頁\總數(shù)九十一頁\編于二十一點求解高斯消去法(逐次消去法)及消去法和矩陣三角分解之間的關(guān)系：§5.2高斯消去法目前十九頁\總數(shù)九十一頁\編于二十一點例2

用消去法解方程組解第1步.將方程(2.2)乘上-2加到方程(2.4)上,消去未知數(shù)x1,得到(2.2)(2.3)(2.4)第2步.將方程(2.3)加到方程(2.5)上去,消去方程(2.5)中的x2,得到與原方程組等價的三角形方程組解為:首先舉一個簡單的例子來說明消去法的基本思想.目前二十頁\總數(shù)九十一頁\編于二十一點上述過程相當(dāng)于思路首先將A化為上三角陣/*upper-triangularmatrix*/，再回代求解/*backwardsubstitution*/。=目前二十一頁\總數(shù)九十一頁\編于二十一點消元記Step1：設(shè)，計算因子將增廣矩陣/*augmentedmatrix*/第i行mi1

第1行，得到其中Stepk：設(shè)，計算因子且計算目前二十二頁\總數(shù)九十一頁\編于二十一點回代注意1：只要A

非奇異，即A1

存在，則可通過逐次消元及行交換，將方程組化為三角形方程組，求出唯一解。共進行？步n

1目前二十三頁\總數(shù)九十一頁\編于二十一點注意2:設(shè)Ax=b,其中A為非奇異矩陣,如果由于A為非奇異矩陣,所以A的第一列一定有元素不等于零.例如目前二十四頁\總數(shù)九十一頁\編于二十一點定理5

設(shè)Ax=b,其中(1)如果則可通過高斯消去法將Ax=b約化為等價的三角形線性方程組(2.10),且計算公式為:①消元計算(k=1,2,…,n-1)②回代計算目前二十五頁\總數(shù)九十一頁\編于二十一點(2)如果A為非奇異矩陣,則可通過高斯消去法(及交換兩行的初等變換)將方程組Ax=b約化為方程組(2.10).定理6

約化的主元素aii(i)

≠0(i=1,2,…,k)的充要條件是矩陣A的所有順序主子式

/*determinantofleadingprincipalsubmatrices*/Di≠0(i=1,2,…,k).即(2.12)推論如果A的順序主子式Dk≠0(k=1,2,…,n-1),則目前二十六頁\總數(shù)九十一頁\編于二十一點§5.2.2三角分解法

/*MatrixFactorization*/高斯消元法的矩陣形式

/*MatrixFormofG.E.*/：Step1:記L1=，則Stepn

1:其中

Lk=目前二十七頁\總數(shù)九十一頁\編于二十一點記為L單位下三角陣/*unitarylower-triangularmatrix*/記

U=目前二十八頁\總數(shù)九十一頁\編于二十一點定理7

若A的所有順序主子式/*determinantofleadingprincipalsubmatrices*/

均不為0，則A

的

LU

分解唯一（其中L

為單位下三角陣）。證明：由§1中定理可知，LU分解存在。下面證明唯一性。若不唯一，則可設(shè)A=L1U1=L2U2

，推出Upper-triangularLower-triangularWithdiagonalentries1注：L

為一般下三角陣而U

為單位上三角陣的分解稱為Crout分解。實際上只要考慮A*的LU

分解，即

，則即是A的Crout分解。目前二十九頁\總數(shù)九十一頁\編于二十一點例3

對于例2,系數(shù)矩陣由高斯消去法,返回主頁目前三十頁\總數(shù)九十一頁\編于二十一點5.3.1列主元消去法:在順序消元過程中,只要即可進行計算,但如果很小,則將導(dǎo)致舍入誤差增長,使解的誤差很大.例4用Gauss消去法求解方程組§5.3高斯主元素消去法目前三十一頁\總數(shù)九十一頁\編于二十一點解:因故方程有唯一解,且精確解為若用Gauss消去法取四位有效數(shù)字計算,可得解比較,誤差很大,若將兩個方程互換為用Gauss消去法取四位有效數(shù)字計算,可得解目前三十二頁\總數(shù)九十一頁\編于二十一點

本例表明通過行交換可避免舍入誤差增長,這就是列主元消去法的基本思想.其計算步驟如下:第1步，在中的第1列選主元，即行為主元,若將

的第i1行與第1行互換，再按消元公式計算得到假定上述過程已進行(k-1)步，得到第k步，在中的第k列選主元，即若則在中將ik行與第k行互換，再按消元公式計算得到對k=1,2,…,n-1,重復(fù)以上過程則得如果某個k出現(xiàn)主元方程沒有唯一解或嚴(yán)重病態(tài)，否則可由(3.2.4)求得解.則表明detA=0，目前三十三頁\總數(shù)九十一頁\編于二十一點它也表明當(dāng)Ａ非奇異時，存在排列矩陣Ｐ(若干初等排列矩陣的乘積),使PA=LU,其中L為單位下三角矩陣,其元素|lij|<=1,U為上三角矩陣.上述每步行交換后再消元相當(dāng)于其中是指標(biāo)為k的初等下三角矩陣,為初等排列矩陣時,表示不換行,經(jīng)過(n-1)步換行與消元,A化為上三角矩陣.即:目前三十四頁\總數(shù)九十一頁\編于二十一點解:例5用列主元消去法解Ａx=b,其中消元消元目前三十五頁\總數(shù)九十一頁\編于二十一點消元結(jié)束．由回代公式求得解此例的精確解為可見結(jié)果精度較高．若不選列主元Ｇauss消去法，求得解，誤差較大．除列主元消去法外，還有一種消去法，是在Ａ的所有元素aij中選主元，稱為全主元消去法．因計算量較大且應(yīng)用列主元已能滿足實際要求，故不再討論．目前很多數(shù)學(xué)軟件庫都有列主元消去法，可直接調(diào)用．目前三十六頁\總數(shù)九十一頁\編于二十一點注:為了減少計算的舍入誤差,使用消去法通常都要選主元.目前最常用的是列主元消去法,也就是每步消元之前選主元,當(dāng)A=(aij)第一步選A中第1列的主元,即max|ai1|=ai1.然后將i1行與第１行互換，再進行消元，以后每步消元做法類似，先選主元，再消元．目前三十七頁\總數(shù)九十一頁\編于二十一點5.3.2高斯若當(dāng)消去法消去對角線上方和下方的元素.假設(shè)已經(jīng)完成k-1步,得到與方程Ax=b等價的方程組返回主頁目前三十八頁\總數(shù)九十一頁\編于二十一點

高斯消去法有很多變形,有的是高斯消去法的改進、改寫，有的是用于某一類特殊性質(zhì)矩陣的高斯消去法的簡化。5.3.1直接三角分解法.

將高斯消去法改寫為緊湊形式,可以直接從A的元素得到計算L,U元素的遞推公式,而不需要任何中間步驟,這就是所謂直接三角分解法.

一旦實現(xiàn)A的LU分解,那么求解Ax=b的問題就等價于求解兩個三角形方程組

(1)Ly=b,求y;(2)Ux=y,求x.§５.4矩陣三角分解法

/*MatrixFactorization*/目前三十九頁\總數(shù)九十一頁\編于二十一點１.不選主元的三角分解法設(shè)A為非奇異矩陣,且有分解式A=LU,其中L為單位下三角,U為上三角即L,U元素可以由n步直接計算定出,其中第r步定出U的第r行和L的第r列元素.由上式有:目前四十頁\總數(shù)九十一頁\編于二十一點故同樣有:

設(shè)已經(jīng)定出U的第1行到第r-1行元素與L的第1列到第r-1列元素.利用矩陣乘法(注意當(dāng)r<k時,lrk=0),有得目前四十一頁\總數(shù)九十一頁\編于二十一點通過比較法直接導(dǎo)出L和

U的計算公式。思路目前四十二頁\總數(shù)九十一頁\編于二十一點固定r:對i=r,r+1,…,n

有l(wèi)ii=1a將r

，i

對換，對r=i,i+1,…,n有b目前四十三頁\總數(shù)九十一頁\編于二十一點結(jié)論:用直接三角分解法解Ax=b(要求A的所有順序主子式都不等于0)的計算公式如下.Step1:u1i=a1i;li1=ai1/u11;(i=2,…,n)計算U的第r行,L的第r列元素(r=2,3,…，n）Step2:求解Ly=b,Ux=y

的計算公式：Step3:Step4:Step5:目前四十四頁\總數(shù)九十一頁\編于二十一點例5

用直接三角分解法解解用分解公式計算得求解目前四十五頁\總數(shù)九十一頁\編于二十一點2.選主元的三角分解法當(dāng)urr=0時計算中斷，或者當(dāng)urr絕對值很小時，按分解公式計算可能引起舍入誤差的積累。但如果A非奇異，可以通過交換A的行實現(xiàn)矩陣A的LU分解，因此可采用與列主元消去法類似的方法，將直接三角分解法修改為（部分）選主元的三角分解法。設(shè)第r-1步分解已完成，這時有目前四十六頁\總數(shù)九十一頁\編于二十一點第r步分解需用到（3.2）及（3.3）式，為避免用小的數(shù)urr做除數(shù)，引進量于是有：取交換A的r行與ir行元素，將調(diào)到(r,r)位置(將（i，j）位置的新元素仍記為lij

及

aii),于是有|lir|<=1(i=r+1,…，n）.由此再進行第r步分解計算。目前四十七頁\總數(shù)九十一頁\編于二十一點5.3.2平方根法當(dāng)A對稱正定時，A的順序主子式故由定理知，A=LU的分解存在且唯一，其中L為單位下三角為了A利用對稱性其中D為對角陣,U0為單位上三角陣,于是又代入到上式,就得到對稱矩陣A的分解式矩陣，U為上三角矩陣，且目前四十八頁\總數(shù)九十一頁\編于二十一點定理9(對稱陣的三角分解定理)設(shè)A為n階對稱陣,且A的順序主子式則A可唯一分解為，其中L為單位下三角矩陣，.

D為對角矩陣定理10(對稱正定矩陣的三角分解或Cholesky分解)

如果A為n階對稱正定矩陣,則存在一個實的非奇異下三角陣L使A=LLT,當(dāng)限定L的對角元素為正時,這種分解是唯一的.將求方程組的解轉(zhuǎn)化為求方程的解LLTx=b.令

，求得方程的解由根據(jù)矩陣乘法，由目前四十九頁\總數(shù)九十一頁\編于二十一點得

i=j有

當(dāng)i＞j,得目前五十頁\總數(shù)九十一頁\編于二十一點例6用平方根法求以下方程組的解.

解先驗證系數(shù)矩陣A對稱正定，對稱顯然，

故A對稱正定，可用Cholesky分解計算，求得

求解

得再求得目前五十一頁\總數(shù)九十一頁\編于二十一點5.3.3追趕法解三對角方程組

/*CroutReductionforTridiagonalLinearSystem*/在一些實際問題中,例如解常微分方程邊值問題,解熱傳導(dǎo)方程以及船體數(shù)學(xué)放樣中建立三次樣條函數(shù)等,都會要求解系數(shù)矩陣為對角占優(yōu)的三對角線方程組.簡記為Ax=f.其中,當(dāng)|i-j|>1時,aij=0,且：目前五十二頁\總數(shù)九十一頁\編于二十一點Step1:對A作Crout分解直接比較等式兩邊的元素，可得到計算公式。Ｌ為下三角，Ｕ為單位上三角目前五十三頁\總數(shù)九十一頁\編于二十一點注意當(dāng)j=1時有對j=2，3，…，n求得L的元素，這就是A的Cholesky分解，然后再解兩個三角方程組，得這就是對稱正定方程組的平方根法

另外，由于

故有這表明分解過程中矩陣L中元素因此平方根法計算是數(shù)值穩(wěn)定的.

的數(shù)量級不增長，目前五十四頁\總數(shù)九十一頁\編于二十一點Step2:追——即解Step3:趕——即解：求解Ａx＝f等價于Step１:計算的遞推公式將計算系數(shù)為追的過程將計算方程組的解為趕的過程目前五十五頁\總數(shù)九十一頁\編于二十一點定理：設(shè)有三對角線方程組Ａx＝f，其中Ａ滿足對角占優(yōu)的條件，則Ａ為非奇異矩陣且追趕法計算公式中的滿足：目前五十六頁\總數(shù)九十一頁\編于二十一點定理

若A

為對角占優(yōu)

/*diagonallydominant*/的三對角陣，且滿足，則追趕法可解以A

為系數(shù)矩陣的方程組。注：

如果A是嚴(yán)格對角占優(yōu)陣，則不要求三對角線上的所有元素非零。

根據(jù)不等式可知：分解過程中，矩陣元素不會過分增大，算法保證穩(wěn)定。

運算量為O(6n)。返回主頁目前五十七頁\總數(shù)九十一頁\編于二十一點5.5.1內(nèi)積與向量范數(shù)為了研究方程組Ax=b解的誤差和迭代法收斂性，需對向量及矩陣的"大小"引進一種度量，就要定義范數(shù)，它是向量"長度"概念的直接推廣，通常用表示n維實向量空間，表示n維復(fù)向量空間.定義2

設(shè)將實數(shù)稱為向量x,y的數(shù)量積.非負實數(shù)稱為向量x的歐氏范數(shù)或2-范數(shù).§5.5向量和矩陣范數(shù)目前五十八頁\總數(shù)九十一頁\編于二十一點定理12設(shè)則內(nèi)積有以下性質(zhì)：(1)(2)(3)(4)(5)(柯西-施瓦茨不等式)等式當(dāng)且僅當(dāng)x與y線性相關(guān)時成立;(6)三角不等式目前五十九頁\總數(shù)九十一頁\編于二十一點定義3(向量范數(shù))如果向量的某個實值函數(shù)滿足條件:則稱目前六十頁\總數(shù)九十一頁\編于二十一點對于由內(nèi)積性質(zhì)可知它滿足定義2的三個條件，故它是一種向量范數(shù).此外還有以下幾種常用的向量范數(shù).容易驗證均滿足定義2的三個條件.更一般的還可定義但只有p=1,2,∞時的三種范數(shù)是常用的向量范數(shù).例如給定則可求出目前六十一頁\總數(shù)九十一頁\編于二十一點定理14設(shè)是則N(x)是向量x的分量上任一種向量范數(shù)，的連續(xù)函數(shù).定理15(向量范數(shù)的等價性)設(shè)是上任意兩種向量范數(shù)，則存在常數(shù)使目前六十二頁\總數(shù)九十一頁\編于二十一點5.5.2矩陣范數(shù)

矩陣可看成n×n維向量，如果直接將向量的2-范數(shù)用于矩陣A，則可定義稱為矩陣A的Frobenius范數(shù)，簡稱F-范數(shù).它顯然滿足向量范數(shù)的三條性質(zhì)，但由于矩陣還有乘法運算，因此矩陣范數(shù)的定義中應(yīng)增加新條件.目前六十三頁\總數(shù)九十一頁\編于二十一點定義4

如果的某個非負實函數(shù)N(A)，記作‖A‖，滿足條件：則稱目前六十四頁\總數(shù)九十一頁\編于二十一點顯然滿足定義中的四個條件,(3)，(4)兩條均可由Cauchy-Schwarz不等式證明，故是一種矩陣范數(shù).

除矩陣自身的運算外，在解方程中矩陣乘向量的運算即Ax，也是必不可少的.因此要求所引進的范數(shù)應(yīng)滿足條件：上式稱為相容性條件.為使引進的矩陣范數(shù)滿足條件(4.5)，我們給出以下定義.(4.5)目前六十五頁\總數(shù)九十一頁\編于二十一點定義6(矩陣的算子范數(shù))設(shè)當(dāng)給定向量范數(shù)時可定義稱為矩陣的算子范數(shù)或從屬范數(shù).(4.6)定理17

設(shè)上的一種向量范數(shù)，則由(4.6)定義的是一種矩陣范數(shù)，且滿足相容性條件目前六十六頁\總數(shù)九十一頁\編于二十一點證明因中有界閉集上的連續(xù)函數(shù)，故在D上有最大值，即使而對故所以從而當(dāng)成立，而x=0時顯然也成立.目前六十七頁\總數(shù)九十一頁\編于二十一點定理17

設(shè)則這里為矩陣的譜半徑.目前六十八頁\總數(shù)九十一頁\編于二十一點例7

已知解從定理可以看出，計算較容易，而計算

時因為要求的特征值，所以較為困難.但當(dāng)A對稱時，有目前六十九頁\總數(shù)九十一頁\編于二十一點定理19定理18

對任何為任一種從屬范數(shù)則反之，對任意ε＞0，至少存在一種從屬范數(shù)使證明:設(shè)為A的特征值，則由得目前七十頁\總數(shù)九十一頁\編于二十一點非奇異，且證明用反證法.假定(I+B)奇異，則齊次方程有非零解而與‖B‖＜1的假設(shè)矛盾，故(I+B)非奇異.

又得取范數(shù)得定理20設(shè)

返回主頁目前七十一頁\總數(shù)九十一頁\編于二十一點矩陣條件數(shù)與擾動方程組誤差界

在解方程組Ax=b時，由于各種原因，A或b往往有誤差，從而使得解也產(chǎn)生誤差.例8

方程組

的準(zhǔn)確解為

,當(dāng)A與b有微小變化時，如變?yōu)榉匠虅t準(zhǔn)確解為

它表明A,b的微小擾動引起方程解x的很大變化，這就是病態(tài)方程.§5.6誤差分析與病態(tài)方程組目前七十二頁\總數(shù)九十一頁\編于二十一點定義7

求解線性方程組Ax=b時，若A或b有微小擾動

解x的誤差

很大,，則稱此方程組為病態(tài)方程組，相應(yīng)的系數(shù)矩陣A稱為病態(tài)矩陣.

反之，若此時

很小,，則稱此方程組為良態(tài)方程組，相應(yīng)的系數(shù)矩陣A稱為良態(tài)矩陣.

注意方程組是否病態(tài)與用什么數(shù)值方法無關(guān)，它是由方程自身性質(zhì)決定的.

在例8中因為行列式

因此出現(xiàn)病態(tài).但有時A從表面上看性質(zhì)很好，也可能是病態(tài)的.目前七十三頁\總數(shù)九十一頁\編于二十一點那么如何判斷A是否病態(tài)?先給出如下定義.例9

方程組Ax=b表示為它的準(zhǔn)確解

A對稱正定且

表面看性質(zhì)"較好"，但若對右端b作微小變化，如方程改為

則解變?yōu)?/p>

這里b的相對誤差大約只有

但解的相對誤差卻很大，故A也是病態(tài)矩陣.目前七十四頁\總數(shù)九十一頁\編于二十一點定義8

設(shè)

非奇異，‖·‖v為矩陣的任一種從屬范數(shù)，則

稱為矩陣A的條件數(shù).

從定義看到矩陣條件數(shù)依賴于范數(shù)的選取，如范數(shù)為2-范數(shù)，

則記為

同理有

等等.目前七十五頁\總數(shù)九十一頁\編于二十一點條件數(shù)有以下性質(zhì)：

(1)(2)(3)U為正交矩陣，則

(4)若

與為A的按模最大與最小特征值，則若A對稱，則

目前七十六頁\總數(shù)九十一頁\編于二十一點下面給出擾動方程組解的誤差分析.先考察b有擾動

則擾動方程為由于Ax=b,故得

于是再由Ax=b,有

即故得目前七十七頁\總數(shù)九十一頁\編于二十一點下面再研究方程Ax=b，當(dāng)A有擾動

時，其解

的誤差分析.

此時擾動方程為

因Ax=b,故有

因存在,若假定

則由定理20可知

非奇異,并有:(5.6)由(5.6)可得

因此(5.7)目前七十八頁\總數(shù)九十一頁\編于二十一點定理22

設(shè)A為非奇異矩陣,Ax=b≠0,且如果則(5.7)式成立.從(5.7)看到，當(dāng)A的條件數(shù)Cond(A)很大時，解的相對誤差

也很大，故方程組為病態(tài).在例9中

而于是條件數(shù)很大，故方程是嚴(yán)重病態(tài)的.

目前七十九頁\總數(shù)九十一頁\編于二十一點例10

Hibert矩陣是一個著名的病態(tài)矩陣，記作

它是一個對稱正定矩陣，當(dāng)n≥3時它是病態(tài)矩陣.例如

故另外還有

等等.因此Hn是嚴(yán)重的病態(tài)矩陣，且n越大

Cond(Hn)越大.目前八十頁\總數(shù)九十一頁\編于二十一點例11

在例10的方程組中可算出A的特征值

故例中實際相對誤差是而根據(jù)(5.6)的誤差估計為這與實際相差不大，即相對誤差放大了將近3000倍.故方程為病態(tài)方程組.目前八十一頁\總數(shù)九十一頁\編于二十一點定理23(事后誤差估計)設(shè)方程組

，則若實際求得解為證明記剩余則它表明如果方程組病態(tài)，即使剩余‖r‖很小，解的相對誤差仍可能很大.目前八十二頁\總數(shù)九十一頁\編于二十一點5.6.2病態(tài)方程組的解法

如果A的條件數(shù)Cond(A)>>1，則Ax=b為病態(tài)方程，但計算Cond(A)時需要求A-1,計算量很大，相當(dāng)于解方程組，在實際中常可通過求解過程直觀地判斷方程組的病態(tài)性質(zhì)，如果解方程時出現(xiàn)下述情況之一，則可能是"病態(tài)"方程組.(1)在列主元消去法中出現(xiàn)小主元；(2)在計算過程中行或列幾乎線性相關(guān)或三角分解中對角元出現(xiàn)近似零的元素；(3)矩陣A的元素數(shù)量級相差很大且無規(guī)律；(4)剩余很小，而解很大，又達不到精度要求.目前八十三頁\總數(shù)九十一頁\編于二十一點

(1)采用高精度運算，減輕病態(tài)影響，例如用雙倍字長運算.對病態(tài)方程組求解可采用以下措施：(2)用預(yù)處理方法改善A的條件數(shù)，即選擇非奇矩陣與Ax=b等價,而(3)平衡方法，當(dāng)A中元素的數(shù)量級相差很大，可采用行均衡或列均衡的方法改善A的條件數(shù).設(shè)非奇異，計算于是求Ax=b等價于求的條件數(shù)可得到改善，這就是行均衡法.目前八十四頁\總數(shù)九十一頁\編于二十一點例12

給定方程組Ax=b為A的條件數(shù)若用行均衡法可取則平衡后的方程用三位有效數(shù)字的列主元消去法求解得目前八十五頁\總數(shù)九十一頁\編于二十一點functionx=threedia(a,b,c,f)N=length(f);x=zeros(1,N);y=zeros(1,N);beta=zeros(1,N);gramma=zeros(1,N);beta(1)=b(1);fori=1:N-1gramma(i)=c(i)/beta(i);beta(i+1)=b(i+1)-a(i+1)*gramma(i);end%追的過程y(1)=f(1)/beta(1);fori=2:Ny(i)=(f(i)-a(i)*y(i-1))/beta(i);end%趕的過程x(N)=y(N);fori=N-1:-1:1x(i)=y(i)-gramma(i)*x(i+1);end

a=[0,-1,-1,-3];>>b=[2,3,2,5];>>c=[-1,-2,-1,0];>>f=[6,1,0,1]';>>x=threedia(a,b,c,f)追趕法求解三對角方程組目前八十六頁\總數(shù)九十一頁\編于二十一點Cholesky方法：

A=[4,-2,4;-2,17,10;4,10,9];

b=[8.7,13.7,-0.7]';

[x,L,D]=Chol_decompose(A,b)L=1.000000-0.50001.000001.00000.75001.0000D=416-4x=-5.1457-3.17275.7344%用Cholesky分解求解%A是對稱矩陣%L是單位下三角陣%D是對角陣%對稱陣A進行三角分解：A=LDL'目前八十七頁\總數(shù)九十一頁\編于二十一點function[x,L,D]=Chol_decompose(A,b)N=length(A);L=zeros(N,N);D=zeros(1,N);fori=1:NL(i,i)=1;endD(1)=A(1,1);fori=2:Nforj=1:i-1ifj==1L(i,j)=A(i,j)/D(j);elsesum1=0;fork=1:j-1sum1=sum1+L(i,k)*D(k)*L(j,k);endL(i,j)=(A(i,j)-sum1)/D(j);endendsum2=0;fork=1:i-1sum2=sum2+L(i,k)^2*D(k);endD(i)=A(i,i)-sum2;end%分別求解線性方程組Ly=b;L'x=y/Dy=zeros(1,N);y(1)=b(1);fori=2:Nsumi=0;fork=1:i-1sumi=sumi+L(i,k)*y(k);endy(i)=b(i)-sumi;endx=zeros(1,N);x(N)=y(N)/D(N);fori=N-1:-1:1sumi=0;fork=i+1:Nsumi=sumi+L(k,i)*x(k);endx(i)=y(i)/D(i)-sumi;end目前八十八頁\總數(shù)九十一頁\編于二十一點用Dollittle三角分解法求解方程組：>>A=[0.001,2,3;-1,3.712,4.623;-2,1.072,5.643];>>b=[1,2,3]';>>[x,L,U]=lu_decompose(A,b)x=-0.4904-0.05100.3675L=1.0e+003*0.001000-1.00000.00100-2.00000.00200.0010U=1.0e+003*0.00000.00200.003002.00373.0046000.0059目前八十九頁\總數(shù)九十一頁\編于二十一點function[x,L,U]=lu_decompose(A,b)%用Dollittle三角分解%