解線性方程組的直接方法

解線性方程組的直接方法這一章討論線性方程組引言與預備知識的數(shù)值解法.第2頁,共155頁，2024年2月25日，星期天在自然科學和工程技術中,很多問題歸結為解線性方程組.例如電學中的網絡問題，船體數(shù)學放樣中建立三次樣條函數(shù)問題，用最小二乘法求實驗數(shù)據(jù)的曲線擬合問題，解非線性方程組問題，用差分法或者有限元方法解常微分方程、偏微分方程邊值問題等都導致求解線性方程組，而這些方程組的系數(shù)矩陣大致分為兩種，一種是低階稠密矩陣(例如，階數(shù)不超過150).另一種是大型稀疏矩陣(即矩陣階數(shù)高且零元素較多).一、引言第3頁,共155頁，2024年2月25日，星期天有的問題的數(shù)學模型中雖不直接表現(xiàn)為含線性方程組,但它的數(shù)值解法中將問題“離散化”或“線性化”為線性方程組.因此線性方程組的求解是數(shù)值分析課程中最基本的內容之一.

關于線性方程組的解法一般有兩大類：第4頁,共155頁，2024年2月25日，星期天但實際計算中由于舍入誤差的存在和影響，這種方法也只能求得線性方程組的近似解.本章將闡述這類算法中最基本的和具有代表性的算法就是高斯消去法,以及它的某些變形和應用.這類方法是解低階稠密矩陣方程組及某些大型稀疏矩陣方程組(例如，大型帶狀方程組)的有效方法.

1.直接法經過有限次的算術運算,可以求得方程組的精確解(假定計算過程沒有舍入誤差).如線性代數(shù)課程中提到的克萊姆算法就是一種直接法.但該法對高階方程組計算量太大,不是一種實用的算法.第5頁,共155頁，2024年2月25日，星期天就是用某種極限過程去逐步逼近方程組精確解的方法.迭代法具有計算機的存儲單元較少、程序設計簡單、原始系數(shù)矩陣在計算過程中始終不變等優(yōu)點，但存在收斂條件和收斂速度問題.迭代法是解大型稀疏矩陣方程組(尤其是由微分方程離散后得到的大型方程組)的重要方法(見第6章).為了討論線性方程組數(shù)值解法，需復習一些基本的矩陣代數(shù)知識.

2.迭代法第6頁,共155頁，2024年2月25日，星期天二、向量和矩陣

基本概念：用Rm×n表示全部m×n實矩陣的向量空間；用Cm×n表示全部m×n復矩陣的向量空間.(由數(shù)排成的矩陣表，稱為m行n列矩陣).第7頁,共155頁，2024年2月25日，星期天其中aj為A的第j列的m維列向量.同理(m維列向量).其中biT為A的第i行的n維行向量.第8頁,共155頁，2024年2月25日，星期天

矩陣的基本運算：

(1)矩陣加法

(2)矩陣與標量的乘法

(3)矩陣與矩陣的乘法

(4)轉置矩陣第9頁,共155頁，2024年2月25日，星期天(5)單位矩陣其中

(6)非奇異矩陣則稱B是A的逆矩陣，記為A-1，且(A-1)T=(AT)-1.如果A-1存在，則A稱為非奇異矩陣.如果A、B均為非奇異矩陣，則有(AB)-1=B-1A-1.第10頁,共155頁，2024年2月25日，星期天

(7)矩陣的行列式設A∈Rn×n，則A的行列式可按任一行(列)展開,其中Aij為aij的代數(shù)余子式，Aij=(-1)i+jMij，為Mij元素aij的余子式.

行列式性質：第11頁,共155頁，2024年2月25日，星期天三、矩陣的特征值與譜半徑

定義1設A=(aij)

Rn×n，若存在數(shù)

(實數(shù)或復數(shù))和非零向量x=(x0,x1,…,xn)T

Rn，使Ax=

x,(0.1)則稱

為A的特征值,x為A對應

的特征向量,A的全體特征值稱為A的譜,記作σ(A),即σ(A)={

1,

2,…,

n}.記稱為A的譜半徑.由(0.1)式

可知可使齊次線性方程組(

I-A)x=0有非零解，故系數(shù)行列式det(

I-A)=0，記第12頁,共155頁，2024年2月25日，星期天p(

)稱為矩陣A的特征多項式，方程(0.3)稱為A的特征方程.因為n次代數(shù)方程p(

)在復數(shù)域中有n個根

1,

2,…,

n，故第13頁,共155頁，2024年2月25日，星期天由(1.3)式中的行列式展開可得故A=(aij)

Rn×n的n個特征值

1,

2,…,

n是它的特征方程(1.3)的n個根，并有稱trA為A的跡.及第14頁,共155頁，2024年2月25日，星期天

A的特征值

和特征向量x還有以下性質：⑴AT與A有相同的特征值

及特征向量x.⑵若A非奇異,則A-1的特征值為

-1,特征向量為x.⑶相似矩陣B=S-1AS有相同的特征多項式.

例1求矩陣的特征值及譜半徑

解矩陣A的特征方程為得A的特征值為

1=

2=2,

3=7，譜半徑

(A)=7.第15頁,共155頁，2024年2月25日，星期天特殊矩陣設A=(aij)∈Rn×n.

(1)對角矩陣如果當i≠j時，aij=0.

(2)三對角矩陣如果當|i-j|>1時，aij=0.

(3)上三角矩陣如果當i>j時，aij=0.

(4)上海森伯格(Hessenberg)矩陣如果當i>j+1時，aij=0.

(5)對稱矩陣如果AT=A.

