第七章常微分方程數(shù)值解法

一、

基本思想7.3龍格-庫塔(Runge-Kutta)方法1.從積分的角度考慮.2.從微分的角度考慮.應用微分中值定理，有，利用所給微分方程，有其中K*=f(ξ,

y(ξ))

稱為區(qū)間[xn,xn+1]上的平均斜率。(7)對差商這樣，只要對平均斜率K*提供一種算法，由式按照這種觀點考察Euler格式，作為平均斜率K*的近似，精度自然很低。便相應地導出一種計算格式。它簡單地取點xn的斜率可將改進的歐拉格式改寫成的算術平均值作為平均斜率。該公式可以看作是用xn和xn+1兩個點處的斜率和由改進型歐拉公式我們可以猜想，如果在[xn,xn+1]內多預測幾個點的斜率，再對他們進行加權平均，可能得到精度更好的平均斜率！再考察改進的Euler格式。二、下面以2階龍格-庫塔方法為例來闡述這種思想考察區(qū)間[xn,xn+1]

上的一點xn+α

=xn+αh,0<α≤1,用xn和xn+α

的斜率K1和K2的加權平均作為平均斜率K*的近似值：即取其中和是待定常數(shù)。若取k1=f(xn,yn)，則問題在于如何確定xn+α

處的斜率k2和常數(shù)λ1和λ2。K*=λ1K1+λ2K2yn+1=yn+h(λ1K1+λ2K2)仿照改進的歐拉方法，用歐拉方法預測y(xn+α)的值，并用它來估計斜率K2：于是得到如下形式的算法：通過適當選取參數(shù)λ1,

λ2和α的值，使得公式具有2階精度??！yn+α

=yn+αhK1K2=f(xn+α,yn+α

)(8)式(8)是為二級二階Runge-Kutta方法的一般形式由泰勒公式展開，要使公式具有2階精度，只需方程組有無窮多解：二級二階R-K方法有無窮多種取?。ㄖ悬c法）（改進的Euler法）三、三級三階R-K方法形式：（9）我們要確定參數(shù)使得（9）式的局部截斷誤差為O(h4).將K2,K3按二元函數(shù)泰勒公式展開，代入(9)中第一式，再將y(xn+1)在xn處作泰勒展開，最后通過比較兩個展開式的系數(shù)，使局部截斷誤差為O(h4).得到如下方程組：此為5個方程，但含有7個未知量的方程組，它有無窮多個解。

庫塔三階方法四、四級四階R-K方法：局部截斷誤差經典龍格-庫塔方法解：例2：用經典的龍格-庫塔方法求解下列初值問題。經典的四階龍格-庫塔公式：0.10.20.30.40.51.09541.18321.26491.34161.41420.60.70.80.91.01.48321.54921.61251.6733

1.7321同保留5位的精確值完全一致：0.10.20.30.40.51.09541.18321.26491.34161.41420.60.70.80.91.01.48321.54921.61251.6733

1.7321

0.11.1000

1.0959

1.0954

0.21.1918

1.1841

1.18320.31.2774

1.2662

1.26490.41.3582

1.3434

1.34160.51.4351

1.4164

1.41420.61.5090

1.4860

1.48320.71.5803

1.5525

1.54920.81.6498

1.6165

1.61250.91.7178

1.6782

1.6733

1.01.7848

1.7379

1.7321

7.4單步法的收斂性和穩(wěn)定性收斂性和穩(wěn)定性從不同角度描述了數(shù)值方法的可靠性.只有既收斂又穩(wěn)定的方法，才能提供比較可靠的計算結果。一、收斂性定義若一種數(shù)值方法對任意固定的xn=x0+nh,當h→0時，有yn→y(xn),則稱該數(shù)值方法是收斂的。例設考察Euler方法求解此問題的收斂性。該問題的精確解為Euler方法計算公式為從而有當h→0時,考察在固定點x*=xn=x0+nh近似值的變化。開始取步長h=h0,x*=x2=x0+2h再取步長h=h0∕2,x*=x4=x0+4h進一步取步長h=h0∕4,x*=x8=x0+8h,

……x0

x1

x2x0

x4

x0

x8

這樣，是固定的。當h→0

