《數值分析》筆記

使用教材:《數值分析》李慶揚,王能超,易大義 編。
大致梳理整本書的內容和重點。

第1章 數值分析與科學計算引論

本章是介紹性質的,主要說了數值分析(乃至計算數學)在學什么,解決什么問題。
內容上,需要掌握誤差的幾個種類,絕對/相對誤差(限)的計算,有效數字。
1.3節(jié)“誤差定性分析與避免誤差危害”:了解數值穩(wěn)定、病態(tài)問題、條件數的概念。尤其是條件數的意思,在后續(xù)數值線性代數中也會用到。避免誤差危害的方法:避免兩相近數相減,避免用絕對值小的數做除數,注意運算次序和減少運算次數。
了解一些古典的數值計算思想:秦九韶算法,迭代法,以直代曲,加權平均的松弛技術(感覺類似外推法?)
本章內容淺顯零散,了解即可。

第2章 插值法

注意區(qū)分插值擬合:插值函數一定確保已有的數據點完全擬合,一般維數(比如多項式次數)較高;擬合則不需要,一般維數較低。

2.1 引言

差值問題的提出:略
多項式擬合:對于n+1個數據點(x,y),用多項式擬合,則使用n次多項式剛好保證存在唯一。(代數方法證明)
雖然結果唯一,但插值多項式系數的求法有許多種,這里介紹拉格朗日插值牛頓插值

2.2 拉格朗日插值

思路:基函數l_k(x)(線性插值基函數):滿足l_j(x_j)=1,l_j(x_i)=0,i\neq j的n+1次多項式(顯然唯一)。插值多項式L(x)=\sum y_kl_k(x)為基函數的線性組合,系數由數據點確定。
基函數方法是將插值問題劃歸于特定條件下容易實現的插值問題,本質上是廣義的坐標系方法。
插值余項:若f^{(n+1)}(x)存在,則R_n(x)=f(x)-L(x)=\frac{f^{(x+1)(\xi)}}{(n+1)!}\prod_{i=0}^n (x-x_i)。
使用中值定理證明。(中值定理和泰勒展開在數值分析中用處很多)。由插值余項的公式得,要保證多項式插值與原函數較為接近,最好原函數的n+1次導數有界。
拉格朗日插值法思路清晰,公式漂亮,但若數據點變化需要重新計算。下面的牛頓插值多項式解決了這個問題。

2.3 均差與牛頓插值多項式

引入“均差”(差商)概念,類比導數。
牛頓法思路:利用均差,依次添加數據點,每添加一個數據點插值多項式P(x)會增加一項。
對數據點的順序沒有要求。
牛頓插值法在添加數據點時計算量較小,但公式沒那么漂亮。

2.4 埃爾米特插值

拉格朗日插值和牛頓插值僅保證數據點f(x_i)=y_i,沒有導數信息。埃爾米特插值用于一些f^{(k)}(x_i)=y_i^k的需求。
計算方法類似牛頓插值,必要時可用待定系數法。
書中僅討論了兩個典型的埃爾米特插值的計算。
注意兩點三次埃爾米特插值三次樣條插值的區(qū)別:兩者都是用三次多項式對兩點做插值,但前者是為了保證f(x_1),f(x_2),f'(x_1),f'(x_2)均與數據相同,后者僅保證f(x_1),f(x_2)與數據相同,導數的條件是為了保證f的光滑性。

2.5 分段低次插值

了解高次插值的病態(tài)性質:龍格現象,高次插值沒有收斂性。
通過分段,降低次數,減少病態(tài)現象。
常用的:分段線性插值,分段三次埃爾米特插值,三次樣條插值。
分段插值具有一致收斂性。

2.6 三次樣條插值

保證數據點的同時插值函數具有二階光滑性,在很多實際問題中有應用。
思路:對于n+1個數據點,分成n段三次多項式,共需4n個參數。限制條件:每個數據點吻合,2n個方程;為保證光滑性,每個連接點兩邊的一階導數、二階導數相同,2(n-1)個方程;還差2個方程,一般為兩個端點的條件(邊界條件),可根據實際問題要求給定。
到此為止已經列出了一個線性方程組,解法看書即可,最終化成三對角矩陣的線性方程組,使用追趕法求解。
性質:插值函數及其一階導數、二階導數均一致收斂與數據的原函數。

第3章 函數逼近與快速傅里葉變換

3.1 函數逼近的基本概念