(6)埃爾米特(Hermit)矩陣設A∈Cn×n，如果當AH=A(AH是A的共軛轉置矩陣).第16頁,共155頁，2024年2月25日，星期天

(7)對稱正定矩陣如果AT=A且對任意非零向量，x∈Rn，(Ax,x)=xTAx>0.

(8)正交矩陣如果A-1=AT.(9)酉矩陣設A∈Cn×n，如果A-1=AH.(10)初等置換陣由單位矩陣I交換第i行與第j行(或第i列與第j列)，得到的矩陣記為Iij，且.

IijA=B(為交換A第i行與第j行得到的矩陣)；

AIij=C(為交換A第i列與第j列得到的矩陣).

(11)置換陣由初等置換陣的乘積得到的矩陣.第17頁,共155頁，2024年2月25日，星期天

定理1設A∈Rn×n，則下述命題等價：(1)對任何b∈Rn，方程組Ax=b有唯一解.(2)齊次方程組Ax=0只有唯一解零解x=0.(3)det(A)≠0.(4)A-1存在.(5)A的秩rank(A)=n.第18頁,共155頁，2024年2月25日，星期天

定理2設A∈Rn×n為對稱正定矩陣，則(1)A為非奇異矩陣，且A-1亦是對稱正定矩陣.(2)設Ak為A的順序主子陣，則Ak(k=1,2,

,n)亦是對稱正定矩陣，其中(3)A的特征值

i>0

(i=1,2,

,n).(4)A的順序主子式都大于零，即Ak的行列式det(Ak)>0(k=1,2,

,n).第19頁,共155頁，2024年2月25日，星期天

定理3設A∈Rn×n為對稱矩陣，如果det(Ak)>0

(k=1,2,

,n)，或A的特征值

i>0

(i=1,2,

,n)，則A為對稱正定矩陣.有重特征值的矩陣不一定相似于對角矩陣，那么一般n階矩陣A在相似變換下能簡化到什么形狀.

定理4(Jordan標準型)設A為n階矩陣，則存在一個非奇異n階矩陣P使得第20頁,共155頁，2024年2月25日，星期天為若當(Jordan)塊.其中第21頁,共155頁，2024年2月25日，星期天(1)當A的若當標準型中所有若當塊Ji均為一階時，此標準型變?yōu)閷蔷仃?(2)如果A的特征值各不相同，則其若當標準型必為對角矩陣diag(

1,

2,

,

n).第22頁,共155頁，2024年2月25日，星期天6.1高斯消去法

本節(jié)介紹高斯消去法(逐次消去法)及消去法和矩陣三角分解之間的關系.雖然高斯消去法是一種古老的求解線性方程組的方法(早在公元前250年我國就掌握了解方程組的消去法)，但由它改進、變形得到的選主元素消去法、三角分解法仍然是目前計算機上常用的有效方法.我們在中學學過消去法,高斯消去法就是它的標準化的、適合在計算機上自動計算的一種方法.第23頁,共155頁，2024年2月25日，星期天6.1.1高斯消去法設有線性方程組或寫成矩陣形式簡記為Ax=b.第24頁,共155頁，2024年2月25日，星期天解為當m=n時，對三角形方程組的解非常容易的.如:上三角矩陣所對應的線性方程組第25頁,共155頁，2024年2月25日，星期天下三角矩陣所對應的線性方程組計算量（乘除法的主要部分）都為n2/2.解為因此，我們將一般的線性方程組化成等價的三角形方程組來求解.第26頁,共155頁，2024年2月25日，星期天首先舉一個簡單的例子來說明消去法的基本思想.

例2用消去法解線性方程組

解第1步，將方程(1.2)乘上-2加到方程(1.4)上去，消去(1.4)中的未知數(shù)x1，得到第2步，將方程(1.3)加到方程(1.5)上去，消去(1.5)中的未知數(shù)x2，得到與原方程組等價的三角形方程組第27頁,共155頁，2024年2月25日，星期天顯然，方程組是(1.6)是容易求解的，解為上述過程相當于對方程的增廣陣做初等行變換其中ri用表示矩陣的第i行.第28頁,共155頁，2024年2月25日，星期天由此看出，用消去法解方程組的基本思想是用逐次消去未知數(shù)的方法把原方程組Ax=b化為與其等價的三角形方程組，而求解三角形方程組可用回代的方法求解.換句話說，上述過程就是用初等行變換將原方程組系數(shù)矩陣化為簡單形式(上三角矩陣)，從而求解原方程組(1.1)的問題轉化為求解簡單方程組的問題.或者說，對系數(shù)矩陣A施行一些行變換(用一些簡單矩陣左乘A)將其約化為上三角矩陣.這就是高斯消去法.第29頁,共155頁，2024年2月25日，星期天下面討論求解一般線性方程組的高斯消去法.由第30頁,共155頁，2024年2月25日，星期天將(1.1)記為A(1)x=b(1)，其中第31頁,共155頁，2024年2月25日，星期天(1)第1步(k=1).

設a11(1)0，首先計算乘數(shù)

mi1=ai1(1)/a11(1)(i=2,3,…,m).用-mi1乘(1.1)的第一個方程，加到第i個(i=2,3,

,m)方程上，消去(1.1)的從第二個方程到第m個方程中的未知數(shù)x1，得到與(1.1)等價的方程組(1.7)第32頁,共155頁，2024年2月25日，星期天簡記為A(2)x=b(2)，其中A(2),b(2)的元素計算公式為(2)第k次消元(k=1,2,

,s,s=min(m-1,n)).

設上述第1步，

，第k-1步消元過程計算已經完成，即已計算好與(1.1)等價的方程組，簡記為A(k)x=b(k)，其中第33頁,共155頁，2024年2月25日，星期天