(即n→∞)時對于此例，Euler方法是收斂的。再考察隱式Euler方法。此問題的隱式Euler格式為這時有從而數(shù)值解這時當h→0時仍然有yn→y(xn).一種方法的收斂性不能只針對某一個具體實例的成立與否，而是必須考慮對所有初值問題都收斂。下面討論一般單步法的收斂性。首先，關于求解初值問題的所有顯式單步方法都可以寫成如下形式其中h為步長，稱為增量函數(shù),它依賴于方程右端函數(shù)f,且不同的單步方法，φ的形式不同。(1)(2)Euler方法：改進的Euler方法：單步法式（2）在節(jié)點xn+1的局部截斷誤差定義為注意，局部截斷誤差en+1是在yn=y(xn)精確成立的假定下，從xn出發(fā)，用式（2）向前推進一步到xn+1而產生的誤差，并不是在點xn+1處，由于使用式（2）求解所引進的全部誤差。(3)事實上，用式（2）計算，每推進一步都會產生誤差，因此近似值yn包含有前n步的累計誤差。把從x0處的初值y0出發(fā)，用單步法（2）推進n+1步，到點xn+1所得到的近似值yn+1與其精確值y(xn+1)的偏差En+1=y(xn+1)-yn+1稱為整體截斷誤差。由上述定義可知，一種單步法是否收斂，就是看其整體截斷誤差En=

y(xn)-yn在h→0

(n→∞)時是否有En→0定義(4)為了估計單步法式（2）的整體截斷誤差En，我們需要下面的引理。引理設α>0,

β≥0為實數(shù)，且實數(shù)序列{ηn}滿足遞推關系則有證利用不等式eα>1+α,用數(shù)學歸納法（略）。定理若初值問題（1）的一個單步方法（2）的局部截斷誤差為O(hp+1)(p>=1)，且（2）式中的增量函數(shù)φ(x,y,h)滿足對y的Lipschitz條件，即存在L>0使得對一切有則單步法（2）的整體截斷誤差為En=O(h

p

).從而單步法收斂。(5)證依據局部截斷誤差的定義，有用上式減去（2）式，得根據定理中的假設條件，有再利用引理的結論可得到（C為正常數(shù)）假定初值是準確的，即E0=

y(a)-y0=0則有，從而有當h→0時，En→0.依據這個定理，判斷一個單步法（2）的收斂性，就歸結為驗證其增量函數(shù)φ(x,y,h)是否對y滿足Lipschitz條件。推論1Euler方法是收斂的。推論2

改進的Euler方法是收斂的。推論3各級Runge-Kutta方法是收斂的。在以上討論中，(2)是p階的方法，即en+1=O(hp+1),一般都是指p至少等于1，即p+1>=2.若將en+1按變量h在h=0作Taylor展開，得到很明顯，en+1=O(hp+1)，而p>=1的充要條件是二、相容性(6)定義若方法（2）增量函數(shù)滿足則稱單步法（2）與初值問題（1）相容。據以上討論，與初值問題(1)相容的方法至少是一階的。

對于與初值問題相容的單步方法(2)，若φ(x,y,h)

滿足對y的Lipschitz條件，則由前面的定理，此方法是收斂的。(7)若把（2）寫成(8)并令x

=

xn

(=

x0+nh)固定，h→0.因為`yn→y(xn),(8)式可寫成即(9)這意味著h→0時，數(shù)值計算格式(8)趨于微分方程(9).這就是“相容”一詞的意義。不難驗證，已討論過的各種方法都是與所對應的初值問題相容的三、穩(wěn)定性定義

設用某一數(shù)值方法計算yn時，所得到的實際計算結果為?n,且由誤差δn=yn-?n引起以后各點處ym(m>n)的誤差為δm,如果總有|δm|≤|δn|，則稱該數(shù)值方法是絕對穩(wěn)定的。對于微分方程若Reλ>0，它是不穩(wěn)定的（這里指微分方程本身的性質）若Reλ<0，它是穩(wěn)定的.但用不同的數(shù)值方法來求解時，還有數(shù)值穩(wěn)定與數(shù)值不穩(wěn)定的兩種可能。所以我們考慮Reλ<0的情形，對于實數(shù)λ，考慮λ<0。當λ為實數(shù)時，Euler方法穩(wěn)定條件：-2<λh<0改進的Euler方法