函數逼近問題就是:在區(qū)間[a,b]上用簡單函數逼近已知復雜函數的問題,簡單函數一般是n次多項式、有理函數、分段低次多項式等。(可以理解為,從任意函數構成的空間到多項式函數空間的投影?)
維爾斯特拉斯定理:閉區(qū)間上連續(xù)函數,均可用有限維的多項式任意逼近(一致范數)。證明方法用到伯恩斯坦多項式,這是一個理論上對[0,1]上連續(xù)函數f任意逼近的多項式,但實際上收斂速度慢。

介紹范數與賦范線性空間,用來描述逼近的程度。引入內積、正交、施瓦茨不等式、權函數、函數內積、函數范數等概念。

最佳逼近:對于某個函數空間H,某個范數,找出H中與目標函數f的差的范數最小的函數P*,為最佳逼近函數。
||f(x)-P*(x)||=min_{P\in H_n}||f(x)-P(x)||

若取P為多項式,范數為一致范數,則為最佳一致逼近多項式;若范數取二范數,則為最佳平方逼近多項式;若f(x)為列表函數(離散情況),范數取二范數,則為最小二乘擬合。

3.2 正交多項式

正交多項式是函數逼近的重要工具,在數值積分中也有重要應用。
正交多項式可以理解為多項式內積空間中的正交基?

概念:帶權\rho (x)正交,正交函數族,標準正交函數族,施密特正交化,正交多項式的遞推關系式,正交多項式零點的定理。

勒讓德多項式P_n(x):區(qū)間為[-1,1],權函數\rho \equiv1時,由\{1,x, \cdots ,x^n,\cdots \}\正交化得到的多項式。
性質:正交性,奇偶性,遞推關系,n個不同零點。

切比雪夫多項式T_n(x):區(qū)間為[-1,1],權函數\rho = \frac{1}{\sqrt{1-x^2}}時,由\{1,x, \cdots ,x^n,\cdots \}\正交化得到的多項式。另一個比較直觀的表示:若令x=cos\theta,則T_n(x)=cos(n\theta),0 \le \theta\ \le \pi
性質:遞推關系,(帶權)正交性,奇偶性,n個零點,首項系數,最佳一致逼近多項式,切比雪夫多項式零點插值(通過改變插值點的選位,從等距點改為切比雪夫多項式零點,可避免龍格現象)

一些其他有用的正交多項式:第二類切比雪夫多項式,蓋拉爾多項式,埃爾米特多項式。

正交函數族的優(yōu)勢在下面一節(jié)中就有體現。

3.3 最佳平方逼近

數學表達式:
||f(x)-S*(x)||_2^2=min_{S(x)\in \varphi} ||f(x)-S(x)||_2^2
連續(xù)形式:min_{S(x)\in \varphi} \int_a^b \rho(x)[f(x)-S(x)]^2dx
對多項式空間,等價于求多元函數的最小值:I(a_0,a_1,\cdots,a_n) = \int_a^b \rho(x) \left [ \sum_{j=0}^n a_j \varphi_j(x) - f(x) \right ]^2 dx
這是一個分析性質比較好的函數,求最小值只要求某個點對所有變量的偏導數等于0,剛好是一個線性方程組,稱為法方程。這樣看起來是找到了解析解,只要用數值方法求解線性方程組就可以,但直接取\varphi_k(x)=x^k時,系數矩陣為希爾伯特矩陣,是一個條件數大的矩陣,對它解線性方程組得到的誤差較大。

為了解決上面的問題,可以用正交函數族作最佳平方逼近,參考傅里葉級數,由正交性質,系數矩陣是對角陣,因此在計算傅里葉級數的系數時不需要解線性方程組。在上面的例子中將x^k換成勒讓德多項式P_k(x)就能解決問題。這就體現出正交函數族對于數值逼近的優(yōu)勢。

切比雪夫級數:將函數f按切比雪夫多項式展開,得到切比雪夫級數,可作為f[-1,1]上的近似最佳一致逼近多項式。

3.4 曲線擬合的最小二乘法

類似于最佳平方逼近的離散形式,在科學實驗中對實驗數據的處理(曲線擬合)經常用到。
數學表達式:min_{S(x)\in \varphi } \sum_{i=0}^m [S(x_i) - y_i]^2
方法與最佳平方逼近類似,法方程,解線性方程組,用正交多項式解決病態(tài)系數矩陣問題。略。

3.5 有理逼近