設akk(k)0，計算乘數(shù)

mik=aik(k)/akk(k)

(i=k+1,

,m).用-mik乘(1.8)的第k個方程加到第i個(i=k+1,

,m)方程上，消去從第k+1個方程到第m個方程中的未知數(shù)xk，得到與(1.1)等價的方程組A(k+1)x=b(k+1)，第34頁,共155頁，2024年2月25日，星期天其中A(k+1),b(k+1)的元素計算公式為，顯然A(k+1)中從第1行到第k行與A(k)相同.

(3)繼續(xù)上述過程，且設aii(i)0(i=1,2,

,s)，直到完成第s步消元計算.最后得到與原方程組等價的簡單方程組A(s+1)x=b(s+1)，其中A(s+1)為上階梯.第35頁,共155頁，2024年2月25日，星期天

特別當m=n時，與原方程組等價的簡單方程組為A(n)x=b(n)，即由(1.1)約化為(1.10)的過程稱為消元過程.

如果A是非奇異矩陣，且akk(k)0(k=1,2,

,n-1)，求解三角形方程組(1.10)，得到求解公式(1.10)的求解過程(1.11)稱為回代過程.第36頁,共155頁，2024年2月25日，星期天

注意：設Ax=b，其中A∈Rn×n為非奇異矩陣，如果a11(1)=0，由于A為非奇異矩陣，所以A的第1列一定有元素不等于零，例如al10，于是可交換兩行元素(即r1?rl)，將al1

調到(1,1)位置，然后進行消元計算，這時A(2)右下角矩陣為n-1階非奇異矩陣，繼續(xù)這過程，高斯消去法照樣可進行計算.總結上述討論即有第37頁,共155頁，2024年2月25日，星期天

定理5設Ax=b，其中A∈Rn×n.

(1)如果akk(k)0(k=1,2,…,n)，則可通過高斯消去法將Ax=b約化為等價的三角形方程組(1.10)，且計算公式為

①消元計算(k=1,2,

,n-1)第38頁,共155頁，2024年2月25日，星期天

②回代計算

(2)如果A為非奇異矩陣，則可通過高斯消去法(及交換兩行的初等變換)將方程組Ax=b約化為等價的三角形方程組(1.10).第39頁,共155頁，2024年2月25日，星期天

算法1(高斯算法)設A∈Rm×n(m>1),s=min(m-1,n)，如果akk(k)0(k=1,2,

,s)，本算法用高斯方法將A約化為上三角形矩陣，且A(k)覆蓋A，乘數(shù)mik覆蓋aik.

對于k=1,2,

,s

(1)如果akk=0，則計算停止；

(2)對于i=k+1,

,m

(a)aik←mik=aik/akk

(b)對于j=k+1,

,n

aij←aij-mik·akj

.第40頁,共155頁，2024年2月25日，星期天

顯然，算法1第k步需要m-k次除法，(m-k)(n-k)次乘法，因此，本算法(從第1步到第s步消元計算總的計算量)大約需要s3/3-(m-k)s2/2+mns次乘法運算(對相當大的s).當m=n時,總共大約需要n3/3次乘法運算.

數(shù)akk(k)在高斯消去法中有著突出的作用，稱為約化的主元素(簡稱主元).第41頁,共155頁，2024年2月25日，星期天

算法2(回代算法)設Ux=b，其中U=(uij)∈Rn×n為非奇異上三角矩陣，本算法計算Ux=b的解.

對于i=n,n-1,

,1

(1)xi←bi

(2)對于j=i+1,

,n

xi←xi-uij·xj(i>n不算)

(3)xi←xi/uii

這個算法需要n(n-1)/2次乘除法運算.第42頁,共155頁，2024年2月25日，星期天例子解方程組解：消元回代得第43頁,共155頁，2024年2月25日，星期天

高斯消去法對于某些簡單的矩陣可能會失敗，例如

由此，需要對算法1進行修改，首先研究原來矩陣A在什么條件下才能保證akk(k)0(k=1,2,

,s).下面的定理給出了這個條件.第44頁,共155頁，2024年2月25日，星期天

定理6約化的主元素aii(i)0(i=1,2,

,k)的充要條件是方陣A的順序主子式Di0(i=1,2,

,k)，即

證明首先利用歸納法證明充分性.顯然,當k=1時，結論成立，假設結論對k-1是成立的，對k由條件設Di0(i=1,2,

,k)，于是由歸納法假設我們有aii(i)0

(i=1,2,

,k-1)，可用高斯消去法將A(1)約化到A(k)，即第45頁,共155頁，2024年2月25日，星期天利用行列式的性質，我們有第46頁,共155頁，2024年2月25日，星期天

由條件有Di0(i=1,2,

,k)，利用(1.13)式，則有aii(i)0(i=1,2,

,k)，定理6的對k成立.

必要性，由條件aii(i)0(i=1,2,

,k)，利用(1.13)式亦可推出Di0(i=1,2,

,k).

定理6得證.

推論如果A的順序主子式Di0(i=1,2,

,n-1),則第47頁,共155頁，2024年2月25日，星期天6.1.2矩陣的三角分解

下面我們借助矩陣理論進一步對消去法做些分析，從而建立高斯消去法與矩陣因式分解的關系.

設(1.1)的系數(shù)矩陣A∈Rn×n的各順序主子式均不為零，對于A施行初等行變換相當于用初等矩陣左乘A，于是對(1.1)施行第1次消元后化為(1.7)，這時A(1)化為A(2)，b(1)化為b(2)，即L1A(1)=A(2),L1b(1)=b(2)，第48頁,共155頁，2024年2月25日，星期天其中

