科學(xué)計算方法6(高斯消元法)

1、設(shè)3種食物每100克中蛋白質(zhì)、碳水化合物和脂肪的含量 如下表; 表中還給出了美國流行的劍橋大學(xué)醫(yī)學(xué)院的簡捷營養(yǎng)處方。 現(xiàn)在的問題是:如果用這3種食物作為每天的主要食物,那么它們的 用量應(yīng)各取多少才能全面準(zhǔn)確地實現(xiàn)這個營養(yǎng)。,大規(guī)模代數(shù)系統(tǒng)的求解是科學(xué)與工程計算的基礎(chǔ)和關(guān)鍵,其計算量也占有非常大的比重, 有的甚至占整個計算的80%以上！,Gene Golub (19322007) 美國科學(xué)院,工程院、藝術(shù)與科學(xué)學(xué)院院士 以及瑞典皇家工程科學(xué)院院士,求解,線性方程組的矩陣形式,A x = b,增廣矩陣:,存在性,如果rank(A,b)=rank(A), 則線性方程組的有解。,秩是定量刻畫矩陣性質(zhì)的

2、指標(biāo), 相對于矩陣行列式更加精細(xì)。,唯一性,如果系數(shù)矩陣非奇異即|A| 0, 則線性方程組的解唯一。,如果系數(shù)矩陣奇異即|A| = 0,則線性方程組的有無窮多解或無解。例如, 2x + y = 3和4x + 2y = 6有無窮多解; 而2x + y = 3和4x + 2y = 0則無解。,求解如下下三角線性方程組Lx=c:,從第一個方程開始, 從前往后每次求解一個方程, 則下三角方程組的解可以表示如下,類似的, 求解如下上三角線性方程組Ux=c:,從最后一個方程開始, 從后往前每次求解一個方程, 則上三角方程組的解可以表示如下,例,求解方程組Ax = b, 其中,假設(shè)系數(shù)矩陣A存在如下LU分解

3、:,解:,首先求解下三角線性方陣組Ly = b:,然后求解上三角方程組Ux = y :,Erds 相信上帝手中有一本包含世間所有精妙證明的天書。上帝相信這本書在 Gauss 手上。,高斯消元法,高斯消元法是求解線性方程組的經(jīng)典方法之一。高斯消元法由兩個部分構(gòu)成: 消元階段和回代階段 。消元階段的目的是將一般線性方程組轉(zhuǎn)化成容易求解的上三角線性方程組。,例如,如何將線性方程組化為更容易求解的線性方程組且不改變線性方程組的解？,交換兩個方程; 某一方程乘以非零常數(shù)c; 把某一方程的c倍加于另一方程。,消元階段:,經(jīng)過第一輪消元:,注釋: 8是主元(pivot), -0.5和0.5是乘子(multi

4、plier)。,經(jīng)過第二輪消元:,注釋: 8是主元(pivot), -0.5是乘子(multiplier)。,更加簡潔的表示:,Guass消元法和LU分解的關(guān)系？,回代階段:,從最后一個方程開始, 從后往前每次求解一個方程, 則三角方程組的解可以表示如下,代碼:,%elimination phase for k = 1:n-1 %pivot row for i= k+1:n %rows to be transformed if A(i,k) = 0 %multiplier lambda = A(i,k)/A(k,k); A(i,k+1:n) = A(i,k+1:n) - lambda*A(k,

5、k+1:n); b(i)= b(i) - lambda*b(k); end end end,%back-substitution phase for k = n:-1:1 b(k) = (b(k) - A(k,k+1:n)*x(k+1:n)/A(k,k); end,消元階段: 除法次數(shù)(n-1)+ (n-2)+1= n(n-1)/2, 乘法次數(shù)(n-1)*n+ (n-2)* (n-1) +1*2= (n3-n)/3, 加(減)法次數(shù)(n-1)*n+ (n-2)* (n-1) +1*2= (n3-n)/3。 共計： 2n3/3+n2/2-7n/6,回代階段: 除法次數(shù) n次, 乘法次數(shù) n(n-

6、1)/2次, 加(減)法次數(shù) n(n-1)/2次 共計： n2,高斯消元法由不相等的兩部分組成:工作量相對多的消元階段和工作量相對少的回代階段。對于大的n, 冪次較低的n相比而言可以忽略不計,最高次冪決定了當(dāng)n趨近于無窮時的極限形態(tài)。換而言之, 對大的n, 低階項對算法的執(zhí)行時間的估計沒有太大的影響, 當(dāng)僅需近似估計執(zhí)行時間時可以忽略不計。,我們經(jīng)常使用記號O表示算法是“是多少階的”。因此高斯消元法為O(n3/3)的算法。復(fù)雜度(complexity)是衡量算法性能優(yōu)劣的主要指標(biāo)。,例1. 在一臺特定的計算機上, 估計用高斯消元法求解5000*5000的線性方程組所需時間。,高斯消元法求解50

7、*50的線性方程組所需時間為t1。,回顧:,高斯消元法是時間復(fù)雜度為O(n3)的算法, 最終在有限步之內(nèi)給出一個解。因此高斯消元法叫做求解線性方程組的直接法。,理論上直接法在有限步之內(nèi)給出精確解(當(dāng)然在一有限精度的計算機上執(zhí)行時得到的結(jié)果將只是近似的)。,O(n-1)n!) (Cramer),O(n3) (Gauss),O(n2.373),O(n2+eps)?,Ref: The smart moneys on numerical analysis, SIAM News, 2012,O(n2.81) (Strassen),My own faith has been strong for years. Back in 1985, I made a $100 bet with Peter Alfeld that a fast matrix inverse would be found by 1995. It wasnt, so I paid him $100 and renewed the bet for another decade. It