前面的多項式逼近對于有界連續(xù)函數逼近效果較好,現在對于在某點附近無界的情況,用有理函數逼近效果較好。
有理函數:R_{nm}(x)=\frac{P_n(x)}{Q_m(x)} = \frac{\sum_{k=0}^n a_kx^k}{\sum_{k=0}^m b_kx^k},即兩個多項式的比。這樣可以得到比如\frac{1}{x}這樣的可以趨于無窮的函數逼近。
計算方法是利用連分式。

帕德逼近:如果有理函數R_{nm}(x)逼近f(x)R_{nm}^{(k)}(0)=f^{(k)}(0),k=0,1,\cdots,(n+m),則稱R_{nm}(x)f(x)x=0處的(n,m)階帕德逼近。這個逼近的好處個人沒用過不太清楚,略。

3.6 三角多項式逼近與快速傅里葉變換

在模型數據具有周期性時,用三角函數特別是正弦函數和余弦函數作為基函數是合適的。
最佳平方三角逼近:傅里葉級數
快速傅里葉變換:略。

評注

正交多項式中的勒讓德多項式和切比雪夫多項式是兩個十分重要且經常使用的正交多項式,應引起高度關注。
當模型為多項式時法方程是病態(tài)的,為此推薦用正交化方法避免解法方程。
如果數據是周期的,使用三角最小二乘或三角插值是合適的,計算用快速傅里葉變換(FFT),它是節(jié)省計算量的一個范例。

第4章 數值積分與數值微分

4.1 數值積分概論

一般的數值積分:\int_a^b f(x) dx \approx \sum_{k=0}^n A_k f(x_k)
其中x_k為求積節(jié)點,A_k為求積系數,即節(jié)點x_k的權重。
代數精度:描述求積公式的誤差的階數,計算方法為計算f(x) = x^k的情況的誤差。
插值型求積公式:先求f(x)的插值多項式L_n(x),再用I_n := \int_a^b L_n(x) dx作為I := \int_a^b f(x) dx的近似值,其誤差可用插值多項式的誤差計算:R[f] = \int_a^b [ f(x) - L_n(x) ] dx = \int_a^b R_n(x) dx。
求積公式是插值型的,當且僅當其至少有n次代數精度。

收斂性:只要插值點充分密集,數值積分就會充分接近真實值。
穩(wěn)定性:插值點數據(x_k,f(x_k))可能存在誤差,穩(wěn)定性是指,只要誤差充分小,數值積分就會充分接近真實值。反之,稱為對誤差敏感。

4.2 牛頓-柯特斯公式

先默認選取等距節(jié)點,因為簡單,此時的公式I_n = (b-a) \sum_{k=0}^n C_k^{(n)} f(x_k)稱為牛頓-柯特斯公式,C_k^{(n)}稱為柯特斯系數。

對于[a,b],選取a、b兩點插值(n=1),稱為梯形公式;選取a、b、(a+b)/2三點(n=3),稱為辛普森公式,5個點(n=4)稱為柯特斯公式。9點以上的牛頓柯特斯公式穩(wěn)定性差,不用。

當n為偶數,牛頓柯特斯公式至少具有n+1階代數精度。

4.3 復合求積公式

類似分段低階插值,先分段再用低階求積公式,稱為復合求積公式。常用的,復合梯形公式(折線段),復合辛普森公式。

4.4 龍貝格求積公式

思路:實際計算時若精度不夠可將步長逐次分半,直到精度足夠為止。

理查森外推法:只要真值與近似值的誤差能表示成h的冪級數,就能使用外推法提高精度。
例如,T(h) = I + a_1 h^2 + a_2 h^4 + \cdots,則可增加一些點將h減半得到T(h/2) = I + \frac{1}{4}a_1 h^2 + \frac{1}{16}a_2 h^4 + \cdots,通過線性組合得到S(h) = \frac{4}{3} T(h/2) - \frac{1}{3} T(h) = I + b_1 h^4 + b_2 h^6 + \cdots,誤差階數就從2階提高到了4階!

龍貝格算法:從復合梯形公式開始,不斷使用外推法,(計算有遞推公式),直到精度足夠為止。

4.5 自適應積分方法

回顧以上符合積分公式,對整個積分區(qū)域使用同樣的階數,適用于被積函數變化不大的積分。
對于在不同區(qū)域變化程度不同的函數,我們可以對不同區(qū)間使用不同的插值點數量(根據誤差自動判斷,龍貝格算法),平坦的部分步長大,變化劇烈的部分步長小,既減少運算量又減少誤差。具體過程略。

4.6 高斯求積公式