一般第k步消元，A(k)化為A(k+1)，b(k)化為b(k+1)，相當于LkA(k)=A(k+1),Lkb(k)=b(k+1)，第49頁,共155頁，2024年2月25日，星期天其中重復以上過程，最后得到(1.14)第50頁,共155頁，2024年2月25日，星期天將上三角矩陣A(n)記為U，由(1.14)得到其中為單位下三角矩陣.第51頁,共155頁，2024年2月25日，星期天

這就是說，高斯消去法實質上產生了一個將A分解為兩個三角形矩陣相乘的因式分解，于是我們得到如下重要定理，它在解方程組的直接法中起著重要作用.

定理7(矩陣的LU分解)設A為n階矩陣，如果A的順序主子式Di0(i=1,2,

,n-1)，則A可分解為一個單位下三角矩陣L和一個上三角矩陣U的乘積，即A=LU,且這種分解是唯一的.第52頁,共155頁，2024年2月25日，星期天

證明根據(jù)以上高斯消去法的矩陣分析，A=LU的存在性已經得到證明，現(xiàn)僅在A為非奇異矩陣的假定下來證明唯一性，當A為奇異矩陣的情況留作練習，設A有兩種分解為A=LU=L1U1，其中L,L1為單位下三角矩陣，U,U1為上三角矩陣.由于L-1,U1-1存在，故有L-1L1=UU1-1.上式右邊為上三角矩陣，左邊為單位下三角矩陣，從而上式兩邊都必須等于單位矩陣I，故有

L=L1，U=U1.證畢.第53頁,共155頁，2024年2月25日，星期天

例3對于例2，系數(shù)矩陣由高斯消去法m21=0，m31=2，m32=-1，故從而得到求矩陣行列式的計算公式為第54頁,共155頁，2024年2月25日，星期天6.1.3列主元消去法

由高斯消去法知道,在消元過程中可能有akk(k)＝0的情況，這時消去法將無法進行；即使主元素akk(k)0但很小時，用其作除數(shù)，會導致其它元素數(shù)量級的嚴重增長和舍入誤差的擴散，最后也使得計算解不可靠.先看一個例子.

高斯消去法也稱主元素消去法(條件detA0)

即當akk(k)=0時，高斯消元法無法進行；或

|akk(k)|<<1時，帶來舍入誤差的擴散。（小主元）第55頁,共155頁，2024年2月25日，星期天解

(方法1)高斯消去法(用4位有效數(shù)字)例4(小主元產生了大誤差)第56頁,共155頁，2024年2月25日，星期天顯然計算得到的解x1=0是一個很壞的結果，不能作為方程組的近似解.其原因是我們在消元計算時用了小主元0.00001。使得約化后的方程組元素數(shù)量級大大增長，經再舍入使得在計算(0.3),(1.2),(1.3)元素時發(fā)生了嚴重的相消情況，因此經消元后得到的三角形方程組就相當不準確了，所以產生了很大的誤差，結果是不能用的.第57頁,共155頁，2024年2月25日，星期天

(方法2)列主元消元法，先交換行第58頁,共155頁，2024年2月25日，星期天

這個例子告訴我們，在采用高斯消去法解方程組時，小主元可能產生麻煩，故應避免絕對值小的主元素akk(k)，對一般矩陣來說，最好每一步選取系數(shù)矩陣(或消元后的低階矩陣)中絕對值最大的元素作為主元素，以使高斯消去法具有較好的數(shù)值穩(wěn)定性，這就是全主元素消去法，在選主元時要花費較多機器時間，目前主要使用的是列主元素消去法，下面介紹列主元素消去法，并假定線性方程組(1.1)的系數(shù)矩陣A∈Rn×n為非奇異的.第59頁,共155頁，2024年2月25日，星期天設方程組(2.1)的增廣矩陣為首先在A的第1列中選取絕對值最大的元素作為主元素，例如第60頁,共155頁，2024年2月25日，星期天然后交換B的第1行與第i1行，經第1次消元計算得(A|b)→(A(2)|b(2))重復上述過程，設已完成第k-1步的選主元素，交換兩行及消元計算，(A|b)約化為其中A(k)的元素仍記為aij，b(k)的元素仍記為bi.第61頁,共155頁，2024年2月25日，星期天第k步的選主元素(在A(k)右下角方陣的第1列內選)，即確定ik，使交換(A(k)|b(k))第k行與ik行的元素，再進行消元計算，最后將原方程組化為

(k=1,2,

,n-1)第62頁,共155頁，2024年2月25日，星期天

回代求解

算法1(列主元素消去法)設Ax=b.本算法用A的具有行交換的列主元素消去法，消元結果沖掉A，乘數(shù)mij沖掉aij，計算解x沖掉常數(shù)項b，行列式存放在det中.

1.det←1

2.對k=1,2,

,n-1

(1)按列選主元第63頁,共155頁，2024年2月25日，星期天

(2)如果aik,k=0，則計算停止(det(A)=0)

(4)消元計算對于i=k+1,

,n

①aik←mik=aik/akk

②對于j=k+1,

,n

aij←aij-mik·akj

(3)如果ik=k，則轉(4)

換行akj?aik,j(j=k,k+1,

,n)

bk?bik

det←-det

③bi←bi-mik·bk

(5)det←akk·

det第64頁,共155頁，2024年2月25日，星期天

3.如果ann=0，則計算停止(det(A)=0)

4.回代求解

(1)bn←bn/ann

(2)對于i=n-1,

,2,1

5.det←ann·

det

例4的方法2用的就是列主元素消去法.第65頁,共155頁，2024年2月25日，星期天

