龍格-庫塔方法

1、第二節(jié) 龍格-庫塔方法基本思想：基本思想：利用利用 在某些特殊點(diǎn)上的在某些特殊點(diǎn)上的函數(shù)值的函數(shù)值的線性線性組合組合來構(gòu)造高階單步法的來構(gòu)造高階單步法的。( , )f x y第二節(jié)第二節(jié) 龍格龍格- -庫塔法庫塔法什么叫平均斜率？什么叫平均斜率？對差商對差商 應(yīng)用微分中值定理，有，應(yīng)用微分中值定理，有，1()()iiy xy xh 1()()()iiiy xy xhy xh 利用微分方程利用微分方程 ，有，有( , )yf x y 1()()(, ()iiiiy xy xhf xh y xh 這里的這里的 稱為。稱為。(, ()iif xh y xh 第二節(jié) 龍格-庫塔方法可將改進(jìn)的歐拉格式改

2、寫成可將改進(jìn)的歐拉格式改寫成112121112(),(,),(,).iiiiiiyKKKyhfxyKyhfxK 的算術(shù)平均值作為平均斜率。的算術(shù)平均值作為平均斜率。該公式可以看作是用該公式可以看作是用 和和 兩個點(diǎn)處的斜率兩個點(diǎn)處的斜率 和和1 iixx 12 KK由改進(jìn)型歐拉公式我們可以猜想，如果在由改進(jìn)型歐拉公式我們可以猜想，如果在1,iix x 內(nèi)多預(yù)測幾個點(diǎn)的斜率，再對他們進(jìn)行加權(quán)平均，內(nèi)多預(yù)測幾個點(diǎn)的斜率，再對他們進(jìn)行加權(quán)平均，可能得到精度更好的平均斜率！可能得到精度更好的平均斜率！第二節(jié) 龍格-庫塔方法下面以下面以2階龍格階龍格- -庫塔方法庫塔方法為例來闡述這種思想為例來闡述這種

3、思想考察區(qū)間考察區(qū)間 上的一點(diǎn)上的一點(diǎn) ，1,iix x 01, ipixxphp 用用 和和 的斜率的斜率 和和 的加權(quán)平均作為平均的加權(quán)平均作為平均 iipxx 12 KK斜率斜率 的近似值：的近似值：*K1122*KKK 即取即取11122()iiyyhKK 其中其中 和和 是待定常數(shù)。若取是待定常數(shù)。若取 ，則，則12 1(,)iiKf xy 問題在于如何確定問題在于如何確定 處的斜率處的斜率 和常數(shù)和常數(shù) 和和 。ipx 2K12 第二節(jié) 龍格-庫塔方法仿照改進(jìn)的歐拉方法，用歐拉方法預(yù)測仿照改進(jìn)的歐拉方法，用歐拉方法預(yù)測 的值，的值，()ipy x 1ipiyyphK 并用它來估計(jì)斜

4、率并用它來估計(jì)斜率 ：2K2(,)ipipKf xy 于是得到如下形式的算法：于是得到如下形式的算法：111221211(),(,),(,).iiiiiiyyhKKKf xyKf xyphK 通過適當(dāng)選取參數(shù)通過適當(dāng)選取參數(shù) 和和 的值，使得公式具有的值，使得公式具有12, p 2階精度！階精度！第二節(jié) 龍格-庫塔方法由泰勒公式展開，要使公式具有由泰勒公式展開，要使公式具有 2 階精度階精度，只需，只需12211 2,p方程組有方程組有無窮無窮多解：多解：二級二級方法有無窮多種方法有無窮多種常見的常見的3種二級方法：種二級方法： 中點(diǎn)法（修正的中點(diǎn)法（修正的Euler法法）1221012,cc

5、a 取取122(,(,)nnnnnnhhyyhf xyf xy 二階二階龍格庫塔方法龍格庫塔方法122112,cca 取取12 (,)(,(,)nnnnnnnnhyyf xyf xh yhf xy 3()O h第二節(jié) 龍格-庫塔方法三級方法：三級方法：N = 3 類似于類似于N = 2的推導(dǎo)方法，可得到的推導(dǎo)方法，可得到1231;ccc 223312;c ac a 22223313;c ac a 323216c a b 4()O h常見的常見的2種三階方法：種三階方法： 庫塔庫塔三階方法三階方法112346()nnhyykkk 1(,);nnkf xy 2122(,)nnhhkf xyk312

6、2(,)nnkf xh yhkhk 第二節(jié) 龍格-庫塔方法 四級方法：四級方法：N = 45()O h局部截?cái)嗾`差局部截?cái)嗾`差常見的常見的2種四階方法：種四階方法：經(jīng)典經(jīng)典龍格龍格- -庫塔庫塔方法方法11234226()nnhyykkkk 1(,)nnkf xy 2122(,)nnhhkf xyk3222(,)nnhhkf xyk43(,)nnkf xh yhk第二節(jié) 龍格-庫塔方法解：解：201( )dyxydxyy 01( , )x 例例2：用用經(jīng)典的經(jīng)典的龍格龍格- -庫塔庫塔方法方法求解下列初值問題求解下列初值問題 。0 1 .h 經(jīng)典的四階經(jīng)典的四階龍格龍格- -庫塔公式：庫塔公式

7、：11234226()nnhyykkkk 12;nnnxkyy 2112222()nnnhxhkykhyk 4332()nnnxhkyhkyhk 3222222();nnnhxhkykhyk 第二節(jié) 龍格-庫塔方法 0.1 0.2 0.3 0.4 0.5 1.0954 1.1832 1.2649 1.3416 1.4142nxny 0.6 0.7 0.8 0.9 1.0 1.4832 1.5492 1.6125 1.6733 1.7321nxny同保留同保留5位的位的精確值精確值完全一致：完全一致： 0.1 0.2 0.3 0.4 0.5 1.0954 1.1832 1.2649 1.3416

8、 1.4142nxny 0.6 0.7 0.8 0.9 1.0 1.4832 1.5492 1.6125 1.6733 1.7321nxny21yx第二節(jié) 龍格-庫塔方法第二節(jié) 龍格-庫塔方法二、二、高階高階和隱式和隱式Runge-Kutta方法方法注注：對于顯式對于顯式N級級R-K方法，最多只能得到方法，最多只能得到N級方法；級方法； N 1,2,3,4 5,6,7 8,9 10,11, N N-1 N-2()p N已經(jīng)證明已經(jīng)證明N級級R-K方法的方法的階階 具有下列關(guān)系：具有下列關(guān)系：()p N2N 若要得到若要得到N階以上方法，則使用階以上方法，則使用N級隱式級隱式R-K方法方法 N級

9、隱式級隱式R-K方法的一般形式：方法的一般形式：11Nnniiiyyhc k 11(,);,.,Nininijjjkf xa h yhb kiN N級隱式級隱式R-K法法可以達(dá)到可以達(dá)到2N階階第二節(jié) 龍格-庫塔方法(1)(1)一一級級二階二階的隱式的隱式中點(diǎn)中點(diǎn)方法：方法：11nnyyhk 1122(,)nnhkhkf xy (2)(2)二二級級四階四階的隱式的隱式R-K方法：方法：1212()nnh kkyy 112131326446() ,()nnhkkf xh yhk221131326464() ,()nnhkkf xh yhk 第二節(jié) 龍格-庫塔方法三、三、變步長變步長方法方法基本基

10、本思想思想：根據(jù)精度：根據(jù)精度自動自動地選擇地選擇步長步長對于對于經(jīng)典經(jīng)典Runge-Kutta方法：方法：0 1 2, , ,n Step1:設(shè)從設(shè)從 出發(fā)，以出發(fā)，以 為步長，經(jīng)過為步長，經(jīng)過一步一步計(jì)算得到計(jì)算得到nxh511( )()hnny xyChStep2:取取 為步長，再從為步長，再從 出發(fā)，經(jīng)過出發(fā)，經(jīng)過兩步兩步計(jì)算得到計(jì)算得到nx2h521122()()( )hnnhy xyC第二節(jié) 龍格-庫塔方法21111116( / )( )()()hnnhnny xyy xy 2111116( / )( )( ()()hhnnnny xyy xy 221111115( / )( /

11、)( )()hhhnnnny xyyy211( / )( )|hhnnyy 記記如果如果 ，則將步長，則將步長折半折半進(jìn)行計(jì)算，直到進(jìn)行計(jì)算，直到 為止為止 此時(shí)取此時(shí)取 為最終結(jié)果；為最終結(jié)果；21( / )hny 如果如果 ，則將步長，則將步長加倍加倍進(jìn)行計(jì)算，直到進(jìn)行計(jì)算，直到 為止為止 此時(shí)將步長此時(shí)將步長折半折半一次計(jì)算，得到的為最終結(jié)果。一次計(jì)算，得到的為最終結(jié)果。第二節(jié) 龍格-庫塔方法一、一、收斂性收斂性 /* *Convergence* */3 單步法的單步法的收斂性收斂性、相容性相容性和和絕對穩(wěn)定性絕對穩(wěn)定性1(, )nnnnyyhxy h 對于初值問題對于初值問題 的一的一

12、種種000( , )();dyf x ydxy xyxx ( ) 1Def單步法單步法 產(chǎn)生的近似解，如果產(chǎn)生的近似解，如果 對于任一對于任一固定固定的的 ，均有，均有 ，0nxxnh 則稱該單步法是則稱該單步法是收斂收斂的。的。0lim()nnhyy x 類似地可以定義類似地可以定義隱式隱式單步法、多步法（單步法、多步法（4）的）的收斂性收斂性第二節(jié) 龍格-庫塔方法3 1 .Th設(shè)初值問題（設(shè)初值問題（*）對應(yīng)的下列）對應(yīng)的下列單步法單步法是是 階的，階的，1(, )nnnnyyhxy h p且函數(shù)且函數(shù) 滿足對滿足對 的的Lipschitz條件，即存在常數(shù)條件，即存在常數(shù)y 0L 1212

13、12| ( , )( , )|,x y hx y hL yyy y 則該則該單步法單步法是收斂的，且是收斂的，且()()pnny xyO h 證明：證明：()nnney xy 記記由由截?cái)嘟財(cái)嗾`差的定義誤差的定義11()()(, (), )nnnnny xy xhxy xhT 11 (, (), )(, )nnnnnnneehxy xhxyhT 第二節(jié) 龍格-庫塔方法因?yàn)橐驗(yàn)閱尾椒▎尾椒ㄊ鞘?階的：階的：p000,hhh 滿足滿足11|pnTCh 11| |pnnneehL eCh |ne 其中其中11,phLCh 212|nnneee 3231|()ne 2101|(.)nnnee 1000

14、1| exp ()|exp ()pnnneL xxeCh LL xx 101exp ()pnCh LL xx ()pneO h 00()h第二節(jié) 龍格-庫塔方法二、二、相容性相容性 /* *Consistency* */()( )( , ( ), )y xhy xhx y x h 0 ( )( ).( ) ( , ( ), ).y xhy xy xhx y x 0( )( , ( ), ).h y xx y x 1()pO h 100( )( , ( ), )py xx y x 對于對于 階方法：階方法：1(, )nnnnyyhxy h p() 若方法（若方法（*）的）的增量增量函數(shù)滿足：函數(shù)

15、滿足：2Def0( , , )( , )x yf x y 則稱該方法與初值問題（則稱該方法與初值問題（*）相容相容。第二節(jié) 龍格-庫塔方法設(shè)方法（設(shè)方法（*）與初值問題（）與初值問題（*）相容相容，且，且 滿足滿足L-條件，條件， 則該方法（則該方法（*）是）是收斂收斂的，即當(dāng)?shù)模串?dāng) 固定，固定， 時(shí)時(shí)nxx 0h()nnyy x1(, )nnnnyyxyhh 0(, (), )()nnnxy xy x ()(, ()nnny xf xy x 再由再由相容性相容性得：得：上式說明：當(dāng)上式說明：當(dāng) 時(shí)，方法（時(shí)，方法（*）趨于）趨于原微分方程原微分方程0h 本章討論的數(shù)值方法都是與原初值問題本

16、章討論的數(shù)值方法都是與原初值問題相容相容的的 第二節(jié) 龍格-庫塔方法三、三、絕對穩(wěn)定性絕對穩(wěn)定性 /* *Absolute Stibility* */計(jì)算過程中產(chǎn)生的計(jì)算過程中產(chǎn)生的舍入誤差舍入誤差對計(jì)算結(jié)果的影響對計(jì)算結(jié)果的影響首先以首先以Euler公式為例，來討論一下公式為例，來討論一下舍入誤差舍入誤差的傳播的傳播:1(,)nnnnyyhf xy 設(shè)設(shè)實(shí)際實(shí)際計(jì)算得到的點(diǎn)計(jì)算得到的點(diǎn) 的的近似近似函數(shù)值為函數(shù)值為 ，nnnyy nxny其中其中 為為精確值精確值， 為誤差為誤差n 1(,)nnnnyyhf xy 111nnnyy 11 (,)(,)(, )nnnnnnynnh f xyf

17、xyhfx 如果如果 ，則誤差是，則誤差是不增不增的，故可認(rèn)為是的，故可認(rèn)為是穩(wěn)定穩(wěn)定的的11|yhf 第二節(jié) 龍格-庫塔方法例如：例如：對于初值問題對于初值問題0()yyy xa 精確解精確解為為0 x xyae 而而實(shí)際求解實(shí)際求解的初值問題為的初值問題為0()yyy xaa 精確解精確解為為0()x xyaa e 在在 處的誤差為處的誤差為nx0nxxae 可見誤差隨著可見誤差隨著 的增加呈的增加呈指數(shù)函數(shù)指數(shù)函數(shù)增長增長nx如果初值問題為如果初值問題為0()yyy xa 精確解精確解為為0 xxyae 第二節(jié) 龍格-庫塔方法實(shí)際求解實(shí)際求解的初值問題為的初值問題為0()yyy xaa

18、精確解精確解為為0()xxyaa e 在在 處的誤差為處的誤差為nx0nxxae 可見誤差隨著可見誤差隨著 的增加呈的增加呈指數(shù)函數(shù)指數(shù)函數(shù)遞減遞減nx當(dāng)當(dāng) 時(shí)，微分方程是時(shí)，微分方程是不穩(wěn)定不穩(wěn)定的；的；0yf 而而 時(shí)，微分方程是時(shí)，微分方程是穩(wěn)定穩(wěn)定的。的。0yf 上面討論的上面討論的穩(wěn)定性穩(wěn)定性，與，與數(shù)值方法數(shù)值方法和方程中和方程中 有關(guān)有關(guān)f第二節(jié) 龍格-庫塔方法實(shí)驗(yàn)實(shí)驗(yàn)方程：方程：0,Re( )yyC 1(, )nnnnyyhxy h 對單步法對單步法 應(yīng)用應(yīng)用實(shí)驗(yàn)實(shí)驗(yàn)方程，方程，1()nnyEh y 3Defh 如果如果 ，當(dāng)，當(dāng) 時(shí)，則稱該時(shí)，則稱該1()Eh 單步法是單步法是絕對穩(wěn)定絕對穩(wěn)定的，在復(fù)平面上復(fù)變量的，在復(fù)平面上復(fù)變量 滿足滿足1()Eh 的區(qū)域，稱為該單步法的絕對穩(wěn)定的區(qū)域，稱為該單步法的絕對穩(wěn)定域域，它與它與實(shí)軸實(shí)軸的的交集交集稱為絕對穩(wěn)定稱為絕對穩(wěn)定區(qū)間區(qū)間。1111()()pnnnTy xyO h 若單步法是若單步法是 階的，則階的，則p由由實(shí)驗(yàn)實(shí)驗(yàn)方程可得：方程可得：1()exp()nny xyh 1exp()()()pnnyhEh yO h ()exp()Ehh 第二節(jié) 龍格-庫塔方法11