以上我們都默認等距節(jié)點,現在我們改變節(jié)點的位置,希望在不改變節(jié)點數量的情況下得到更高階的誤差,這就是高斯求積公式。
如果求積公式(n+1個節(jié)點)具有2n+1次代數精度(n+1個節(jié)點的理論最高精度),則稱其節(jié)點為高斯點,公式為高斯型求積公式。

如何取高斯點?
由定理證明。。。高斯點是正交多項式的零點。
性質:高斯求積公式收斂,穩(wěn)定。

高斯-勒讓德求積公式:使用勒讓德多項式作為高斯求積公式需要的正交多項式。
對于區(qū)間不是[-1,1]的情況,做變換x=\frac{b-a}{2}t + \frac{a+b}{2}即可。

高斯-切比雪夫求積公式:與上面類似,應對不同情況。
還有其他正交多項式。

注:高斯求積公式不一定要用到區(qū)間的端點,對一些在端點出趨于無窮的函數積分(反常積分)近似較好。

4.7 多重積分

寫成累次積分,化成一重積分的情況。略。

4.8 數值微分

思路1:用差商近似導數,比如中點公式f'(a) \approx \frac{f(a+h) - f(a-h)}{2h},誤差用泰勒展開計算。

思路2:插值型求導公式:先求f的插值函數P_n(x),再用P'_n(x)作為導數的近似值。

思路3:三次樣條求導。利用三次樣條函數導數值也收斂的性質,可作為數值微分的方法。

思路4:數值微分的外推算法,對思路1的改進,略。

第5章 解線性方程組的直接方法

5.1 引言與預備知識

線性方程組的系數矩陣可大致分為:低階稠密矩陣,大型稀疏矩陣(零元素多)。
方法大致分兩種:直接法(對低階稠密矩陣有效),迭代法(對大型稀疏矩陣有效)。(本章只講前者)
特殊形式的矩陣(三角,三對角等)有特殊方法。

概念:矩陣特征值,譜半徑,特殊矩陣,若當型矩陣

5.2 高斯消去法

分為前代法和回代法兩步,分別求解三角矩陣方程,LU分解(A=LU)。具體可參考《數值線性代數》。

約化的主元素都不為零(普通高斯消去法能進行下去)的充要條件是,矩陣A的順序主子式都不為零。

對于約會的主元有零的情況:列主元消去法,每次循環(huán)時增加一步行變換,可避免主元是零或接近零,減少誤差。只要矩陣可逆就能計算下去,實用。列主元的LU分解:PA=LU

5.3 矩陣三角分解法

普通的LU分解和列主元LU分解:上面提到了,具體步驟略。

平方根法:對于對稱正定矩陣的專用方法。A=LDL^T=L_1L_1^T楚列斯基分解)。具體步驟略。

追趕法:對于三對角矩陣的專用方法。具體步驟略。

特殊矩陣的專用方法比平凡的方法效果更好(時間、空間復雜度低),但適用范圍小。

5.4 向量和矩陣范數

度量,為矩陣/向量的誤差分析做鋪墊。

矩陣的算子范數/誘導范數:由向量范數擴展到矩陣,具有相容性。也有不是算子范數的矩陣范數,比如F范數(弗羅貝尼烏斯范數)。

有一些定理需要了解。

5.5 誤差分析

病態(tài)方程組、病態(tài)矩陣:對矩陣A或常數項b的變化敏感。條件數刻畫了病態(tài)程度。希爾伯特矩陣是病態(tài)矩陣。若條件數計算不方便,也有幾條判斷病態(tài)的簡單方法。

面對病態(tài)問題,首先是列主元方法,若解決不了還可以采用高精度算術運算或采用預處理方法

第6章 解線性方程組的迭代法

6.1 迭代法的基本概念

迭代法的格式:x^{(k+1)} = Bx^{(k)} + f。
先(任意)取一個初值x^{(0)},再有Ax=b的A和b確定迭代格式,不斷迭代,直到達到一定次數或誤差達到要求為止。
如果\lim _{k \rightarrow \infty } x^{(k)}存在,則稱迭代法收斂。

定義向量序列的極限和矩陣序列的極限。
\lim _{k \rightarrow \infty }B^k=0等價于\rho(B)<1等價于B的某個范數||B||<1。

對于任意初值x^{(0)},迭代法x^{(k+1)} = Bx^{(k)} + f收斂的充要條件是\lim _{k \rightarrow \infty }B^k=0

迭代法收斂速度的度量:平均收斂速度,漸進收斂速度。