下面用矩陣運算來描述解(1.1)的列主元素消去法.列主元素消去法為其中Lk的元素滿足|mik|≤1

(k=1,2,

,n-1)，是初等置換矩陣.

利用(1.15)得到(1.15)第66頁,共155頁，2024年2月25日，星期天簡記為其中下面就n=4來考察一下矩陣其中可知亦為單位下三角矩陣，其元素的絕對值不超過1.記(1.16)第67頁,共155頁，2024年2月25日，星期天由(1.16)得到PA=LU.其中P為排列矩陣，L為單位下三角矩陣，U為上三角矩陣.這說明對(1.1)應用列主元消去法相當于對(A|b)先進行一系列行交換后對PAx=Pb再應用高斯消去法.在實際計算中我們只能在計算過程中做行的交換.總結以上的討論我們有第68頁,共155頁，2024年2月25日，星期天

定理8(列主元素的三角分解定理)如果A為非奇異矩陣，則存在排列矩陣P使PA=LU，其中L為單位下三角矩陣，U為上三角矩陣.在編程過程中，L元素存放在數(shù)組A的下三角部分，U元素存放在數(shù)組A的上三角部分，由紀錄主行的整型數(shù)組Ip(n)可知P的情況.第69頁,共155頁，2024年2月25日，星期天消元法的應用1.求行列式高斯消元法2.求逆矩陣(用高斯-若當消去法)第70頁,共155頁，2024年2月25日，星期天例子分別用高斯消元法和列主元消元法求矩陣的行列式.解：高斯消元法大數(shù)“吃掉”了小數(shù)！第71頁,共155頁，2024年2月25日，星期天列主元消元法第72頁,共155頁，2024年2月25日，星期天

全主元素消去法介紹在消元計算過程中，若每次選主元素不局限在方框列中，而在整個主子矩陣中選取,便稱為全主元高斯消去法,此時增加的步驟為:

①，確定ik

,jk,若=0，給出|Ａ|=0的信息，退出計算.第73頁,共155頁，2024年2月25日，星期天

行交換:(kjn)

列交換:(kin)②作如下行交換和列交換值得注意的是,在全主元的消去法中，由于進行了列交換，x各分量的順序已被打亂.因此必須在每次列交換的同時，讓機器“記住”作了一次怎樣的交換，在回代解出后將x各分量換回原來相應的位置，這樣增加了程序設計的復雜性.此外作①步比較大小時，全主元消去法將耗用更多的機時.第74頁,共155頁，2024年2月25日，星期天但全主元消去法比列主元消去法數(shù)值穩(wěn)定性更好一些.實際應用中，這兩種選主元技術都在使用.選主元素的高斯消元法是一種實用的算法，它可以應用于任意的方程組Ax=b，只要|A|≠0.第75頁,共155頁，2024年2月25日，星期天6.2.追趕法在一些實際問題中,例如解常微分方程邊值問題,解熱傳導方程以及船體數(shù)學放樣中建立三次樣條插值函數(shù)，都會要求解系數(shù)矩陣為對角占優(yōu)的三對角線方程組Ax=f，即:(2.1)第76頁,共155頁，2024年2月25日，星期天其中,當|i-j|>1時,aij=0,且滿足如下的對角占優(yōu)條件:①|b1|>|c1|>0;②|bi|≥|ai|+|ci|,ai,ci≠0,(i=2,3,…,n-1).③|bn|>|an|>0.我們利用矩陣的直接三角分解法來推導解三對角線方程組(2.1)的計算公式.由系數(shù)陣A的特點，可以將A分解為兩個三角矩陣的乘積，即A=LU,其中取L下三角陣,取U為單位上三角陣,這樣求解方程組Ax=f的方法稱為追趕法.設第77頁,共155頁，2024年2月25日，星期天其中為待定系數(shù)，比較(2.2)兩邊即得(2.2)第78頁,共155頁，2024年2月25日，星期天下面我們用歸納法證明(2.4)對i=1是成立的.現(xiàn)設(2.4)對i-1成立，求證對i亦成立.(2.3)(2.3)(2.4)第79頁,共155頁，2024年2月25日，星期天由歸納法假設0<|

i-1|<1，又由(2.4)及A的假設條件有(2.4)第80頁,共155頁，2024年2月25日，星期天簡記為Ax=f,A滿足條件（1）（2）（3）小結：第81頁,共155頁，2024年2月25日，星期天用歸納法可以證明，滿足上述條件的三對角線性方程組的系數(shù)矩陣A非奇異，所以可以利用矩陣的直接三角分解法來推導解三對角線性方程組的計算公式，用克洛特分解法，將A分解成兩個三角陣的乘積，設A=LU：

按乘法展開

則可計算

可依次計算當，由上式可惟一確定L和U第82頁,共155頁，2024年2月25日，星期天例6.9用追趕法求解三對角方程組

解由Ly=f

解出y又由Ux=y解出x第83頁,共155頁，2024年2月25日，星期天求解Ax=f等價于求解兩個三角形方程組.①Ly=f,求y;②Ux=y,求x.從而得到解三對角線方程組的追趕法公式.(1).計算{

i},{

i}的遞推公式小結：第84頁,共155頁，2024年2月25日，星期天我們將計算系數(shù)及的過程稱為追的過程,將計算方程組的解的過程稱為趕的過程.合起來就是追趕法.追趕法公式實際上就是把高斯消去法用到求解三對角線方程組上去的結果.這時由于A特別簡單,因此使得求解的計算公式非常簡單,而且計算量僅為5n-4次乘除法,而另外增加一個方程組Ax=f2僅增加3n-2次乘除法運算,易見追趕法的計算量是比較小的.第85頁,共155頁，2024年2月25日，星期天總結上述討論有下面定理

