武漢大學(xué)高等工程數(shù)學(xué)課件演示文檔

《武漢大學(xué)高等工程數(shù)學(xué)課件演示文檔》由會(huì)員分享，可在線閱讀，更多相關(guān)《武漢大學(xué)高等工程數(shù)學(xué)課件演示文檔（362頁(yè)珍藏版）》請(qǐng)?jiān)谘b配圖網(wǎng)上搜索。

.,,第 2 章 解線性代數(shù)方程組的直接法與迭代法,.,引言,在自然科學(xué)和工程技術(shù)中，很多問(wèn)題的解決常常歸結(jié)為解線性代數(shù)方程組。例如: 船體數(shù)學(xué)放樣中建立三次樣條函數(shù)問(wèn)題， 用最小二乘法求實(shí)驗(yàn)數(shù)據(jù)的曲線擬合問(wèn)題， 非線性方程組求解的問(wèn)題， 用差分法或者有限元法解常微分方程、偏微分方程邊值問(wèn)題等 都導(dǎo)致求解線性方程組，且常常為求解大型線性方程組。,.,引言,與上述這些問(wèn)題有關(guān)的計(jì)算方法包括： 求解線性方程組的數(shù)值方法與計(jì)算矩陣的特征值及特征向量的數(shù)值方法。 在本章中，我們將注意力集中在線性方程組求解的直接方法和迭代方法。 也就是我們要學(xué)習(xí)的求解線性方程組的兩類(lèi)方法。,.,引言,線性方程組的這兩類(lèi)數(shù)值解法有以下特點(diǎn)。 直接法：經(jīng)過(guò)有限步算術(shù)運(yùn)算，可求得方程組的精確解的方法（若在計(jì)算過(guò)程中沒(méi)有舍入誤差）； 迭代法：用某種極限過(guò)程去逐步逼近線性方程組精確解的方法； 迭代法具有占用存儲(chǔ)單元少，程序設(shè)計(jì)簡(jiǎn)單，原始系數(shù)矩陣在迭代過(guò)程中不變等優(yōu)點(diǎn)，但存在收斂性及收斂速度等問(wèn)題。 下面由我們從小熟知的消元法開(kāi)始,.,2.1 高斯消元法,設(shè)線性方程組,,簡(jiǎn)記為 AX=b（矩陣形式/向量形式）,.,高斯消元法,,如果AX=b是n階線性代數(shù)方程組，則,.,高斯消元法,克萊姆法則在理論上有著重大意義，但在實(shí)際應(yīng)用中存在很大的困難，在線性代數(shù)中，為解決這一困難給出了高斯消元法。 我們解釋一下為何存在很大的困難,.,看來(lái)我們確實(shí)有必要學(xué)習(xí)方程組的求解方法， 從一個(gè)簡(jiǎn)單的例子開(kāi)始,Gauss消去法的計(jì)算量—這是我們介紹完Gauss消去法之后，對(duì)其計(jì)算量的分析，大家先看看 高斯法解n階線性方程組需要的總乘除法次數(shù)為 注1：當(dāng)n = 20時(shí)約為2670次，比克萊姆法則9.7?1020次大大減少。 注2：世界上最快的超級(jí)計(jì)算機(jī)每秒千萬(wàn)億次的浮點(diǎn)運(yùn)算能力,即：1015。 注3：目前最好的臺(tái)式或筆記本電腦每秒萬(wàn)億次的浮點(diǎn)運(yùn)算能力，即：1012 。 請(qǐng)同學(xué)們估算一下，用克萊姆法則解20階方程組各需要耗時(shí)多少？?。?！ 超級(jí)計(jì)算機(jī)9.7?1020 / 1015 /3600 =2.8小時(shí) 臺(tái)式或筆記本電腦9.7?1020 / 1012 /3600 =2800小時(shí),另外，還有存儲(chǔ)空間相關(guān)的空間復(fù)雜性問(wèn)題,.,例題,例1.用消元法解方程組,.,例題,第一步：-2 x（1）+(3)得,.,例題,第二步：1 x（2）+（4）,,解方程組（1）（2）（5）得：x=[1,2,3]T,.,,這樣，一般方程組的求解問(wèn)題就轉(zhuǎn)化為一個(gè)上三角形方程組的求解問(wèn)題。當(dāng)然，如果你愿意，同樣也可以轉(zhuǎn)化為一個(gè)下三角形方程組的求解問(wèn)題。,總結(jié),下面我們就來(lái)看看三角形方程組的求解方法！,.,2.1.1 高斯順序消元法,對(duì)下三角形方程的求解 設(shè) （1）,,.,高斯順序消元法,由（1）得,,,.,,對(duì)上三角方程組 設(shè),,.,,由（2）式回代得,,,.,一般線性代數(shù)方程組解法的思考,通過(guò)例子，我們得出了一般線性代數(shù)方程組的求解可以轉(zhuǎn)化為上或下三角形方程組來(lái)求解的結(jié)論，上面我們又研究了上或下三角形方程組的求解方法，下面是討論一般線性代數(shù)方程組求解的高斯消去法的時(shí)候了！。。。。。。,.,高斯順序消去法,對(duì)一般線性代數(shù)方程組設(shè) Ax=b.,,,記A(1)=A b(1)=b,1、第一次消元。設(shè),.,高斯順序消去法,,,.,高斯順序消去法,設(shè)第k-1次消元得A(k)x=b(k) 其中,,下面進(jìn)行第k次消元………,.,高斯順序消去法,第k次消元,,.,高斯順序消去法,最后,在完成n-1次消元后得到如下上三角形矩陣,這樣，一般方程組就變成了所謂的三角形方程組，問(wèn)題也就轉(zhuǎn)化為三角形方程組的求解問(wèn)題。,.,高斯順序消去法,上述過(guò)程，也就是對(duì)于方程組AX=b系數(shù)矩陣做如下變換：,,得到上三角方程組,.,高斯順序消去法,,,.,高斯順序消去法,,,.,高斯順序消去法,注意：算法對(duì)系數(shù)矩陣對(duì)角元素的非零要求,.,高斯順序消去法,.,高斯順序消去法算法框圖,,高斯順序消去法框圖,.,高斯消去法的計(jì)算量,,,.,Gauss消去法的計(jì)算量 以乘除法的次數(shù)為主 (1) 消元過(guò)程： 第k步時(shí)(n – k) + (n – k) (n – k + 1) = (n – k) (n – k + 2) 共有,.,Gauss消去法的計(jì)算量 (2) 回代過(guò)程： 求xi中，乘n–i次，除1次，共n–i+1次（i = 1，…，n–1） 共有,.,Gauss消去法的計(jì)算量 (3) 高斯法解n階線性方程組需要的總乘除法次數(shù)為 注1：當(dāng)n = 20時(shí)約為2670次，比克萊姆法則9.7?1020次大大減少。 注2：世界上最快的超級(jí)計(jì)算機(jī)每秒萬(wàn)萬(wàn)億次的浮點(diǎn)運(yùn)算能力,即：1016。 注3：目前最好的臺(tái)式或筆記本電腦每秒萬(wàn)億次的浮點(diǎn)運(yùn)算能力，即：1012 。 請(qǐng)同學(xué)們估算一下，用克萊姆法則解20階方程組各需要耗時(shí)多少？?。?！,.,高斯順序消去法條件,,,.,進(jìn)一步的考慮—解決辦法 (1) 若消元過(guò)程中出現(xiàn)akk(k) = 0，則無(wú)法繼續(xù) (2) 若akk(k) ≠ 0，但較小，則小主元做除數(shù)將產(chǎn)生大誤差 (3) 每次消元都選擇絕對(duì)值最大者作主元，稱(chēng)為高斯主元消去法 (4) 通常第k步選擇 ——第k列主對(duì)角元以下元素絕對(duì)值最大者作主元（該行與第k行對(duì)調(diào)），稱(chēng)為列主元消去法。,.,2.1.2 高斯主元素消去法,Gauss列主元消元法 從第一列中選出絕對(duì)值最大的元素，如果 ，則 交換第一行 和第i行，即,,,,,,交換,.,高斯列主元消去法算法,,.,高斯列主元消去法,,,,第k步 從 的第k列 ， ， 中選取絕對(duì)值最大項(xiàng)，記錄所在行，即 若 交換第k行與l行的所有對(duì)應(yīng)元素，再進(jìn)行順序消元。 下面，詳細(xì)描述一下列選主元的高斯消去法,,,.,高斯列主元消去法,.,高斯列主元消去法框圖,.,,注意到算法選主元素的第2步，是非正常結(jié)束，即 出現(xiàn)上述問(wèn)題的原因是在某步消元過(guò)程之前的列選主元失敗，無(wú)法選到合適的非零的適當(dāng)大的主元素。 解決辦法，在更大的范圍，即余下的n-k階子陣中來(lái)選擇主元素。 在第k到n列和第k到n行的子陣中選元素絕對(duì)值最大者作主元，稱(chēng)為全主元消去法。,.,2. 全主元消去法,先看一個(gè)例子.求解方程組,,.,全主元消去法,,.,全主元消去法,,.,全主元消去法,,很明顯，用列選主元法所得到 的解要比順序高斯消去法好， 但仍然有較大的誤差，下面介 紹所謂的全選主元高斯消去法：,.,全主元消去法,.,全主元消去法,.,全主元消去法,.,全主元消去法,,.,,從此可以看出，全主元法與列主元相比，其最大 的差別是需要根據(jù)列的交換情況，按照相反的順 序?qū)⒔獾姆至康捻樞蛑匦抡{(diào)整！,全主元消去法,.,全主元消去法,全主元消去法的嚴(yán)格 算法描述如下：,.,Gauss全主元消元算法,.,Gauss全主元消元算法,.,Gauss全主元消元算法,.,Gauss全主元消元算法,.,小結(jié),比較而言,Gauss順序消去法條件苛刻，且數(shù)值不穩(wěn)定；Gauss全主元消去法工作量偏大，需要比較 個(gè)元素及行列交換工作，算法復(fù)雜。因此從算法優(yōu)化的角度考慮， Gauss列主元消去法比較好。,.,思考,在上述討論中，Gauss順序消去法、Gauss列主元消去法、 Gauss全主元消去法的表示方法或算法描述過(guò)程中主要采用的是分量的形式。,但本章一開(kāi)始，我們就漂亮地將分量形式的方程組表示為了矩陣的形式“AX=b”，而從上述過(guò)程來(lái)看，無(wú)論是消元過(guò)程，還是回代求解，都可以同樣漂亮地表示為矩陣的形式。,也就是說(shuō)，至少在形式上方程組的求解問(wèn)題可以借助矩陣運(yùn)算以及相關(guān)的理論、方法。下面就來(lái)探討方程組求解的矩陣三角分解方法。首先，看看矩陣的三角分解！！,.,2.2 矩陣的三角分解法,我們知道，對(duì)矩陣進(jìn)行一次行的初等變換，就相當(dāng)于用相應(yīng)的初等矩陣去左乘原來(lái)的矩陣。因此我們用矩陣乘法來(lái)考察Gauss消元法。,,即可得到求解線性方程組的另一種直接法：矩陣的三角分解法。,.,三角分解法解方程的原理 若A = LU，則Ax = b ? LUx = b ? 其中L、U為三角形矩陣，則方程組易解 定義： (1) 稱(chēng) 為單位下三角陣 (2) 設(shè)L為單位下三角陣，U為上三角陣，若A = LU，則稱(chēng)A可三角（LU）分解，并稱(chēng)A = LU為A的三角分解（杜利特爾分解） (3) A也可以分解為一個(gè)下三角陣L和一個(gè)單位上三角陣U，此時(shí)， A的三角分解（LU）被稱(chēng)為克勞特分解。,.,三角分解可順利進(jìn)行的條件 定理：(1) A = (aij)n?n有唯一LU分解 ? A的順序主子式?k ≠ 0，k = 1，2，...，n – 1 (2) 若A = (aij)n?n可逆，則存在置換陣P，使PA的n個(gè)順序主子式全不等于0 注：由Ax = b ? PAx = Pb知，也就是說(shuō)，經(jīng)適當(dāng)行交換后，A總可三角分解,.,矩陣三角分解過(guò)程 設(shè)A的各階順序主子式不為0，且A = LU，即 (1) L陣主對(duì)角線（含）上邊：當(dāng)列標(biāo)m > 行標(biāo)k時(shí)，lkm = 0 j = k，k+1，...，n；k = 1，2，...，n (2)U陣主對(duì)角線下邊：當(dāng)行標(biāo)m > 列標(biāo)k時(shí)，ukm = 0 i = k+1，k+2，...，n；k = 1，2，...，n–1,.,矩陣三角分解計(jì)算公式,設(shè)A的各階順序主子式不為0，且A = LU，即 (1) 主對(duì)角線（含）上邊：當(dāng)列標(biāo)m > 行標(biāo)k時(shí)，lkm = 0,,,(2) 主對(duì)角線下邊：當(dāng)行標(biāo)m > 列標(biāo)k時(shí)，ukm = 0,欲求lik和ukj,.,矩陣三角分解計(jì)算公式----k=1時(shí) (1) 主對(duì)角線（含）上邊：當(dāng)列標(biāo)m > 行標(biāo)k時(shí)，lkm = 0 j=k,k+1, ...,n；k=1,2,...,n (2) 主對(duì)角線下邊：當(dāng)行標(biāo)m > 列標(biāo)k時(shí)，ukm = 0 i=k+1,k+2,...,n；k=1,2,...,n–1 欲求lik和ukj 設(shè)k = 1，那么，就可以計(jì)算L的第一行和U的第一列，有：a1j = l11u1j ? u1j = a1j j = 1，2，...，n（ U的第一行） ai1 = li1u11 ? li1 = ai1/u11 i = 2，3，...，n（ L的第一列）,,,.,矩陣三角分解計(jì)算過(guò)程----逐行（列）計(jì)算L、U 一般地，逐行（列）計(jì)算L、U 最后：,.,三角分解后求解兩個(gè)三角形方程組 下三角 Ly = b： 上三角 Rx = y： 緊湊格式：,.,事實(shí)上，上述過(guò)程就是高斯消去法的矩陣形式,.,第2步,,.,第n-1步完成后，就可以將A變換為上三角陣U,.,所以，A就有三角分解LU,.,一個(gè)例子，n=3的情況,.,Doolittle分解的一個(gè)例子，n=3的情況,.,Doolittle分解，一般情況,,.,Doolittle分解（L第一列，U第一行）,.,Doolittle分解（L第2列，U第2行）,.,Doolittle分解（L第k列，U第k行）,.,Doolittle分解（L第k列，U第k行）,,.,Doolittle分解,一般計(jì)算公式,.,Doolittle分解后求解兩個(gè)三角形方程組,.,例題,,.,例題,.,例題,.,例題,.,Doolittle分解算法,.,Doolittle分解算法,.,選主元的三角分解法,針對(duì)方程組Ax=b，矩陣A的三角分解能順利進(jìn)行的條件是A的各階順序主子式不等于零。同時(shí)，我們也已經(jīng)指出，三角分解法實(shí)質(zhì)上是高斯消去法的矩陣表示。 因此，高斯消去法求解方程組時(shí)存在的數(shù)值穩(wěn)定性問(wèn)題（除數(shù)為零或很?。?，三角分解法同樣存在。 為此，有必要與高斯消去法一樣進(jìn)行主元的選擇，以便能保證三角分解的順利進(jìn)行，使三角分解法具有好的數(shù)值穩(wěn)定性。 因而就有了選主元的三角分解法，只需要在直接三角分解法中引入選主元的算法，就可以實(shí)現(xiàn)選主元的三角分解法。這里就不詳細(xì)討論了！,下面針對(duì)一類(lèi)特殊的方程組，介紹相應(yīng)的三角分解法……,.,2.2.3 平方根法和改進(jìn)的平方根法--對(duì)稱(chēng)矩陣的Cholesky（喬列斯基）分解,在科學(xué)與工程計(jì)算中，線性方程組大多數(shù)的系數(shù)矩陣具有對(duì)稱(chēng)正定的性質(zhì)，因此基于對(duì)稱(chēng)正定矩陣這一特性的三角分解方法是求解對(duì)稱(chēng)正定方程組的一種有效方法。 由于系數(shù)矩陣A的良好性質(zhì)，分解過(guò)程無(wú)需選主元，具有良好的數(shù)值穩(wěn)定性。 首先，我們討論對(duì)稱(chēng)正定矩陣的三角分解，即Cholesky（喬列斯基）分解。,.,對(duì)稱(chēng)矩陣的Cholesky （喬列斯基）分解,A對(duì)稱(chēng)：AT=A A正定：A的各階順序主子式均大于零。即,.,對(duì)稱(chēng)矩陣的Cholesky分解,由Doolittle分解，對(duì)稱(chēng)正定矩陣A有唯一分解,,.,對(duì)稱(chēng)矩陣的Cholesky分解,定理2.2.4 設(shè)A為對(duì)稱(chēng)正定矩陣，則存在唯一分解A=LDLT，其中L為單位下三角陣，D=diag(d1,d2,…,dn)且di>0(i=1,…,n),.,對(duì)稱(chēng)矩陣的Cholesky分解,證明：,.,對(duì)稱(chēng)矩陣的Cholesky分解,所以，非奇異的上三角陣U1可化為,.,對(duì)稱(chēng)矩陣的Cholesky分解,.,對(duì)稱(chēng)矩陣的Cholesky分解,.,對(duì)稱(chēng)矩陣的Cholesky分解,推論：設(shè)A為對(duì)稱(chēng)正定矩陣，則存在唯一分解 其中L為具有主對(duì)角元素為正數(shù)的下三角矩陣。LLT稱(chēng)為對(duì)稱(chēng)正定矩陣A的Cholesky分解。,.,對(duì)稱(chēng)矩陣的Cholesky分解,證明：,,推論：設(shè)A為對(duì)稱(chēng)正定矩陣，則存在唯一分解 其中L為具有主對(duì)角元素為正數(shù)的下三角矩陣。,.,Cholesky分解的求法,,.,Cholesky分解的求法,.,Cholesky分解的求法,.,Cholesky分解法,,Cholesky分解法缺點(diǎn)及優(yōu)點(diǎn) 優(yōu)點(diǎn)：可以減少存儲(chǔ)單元。 缺點(diǎn)：存在開(kāi)方運(yùn)算，可能會(huì)出現(xiàn)根號(hào)下負(fù)數(shù)，同時(shí)計(jì)算量也大增。 怎么改進(jìn)？！,.,記得我們?cè)谟蒐DLT獲得LLT分解時(shí)，人為地將D分為D1/2* D1/2，這是導(dǎo)致平方根法計(jì)算量大等問(wèn)題的根源。因而，改進(jìn)的想法首先就是避免這種分割。具體地：,Cholesky分解法，怎么改進(jìn)？！,.,改進(jìn)Cholesky分解法,基于分解A=LDLT，其中,.,改進(jìn)的cholesky分解,.,改進(jìn)的cholesky分解,.,,.,例題,.,例題,.,例題,.,例題,.,改進(jìn)平方根法,A=LDLT 分解（對(duì)應(yīng)改進(jìn)平方根法） ，既適合于解對(duì)稱(chēng)正定方程組，也適合求解A為對(duì)稱(chēng)，而各階順序主子式不為零的方程組 而對(duì)A=LLT分解（對(duì)應(yīng)平方根法）只適合于對(duì)稱(chēng)正定方程組 上面，就是求解系數(shù)矩陣對(duì)稱(chēng)正定的方程組的平方根法和改進(jìn)平方根法。下面針對(duì)另外一類(lèi)特殊方程組求解的方法：追趕法！,.,2.2.4 三對(duì)角方程組求解的追趕法,.,三對(duì)角方程組求解的追趕法,.,三對(duì)角方程組求解的追趕法,.,三對(duì)角方程組求解的追趕法,.,三對(duì)角方程組求解的追趕法,其計(jì)算工作量為5n-4次乘除法（包括三角分解2n-2， 求解上三角方程n-1,回代2n-1）。工作量小. 其實(shí)現(xiàn)的條件為qi不為零。有以下定理可得三對(duì)角方程組求解的充分性條件。,.,三對(duì)角方程組求解的追趕法,.,解三對(duì)角矩陣線性方程組的追趕法程序框圖,.,從三對(duì)角矩陣帶狀矩陣：思考,如果方程組的系數(shù)矩陣的上半帶寬為s、下半帶寬為r,即,有如下三角分解的結(jié)果,.,從三對(duì)角矩陣帶狀矩陣：思考,有興趣的同學(xué)請(qǐng)自己針對(duì)你在科學(xué)研究中碰到的實(shí)際問(wèn)題，分析設(shè)計(jì)相應(yīng)的三角分解算法,.,直接法時(shí)間復(fù)雜性,.,解線性代數(shù)方程組的小結(jié),高斯消去法 順序高斯消去法 列選主元高斯消去法 全選主元高斯消去法 三角分解法 杜利特爾、克勞特三角分解法 平方根法、改進(jìn)的平方根法 追趕法 進(jìn)一步需要考慮的問(wèn)題：上述方法求解方程組可能存在的問(wèn)題。我們從研究方程組性態(tài)的工具----向量范數(shù)、矩陣范數(shù)開(kāi)始,.,2.3 向量和矩陣的范數(shù),為了研究線性方程組近似解的誤差估計(jì)和下一章要講的迭代法的收斂性，我們需要對(duì)Rn(n維向量空間)中的向量或Rnxn中矩陣的“大小”引入一種度量： ——向量和矩陣的范數(shù)。,.,向量和矩陣的范數(shù),在一維數(shù)軸上，實(shí)軸上任意一點(diǎn)x到原點(diǎn)的距離用|x|表示。而任意兩點(diǎn)x1，x2之間距離用| x1-x2 |表示。,,,.,向量和矩陣的范數(shù),而在二維平面上，平面上任意一點(diǎn)P(x,y)到原點(diǎn)的距離用 表示。而平面上任意兩點(diǎn)P1(x1,y1),P2(x2,y2)的距離用 表示。,推廣到n維空間，則稱(chēng)為向量范數(shù)。首先引入 向量范數(shù)的定義,.,向量范數(shù),,,.,常見(jiàn)的向量范數(shù),可以驗(yàn)證上述三種定義的向量空間的常用范數(shù)都滿足向量范數(shù)的定義中的三個(gè)條件，均為向量范數(shù)。分別稱(chēng)為1-范數(shù)，2-范數(shù)，∞-范數(shù)。,.,常見(jiàn)的向量范數(shù),參見(jiàn)教材61頁(yè)，驗(yàn)證2-范數(shù)滿足向量范數(shù)的三個(gè)條件,首先是非負(fù)性，對(duì)x≠0,有2-范數(shù)>0,對(duì)x=0,有2-范數(shù)=0；,其次是齊次性，對(duì)任意c,有cx的2-范數(shù)等于|c|乘x的2-范數(shù)；,然后是三角不等式，這里要用到Cauchy-Schwarz不等式，證明如下,.,常見(jiàn)的向量范數(shù),因而，2-范數(shù)是向量范數(shù)。,.,向量范數(shù)性質(zhì),向量范數(shù)是等價(jià)的。,.,向量范數(shù)性質(zhì),等價(jià)性質(zhì)：,,下面證明1）,還有其它形式的向量范數(shù)之間的等價(jià)關(guān)系，。,.,向量范數(shù)性質(zhì),等價(jià)性質(zhì)1）的證明：,,.,向量范數(shù)性質(zhì),等價(jià)性質(zhì)證明：,,.,2.3.2 矩陣范數(shù),.,相容范數(shù),,.,相容范數(shù),,.,算子范數(shù)是矩陣范數(shù),下面證明算子范數(shù)是矩陣范數(shù)！,.,算子范數(shù)是矩陣范數(shù),.,算子范數(shù)是矩陣范數(shù),.,算子范數(shù)是矩陣范數(shù),.,常見(jiàn)的矩陣范數(shù),.,常見(jiàn)的矩陣范數(shù),.,常見(jiàn)的矩陣范數(shù),.,常見(jiàn)的矩陣范數(shù),.,例題,.,2.3.3 矩陣的譜半徑和矩陣序列收斂性,.,譜半徑和矩陣序列的收斂性,.,例題,.,例題,.,,.,2.4 病態(tài)方程組與矩陣的條件數(shù),.,2.4.1 病態(tài)方程組與擾動(dòng)方程組的誤差分析,.,病態(tài)方程組與擾動(dòng)方程組的誤差分析,.,病態(tài)方程組與擾動(dòng)方程組的誤差分析,.,病態(tài)方程組與擾動(dòng)方程組的誤差分析,.,病態(tài)方程組與擾動(dòng)方程組的誤差分析,.,病態(tài)方程組,擾動(dòng)方程 由于計(jì)算機(jī)字長(zhǎng)限制，在解AX=b時(shí)，舍入誤差是不可避免的。因此我們只能得出方程的近似解 。 是方程組（A+△A）x=b+△b (1) 這就是所謂的擾動(dòng)方程。,,.,病態(tài)方程組,稱(chēng)（A+△A）x=b+△b為方程Ax=b的擾動(dòng)方程。其中△A， △b為由舍入誤差所產(chǎn)生的擾動(dòng)矩陣和擾動(dòng)向量。 當(dāng)△A， △b有微小擾動(dòng)，解得擾動(dòng)方程（1）的解與Ax=b的解x的相對(duì)誤差不大時(shí)，稱(chēng)為良態(tài)方程，否則為病態(tài)方程。 下面運(yùn)用前面介紹的向量范數(shù)、矩陣范數(shù)來(lái)研究方程組解得誤差傳播規(guī)律。,.,擾動(dòng)方程組的誤差界----右端項(xiàng)擾動(dòng)△b,.,擾動(dòng)方程組的誤差界----系數(shù)矩陣擾動(dòng)△A,.,擾動(dòng)方程組的誤差界----系數(shù)矩陣擾動(dòng)△A 、右端項(xiàng)擾動(dòng)△b,.,總結(jié)將上述三種情況，我們就能得到在系數(shù)矩陣有擾動(dòng)△A 和/或右端項(xiàng)有擾動(dòng)△b時(shí)，方程組的解的誤差限！在上述三種情況中都出現(xiàn)了一個(gè)決定誤差界大小的量，即||A||||A-1||。也就是這個(gè)數(shù)決定了方程組得性態(tài)，我們將它稱(chēng)之為條件數(shù)。,,,.,2.4.2 矩陣的條件數(shù),.,矩陣的條件數(shù)的性質(zhì),.,擾動(dòng)方程組的誤差界----定理,.,用條件數(shù)來(lái)刻畫(huà)方程組的性態(tài),.,例題，計(jì)算矩陣的條件數(shù),,這里，我們用了條件數(shù)的性質(zhì)來(lái)計(jì)算對(duì)稱(chēng)矩陣的條件數(shù)，也可以用條件數(shù)的定義來(lái)計(jì)算。即先求A的逆矩陣，再求A的條件數(shù)！,.,例：Hilbert 陣,注：現(xiàn)在用Matlab數(shù)學(xué)軟件可以很方便求矩陣的條件數(shù)!,對(duì)于嚴(yán)重的病態(tài)線性方程組，即使原始數(shù)據(jù)A和b都沒(méi)有誤差，但如果在求解過(guò)程中有舍入誤差，所得到的解也會(huì)有很大的相對(duì)誤差。,.,什么樣的線性方程組可能是病態(tài)的？,注：一般判斷矩陣是否病態(tài)，并不計(jì)算A?1，而由經(jīng)驗(yàn)得出。 ? 行列式很大或很?。ㄈ缒承┬小⒘薪葡嚓P(guān)）； ? 元素間相差大數(shù)量級(jí)，且無(wú)規(guī)則； ? 主元消去過(guò)程中出現(xiàn)小主元； ? 特征值相差大數(shù)量級(jí)。,.,緩和甚至解決線性方程組病態(tài)的手段,采用高精度的算法（不是總有效）； 通過(guò)行或列的平衡來(lái)降低矩陣的條件數(shù)（平衡法）； 殘差校正法； 通過(guò)矩陣的預(yù)條件處理來(lái)降低矩陣的條件數(shù)。,.,殘差校正法步驟：,求解Ax＝b，得到x的近似解x1； 校正x1，令r1＝b－Ax1，解AΔ x1 ＝ r1； x2＝x1＋ Δ x1 ； 令r2＝b－Ax2，解AΔ x2 ＝ r2； x3＝x2＋ Δ x2 ； 。。。 令rk＝b－Axk，解AΔ xk ＝ rk； Xk＋1＝xk＋ Δ xk ；。。。,.,預(yù)條件技術(shù),80年代提出，純代數(shù)觀點(diǎn) M 近似 A, 求解 　　　　　　　　　　(左預(yù)條件) M y = c 易解 相當(dāng)于化學(xué)反應(yīng)中尋找高效、廉價(jià)的催化劑，能極大地提高迭代方法的速度。,.,預(yù)條件技術(shù)(續(xù)),代數(shù)預(yù)條件技術(shù) ILU、SPAI、SOR、多項(xiàng)式 幾何預(yù)條件技術(shù) 某方向壓縮粗化、網(wǎng)格均勻化、區(qū)域規(guī)則化近似 分析預(yù)條件技術(shù) 變系數(shù)常數(shù)化、Green函數(shù)稀疏近似、強(qiáng)化橢圓型 物理預(yù)條件 量綱平衡、物理參數(shù)逼近、狀態(tài)方程近似、某些物理項(xiàng)簡(jiǎn)化、算子分裂 高級(jí)預(yù)條件 MG、DDM、快速變換（FFT）、基底變換法,.,解線性代數(shù)方程組的迭代解法,.,In Scientific Computing ↓ Large Linear Systems Ax=b as sub-problems/ as intermediate steps,Conjugate Gradient method for symmetric systems,,.,173,直接法和迭代法的比較,直接法: 經(jīng)過(guò)有限次運(yùn)算后可求得方程組精確解的方法(不計(jì)舍入誤差!),迭代法：從解的某個(gè)近似值出發(fā)，通過(guò)構(gòu)造一個(gè)無(wú)窮列去逼近精確解的方法。（一般有限步內(nèi)得不到精確解）,直接法比較適用于中小型方程組。對(duì)高階方程組，既使系數(shù)矩陣是稀疏的，但在運(yùn)算中很難保持稀疏性，因而有存儲(chǔ)量大，程序復(fù)雜等不足。 迭代法則能保持矩陣的稀疏性，具有計(jì)算簡(jiǎn)單，編制程序容易的優(yōu)點(diǎn)，并在許多情況下收斂較快。故能有效地解一些高階方程組。,.,解線性代數(shù)方程組的迭代法,引進(jìn)迭代法的目的： 解大型方程組； 存儲(chǔ)的考慮； 解病態(tài)方程 。,.,2.5 解線性代數(shù)方程組的迭代法,迭代發(fā)的一般理論,Jacobi迭代法、G-S迭代法,松弛迭代法,迭代法的收斂條件,迭代法收斂的其它判定方法,解線性方程組的迭代法,解線性代數(shù)方程組的極小化方法,.,2.5 解線性代數(shù)方程組的迭代法,迭代法的一般理論,Jacobi迭代法、G-S迭代法,松弛迭代法,解線性方程組的迭代法,2.5.1 一般迭代法,.,引例,解,迭代法的基本思想,.,或?qū)憺? ，,其中,迭代法的基本思想,.,迭代法的基本思想,.,迭代法的基本思想,.,我們?nèi)稳〕跏贾?，例? ，將這些值代入（2）式右邊（若（2）式為等式即求得方程組的解，但一般不滿足），得到新的值,得到向量序列和一般的計(jì)算公式（迭代公式）,迭代法的基本思想,.,其中k表示迭代次數(shù)（k=0,1,2, …）.,,迭代到第10次有,迭代法的基本思想,.,即建立一種從已知近似解到新的近似解的計(jì)算公式。,如何建立這樣的迭代計(jì)算公式呢？,迭代法的基本思想,迭代法求方程組的解的基本思想即是：構(gòu)造 一串收斂到解的序列，,.,思路,基本迭代法,.,基本迭代法,.,將線性方程組 Ax=b 化為 x=Mx+g,再由此構(gòu)造向量序列{x (k)}: x(k+1)=Mx (k)+g,若{x (k)}收斂至某個(gè)向量x *,則可得向量x *就是所求方程組 Ax=b 的準(zhǔn)確解.,迭代法的基本思想,我們已經(jīng)回答了 “如何建立這樣的迭代計(jì)算公式呢？”的問(wèn)題。方法是：,下面就來(lái)構(gòu)造我們這一章中將要介紹的第一種迭代法,.,Jacobi迭代法,.,Jacobi迭代法,.,Jacobi迭代法,.,Jacobi迭代法,.,求解一般線性代數(shù)方程組的Jacobi迭代法的計(jì)算公式,Jacobi迭代法,.,舉 例,按照jacobi法的一般計(jì)算公式寫(xiě)出本例的jacobi迭代式,.,Jacobi迭代法的特點(diǎn),Jacobi迭代法的特點(diǎn),.,Jacobi迭代法,.,Gauss-Seidel迭代法,.,Gauss-Seidel迭代法,.,Gauss-Seidel迭代法,.,即：,Gauss-Seidel迭代法,.,用高斯-賽德?tīng)柕ń夥匠探M,例4.1,Gauss-Seidel迭代法舉例,.,如此繼續(xù)下去…………,Gauss-Seidel迭代法的特點(diǎn),.,簡(jiǎn)單，每迭代一次只需進(jìn)行一次矩陣和向量乘法， 不改變系數(shù)矩陣的稀疏性， 僅需要一組存儲(chǔ)單元 。 Gauss-Seidel迭代法一般比Jacobi迭代法快，但也有比 Jacobi迭代法慢的，甚至有Jacobi迭代法收斂而Gauss- Seidel迭代法不收斂的。,Gauss-Seidel迭代法特點(diǎn),.,舉 例,請(qǐng)按照下述Gauss-Seidel法的一般計(jì)算公式寫(xiě)出該例的Gauss-Seidel迭代式,.,上面我們由介紹了一般迭代法的構(gòu)造方法，又由一般迭代方法具體引出了Jacobi迭代法，進(jìn)一步考慮迭代過(guò)程中新信息的利用，得到了Gauss-seidel迭代法。 當(dāng)然，這兩種特殊的迭代法同樣可以直接由系數(shù)矩陣的和、差分解導(dǎo)出。 下面我們?cè)诜治鲞@兩者方法存在誤差的基礎(chǔ)上，提出一種新的方法，即所謂松弛迭代法。首先來(lái)看松弛迭代法的基本思想。,松 弛 法,.,松弛法的基本思想,我們從為Gauss-Seidel 迭代法加速的角度來(lái)考慮提出松弛法。,Gauss-Seidel 迭代公式得到，于是有,這樣我們就可以得到這兩次近似解的差別△x,即,.,松弛法的基本思想,具體地,.,可以把 看作Gauss-Seidel 迭代的修正項(xiàng)，即第k次,近似解 以此項(xiàng)修正后得到新的近似解。,我們有理由相信，新的近似解，應(yīng)有更好的精度，誤 差更小。當(dāng)然，考慮到△x也是利用兩次近似估計(jì)得來(lái) 的，并不是真正的第k次近似解的誤差，所以在實(shí)際修 正時(shí)，更好的做法是對(duì)修正項(xiàng)△x加權(quán)！具體地有：,松弛法的基本思想,具體地，我們,.,得到新的近似解，具體公式為,松弛法是將 乘上一個(gè)參數(shù)因子 作為修正項(xiàng)而,松弛法的基本思想,所以，松弛法就可以進(jìn)一步表示為：,.,,稱(chēng)為松弛因子,當(dāng) 時(shí)稱(chēng)為低松弛；,時(shí)稱(chēng)為超松弛法。,松弛法的基本思想,上面我們運(yùn)用誤差估計(jì)—校正法，得到了松弛法的分 量形式，這已經(jīng)可以方便地讓我們編程計(jì)算了，但分 量形式不利于分析推理，下面將松弛法表為矩陣形式：,.,選取分裂矩陣M為帶參數(shù)下三角陣,其中 為可選擇的松弛因子.,松弛法的矩陣表示,.,松弛法的矩陣表示,通過(guò)簡(jiǎn)單的推導(dǎo)，得到,.,松弛法的矩陣表示,對(duì)照前面給出的分量形式！,.,SOR方法的 計(jì)算公式,SOR 方 法,對(duì)照前面給出的分量形式！,.,例3,取 用松弛法解方程組,舉 例,解：由,.,1. 1. 1. 1. 1. 1.56 1. 1.392 1.532 1.2744 1.3528 1.602 1.1372 1.40376 1.59164 1.22775 1.39396 1.60108 1.18467 1.40106 1.59889 1.20687 1.39917 1.60024 1.19667 1.40021 1.59984,1.20148 1.39988 1.60004 1.19932 1.40004 1.59998 1.2003 1.39998 1.60001 1.19987 1.40001 1.6 1.20006 1.4 1.6 1.19998 1.4 1.6 1.20001 1.4 1.6 1.2 1.4 1.6,計(jì)算結(jié)果,松弛法的程序,.,驗(yàn)證： 松弛因子與迭代次數(shù)的關(guān)系,例4,舉 例,.,舉 例,.,舉 例,.,舉 例,請(qǐng)按照下述松弛法的一般計(jì)算公式寫(xiě)出該例的松弛迭代式,.,2.5.2 迭代法的收斂條件,1、向量序列收斂與矩陣序列收斂,2、迭代法收斂概念,3、迭代法收斂定理,4、舉例,迭代法的收斂條件,.,不一定!!!,定義4.1,關(guān)注的問(wèn)題及幾個(gè)相關(guān)概念,.,,,,,關(guān)注的問(wèn)題及幾個(gè)相關(guān)概念,.,證明：,上述定理表明，向量序列和矩陣序列的收斂可歸結(jié)為對(duì)應(yīng)分量或?qū)?yīng)元素序列的收斂。,矩陣有類(lèi)似結(jié)果,關(guān)注的問(wèn)題及幾個(gè)相關(guān)概念,.,定義4.2,迭代矩陣,.,定義4.3,矩陣的譜半徑,回顧我們上一章中學(xué)習(xí)過(guò)的矩陣的譜半徑的定義,譜半徑是我們討論迭代法收斂性的重要工具！它 具有以下性質(zhì)：,.,矩陣的譜半徑,.,零矩陣的充要條件,有了譜半徑，就可以給出前面 定義的向量序列、矩陣序列收 斂的條件了！ 當(dāng)然也就可以討論迭代法的收 斂性問(wèn)題了。,.,定理4.1,證明：（必要性）,零矩陣的充要條件,.,（充分性）,零矩陣的充要條件,.,零矩陣的充要條件,有了上述向量序列、矩陣序列收 斂的條件，也就等價(jià)地有了---- 一般迭代法收斂條件了。我們有迭代法的如下基本定理：,.,定理4.2,,證明：,迭代法的收斂定理,.,零矩陣的充要條件,有了一般迭代法收斂的充要條件，從理論上解決了迭代法收斂性的判定問(wèn)題。但收斂的速度、誤差傳播的規(guī)律或者誤差估計(jì)的問(wèn)題仍然沒(méi)有解決。,還記得我們有 的結(jié)果嗎？如果我們能將迭代矩陣的條件加強(qiáng)為 ，我們應(yīng)該有希望得到我們想要的誤差估計(jì)：,.,和,誤 差 估 計(jì),.,誤 差 估 計(jì),.,有,由此可得,誤 差 估 計(jì),.,所以,誤 差 估 計(jì),由于,.,取范數(shù),誤 差 估 計(jì),證畢,.,誤差估計(jì)式,給出了停止迭代的一個(gè)判別準(zhǔn)則，即用相鄰兩次迭代結(jié)果的差的范數(shù)來(lái)判別，它越小，迭代序列與精確解越接近。,由誤差估計(jì)式,誤 差 估 計(jì),.,例：方程組,解：用雅可比迭代法,誤 差 估 計(jì) 舉 例,.,誤 差 估 計(jì) 舉 例,.,迭代法收斂性的判定,2、若A為嚴(yán)格對(duì)角占優(yōu)陣或不可約弱對(duì)角占優(yōu)陣，則Jacobi迭代法和Gauss-Seidel迭代法均收斂。,小 結(jié),.,三種迭代法收斂性判定的其它定理,前面，我們已經(jīng)給出了一個(gè)一般迭代法收斂判定的充要條件和一個(gè)迭代矩陣的某種范數(shù)小于1的特定迭代法收斂的充分條件。下面，我們針對(duì)一些更特別的方程組的迭代法來(lái)給出幾個(gè)收斂判定定理。,.,若n階方陣 滿足,定義1,定義2,如果矩陣A不能通過(guò)行的互換和相應(yīng)的列互換成,為形式 ，其中 為方陣，則稱(chēng)A為,不可約。,且至少有一個(gè)i值,使得上式中不等號(hào)嚴(yán)格成立，則稱(chēng)A,為弱對(duì)角占優(yōu)陣。若對(duì)所有i,上式不等號(hào)均嚴(yán)格成立 ，,則稱(chēng)A為嚴(yán)格對(duì)角占優(yōu)陣。,定 義,例,.,定理,定 理,設(shè)有線性方程組 ，下列結(jié)論成立,1. 若A為嚴(yán)格對(duì)角占優(yōu)陣或不可約弱對(duì)角占優(yōu)陣，則,Jacobi 迭代法和Gass-Seidel迭代法均收斂。,2.若A為嚴(yán)格對(duì)角占優(yōu)陣， 則松弛法收斂。,3.若A為對(duì)稱(chēng)正定陣，則松弛法收斂的充要條件為,上兩例中：,A為嚴(yán)格對(duì)角占優(yōu)陣， Jacobi 迭代法和Gass-Seidel迭代均收斂。B為非嚴(yán)格對(duì)角占優(yōu)陣， 故松弛法收斂。,.,例5,解：因A為對(duì)稱(chēng)且各階主子式皆大于零，故A為對(duì)稱(chēng)正定,均收斂。A不是弱對(duì)角占優(yōu)陣，故不能用條件1判斷。,陣。由判別條件3，Gauss-Seidel迭代法與松弛法,Jacobi迭代法迭代矩陣,舉 例,.,其特征方程,舉 例,.,Jacobi與Gauss-Seidel迭代的迭代矩陣分別為,改變方程組中方程的次序，即將系數(shù)矩陣作行,交換，會(huì)改變迭代法的收斂性,例,舉 例,.,若交換方程的次序，得 Ax=b 的同解方程組,為嚴(yán)格對(duì)角占優(yōu)矩陣，因而對(duì)方程組 用,Jacobi 與Gass-Seidel 迭代求解均收斂。,舉 例,.,對(duì)稱(chēng)超松弛迭代法（SSOR法）,.,,對(duì)稱(chēng)超松弛迭代法（SSOR法）,.,SSOR迭代法的矩陣形式：,注： （1）關(guān)于 的收斂條件和準(zhǔn)則與SOR方法相同； （2）收斂快慢對(duì) 的選取不敏感。,.,塊超松弛迭代法（BSOR法）,.,案 例 1 (熱傳導(dǎo)問(wèn)題),設(shè)有一維熱傳導(dǎo)方程的初邊值問(wèn)題,試用數(shù)值方法求出t=0.2時(shí)刻金屬桿的溫度分布.,應(yīng)用實(shí)例,.,,解:(1)對(duì)空間進(jìn)行離散.,(2)對(duì)微分算子進(jìn)行離散.,,,,,,,,,,,,,,t 0.02 0.01 0,0 0.1 0.2 …. 0.9 1 x,,.,,采用無(wú)條件穩(wěn)定的Crank-Nicholson格式,則有,.,或,加上邊界條件后有,.,加上邊界條件后有,其矩陣形式為 ATn+1=dn,.,案 例 2,數(shù)值求解正方形域上的Poisson方程邊值問(wèn)題,.,,解:(1)剖分求解域.,(2)對(duì)微分算子進(jìn)行離散.,,,,,,,,,,,,,,0 1 2 …. N N+1 X,,Y N+1 N ： 2 1 0,.,在每個(gè)點(diǎn)(xi,yj)上的有限差分方程為,在邊界上,又稱(chēng)為五點(diǎn)差分格式,.,對(duì)非邊界點(diǎn)進(jìn)行編號(hào): 順序?yàn)?----從下往上,從左往右,相應(yīng)的解向量和右端向量分別為,.,.,.,差分方程組的矩陣形式為 Au=f,如果把每一條線上的網(wǎng)點(diǎn)看作一個(gè)組,如,.,其中,可用塊迭代法 求解Au=f,.,.,2.5.3解線性代數(shù)方程組的極小化方法,一、與線性方程組等價(jià)的變分問(wèn)題,三、共軛斜量法(共軛梯度法),四、預(yù)條件共軛斜量法,二、最速下降法,.,共軛梯度法最早是由計(jì)算數(shù)學(xué)家Hestenes 和幾何學(xué)家Stiefel 在20世紀(jì)50年代初為求解線性方程組而各自獨(dú)立提出的。他們合作的文章被公認(rèn)為共軛梯度法的奠基之作，詳細(xì)討論了求解線性方程組的共軛梯度法的性質(zhì)以及它和其他方法的關(guān)系。 在A為對(duì)稱(chēng)正定陣時(shí)，上述線性方程組等價(jià)于最優(yōu)化問(wèn)題。,,解線性代數(shù)方程組的極小化方法,.,由此，Hestenes和Stiefel的共軛梯度方法也可視為求二次函數(shù)極小值的共軛梯度法。 1964年Fletcher和Reeves 將此方法推廣到非線性?xún)?yōu)化，得到求一般函數(shù)極小值的共軛梯度法。,,解線性代數(shù)方程組的極小化方法,.,提到最優(yōu)化問(wèn)題，首先介紹何為最速下降法。 考慮線性方程組Ax=b的求解問(wèn)題，其中A是給定的n階對(duì)稱(chēng)正定矩陣，b是給定的n維向量。 為此我們定義二次泛函,解線性代數(shù)方程組的極小化方法,.,定理：設(shè)A對(duì)稱(chēng)正定，求方程組Ax=b的解，等價(jià)于求二次泛函 的極小點(diǎn) 由定理，求解線性方程組的問(wèn)題就轉(zhuǎn)化為求二次泛函的極小點(diǎn)的問(wèn)題，也可表為與線性方程組等價(jià)的變分問(wèn)題：,解線性代數(shù)方程組的極小化方法,.,一、與線性方程組等價(jià)的變分問(wèn)題,,.,一、與線性方程組等價(jià)的變分問(wèn)題,.,一、與線性方程組等價(jià)的變分問(wèn)題,.,一、與線性方程組等價(jià)的變分問(wèn)題,.,一、與線性方程組等價(jià)的變分問(wèn)題,.,,一、與線性方程組等價(jià)的變分問(wèn)題,.,.,一、與線性方程組等價(jià)的變分問(wèn)題,.,一、與線性方程組等價(jià)的變分問(wèn)題,.,一、與線性方程組等價(jià)的變分問(wèn)題,.,二、最速下降法,.,二、最速下降法,.,二、最速下降法,.,,,二、最速下降法,.,三、共軛斜量法 (CG) (共軛梯度法),.,三、共軛斜量法 (CG) (共軛梯度法),.,三、共軛斜量法 (CG) (共軛梯度法),.,三、共軛斜量法 (CG) (共軛梯度法),.,,三、共軛斜量法 (CG) (共軛梯度法),.,三、共軛斜量法 (CG) (共軛梯度法),.,三、共軛斜量法 (CG) (共軛梯度法),.,公式化簡(jiǎn),.,三、共軛斜量法 (CG) (共軛梯度法),.,三、共軛斜量法 (CG) (共軛梯度法),.,三、共軛斜量法 (CG) (共軛梯度法),.,公式化簡(jiǎn),三、共軛斜量法 (CG) (共軛梯度法),.,三、共軛斜量法 (CG) (共軛梯度法),.,.,四、預(yù)條件共軛斜量法(PCG),.,,預(yù)條件共軛斜量法,實(shí)際計(jì)算，可通過(guò)變換，轉(zhuǎn)化成用原方程組的量來(lái)計(jì)算。,.,.,,預(yù)優(yōu)矩陣的選取,下面介紹幾種選取預(yù)優(yōu)矩陣的方案:,.,幾種選取預(yù)優(yōu)矩陣的方案,.,幾種選取預(yù)優(yōu)矩陣的方案,.,幾種選取預(yù)優(yōu)矩陣的方案,.,三種預(yù)條件共軛斜量法在MATLAB中的調(diào)用格式: 1.不用預(yù)優(yōu)矩陣的共軛斜量法 x=pcg(a,b,tol,kmax) 2.用預(yù)優(yōu)矩陣的共軛斜量法 (1)x=pcg(a,b,tol,kmax,m) (2) r=chol(m) x=pcg(a,b,tol,kmax,r’,r,x0) 3.未給定預(yù)優(yōu)矩陣的共軛斜量法 r=cholinc(sa,’0’) x=pcg(a,b,tol,kmax,r’,r,x0),.,CG(Conjugate Gradient)（對(duì)稱(chēng)正定） MRES(Minimal RESidual) (對(duì)稱(chēng)不正定） GMRES(Generilized Minimal RESidual） （不對(duì)稱(chēng)不正定） BiCG(不對(duì)稱(chēng)不正定） BiCGSTAB (不對(duì)稱(chēng)不正定）,幾種常用的Krylov-subspace方法,.,背景,大規(guī)模線性方程組的求解 很多科學(xué)工程計(jì)算問(wèn)題都轉(zhuǎn)化為求解方程組Ax=b.如偏微分方程組的差分格式，有限元方法離散得到剛度矩陣. 大規(guī)模矩陣特征值和特征向量的計(jì)算 工程計(jì)算領(lǐng)域十分常見(jiàn)。如量子物理中的Kohn-Sham方程求解化為哈密頓矩陣某些關(guān)鍵特征值對(duì)的計(jì)算.,Krylov-subspace方法,.,投影方法,線性方程組的投影方法 方程組Ax=b,A是n×n的矩陣. 給定初始x(0),在m維空間K(右子空間）中尋找x的近似解x(1) 滿足殘向量r=b-Ax(1) 與m維空間L(左子空間）正交，即 b-Ax(1) ⊥L 此條件稱(chēng)為Petrov-Galerkin條件 . 當(dāng)空間K=L時(shí)，稱(chēng)相應(yīng)的投影法為正交投影法，否則稱(chēng)為斜交投影法.,,,,,,,,,,K(L),X(0),X(1),K,L,X(0),X(1),.,解方程組的投影法的矩陣表示 　　設(shè)n×m階矩陣V=[v(1), v(2), …v(m)]與W=[w(1), w(2), …w(m)]的列分別構(gòu)成K與L的一組基 .記z=x(1)-x(0),z=Vy,有 　 當(dāng)　　　非奇異時(shí)有 從而得到迭代公式,,.,投影方法的最優(yōu)性 １. (誤差投影)設(shè)A為對(duì)稱(chēng)正定矩陣, x(0)為初始近似解,且K=L,則x(1)為采用投影方法得到的新近似解的充要條件是 其中, 且x為問(wèn)題的精確解.,,.,,２. (殘量投影)設(shè)A為任意方陣, x(0)為初始近似解,且 L=AK, 則x(1)為采用投影方法得到的新近似解的充要條件是 其中,.,,矩陣特征值的投影方法 　 對(duì)于特征值問(wèn)題Ax=λx,其中A是n×n的矩陣，斜交投影法是在m維右子空間K中尋找xi和復(fù)數(shù)λi滿足Axi- λixi ⊥L,其中L為m維左子空間.當(dāng)L=K時(shí)，稱(chēng)此投影方法為正交投影法.,.,Kyrlov子空間,定義為m維Krylov子空間 　μ為v的次數(shù)，即使得q(A)v＝０的非零首一多項(xiàng)式的最低次數(shù),.,Krylov子空間的標(biāo)準(zhǔn)正交基,Arnoldi方法 基于Gram-Schmit正交化方法 首先，選取一個(gè)Euclid范數(shù)為1的向量v(1),對(duì) ,通常可取 在已知v(1),v(2),……v(j)的情況下,不妨設(shè)v(1),v(2),……v(j),Av(j)線性無(wú)關(guān)(否則構(gòu)造完畢),則可求出與每個(gè)都正交的向量,.,,而不難看出　　　　　　　　　　　，再記　　　　　　　　　　　　　　　　　　　　　　　　　　　 　　　　　　　　　　,得到與v(1),v(2),……v(j)都正交的向量 　　　　　　重復(fù)此過(guò)程，即可得到一組標(biāo)準(zhǔn)正交基.若期間某個(gè)j使得hj+1,j=0, 則說(shuō)明v的次數(shù)是j,且Kj是A的不變子空間.,.,,,.,,定理 　 如果記以v(1),v(2),……v(m)為列構(gòu)成的矩陣為Vm,由hij定義的(m+1)×m階上Hessenberg矩陣為 刪除最后一行得到的矩陣為Hm,則 在Arnoldi算法中，可能有較大舍入誤差，改寫(xiě),.,.,Krylov子空間方法解線性方程組,誤差投影型方法 取L=K時(shí)的正交投影法 　 １）非對(duì)稱(chēng)矩陣的FOM方法,.,.,,不難看出，當(dāng)采用上述FOM算法時(shí)，需要存儲(chǔ)所有的vi，(i=1,2,…m)，當(dāng)m增大時(shí)，存儲(chǔ)量以O(shè)(mn)量級(jí)增大.而FOM計(jì)算量是O(m2 n).可見(jiàn)其代價(jià)十分高昂.因此我們考慮重啟的?FOM算法,.,.,Krylov子空間方法解線性方程組,誤差投影型方法 取L=K時(shí)的正交投影法 　 １）非對(duì)稱(chēng)矩陣的FOM方法 2) 非對(duì)稱(chēng)矩陣的IOM方法和DIOM方法,.,,所謂不完全正交化方法(IOM),是指在正交化過(guò)程中,v(j+1)僅與最近k個(gè)v(j-k+1),…v(j-1),v(j)正交,這樣做雖然破壞的正交性,但是降低了計(jì)算量.當(dāng)然k選得越小，對(duì)每個(gè)j對(duì)應(yīng)的計(jì)算量也越小,但可能要選更大的m才能取得滿足精度要求的近似解. IOM算法僅僅是把FOM算法的第三步改為 (iii’)對(duì)i=max(1,j-k+1),…,j, 計(jì)算hij與w(j).,.,,但采用IOM后，仍然需要存儲(chǔ)v(1), v(2), …v(m),因?yàn)樵诘?vi)步　　　　　　　　　中仍然需要這些向量. 解決這個(gè)問(wèn)題可以考慮采用H的LU分解，通過(guò)自身分解的迭代更新以減少每一步的存儲(chǔ)量,.,.,.,Krylov子空間方法解線性方程組,誤差投影型方法 取L=K時(shí)的正交投影法 　 １）非對(duì)稱(chēng)矩陣的FOM方法 2 ) 非對(duì)稱(chēng)矩陣的IOM方法和DIOM方法 ３）對(duì)稱(chēng)矩陣Lanczos算法,.,,A是非對(duì)稱(chēng)陣時(shí),H是Hessenberg陣，則當(dāng)A是對(duì)稱(chēng)陣時(shí)，H也為對(duì)稱(chēng)陣，故應(yīng)為三對(duì)角陣，記為 對(duì)它采用修正Arnoldi正交化過(guò)程，就得到Lanczos方法,.,.,Krylov子空間方法解線性方程組,誤差投影型方法 取L=K時(shí)的正交投影法 　 １）非對(duì)稱(chēng)矩陣的FOM方法 2 ) 非對(duì)稱(chēng)矩陣的IOM方法和DIOM方法 ３）對(duì)稱(chēng)矩陣Lanczos算法 4 ）正定對(duì)稱(chēng)陣的CG法,.,,CG法是解正定對(duì)稱(chēng)線性方程組最有效的方法之一 該方法充分利用了矩陣A的稀疏性，每次迭代的主要計(jì)算量是向量乘法．而實(shí)際上，CG方法也是一種Krylov子空間方法 后面將給出一個(gè)數(shù)值例子,.,.,Krylov子空間方法解線性方程組,殘量投影型方法 取L=AK時(shí)的斜交投影法 　 １）GMERS方法,.,,GMERS方法把方程組的求解化為一個(gè)極小問(wèn)題的求解，回顧之前介紹投影類(lèi)方法,.,回顧,２. (殘量投影)設(shè)A為任意方陣, x(0)為初始近似解,且 L=AK, 則x(1)為采用投影方法得到的新近似解的充要條件是 其中,.,,GMERS方法在Arnoldi正交化過(guò)程之后把方程組的求解化為一個(gè)極小問(wèn)題的求解，回顧之前介紹投影類(lèi)方法 不去直接求x(1),轉(zhuǎn)而去解此極小問(wèn)題，這是GMERS方法的獨(dú)到之處.再由之前的定理,.,.,.,Krylov子空間方法解線性方程組,殘量投影型方法 取L=AK時(shí)的斜交投影法 　 １）GMERS方法 2 ) 重啟型GMERS方法、QGMERS、DGMERS,.,,GMERS算法具有很好的收斂性，在不考慮舍入誤差的情況下，該算法能在在有限步得到精確解 .但是，和FOM算法類(lèi)似，它需要保存大量正交基，因此我們可以采取類(lèi)似重啟型FOM的算法實(shí)現(xiàn)重啟型GMERS算法,.,.,,同樣可以采取不完全正交的方法，在正交化過(guò)程中,v(j+1)僅與最近k個(gè)v(j-k+1),…v(j-1),v(j)正交就能得到QGMERS方法. 為了避免存儲(chǔ)v(1),v(2)…v(m-k)，可以對(duì) 進(jìn)行QR分解.這種算法被稱(chēng)為DQGMERS,.,.,Krylov子空間方法解矩陣特征值問(wèn)題,Arnoldi方法求解非對(duì)稱(chēng)矩陣特征值 Lanczos方法求解對(duì)稱(chēng)矩陣特征值,.,,Arnoldi方法 再次回顧前面的定理 而 所以有以下結(jié)論：,.,,1）若 ,則 是A的不變子空間,從而Hm的所有特征值是A的特征值子集 2）若 ,設(shè)Hm的特征值對(duì)是 ,即 ,由上述定理可得,,.,,以上討論說(shuō)明,可以用Hm的特征值 (又稱(chēng)Ritz值)來(lái)近似A的特征值,用向量 (又稱(chēng)Ritz向量)來(lái)近似相應(yīng)的特征向量.事實(shí)上, Hm為A在Kyrlov子空間上的正交投影. 由于Hm是上Hessenberg陣，它的特征值求解難度遠(yuǎn)小于原問(wèn)題,例如可以采取分解法來(lái)求解.,.,.,,Arnoldi算法構(gòu)造標(biāo)準(zhǔn)正交基 1需要存儲(chǔ)所有的基向量,當(dāng)m很大時(shí),存儲(chǔ)量大 2理論上為了保證收斂速度,m越大越好 矛盾！,.,,Lanczos 方法 A是對(duì)稱(chēng)陣時(shí),Hm是三對(duì)角陣，仍然采用Hm的特征值來(lái)近似A的特征值，相應(yīng)的Ritz向量來(lái)近似A的特征向量. 問(wèn)題轉(zhuǎn)化為三對(duì)角陣特征值的求解問(wèn)題，而它的求解，復(fù)雜度大大減小，有很多成熟的辦法，如分而治之法，QR法等，可以在并行機(jī)上達(dá)到O(logN)的量級(jí).,.,.,,對(duì)稱(chēng)Lanczos算法,在沒(méi)有舍入誤差的情形下,得到的是標(biāo)準(zhǔn)正交基相互正交的,并且至多n步必然終止.但是在出現(xiàn)舍入誤差的情況下，計(jì)算得到的將失去正交性.因此,長(zhǎng)期以來(lái)被人們認(rèn)為是不穩(wěn)定的.直到1971年,Paige通過(guò)仔細(xì)的舍入誤差分析,發(fā)現(xiàn)失去正交性恰與近似特征值的精度提高有關(guān).之后,經(jīng)過(guò)大量的理論分析和數(shù)值實(shí)驗(yàn),人們充分認(rèn)識(shí)到Lanczos方法是求解大型稀疏對(duì)稱(chēng)矩陣特征值問(wèn)題的最有效方法之一. 目前，Lanczos方法是求大型稀疏對(duì)稱(chēng)矩陣部分極端特征值的最有效方法．,