Runge-Kutta法

1、1 7.2 Runge-Kutta法法龍格龍格-庫塔庫塔(Runge-Kutta)方法簡稱方法簡稱R-K法，是一種應用較廣的法，是一種應用較廣的高精度的單步法。高精度的單步法。所謂單步法就是在計算所謂單步法就是在計算yi時只用到前一步信息時只用到前一步信息yi-1的方法。的方法。本節(jié)介紹本節(jié)介紹R-K法的構(gòu)造原理、常用公式。法的構(gòu)造原理、常用公式。2一、Runge-Kutta方法的構(gòu)造原理0)(),(yaybxayxfy對于常微分方程的初值問題的解),(xyy 1,iixx在區(qū)間上使用微分中值定理 有11( )()( )()iiiiy xy xyxx1(,)iixx1( )()( )iiy x

2、y xhy即-(3)1iihxx其中，為步長31( ),iiKyy xxx可以認為是在區(qū)間上的平均斜率引入記號( )Ky , ( )fy1( ),iiy xxxK只要使用適當?shù)姆椒ㄇ蟪鲈趨^(qū)間上平均斜率的近似值就可得到相應的Runge-Kutta方法1nxnxxy)(xyy K4二、低階Runge-Kutta方法1nxnxxy)(xyy 如下圖111.( )( ),iiiy xxy xxxK如果以在處的斜率作為在上的平均斜率即1()iKy x11, ()iif xy x則(3)式化為111(,)iiiiyyhf xy11(,)iif xy即Euler方法Euler方法也稱為一階Runge-Kut

3、ta方法2( )()ie hO h由于-(4)KK51nxnxxy)(xyy11212.( )( ),iiiiy xxxKKy xxx如果以在和 處的斜率和的算術(shù)平均值作為在上的平均斜率11()iKy x11(,)iif xy2( )iKy x , ( )iif x y x( ,)iif x y11( ,)iif x yhK(由(4)式)令221KKK則(3)式化為K1K2K611211211()()222iiihyyhKKyKK111(,)iiKf xy211( ,)iiKf x yhK)(00 xyy -(5)即改進Euler公式,也稱為二階Runge-Kutta法3( )()ie hO

4、h7三、高階Runge-Kutta方法112,iimxxKKK如果上多預估幾個點上的斜率值、1( ),iiy xxx用它們的加權(quán)平均值作為在上的平均斜率1nxnx21hxnxy)(xyy K1K2KmK111221222()( ,)(,)(,)iimmiiiimimimyyhKKKKf x yKf xh yhKf xh yh01,()jjjjjyhy xh其中是的預估值。8RKjjjh上述公式就是公式。式中，都是待定參數(shù)，可以通過是公式的局部截斷誤差關(guān)于 的階數(shù)盡量高來確定。下面以二階R-K 為例，說明R-K公式的構(gòu)造過程。1112212221()( ,)(,)iiiiiiyyhKKKf x

5、yKf xh yhK122RK上述公式就是二階公式。式中，都是待定參數(shù)。-(6)1226為確定， ，需分析式（ ）的局部截斷誤差。91112212221212221222( ,)(,)( ,)( ,)( ,)( ,)( ,)()( ,)( ,)()( ,)( ,)iiiiiiiiiiiiiiiiiiiiiiiiiiiyyhKhKyhf x yhf xh yhKf x yf x yyhf x yhf x yhf x yo hxyf x yf x yyhf x yhf x yxy 3()o h( )iiyy x假定1122322231222()( ,)( ,)( ,)( ,)()( )()( )(

6、 )()( )( )( )( )( )iiiiiiiiiiiiiyhf xf xf xhf xo hxyy xhy xh yy xy xyxo hxxxyy 另一方面231()( )( )( )()2iiiihy xy xhy xyxo h1011231222( )()11 ()( )( )()2niiiie hy xyhy xhy xo h 6式（ ）的局部截斷誤差:1222112 22R K1110121212有無窮多組解，從而可以得到許多具體的二階 - 公式如：，111123(4)6iihyyKKK111(,)iiKf xy2111(,)22iihhKf xyK31121(,(2)iiK

7、f xh yhKK-(7)(00 xyy (7)式稱為三階Runge-Kutta方法用類似的方法也可以構(gòu)造三階R-K公式4( )( )()iiiy xye hO h因而三階R-K方法(7)具有3階精度1211234(22)6iihyyKKKK111(,)iiKf xy2111(,)22iihhKf xyK4113(,)iiKf xh yhK)(00 xyy 類似于(7)式,還可構(gòu)造四階(經(jīng)典)Runge=Kutta方法3112(,)22iihhKf xyK5( )()ie hO h-(8)因而方法(8)有4階精度13例1. 使用高階R-K方法計算初值問題1)0(5 . 002yxyy. 1 .

8、 0h取解:(1) 使用三階R-K方法時1n201yK 12102)21 . 0(KyK1025. 121203)2( 1 . 0(KKyK2555. 1)4(61 . 032101KKKyy1111. 1RK.m14其余結(jié)果如下:(2) 如果使用四階R-K方法時1n201yK 12102)21 . 0(KyK1025. 1 n xn k1 k2 k3 yn 1.0000 0.1000 1.0000 1.1025 1.2555 1.1111 2.0000 0.2000 1.2345 1.3755 1.5945 1.2499 3.0000 0.3000 1.5624 1.7637 2.0922

9、1.4284 4.0000 0.4000 2.0404 2.3423 2.8658 1.6664 5.0000 0.5000 2.7768 3.2587 4.1634 1.9993152203)21 . 0(KyK1133. 12304)1 . 0(KyK2351. 1)22(61 . 0432101KKKKyy1111. 1其余結(jié)果如下: i x i K1 K2 K3 K4 yi 1.0000 0.1000 1.0000 1.1025 1.1133 1.2351 1.1111 2.0000 0.2000 1.2346 1.3756 1.3921 1.5633 1.2500 3.0000 0.

10、3000 1.5625 1.7639 1.7908 2.0423 1.4286 4.0000 0.4000 2.0408 2.3428 2.3892 2.7805 1.6667 5.0000 0.5000 2.7777 3.2600 3.3476 4.0057 2.00001600.050.10.150.20.250.30.350.40.450.511.21.41.61.822.2精確 解 00.050.10.150.20.250.30.350.40.450.511.21.41.61.822.2精確 解 Euler法 00.050.10.150.20.250.30.350.40.450.511.21.41.61.822.2精確 解 Euler法 預測 -校 正 系 統(tǒng) 00.050.10.150.20.250.30.350.40.450.511.21.41.61.822.2精確 解 Euler法 預測 -校 正 系 統(tǒng) 改 進 Eul er法 00.050.10.150.20.250.30.350.40.450.511.21.41.61.822.2精確 解 Euler法 預測 -校 正 系 統(tǒng) 改 進 Eul er法 3階RK法 00.050.10.150.20.250.30.350.40.450.5