定理11設有三對角線方程組Ax=f,其中A滿足條件(1),(2),(3),則A為非奇異矩陣且追趕法計算公式中的{

i},{

i}滿足由定理11的(1),(2)說明追趕法計算公式中不會出現(xiàn)中間結果數(shù)量級的巨大增長和舍入誤差的嚴重累積.第86頁,共155頁，2024年2月25日，星期天補充：矩陣三角分解法高斯消去法有很多變形，有的是高斯消去法的改進、改寫，有的是用于某一類特殊矩陣的高斯消去法的簡化.下面我們將介紹矩陣的直接三角分解法，解特殊方程組用的平方根法及追趕法.

定義

如果L為單位下三角陣，U為上三角陣，則稱A=LU為杜里特爾(Doolittle)分解；如果L為下三角陣，U為單位上三角陣，則稱A=LU為克勞特(Crout)分解.第87頁,共155頁，2024年2月25日，星期天1、直接三角分解法(LU分解)

在6.1已經通過高斯消去法得到一個將A分解為一個單位下三角矩陣L和一個上三角矩陣U的乘積，A=LU，其中并由定理7得到這種分解是唯一的.第88頁,共155頁，2024年2月25日，星期天

將高斯消去法改寫為緊湊形式，可以直接從矩陣A的元素得到計算L、U元素的遞推公式，而不需要任何中間步驟，這就是所謂直接三角分解法.一旦實現(xiàn)了矩陣A的LU分解，那么求解Ax=b的問題就等價于求兩個三角形方程組.①Ly=b，求y；②Ux=y，求x.

1.不選主元的三角分解法設A為非奇異矩陣，且有分解式A=LU，其中L為單位下三角矩陣，U為上三角矩陣，即第89頁,共155頁，2024年2月25日，星期天其中第90頁,共155頁，2024年2月25日，星期天a11a12

a1k

a1n

u11u12

u1k

u1n

第1框a21a22

a2k

a2n

l21u22

u2k

u2n第2框

ak1ak2

akk

akn

lk1lk2

ukk

ukn第k框

an1an2

ank

ann

ln1ln2

lnk

unn第n框比較式A=LU

兩端的元素,按下圖所示順序逐框進行，先求ukj，后求lik.由第一框可得

第91頁,共155頁，2024年2月25日，星期天得假設前k-1框元素已求出，則由第92頁,共155頁，2024年2月25日，星期天有了矩陣A

的LU分解計算公式,解線性方程組Ax=b

就轉化為依次解下三角方程組

Ly=b

與上三角方程組Ux=y

其計算公式如下：

Ly=b

Ux=y第93頁,共155頁，2024年2月25日，星期天例6.10求矩陣解用緊湊形式計算的LU(Doolittle)分解.

u11=2

u12=1u13=4

第94頁,共155頁，2024年2月25日，星期天所以u11=2u12=1u13=4

得行列式第95頁,共155頁，2024年2月25日，星期天由于在計算機實現(xiàn)時當uri計算好后ari就不用了，因此計算好L,U的元素后就存放在A的相應位置.例如最后在存放A的數(shù)組中得到L,U的元素.由直接三角分解計算公式，需要計算形如∑aibi的式子，可采用“雙精度累加”，以提高精度.第96頁,共155頁，2024年2月25日，星期天如果已經實現(xiàn)了A=LU的分解計算，且L,U保存在A的相應位置，則用直接三角分解法解具有相同系數(shù)的方程組Ax=(b1,b2,…,bm)是相當方便的，每解一個方程組Ax=bj僅需要增加次乘除法運算.矩陣A的分解公式(3.2),(3.3)又稱為杜里特爾(Doolittle)分解.直接分解法大約需要n3/3次乘除法，和高斯消去法計算量基本相同.第97頁,共155頁，2024年2月25日，星期天

2.選主元的三角分解法從直接三角分解公式了看出當urr=0時計算將中斷，或者當urr絕對值很小時，按分解公式直接計算可能引起舍入誤差的積累.但如果A非奇異，我們可通過交換A的行實現(xiàn)矩陣PA的LU分解，因此可采用與列主元消去法類似的方法(可以證明下述方法與列主元消去法等價)，將直接三角分解法修改為(部分)選主元的三角分解法.第98頁,共155頁，2024年2月25日，星期天設第r-1步分解已完成，這時有設r步分解需用到(3.2)式及(3.3)式，為了避免用小的數(shù)urr作除數(shù)，引進量第99頁,共155頁，2024年2月25日，星期天于是有取交換A的r行與ir行元素，將調到(r,r)位置(將(i,j)位置的新元素仍記為lij及aij

),于是有|lir|

1(i=r+1,…,n).由此再進行第r步分解計算.第100頁,共155頁，2024年2月25日，星期天算法2(選主元的三角分解法)設Ax=b，其中A為非奇異矩陣.本算法采用選主元的三角分解法，用PA的三角分解沖掉A，用整型數(shù)組Ip(n)記錄主行，解x存放在b內.

1.對于r=1,2,…,n

(1)計算si

(2)選主元第101頁,共155頁，2024年2月25日，星期天上述計算過程完成后就實現(xiàn)了PA的LU分解，且U保存在A的上三角部分，L保存在A的下三角部分，排列矩陣P由Ip(n)最后記錄可知.(這時有|lir|

1)

(3)交換A的r行與ir行元素

(4)計算U的第r行元素，L的第r列元素第102頁,共155頁，2024年2月25日，星期天

求解方程組Ly=Pb及Ux=y.

2.對于i=1,2,

,n-1

(1)t←Ip(i)

(2)如果i=t則轉(3)