6.2 雅克比迭代法與高斯-賽德爾迭代法

對于Ax=b,令A=D-L-U,D為對角矩陣,L為下三交矩陣,U為上三角矩陣。
則有Dx = (L+U)x+b,即x = D^{-1}(L+U)x+D^{-1}b,此為雅克比迭代法的迭代格式。
用方程組表示,就是給定x^{(k)},對每個元素x^{(k+1)}(i)迭代時都使用x^{(k)}的元素,全部算完再更新為x^{(k+1)}

類似的,有(D-L)x = Ux+b,即x = (D-L)^{-1}Ux+(D-L)^{-1}b,此為高斯-賽德爾迭代法的迭代格式。
用方程組表示,就是計算每個元素x^{(k+1)}(i)盡可能使用最新的結果,即x^{k+1}(j)x^{k}(l),\forall j<i, \forall l >i。

收斂性:(可用迭代法收斂條件判斷)
概念:對角占優(yōu)矩陣,嚴格對角占優(yōu)矩陣,可約矩陣,不可約矩陣。
對角占優(yōu)定理:如果A為嚴格對角矩陣或不可約對角占優(yōu)矩陣,則A為非奇異矩陣。
如果A為嚴格對角矩陣或不可約對角占優(yōu)矩陣,雅克比迭代法和高斯賽德爾迭代法均收斂。

6.3 超松弛迭代法

高斯賽德爾迭代法的改進,SOR超松弛迭代法。略。

6.4 共軛梯度法

很重要??!需要仔細看書,在此就不細說了,只寫需要掌握的部分。

對于對稱正定矩陣A,考慮二次函數\varphi:R^n \rightarrow R
\varphi(x) = \frac{1}{2}(Ax,x) - (b,x) = \frac{1}{2}\sum_{i=1}^n \sum_{j=1}^n a_{ij}x_ix_j - \sum_{j=1}^n b_j x_j
\varphi(x)是可導的凸函數,其去最值的時候等價于梯度等于零,即\bigtriangledown \varphi(x) = Ax-b=0
因此求解線性方程組的問題等價于求二次函數最小值的問題,求解方法是構造一個向量序列\{ x^{(k)} \}\使二次函數趨向于最小值。

最速下降法:每次迭代方向為\varphi(x)x^{(k)}處的梯度方向,確定方向后取**該方向上\varphi(x)的最小值點為迭代結果x^{(k+1)}。即以梯度方向為迭代方向,再用線搜索確定步長。
計算:梯度:先求出梯度函數即可。線搜索:令導數等于零可求得步長。
性質:兩個相鄰的搜索方向正交,因此快接近最優(yōu)值時會出現鋸齒狀的路線,速度變慢。收斂速度與矩陣A的最大最小特征值的差有關,特征值越接近收斂速度越快。
速度并不快的方法,但理論簡單,后續(xù)的方法都是在其上的改進。

共軛梯度法(CG方法):求解大型稀疏對稱正定方程組十分有效的方法。
思路:搜索方向從單純的梯度方向改進為,x^{(k)}處梯度方向和上一步的收縮方向組成的平面,這樣需要兩個參數來確定步長。
性質:誤差r^{(k)}:=Ax^{(k)}-b序列是正交向量組,搜索方向p^{(k)}序列是A-共軛向量組。
r^{(k)}正交性,對n維線性方程組,理論上最多n步便可得到精確解,可以說CG方法是直接法;但由于舍入誤差,很難保證r^{(k)}的正交性,一般還是以誤差到達精度要求時終止,也可以說CG方法是迭代法。

共軛梯度法的細節(jié)可參考《數值線性代數》,它還有一些其他的性質。
目前討論的共軛梯度法是關于對稱正定矩陣A的,對于一般矩陣也可以稍加修改后使用,具體參考《最優(yōu)化理論與方法》。

第7章 非線性方程與方程組的數值解法

7.1 方程求根與二分法

二分法:最樸素的方法,在函數連續(xù)且端點函數值異號的情況下,不斷壓縮區(qū)間長度直到達到精度要求為止。

7.2 不動點迭代法及其收斂性。

x = \varphi(x),則稱x\varphi(x)的不動點。若迭代法收斂,則最終收斂到其不動點。

不動點的存在性與迭代法收斂性:不動點存在唯一性定理
\varphi(x)[a,b]連續(xù),且a \le \varphi(x) \le b,且存在正常數L<1使\forall x,y \in [a,b]. |\varphi(x) - \varphi(y)| \le L |x-y|,則\varphi(x)[a,b]上存在唯一的不動點x^*