bi

←→bt

(3)(繼續(xù)循環(huán))

3.

4.第103頁,共155頁，2024年2月25日，星期天

利用算法2的結果(實現(xiàn)PA=LU三角分解)，則可以計算A的逆矩陣A-1=U-1L-1P.

利用PA的三角分解計算A-1步驟：

(1)計算上三角矩陣的逆陣U-1；

(2)計算U-1L-1；

(3)交換U-1L-1列(利用Ip(n)最后記錄).

上述方法求A-1大約需要n3次乘除運算.第104頁,共155頁，2024年2月25日，星期天6.3平方根法實際問題中Ax=b，系數(shù)矩陣A大多是對稱正定矩陣，所謂平方根法，就是利用對稱正定矩陣的三角分解而得到求解對稱正定方程組的一種有效方法，目前在計算機上廣泛應用平方根法解此類方程組.

對稱矩陣的LDLT

分解法當A為對稱矩陣,且其順序主子式均不為0時,由定理7可知A有唯一的LU分解為(3.1)式.為了利用A的對稱性，將U再分解為第105頁,共155頁，2024年2月25日，星期天其中D為對角矩陣，U0為上三角矩陣，于是

A=LU=LDU0.(3.6)又A=AT=U0TDLT,由分解的唯一性即得U0T=L.代入(3.6)得到對稱矩陣A的分解式A=LDLT.總結上述討論有下面定理.第106頁,共155頁，2024年2月25日，星期天

定理9(對稱矩陣的三角分解定理)設A為n階對稱矩陣，且A的各階順序主子式均不為零，則A可以唯一分解為A=LDLT,其中L是單位下三角矩陣，D是對角矩陣.如果A為對稱正定矩陣，則A的分解式A=LDLT中D的對角元素di均為正數(shù).事實上,由A的對稱正定性,定理6的推論成立,即第107頁,共155頁，2024年2月25日，星期天于是有由定理6得到其中L1=LD1/2為下三角矩陣.第108頁,共155頁，2024年2月25日，星期天

定理10(對稱正定矩陣的三角分解或喬雷斯基(Cholesky)分解)如果A為n階對稱正定矩陣，則存在一個實的非奇異下三角矩陣L使A=LLT，當限定L的對角元素為正時，這種分解是唯一的.下面我們用直接分解方法來確定計算L元素的遞推公式，因為第109頁,共155頁，2024年2月25日，星期天其中l(wèi)ii>0(i=1,2,

,n).由矩陣乘法及l(fā)jk=0(當j<k時),得于是得到解對稱正定方程組Ax=b的平方根法計算公式對于j=1,2,

,n第110頁,共155頁，2024年2月25日，星期天求解Ax=b，即求解兩個三角形方程組①Ly=b，求y；②LTx=y，求x.由計算公式(3.7)知所以有第111頁,共155頁，2024年2月25日，星期天于是上面分析說明，分析過程中元素ljk的數(shù)量級不會增長且對角元素ljj恒為正數(shù).于是有結論不選主元素的平方根法是一個數(shù)值穩(wěn)定的方法.當求出L的第j列元素時，LT的第j行元素亦算出.所以平方根約需n3/6次乘除法，大約為一般直接LU分解法計算量的一半.第112頁,共155頁，2024年2月25日，星期天例6.11用平方根法求解對稱正定方程組解首先對A進行Cholesky分解求解Ly=b，得y1=2,y2=3.5,y3=1.求解LTx=y，得x1=1,x2=1,x3=1.第113頁,共155頁，2024年2月25日，星期天例6.12平方根法求解方程組

解:因方程組系數(shù)矩陣對稱正定,設A=,即：由Ly=b解得由解得

由此例可以看出，平方根法解正定方程組的缺點是需要進行開方運算。為避免開方運算，我們改用單位三角陣作為分解陣，即把對稱正定矩陣A分解成

的形式，其中

第114頁,共155頁，2024年2月25日，星期天為對角陣，而

是單位下三角陣,這里分解公式為第115頁,共155頁，2024年2月25日，星期天據(jù)此可逐行計算

運用這種矩陣分解方法,方程組Ax=b即可歸結為求解兩個上三角方程組

和其計算公式分別為

和求解方程組的上述算法稱為改進的平方根法。這種方法總的計算量約為，即僅為高斯消去法計算量的一半。第116頁,共155頁，2024年2月25日，星期天總之，用平方根法解對稱正定方程組時，計算L的元素ljj需要用到開方運算.為了避免開方，采用定理10的分解式A=LDLT.即由矩陣乘法，并注意ljj=1,ljk=0(j<k)，得第117頁,共155頁，2024年2月25日，星期天于是得到L的元素及D的對角元素公式：為了避免重復計算，我們引進tij=lijdj

.由(3.9)得到按行計算L,T(T=LD)元素的公式：d1=a11

對于i=1,2,…,n.第118頁,共155頁，2024年2月25日，星期天求解Ly=b及DLTx=y的計算公式為計算公式(3.10),(3.11)稱為改進的平方根法.對于i=1,2,…,n.第119頁,共155頁，2024年2月25日，星期天6.4誤差分析6.4.1矩陣的條件數(shù)由于A(或b)元素是測量得到的,或者是計算的結果,在第一種情況A(或b)常帶有某些觀測誤差,在后一種情況A(或b)又包含有舍入誤差.因此我們處理的實際矩陣是A+A(或b+b),下面我們來研究數(shù)據(jù)A(或b)的微小誤差對解的影響.即考慮估計x-y,其中y是(A+A)y=b的解.考慮線性方程組Ax=b,其中設A為非奇異矩陣,x為方程組的精確解.第120頁,共155頁，2024年2月25日，星期天首先考察一個例子.

例8設有方程組記為Ax=b,它的精確解為x=(2,0)T

.現(xiàn)在考慮常數(shù)項的微小變化對方程組解的影響,即考察方程組第121頁,共155頁，2024年2月25日，星期天也可表示為A(x+x)=b+b,其中b=(0,0.0001)T,y=x+x,x為(6.1)的解.顯然方程組(5.2)的解為x+x=(1,1)T.我們看到(5.1)的常數(shù)項b的第2個分量只有1/1000的微小變化,方程組的解卻變化很大.這樣的方程組稱為病態(tài)方程組.

定義7如果矩陣A或常數(shù)項b

的微小變化(小擾動),引起方程組Ax=b解的巨大變化,則稱此方程組為“病態(tài)”方程組,其系數(shù)矩陣A稱為“病態(tài)”矩陣(相對于方程組而言),否則稱方程組為“良態(tài)”方程組,A稱為“良態(tài)”矩陣.第122頁,共155頁，2024年2月25日，星期天應該注意,矩陣的“病態(tài)”性質是矩陣本身的特性,下面我們希望找出刻畫矩陣“病態(tài)”性質的量.設有方程組

Ax=b,(5.3)其中A為非奇異矩陣,x為(5.3)的精確解.以下我們研究方程組的系數(shù)矩陣A(或b)的微小誤差(小擾動)時對解的影響.第123頁,共155頁，2024年2月25日，星期天

(1)現(xiàn)設A是精確的,x為Ax=b的精確解，當方程組右端有誤差b,受擾解為

x+x,則

A(x+x)=b+b,Ax=b,x=A-1b,

||x||≤||A-1||||b||.(5.4)由(5.3)有||b||≤||A||||x||.于是由(5.4)及(5.5),得第124頁,共155頁，2024年2月25日，星期天定理21設A是非奇異矩陣,Ax=b≠0

,且

A(x+x)=b+b,則上式給出了解x的相對誤差的上界,常數(shù)項b的相對誤差在解中放大||A-1||||A||倍.

(2)現(xiàn)設b是精確的,x為Ax=b的精確解，當A有微小誤差(小擾動)A,受擾解為

x+x,則(A+A)(x+x)=b,(A+A)x=-(A)x.(5.6)第125頁,共155頁，2024年2月25日，星期天如果A不受限制的話,可能A+A奇異,而(A+A)=A(I+A-1A).由定理20知,||A-1A||<1時,(I+A-1A)-1存在.由(5.6)式x=-(I+A-1A)-1A-1(A)x.因此設||A-1||||A||<1,即得第126頁,共155頁，2024年2月25日，星期天定理22設A是非奇異矩陣,Ax=b≠0

,且

(A+A)(x+x)=b.如果||A-1||||A||<1,則(5.7)式成立.如果A充分小,且在條件||A-1||||A||<1下,那么(5.7)式說明矩陣A的相對誤差||A||/||A||在解中可能放大||A-1||||A||倍.總之,量||A-1||||A||越小,由A(或b)的相對誤差引起的解的相對誤差就越小;量||A-1||||A||越大,解的相對誤差就可能越大.所以量||A-1||||A||事實上刻畫了解對原始數(shù)據(jù)變化的靈敏程度,即刻畫了方程組的“病態(tài)”程度,于是引進下述定義:第127頁,共155頁，2024年2月25日，星期天

(3)現(xiàn)設x為Ax=b的精確解，當A有微小誤差(小擾動)A,而b同時也有微小誤差b(小擾動)時,受擾解為

x+x,則還可以推出相對誤差估計式為

補充第128頁,共155頁，2024年2月25日，星期天定義8設A是非奇異矩陣,稱數(shù)cond(A)v=||A-1||v||A||v(v=1,2或

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

.由此看出矩陣的條件數(shù)與范數(shù)有關.矩陣的條件數(shù)是一個十分重要的概念.由上面討論知,當A的條件數(shù)相對的大,即cond(A)>>1時,則(5.3)是“病態(tài)”的(即A是“病態(tài)”矩陣,或者說A是壞條件的,相對于方程組),當A的條件數(shù)相對的小,則(5.3)是“良態(tài)”的(或者說A是好條件的).注意,方程組病態(tài)性質是方程組本身的特性.A的條件數(shù)越大,方程組的病態(tài)程度越嚴重,也就越難用一般的計算方法求得比較準確的解.第129頁,共155頁，2024年2月25日，星期天例如對前面例8的矩陣作分析由于條件數(shù)cond(A)1很大,可見矩陣A的病態(tài)程度十分嚴重,故由此方程組的解誤差非常大.計算其條件數(shù)第130頁,共155頁，2024年2月25日，星期天通常使用的條件數(shù),有(1)cond(A)∞

=||A-1||

||A||

;(2)A的譜條件數(shù);當A為對稱矩陣時其中

1,

n為A的絕對值最大和絕對值最小的特征值.第131頁,共155頁，2024年2月25日，星期天

條件數(shù)的性質:(1).對任何非奇異矩陣,都有cond(A)v≥1.事實上

cond(cA)v=cond(A)v.(2).設A為非奇異矩陣且c≠0(常數(shù)),則(3).如果A為正交矩陣,則cond(A)2=1;如果A為非奇異矩陣,R為正交矩陣,則

cond(RA)2=cond(AR)2=cond(A)2.第132頁,共155頁，2024年2月25日，星期天由性質1知,1≤cond(A),即cond(A)總是大于等于1的數(shù).條件數(shù)反映了方程組的“病態(tài)程度”.條件數(shù)越小,方程組的狀態(tài)越好,條