局部收斂性:若varphi(x)有不動點x^*,如果存在不動點附近的鄰域,使得當初值在此鄰域中時迭代法一定收斂到不動點,則稱迭代法局部收斂。
局部收斂的充分條件:\varphi'(x)在不動點附近連續(xù)且絕對值小于1,則局部收斂。

收斂速度:迭代誤差序列(一個無窮小量)的階數,p階收斂。

7.3 迭代收斂的加速方法

已有迭代格式\varphi,通過一些變形得到新的迭代格式\phi[\varphi],使新的迭代格式收斂速度(階數)加快。
埃特金加速收斂方法,斯特芬森迭代法。

7.4 牛頓法

思路:利用導數(梯度)信息,對原函數線性化(泰勒展開)來近似,再求零點。
x_{k+1} = x_k - \frac{f(x_k}{f'(x_k)}
有幾何解釋:切線法。
至少是二階收斂,優(yōu)點是收斂速度快,缺點是每一步計算量稍大(要計算導數值),而且對初值有要求(不是全局收斂,需要初值在不動點附近)
可以擴展到高維情況,為了克服缺點有一些變種(簡化牛頓法,牛頓下山法),重根情形可通過變形加速。

7.5 弦截法與拋物線法

以上是用前一個點迭代后一個點,稱為單點迭代法。
弦截法與拋物線法使用前幾個點計算后一個點,稱為多點迭代法,從而可以不使用導數信息且保持高階收斂速度。

7.6 求根問題的敏感性與多項式的零點

知道敏感性的意思即可。

求多項式的全部零點:以上方法只能求一個零點??梢悦壳笠粋€零點,利用因式分解給原方程降次。推薦使用拋物線法,可求復根。

7.7 非線性方程組的數值解法

在一維的基礎上向高維推廣。使用向量范數代替誤差,使用雅可比矩陣代替一維導數,多變量方程的不動點迭代法,不動點存在唯一性定理(壓縮映射原理),p階收斂,牛頓迭代法等。

第8章 矩陣特征值計算

8.1 特征值性質和估計

特征值的基本性質,略。

格什戈林圓盤定理:設A = (a_{ij})_{n*n},則A的每一個特征值必屬于下述某個圓盤之中:
| \lambda - a_{ii} | \leq r_i = \sum_{j=1,j\not= i}^n | a_{ij}|, \;\;\;\; i=1,2,\cdots,n
或者說,A的特征值都在復平面上n個圓盤的并集中。
如果Am個圓盤組成一個連通的并集S,且S與余下的n-m個圓盤是分離的,則S內恰包含Am個特征值。

圓盤定理可初步估計特征值的范圍。

8.2 冪法及反冪法

冪法是一種計算矩陣主特征值(模最大的特征值)及其特征向量的迭代方法,適用于大型稀疏矩陣。
類似的,反冪法計算模最小的特征值及其特征向量,用于計算海森伯格矩陣或三對角陣的對于一個給定近似特征值的特征向量。

冪法:原理和編程都很簡單,看書即可。前提:模最大的特征值只有一個(不是重特征值)。

加速方法:
1.原點平移法,引進矩陣B=A-pI,B是A的平移,特征值也對應平移且特征向量不變,可以改變模最大的特征值、收斂速度。缺點:p的確定困難。
2.瑞利商加速。

反冪法:計算A的模最小特征值,等價于計算A^{-1}的模最大特征值的倒數。(實際上不需要計算A^{-1},只有解線性方程組即可)
若已知某個特征值的近似值p(由其他方法求出),可以對矩陣平移B=A-pI,則B就有一個接近0的特征值,可使用反冪法求這個值和特征向量,也就是求得原矩陣A的更精確的特征值和對應的特征向量。

8.3 正交變換與矩陣分解

正交變換是計算矩陣特征值的有力工具,一般可以先通過正交變換將矩陣變?yōu)樯先蔷仃嚮蚪粕先堑木仃?,再用迭代法快速求出所有特征值(QR方法)。
本節(jié)介紹豪斯霍爾德(Householder)變換吉文斯(Givens)變換

豪斯霍爾德變換:思路:使用反射,將x反射到e_1。計算量稍大,但一次改變向量的所有分量。

吉文斯變換:思路:使用旋轉,將(\cdots,x_i,\cdots,x_j,\cdots)旋轉到(\cdots,1,\cdots,0,\cdots)。計算量較小,但一次只改變兩個分量。

矩陣的QR分解和舒爾分解
QR分解:設n階方陣A非奇異,則可以通過一系列變換(比如n-1次豪斯霍爾德變換)將其變?yōu)樯先顷?,?img class="math-inline" src="https://math.jianshu.com/math?formula=PA%3DR" alt="PA=R" mathimg="1">,其中P為正交矩陣,R為上三角矩陣。因此令Q=P^{-1},就有A=QR的分解。我們只需再求上三角矩陣R的特征值,就可以得到A的特征值。
性質:對n階非奇異方陣A,QR分解存在且當R對角元為正時唯一。

實舒爾分解:解決實矩陣的復特征值問題,實矩陣可以約化為上海森伯格矩陣(上三角和對角線元素的下面一個元素可以非零,其余是零)。
Q^TAQ = R
R為類上三角矩陣,其對角塊R_{ii}為一階或二階方陣。一階對應實特征值,二階對應兩個共軛復特征值。計算Q的方法也是使用正交變換。

8.4 QR方法

QR方法是一種變換方法,是計算一般矩陣(中小型矩陣)全部特征值問題的最有效按方法之一。
對一般矩陣,先用豪斯霍爾德方法將A化為上海森伯格矩陣B,再用QR方法計算B的全部特征值(,再使用反冪法提高精度并求出對應的特征向量)。
迭代過程:對A_k做QR分解:A_k=Q_kR_k,再令A_{k+1}=R_kQ_K = Q_k^T R_kQ_k。

具體細節(jié),比如迭代法何時終止的判定,略。
還有在其基礎上的改進,比如帶原點位移的QR方法,隱式QR方法,比較復雜,略。

第9章 常微分方程初值問題數值方法

9.1 引言

常微分方程初值問題:
y' = f(x,y),\;\;\;\; x\in [x_0,b] \\ y(x_0) = y_0
稱方程的解y(x)為積分曲線。

李普希茲條件:如果存在實數L>0,使得
| f(x,y_1) - f(x,y_2) | \le L | y_1 -y_2|, \;\;\;\; \forall y_1,y_2 \in R
則稱f關于y滿足李普希茲條件。
這個條件用于常微分方程解的存在唯一性,具體參考《常微分方程》。在本章中也用于數值方法的收斂性證明。

對于常微分方程初值問題的數值方法,就是尋求一系列離散節(jié)點x_1<x_2<\cdots<x_n<\cdots上的近似值y_1<y_2<\cdots<y_n<\cdots
為了簡便,我們取x_i = x_0 + i*h的等距節(jié)點,h稱為步長。

一般使用前面的點(已知的或已求出近似值的)來計算后面的點。使用前一個點的叫做單步法,使用前n個點的叫做多步法,這和7.5節(jié)中的多點迭代法類似。
本章我們需要研究的問題有:計算方法(公式),公式的局部階段誤差和階,誤差估計及收斂性,穩(wěn)定性。

9.2 簡單的數值方法

理論上有y_{n+1} = y_n + \int_{x_n}^{x_{n+1}} y'(x)dx = y_n + \int_{x_n}^{x_{n+1}} f(x,y(x))dx

歐拉法y_{n+1} = y_n + hf(x_n,y_n),即使用左矩形公式計算數值積分,\int_{x_n}^{x_{n+1}} f(x,y(x))dx \approx hf(x_n,y_n)。
這是最簡單的思路最清晰的方法,計算量小,當然誤差也比較大(截斷誤差相當于數值積分的誤差)。
這是一個不需要解方程就可以直接計算的方法,即y_{n+1} = y_n + hf(x_n,y_n)等號右邊沒有未知量,這樣的方法稱為顯式方法。

后退歐拉法(隱式歐拉法)y_{n+1} = y_n + hf(x_{n+1},y_{n+1}),思路與歐拉法相同,使用右矩形公式計算數值積分。此時等號右邊有未知量,需要通過解方程來求得y_{n+1},這樣的方法稱為隱式方法,類似于隱函數和顯函數的關系。
本方法的收斂性證明用到李普希茲條件
一般來說,隱式方法雖然計算量稍大,但穩(wěn)定性好。

梯形方法:類似的,可以使用梯形公式計算數值積分。

改進歐拉公式:每次遞推分成兩部:先用歐拉公式求得“預測值”,再用梯形公式求得“校正值”,將梯形公式的隱式方法改為了兩部顯式計算,減少了計算量。

局部截斷誤差:加上目前的數據是正確的,下一個節(jié)點的數值解和精確值的差稱為局部截斷誤差,以此來度量算法的收斂性和收斂速度。一般局部截斷誤差與步長h相關,即步長越小誤差越小,若與h^{p+1}線性相關則稱為p階精度。

9.3 龍格-庫塔方法

為了提高精度而增加了使用的點,增多的是[x_n,x_{n+1}]內的點,不是x_n之前的節(jié)點,因此還是單步法。類似改進的歐拉法,公式如下;
y_{n+1} = y_n + h\varphi(x_n,y_n,h) \\ \varphi(x_n,y_n,h) = \sum_{i=1}^r c_i K_i \\ K_1 = f(x_n,y_n) \\ K_i = f \left (x_n+\lambda_i h, y_n+h \sum_{j=1}^{i-1}\mu_{ij}K_j \right)
稱為r階顯式龍格庫塔法,簡稱RK法,c_i,\lambda_i,\mu_{ij}都是常參數。其中r=1時即為歐拉法,r=2時的一種為改進的歐拉法。

二階顯式R-K方法:以二階為例說明參數的選擇:希望確定這些參數使得p階數盡量高,可使用待定系數法,使用泰勒公式并令幾個低階的誤差項參數為零,從而得到幾個方程,解得參數。由于方程的解不一定唯一,參數不唯一。
三階、四階類似。四階顯式龍格庫塔方法比較常用。
由于龍格庫塔方法的推導基于泰勒展開,因而它要求解具有較好的光滑性。

變步長的龍格庫塔方法:略。

9.4 單步法的收斂性與穩(wěn)定性

收斂性:步長趨于0時誤差也趨于0,則稱為具有收斂性。

定理:若單步法具有p階精度且增量函數\varphi(x,y,h)關于y滿足李普希茲條件,初值也準確,則y(x_n) - y = O(h^p),整體截斷誤差p階收斂。
因此,一般通過判斷是否符合李普希茲條件來判斷是否收斂。

穩(wěn)定性:若計算出現的微小誤差不會惡性增長以至于“淹沒”真解,則稱具有穩(wěn)定性。
而穩(wěn)定性不僅與方法有關,有時也與h的取值范圍有關(在《偏微分方程數值解》中體現比較明顯,有的方法是全局收斂,有的方法是條件收斂,即必須步長/步長比較小時才收斂),收斂時h的取值范圍叫做絕對穩(wěn)定域/絕對穩(wěn)定區(qū)間

9.5 線性多步法

思想:充分利用前面多步的節(jié)點值信息來預測下一個節(jié)點的值,期望會獲得較高的精度。
方法:基于數值積分的方法或基于泰勒展開的方法。

線性多步法的一般公式:
y_{n+k} = \sum_{i=0}^{k-1} \alpha_i y_{n+i} +h\sum_{i=0}^k \beta_i f_{n+i}
類似龍格庫塔公式,可以使用待定系數法和泰勒展開確定參數。

阿當姆斯顯式與隱式公式:公式如下
y_{n+k} = y_{n+k-1} + h\sum_{i=0}^k\beta_if_{n+i}
其中\beta_k=0時為顯式方法,反之為隱式方法。也稱為阿當姆斯-巴什福斯(AB)公式阿當姆斯-蒙爾頓(AM)公式??赏ㄟ^泰勒展開并求解線性方程組確定其系數。

還有其他的線性多步法公式,比如米爾尼方法,辛普森方法,漢明方法等。

預測校正方法:類似改進的歐拉方法,分為預測(顯式公式)和校正(隱式公式)兩步的多步方法,比如阿當姆斯四階預測校正格式(PECE)。

一般先用四階龍格庫塔方法(單步法里最常用的,效果較好但計算量大)算出開始的幾個點,之后使用多步預測校正方法(計算量?。?/p>

9.6 線性多步法的收斂性與穩(wěn)定性

類似的,分為收斂性和穩(wěn)定性,略。

9.7 一階方程組與剛性方程組

略。

最后編輯于
?著作權歸作者所有,轉載或內容合作請聯系作者
【社區(qū)內容提示】社區(qū)部分內容疑似由AI輔助生成,瀏覽時請結合常識與多方信息審慎甄別。
平臺聲明:文章內容(如有圖片或視頻亦包括在內)由作者上傳并發(fā)布,文章內容僅代表作者本人觀點,簡書系信息發(fā)布平臺,僅提供信息存儲服務。

友情鏈接更多精彩內容