第9章 常微分方程初值問題數(shù)值解法

1、上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)第第9章章 常微分方程初值問題數(shù)值解法常微分方程初值問題數(shù)值解法 9.1 引言引言 9.2 簡(jiǎn)單的數(shù)值方法簡(jiǎn)單的數(shù)值方法 9.3 龍格龍格- -庫(kù)塔方法庫(kù)塔方法 9.4 單步法的收斂性與穩(wěn)定性單步法的收斂性與穩(wěn)定性 9.5 線性多步法線性多步法 9.6 線性多步法的收斂性與穩(wěn)定性線性多步法的收斂性與穩(wěn)定性 9.7 一階方程組與剛性方程組一階方程組與剛性方程組上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)9.1 引引 言言 科學(xué)技術(shù)中很多問題都可用微分方程的科學(xué)技術(shù)中很多問題都可用微分方程的定解問定解問題題來描述，主要有來描述，主要有初值問題初值

2、問題與與邊值問題邊值問題兩大類，兩大類，本章本章只考慮初值問題只考慮初值問題. 常微分方程初值問題中最簡(jiǎn)單的例常微分方程初值問題中最簡(jiǎn)單的例子是人口模型，設(shè)某特定區(qū)域在子是人口模型，設(shè)某特定區(qū)域在t0時(shí)刻人口為時(shí)刻人口為y(t0)=y0已知的，該區(qū)域的人口自然增長(zhǎng)率為已知的，該區(qū)域的人口自然增長(zhǎng)率為 ，人口增長(zhǎng)與，人口增長(zhǎng)與人口總數(shù)成正比，所以人口總數(shù)成正比，所以t時(shí)刻的人口總數(shù)時(shí)刻的人口總數(shù)y(t)滿足以下滿足以下微分方程微分方程 .)(),(00ytytyy 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)很多物理系統(tǒng)與時(shí)間有關(guān)，從衛(wèi)星運(yùn)行軌道到單很多物理系統(tǒng)與時(shí)間有關(guān)，從衛(wèi)星運(yùn)行軌道到單

3、擺運(yùn)動(dòng)，從化學(xué)反應(yīng)到物種競(jìng)爭(zhēng)都是隨時(shí)間的延續(xù)而擺運(yùn)動(dòng)，從化學(xué)反應(yīng)到物種競(jìng)爭(zhēng)都是隨時(shí)間的延續(xù)而不斷變化的不斷變化的. 解常微分方程是描述連續(xù)變化的數(shù)學(xué)語(yǔ)解常微分方程是描述連續(xù)變化的數(shù)學(xué)語(yǔ)言，微分方程的求解就是確定滿足給定方程的可微函言，微分方程的求解就是確定滿足給定方程的可微函數(shù)數(shù)y(t)，研究它的數(shù)值方法是本章的主要目的，研究它的數(shù)值方法是本章的主要目的. 考慮考慮一一階常微分方程的初值問題階常微分方程的初值問題 )2 . 1(.)()1 . 1(,),(000yxybxxyxfy則稱則稱f關(guān)于關(guān)于y滿足滿足利普希茨利普希茨(Lipschitz)條件條件，L稱為稱為y的的利普希茨常數(shù)利普希茨常

4、數(shù)(簡(jiǎn)稱簡(jiǎn)稱Lips.常數(shù)常數(shù)).如果存在實(shí)數(shù)如果存在實(shí)數(shù)L0，使得，使得)3 . 1(.,),(),(212121RyyyyLyxfyxf 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 定理定理1 設(shè)設(shè)f在區(qū)域在區(qū)域D=(x, y)|a x b,y R上連續(xù)上連續(xù), 關(guān)關(guān)于于y滿足利普希茨條件，則對(duì)任意滿足利普希茨條件，則對(duì)任意x0 a, b, y0 R，常，常微分方程初值問題微分方程初值問題(1.1)式和式和(1.2)式當(dāng)式當(dāng)x a, b時(shí)存在時(shí)存在唯一的連續(xù)可微解唯一的連續(xù)可微解y(x) .解的存在唯一性定理是常微分方程理論的基本內(nèi)解的存在唯一性定理是常微分方程理論的基本內(nèi)容，也是數(shù)

5、值方法的出發(fā)點(diǎn)，此外還要考慮方程的解容，也是數(shù)值方法的出發(fā)點(diǎn)，此外還要考慮方程的解對(duì)擾動(dòng)的敏感性，它有以下結(jié)論對(duì)擾動(dòng)的敏感性，它有以下結(jié)論. 定理定理2 設(shè)設(shè)f在區(qū)域在區(qū)域D (如定理如定理1所定義所定義) 上連續(xù)上連續(xù), 且且關(guān)于關(guān)于y滿足利普希茨條件，設(shè)初值問題滿足利普希茨條件，設(shè)初值問題.)(),()(0sxyyxfxy 的解為的解為y(x, s)，則，則.),(),(21210ssesxysxyxxL 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)這個(gè)定理表明解對(duì)初值依賴的敏感性，它與右端這個(gè)定理表明解對(duì)初值依賴的敏感性，它與右端函數(shù)函數(shù)f有關(guān)，當(dāng)有關(guān)，當(dāng)f的的Lips.常數(shù)常數(shù)L比

6、較小時(shí)，解對(duì)初值和比較小時(shí)，解對(duì)初值和右端函數(shù)相對(duì)不敏感，可視為好條件右端函數(shù)相對(duì)不敏感，可視為好條件. 若若L較大則可較大則可認(rèn)為壞條件，即病態(tài)問題認(rèn)為壞條件，即病態(tài)問題.如果右端函數(shù)可導(dǎo)，由中值定理有如果右端函數(shù)可導(dǎo)，由中值定理有.,),(),(),(212121之間之間在在yyyyyxfyxfyxf 若假定若假定 在域在域D內(nèi)有界內(nèi)有界, 設(shè)設(shè) , 則則yyxf ),(Lyyxf ),(.),(),(2121yyLyxfyxf 它表明它表明f滿足利普希茨條件，且滿足利普希茨條件，且L的大小反映了右端函的大小反映了右端函數(shù)數(shù)f關(guān)于關(guān)于y變化的快慢，刻畫了初值問變化的快慢，刻畫了初值問(1.

7、1)式和式和(1.2)式式是否為好條件是否為好條件. 這在數(shù)值求解中也是很重要的這在數(shù)值求解中也是很重要的.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 雖然求解常微分方程有各種各樣的解析方法，但雖然求解常微分方程有各種各樣的解析方法，但解析方法只能用來求解一些特殊類型的方程，實(shí)際問解析方法只能用來求解一些特殊類型的方程，實(shí)際問題中歸結(jié)出來的微分方程主要靠數(shù)值解法題中歸結(jié)出來的微分方程主要靠數(shù)值解法. 所謂所謂數(shù)值解法數(shù)值解法, 就是尋求解就是尋求解y(x)在一系列離散節(jié)點(diǎn)在一系列離散節(jié)點(diǎn) 121nnxxxx上的近似值上的近似值 y1,y2,yn,yn+1,. 相鄰兩個(gè)節(jié)點(diǎn)的間距相鄰兩個(gè)節(jié)

8、點(diǎn)的間距hn=xn+1- -xn稱為稱為步長(zhǎng)步長(zhǎng). 今后如不特別說明，總是假定今后如不特別說明，總是假定 hi=h(i=1,2,)為為常數(shù)常數(shù), 這時(shí)節(jié)點(diǎn)為這時(shí)節(jié)點(diǎn)為xn=x0+nh(i=0,1,2,) (等距節(jié)點(diǎn)等距節(jié)點(diǎn)).上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 初值問題的初值問題的數(shù)值解法數(shù)值解法有個(gè)有個(gè)基本特點(diǎn)基本特點(diǎn)，他們都采，他們都采取取“步進(jìn)式步進(jìn)式”，即求解過程順著節(jié)點(diǎn)排列的次序一步，即求解過程順著節(jié)點(diǎn)排列的次序一步一步地向前推進(jìn)一步地向前推進(jìn). 描述這類算法，只要給出用已知信描述這類算法，只要給出用已知信息息yn,yn- -1,yn- -2,計(jì)算計(jì)算yn+1的遞推公式

9、的遞推公式. 本章首先要對(duì)常微分方程本章首先要對(duì)常微分方程(1.1)離散化，建立求解離散化，建立求解數(shù)值解的遞推公式數(shù)值解的遞推公式. 一類是計(jì)算一類是計(jì)算yn+1時(shí)只用到前一點(diǎn)的時(shí)只用到前一點(diǎn)的值值yn，稱為，稱為單步法單步法. 另一類是用到另一類是用到y(tǒng)n+1前面前面 k 點(diǎn)的值點(diǎn)的值yn,yn-1, yn-k+1，稱為，稱為k步法步法. 其次，要研究公式的其次，要研究公式的局局部截?cái)嗾`差部截?cái)嗾`差和和階階，數(shù)值解，數(shù)值解yn與精確解與精確解y(xn)的的誤差估計(jì)誤差估計(jì)及及收斂性收斂性，還有遞推公式的，還有遞推公式的計(jì)算穩(wěn)定性計(jì)算穩(wěn)定性等問題等問題.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下

10、頁(yè)下頁(yè)下頁(yè)9.2 簡(jiǎn)單的數(shù)值方法簡(jiǎn)單的數(shù)值方法9.2.1 歐拉法與后退歐拉法歐拉法與后退歐拉法 我們知道，在我們知道，在xy平面上，微分方程平面上，微分方程(1.1)式式的解的解y=f(x)稱作它的稱作它的積分曲線積分曲線，積分曲線積分曲線上一點(diǎn)上一點(diǎn)(x, y)的切的切線斜率等于函數(shù)線斜率等于函數(shù)f(x, y)的值的值. 如果按如果按f(x, y)在在xy平面上平面上建立一個(gè)方向場(chǎng)，那么，建立一個(gè)方向場(chǎng)，那么，積分曲線積分曲線上每一點(diǎn)的切線上每一點(diǎn)的切線方向均與方向場(chǎng)在該點(diǎn)的方向相一致方向均與方向場(chǎng)在該點(diǎn)的方向相一致.基于上述幾何解釋，我們從初始點(diǎn)基于上述幾何解釋，我們從初始點(diǎn)P0(x0,

11、y0)出發(fā)出發(fā),先依方向場(chǎng)在該點(diǎn)的方向推進(jìn)到先依方向場(chǎng)在該點(diǎn)的方向推進(jìn)到x=x1上一點(diǎn)上一點(diǎn)P1，然后，然后再?gòu)脑購(gòu)腜1點(diǎn)依方向場(chǎng)在該點(diǎn)的方向推進(jìn)到點(diǎn)依方向場(chǎng)在該點(diǎn)的方向推進(jìn)到 x=x2 上一點(diǎn)上一點(diǎn)P2 , 循環(huán)前進(jìn)做出一條循環(huán)前進(jìn)做出一條折線折線P0 P1 P2.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 一般地，設(shè)已做出該折線的頂點(diǎn)一般地，設(shè)已做出該折線的頂點(diǎn)Pn, ,過過Pn(xn, yn)依方向場(chǎng)的方向再推進(jìn)到依方向場(chǎng)的方向再推進(jìn)到Pn+1(xn+1, yn+1)，顯然兩個(gè)，顯然兩個(gè)頂點(diǎn)頂點(diǎn)Pn, ,Pn+1的坐標(biāo)有關(guān)系的坐標(biāo)有關(guān)系),(,()(),(111nnnnnnnnnn

12、nxyxfxyyxfhyyxxyy ) 1 .2 (),(1nnnnyxhfyy 這就是著名的這就是著名的( (顯式顯式) )歐拉歐拉( (Euler) )公式公式. . 若初值若初值y0已已知，則依公式知，則依公式(2.1)可逐次逐步算出各點(diǎn)數(shù)值解可逐次逐步算出各點(diǎn)數(shù)值解. .即即),(0001yxhfyy ),(1112yxhfyy 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 例例1 用歐拉公式求解初值問題用歐拉公式求解初值問題)2 . 2(. 1)0(),10(2 yxyxyy 解解 取步長(zhǎng)取步長(zhǎng)h=0.1，歐拉公式的具體形式為，歐拉公式的具體形式為)2(1nnnnnyxyhyy

13、其中其中xn=nh=0.1n (n=0,1,10), 已知已知y0 =1, 由此式可得由此式可得191818. 1)1 . 12 . 01 . 1 ( 1 . 01 . 1)2(1 . 11 . 01)2(1111200001 yxyhyyyxyhyy上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)依次計(jì)算下去，依次計(jì)算下去，部分計(jì)算結(jié)果部分計(jì)算結(jié)果見下表見下表. xy21 與準(zhǔn)確解與準(zhǔn)確解 相比，可看出歐拉公式的計(jì)算結(jié)相比，可看出歐拉公式的計(jì)算結(jié)果精度很差果精度很差. xn 歐拉公式數(shù)值解歐拉公式數(shù)值解yn準(zhǔn)確解準(zhǔn)確解y(xn) 誤差誤差 0.2 0.4 0.6 0.8 1.0 1.1918

14、18 1.358213 1.508966 1.649783 1.784770 1.183216 1.341641 1.483240 1.612452 1.732051 0.008602 0.016572 0.025726 0.037331 0.052719上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 歐拉公式具有明顯的幾何意義歐拉公式具有明顯的幾何意義, , 就是就是用折線近似用折線近似代替方程的解曲線代替方程的解曲線，因而常稱公式，因而常稱公式(2.1)為為歐拉折線法歐拉折線法. .( )yy xxynx1nxnp1np1np x 還可以通過幾何直觀來考察歐拉方法的精度還可以通過幾何直觀

15、來考察歐拉方法的精度. .假假設(shè)設(shè)yn=y(xn), ,即頂點(diǎn)即頂點(diǎn)Pn落在積分曲線落在積分曲線y=y(x)上，那么，上，那么，按歐拉方法做出的折線按歐拉方法做出的折線PnPn+1便是便是y=y(x)過點(diǎn)過點(diǎn)Pn的切線的切線. .從圖形上看從圖形上看, ,這這樣定出的頂點(diǎn)樣定出的頂點(diǎn)Pn+1顯著顯著地偏離了原來的積分曲地偏離了原來的積分曲線，可見歐拉方法是線，可見歐拉方法是相相當(dāng)粗糙當(dāng)粗糙的的. .上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 為了分析計(jì)算公式的精度，通?？捎锰├照归_為了分析計(jì)算公式的精度，通?？捎锰├照归_將將y(xn+1)在在xn處展開，則有處展開，則有).,()(2),

16、()()(2)()()()(1221 nnnnnnnnnnnnxxyhyxhfxyyhxyhxyhxyxy 在在yn=y(xn)的前提下，的前提下，f(xn,yn )=f(xn,y(xn)=y ( (xn n) ). .于是于是可得歐拉法可得歐拉法(2.1)的的公式誤差公式誤差為為) 3 . 2(),(2)(2)(2211nnnnxyhyhyxy 稱為此方法的稱為此方法的局部截?cái)嗾`差局部截?cái)嗾`差. .上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 如果對(duì)方程如果對(duì)方程(1.1)從從xn到到xn+1積分，得積分，得) 4 . 2(.)(,()()(11 nnxxnndttytfxyxy右端積分

17、用右端積分用左矩形公式左矩形公式hf(xn,y(xn)近似，再以近似，再以yn代替代替y(xn)，yn+1代替代替y(xn+1)也得到歐拉公式也得到歐拉公式(2.1)，局部截，局部截?cái)嗾`差也是斷誤差也是(2.3).稱為稱為( (隱式隱式) )后退的歐拉公式后退的歐拉公式. .它也可以通過利用均差它也可以通過利用均差近似導(dǎo)數(shù)近似導(dǎo)數(shù)y (xn+1)，即，即 如果右端積分用如果右端積分用右矩形公式右矩形公式hf(xn+1,y(xn+1)近似，近似，則得到另一個(gè)公式則得到另一個(gè)公式) 5 . 2 (),(111 nnnnyxhfyy直接得到直接得到. .),(,()()()(1111 nnnnnnn

18、xyxfxyxxxyxy上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 后退的歐拉公式與歐拉公式有著本質(zhì)的區(qū)別后退的歐拉公式與歐拉公式有著本質(zhì)的區(qū)別, 后后者是關(guān)于者是關(guān)于yn+1的一個(gè)直接計(jì)算公式，這類公式稱作是的一個(gè)直接計(jì)算公式，這類公式稱作是顯式的顯式的；前者公式的；前者公式的右端含有未知的右端含有未知的yn+1，它實(shí)際上，它實(shí)際上是關(guān)于是關(guān)于yn+1的一個(gè)函數(shù)方程的一個(gè)函數(shù)方程, ,這類方程稱作是這類方程稱作是隱式的隱式的. . 顯式顯式與與隱式隱式兩類方法各有特點(diǎn)，考慮到數(shù)值穩(wěn)兩類方法各有特點(diǎn)，考慮到數(shù)值穩(wěn)定性等其他因素，人們有時(shí)需要選用定性等其他因素，人們有時(shí)需要選用隱式隱式方法

19、，但方法，但使用使用顯式顯式算法遠(yuǎn)比算法遠(yuǎn)比隱式隱式方便方便. . 隱式方程隱式方程(2.5)通常用迭代法求解，而迭代過程通常用迭代法求解，而迭代過程的實(shí)質(zhì)是的實(shí)質(zhì)是逐步逐步顯式化顯式化. .上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 設(shè)用歐拉公式設(shè)用歐拉公式),() 0(1nnnnyxhfyy 給出迭代初值給出迭代初值 ，用它代入，用它代入(2.5)式的式的右端，使之轉(zhuǎn)右端，使之轉(zhuǎn)化為顯式，直接計(jì)算得化為顯式，直接計(jì)算得) 0(1 ny),() 0(11) 1 (1 nnnnyxhfyy然后再用然后再用 代入代入(2.5)式，又有式，又有) 1 (1 ny).,() 1 (11) 2(

20、1 nnnnyxhfyy如此反復(fù)進(jìn)行，得如此反復(fù)進(jìn)行，得) 6 . 2 ()., 1 , 0(),()(11) 1(1 kyxhfyyknnnkn上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)由于由于f(x, y)對(duì)對(duì)y滿足滿足Lipschitz條件條件(1.3). 由由(2.6)減減(2.5)得得),(),(11)(111) 1(1 nnknnnknyxfyxfhyy.1)(1 nknyyhL由此可知，只要由此可知，只要hL , ,我們反復(fù)將步長(zhǎng)我們反復(fù)將步長(zhǎng)折半計(jì)算折半計(jì)算, ,直至直至 為為止止, ,這時(shí)取最終得到的這時(shí)取最終得到的 作為作為結(jié)果；結(jié)果； 21hny(2) 如果如果 為

21、止，這時(shí)再將步長(zhǎng)折半計(jì)算一次，就得到所為止，這時(shí)再將步長(zhǎng)折半計(jì)算一次，就得到所要的結(jié)果要的結(jié)果. .變步長(zhǎng)方法還可利用變步長(zhǎng)方法還可利用p階與階與p+1階公式的局部截?cái)嚯A公式的局部截?cái)嗾`差得到誤差控制與變步長(zhǎng)的具體方法誤差得到誤差控制與變步長(zhǎng)的具體方法, ,參見文獻(xiàn)參見文獻(xiàn)8.8.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)9.4 單步法的收斂性與穩(wěn)定性單步法的收斂性與穩(wěn)定性9.4.1 收斂性與相容性收斂性與相容性 數(shù)值解法的基本思想是，通過某種離散化手段數(shù)值解法的基本思想是，通過某種離散化手段將微分方程將微分方程(1.1)轉(zhuǎn)化為差分方程轉(zhuǎn)化為差分方程,如單步法如單步法(2.11),即即

22、) 1 . 4(. ),(1hyxhyynnnn 它在點(diǎn)它在點(diǎn)xn處的解為處的解為yn，而初值問題，而初值問題(1.1), (1.2) 在點(diǎn)在點(diǎn)xn處的精確解為處的精確解為y( (xn) )，記，記en=y(xn)- -yn稱為稱為整體截?cái)嗾`整體截?cái)嗾`差差. 收斂性收斂性就是討論當(dāng)就是討論當(dāng) x=xn 固定且固定且 時(shí)時(shí)en0的問題的問題.00 nxxhn上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 定義定義3 若一種數(shù)值方法對(duì)于固定的若一種數(shù)值方法對(duì)于固定的xn=x0+nh, 當(dāng)當(dāng)h0時(shí)有時(shí)有yny(xn)，其中，其中y(x)是是(1.1),(1.2)的準(zhǔn)確解，的準(zhǔn)確解，則稱該方法是則稱

23、該方法是收斂收斂的的. 顯然數(shù)值方法收斂是指顯然數(shù)值方法收斂是指en=y(xn)- -yn0，對(duì)單步，對(duì)單步法法(4.1)有下述收斂性定理：有下述收斂性定理： 定理定理3 假設(shè)單步法假設(shè)單步法(4.1)具有具有p階精度，且增量函階精度，且增量函數(shù)數(shù) (x, y ,h)關(guān)于關(guān)于y滿足利普希次條件滿足利普希次條件又設(shè)初值又設(shè)初值y0是準(zhǔn)確的是準(zhǔn)確的, 即即y0=f(x0), 則其則其整體截?cái)嗾`差整體截?cái)嗾`差)2 . 4(.),(),(yyLhyxhyx ) 3 . 4().()(pnnhOyxy 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 證明證明 設(shè)以設(shè)以yn+1表示取表示取yn=y(xn

24、)用公式用公式(4.1)求得的求得的結(jié)果，即結(jié)果，即) 4 . 4(),),(,()(1hxyxhxyynnnn 則則y(xn)- -yn+1為局部截?cái)嗾`差，由于所給方法具為局部截?cái)嗾`差，由于所給方法具有有p階精度，按定義階精度，按定義2，存在定數(shù)，存在定數(shù)C ，使，使.)(111 pnnChyxy又由式又由式(4.1)與與(4.4)，得，得. ),(),(,()(11hyxhxyxhyxyyynnnnnnnn 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)利用利普希次條件利用利普希次條件(4.2)，有，有.)()1 (11nnnnyxyhLyy 從而有從而有.)()1 ()()(11111

25、11 pnnnnnnnnChyxyhLyxyyyyxy 即對(duì)整體截?cái)嗾`差即對(duì)整體截?cái)嗾`差en=y(xn)- -yn成立下列成立下列遞推關(guān)系式遞推關(guān)系式) 5 . 4(.)1 (11 pnnChehLe 據(jù)此不等式反復(fù)遞推，可得據(jù)此不等式反復(fù)遞推，可得) 6 . 4(.1)1()1 (0 npnnhLLChehLe 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)由此可以斷定，如果初值是準(zhǔn)確的，即由此可以斷定，如果初值是準(zhǔn)確的，即e0=0，則，則(4.3)式成立式成立. 定理證畢定理證畢.再注意到當(dāng)再注意到當(dāng)x=x0+nh T時(shí)時(shí) ,)1 ( TLnhLneehL 最終得下列估計(jì)式最終得下列估計(jì)

26、式) 7 . 4().1(0 TLpTLneLCheee依據(jù)這一定理，判斷單步法依據(jù)這一定理，判斷單步法(4.1)的收斂性，歸結(jié)的收斂性，歸結(jié)為驗(yàn)證增量函數(shù)為驗(yàn)證增量函數(shù) 能否滿足利普希次條件能否滿足利普希次條件(4.2).對(duì)于歐拉方法，由于其增量函數(shù)對(duì)于歐拉方法，由于其增量函數(shù) 就是就是f(x, y), 故當(dāng)故當(dāng)f(x, y)關(guān)于關(guān)于y滿足利普希次條件時(shí)它是收斂的滿足利普希次條件時(shí)它是收斂的.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)再考察改進(jìn)的歐拉方法，其增量函數(shù)再考察改進(jìn)的歐拉方法，其增量函數(shù) 已由已由(3.2)式給出，假定式給出，假定f(x, y)關(guān)于關(guān)于y滿足利普希次條件，即滿

27、足利普希次條件，即 ),(),(hyxhyx 這時(shí)有這時(shí)有.),(),(yyLyxfyxf ),(),(21yxfyxf .2121),(),(212121),(),(2121),(,(),(,(212yyLhLyyhLyyLyxfyxfhLyyLyyLyxhfyyxhfyLyyLyxhfyhxfyxhfyhxf 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)即即.21),(),(yyLhLhyxhyx 設(shè)限定步長(zhǎng)設(shè)限定步長(zhǎng)hh0(h0為定數(shù)為定數(shù))，上式表明，上式表明關(guān)于關(guān)于y的利普的利普希次常數(shù)為希次常數(shù)為.21 LhLL 因此改進(jìn)的歐拉方法也是收斂的因此改進(jìn)的歐拉方法也是收斂的.類似地

28、類似地, 不難驗(yàn)證其它龍格不難驗(yàn)證其它龍格- -庫(kù)塔方法的收斂性庫(kù)塔方法的收斂性.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 定理定理3表明表明p 1時(shí)單步法收斂時(shí)單步法收斂, 并且當(dāng)并且當(dāng)y(x)是初值是初值問題問題(1.1),(1.2)的解的解, (4.1)具有具有p階精度時(shí)階精度時(shí), 則有展開式則有展開式),(,()()(1hxyxhxyhxyTn ) 0),(,() 0),(,(2)()(2 hxyxxyxhhxyhxyx ).()0),(,()(2hOxyxxyh 所以所以p 1的充分必要條件是的充分必要條件是 ，而而 ，于是可給出如下定義：，于是可給出如下定義：0) 0),(

29、,()( xyxxy )(,()(xyxfxy 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 定義定義4 若單步法若單步法(4.1)的增量函數(shù)的增量函數(shù) 滿足滿足 以上討論表明以上討論表明p階方法階方法(4.1)當(dāng)當(dāng)p 1時(shí)與時(shí)與(1.1), (1.2)相容，反之相容，反之相容方法至少是相容方法至少是1階的階的.),()0 ,(yxfyx 則稱單步法則稱單步法(4.1)與初值問題與初值問題(1.1),(1.2)相容相容. 相容性是指方法逼近微分方程相容性是指方法逼近微分方程(1.1),即即(1.1)離散離散化得到的數(shù)值方法化得到的數(shù)值方法, 當(dāng)當(dāng)h0時(shí)可得到時(shí)可得到y(tǒng) (x)=f(x, y

30、) .于是有下面定理于是有下面定理. 定理定理4 p階方法階方法(4.1)與初值問題與初值問題(1.1),(1.2)相容的相容的充分必要條件是充分必要條件是p 1. 于是由定理于是由定理3可知方法可知方法(4.1)收斂的收斂的充分必要條件充分必要條件是是此方法是相容此方法是相容的的.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)9.4.2 絕對(duì)穩(wěn)定性與絕對(duì)穩(wěn)定域絕對(duì)穩(wěn)定性與絕對(duì)穩(wěn)定域 前面關(guān)于收斂性的討論有個(gè)前提，必須假定數(shù)前面關(guān)于收斂性的討論有個(gè)前提，必須假定數(shù)值方法本身的計(jì)算是準(zhǔn)確的值方法本身的計(jì)算是準(zhǔn)確的. 實(shí)際情形并不是這樣，實(shí)際情形并不是這樣，差分方程的求解還會(huì)有差分方程的求解還會(huì)

31、有計(jì)算誤差計(jì)算誤差. 譬如由于數(shù)字舍入譬如由于數(shù)字舍入而引起的而引起的小擾動(dòng)小擾動(dòng). 這類小擾動(dòng)在傳播過程中會(huì)不會(huì)惡這類小擾動(dòng)在傳播過程中會(huì)不會(huì)惡性增長(zhǎng)，以至于性增長(zhǎng)，以至于“淹沒淹沒”了差分方程的了差分方程的“真解真解”呢？呢？這就是這就是差分方程的穩(wěn)定性問題差分方程的穩(wěn)定性問題. 在實(shí)際計(jì)算時(shí)，我們?cè)趯?shí)際計(jì)算時(shí)，我們希望某一步產(chǎn)生的擾動(dòng)值，在后面的計(jì)算中希望某一步產(chǎn)生的擾動(dòng)值，在后面的計(jì)算中能夠被能夠被控制控制，甚至是，甚至是逐步衰減逐步衰減的的.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 定義定義5 若一種數(shù)值方法在節(jié)點(diǎn)值若一種數(shù)值方法在節(jié)點(diǎn)值yn上大小為上大小為 擾擾動(dòng)，于以后各

32、節(jié)點(diǎn)值動(dòng)，于以后各節(jié)點(diǎn)值ym(mn)上產(chǎn)生的偏差均不超過上產(chǎn)生的偏差均不超過 ，則稱該方法是，則稱該方法是穩(wěn)定穩(wěn)定的的. 下面以歐拉法為例考察計(jì)算穩(wěn)定性下面以歐拉法為例考察計(jì)算穩(wěn)定性. 例例4 用歐拉公式求解初值問題用歐拉公式求解初值問題 . 1)0(,100yyy 解解 用歐拉法解方程用歐拉法解方程y=- -100y 得得其準(zhǔn)確解其準(zhǔn)確解 是一個(gè)按指數(shù)曲線衰減很快的是一個(gè)按指數(shù)曲線衰減很快的函數(shù)函數(shù).(見書見書p294的圖的圖9-3)xexy100)( .)1001(1nnyhy 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)若取步長(zhǎng)若取步長(zhǎng)h=0.025，則歐拉公式的具體形式為，則歐拉公

33、式的具體形式為.5 . 11nnyy 節(jié)點(diǎn)節(jié)點(diǎn)xn歐拉方法歐拉方法yn后退歐拉方法后退歐拉方法yn0.0250.0500.0750.100- -1.5 2.25 - -3.375 5.06250.28570.08160.02330.0067計(jì)算結(jié)果見表計(jì)算結(jié)果見表, , 明顯計(jì)算過程不穩(wěn)定明顯計(jì)算過程不穩(wěn)定, , 但取但取h=0.005, , yn+1=0 0.5yn, , 則計(jì)算過程穩(wěn)定則計(jì)算過程穩(wěn)定. .對(duì)后退的歐拉公式，取對(duì)后退的歐拉公式，取h=0.025時(shí)，則計(jì)算公式時(shí)，則計(jì)算公式為為yn+1=(1/3.5)yn . .計(jì)算結(jié)果見表計(jì)算結(jié)果見表, , 這時(shí)計(jì)算過程是穩(wěn)這時(shí)計(jì)算過程是穩(wěn)定

34、的定的. .上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 例題表明穩(wěn)定性不但例題表明穩(wěn)定性不但與方法有關(guān)與方法有關(guān)，也，也與步長(zhǎng)與步長(zhǎng)h有有關(guān)，當(dāng)然關(guān)，當(dāng)然與方程中的與方程中的f(x, y)有關(guān)有關(guān). 為了只考察數(shù)值方為了只考察數(shù)值方法本身，通常只檢驗(yàn)數(shù)值方法用于解法本身，通常只檢驗(yàn)數(shù)值方法用于解模型方程模型方程的穩(wěn)的穩(wěn)定性，定性，模型方程模型方程為為)8 . 4(, yy 其中其中 為復(fù)數(shù)，這個(gè)方程分析較簡(jiǎn)單，對(duì)一般方程可為復(fù)數(shù)，這個(gè)方程分析較簡(jiǎn)單，對(duì)一般方程可以通過局部線性化化為這種形式，例如在以通過局部線性化化為這種形式，例如在(x, y)的鄰的鄰域，可展開為域，可展開為,)(,()

35、(,(),(),( yyyxfxxyxfyxfyxfyyx上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)略去高階項(xiàng)，再做變換即可得到略去高階項(xiàng)，再做變換即可得到u = u的形式的形式. 對(duì)于對(duì)于m個(gè)方程的常微分方程組個(gè)方程的常微分方程組, 可線性化為可線性化為y =Ay, 這里矩這里矩陣陣A為為mm雅可比矩陣雅可比矩陣(fi/yj)，若矩陣，若矩陣A有有m個(gè)特征個(gè)特征值值 1, 2, m，其中，其中 i可能是復(fù)數(shù)，所以，為了使?？赡苁菑?fù)數(shù)，所以，為了使模型方程結(jié)果能推廣到常微分方程組，方程型方程結(jié)果能推廣到常微分方程組，方程(4.8)式中式中 為復(fù)數(shù)為復(fù)數(shù). 為保證微分方程本身的穩(wěn)定性，還應(yīng)

36、假定為保證微分方程本身的穩(wěn)定性，還應(yīng)假定Re( )0.下面先研究歐拉方法的穩(wěn)定性下面先研究歐拉方法的穩(wěn)定性. 模型方程模型方程y = y的歐拉公式為的歐拉公式為)9 . 4(.)1(1nnyhy 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)設(shè)在節(jié)點(diǎn)設(shè)在節(jié)點(diǎn) yn 上有一擾動(dòng)值上有一擾動(dòng)值 n，它的傳播使節(jié)點(diǎn)值，它的傳播使節(jié)點(diǎn)值yn+1產(chǎn)生大小為的擾動(dòng)值產(chǎn)生大小為的擾動(dòng)值 n+1，假設(shè)用，假設(shè)用y*n=yn+ n，按歐拉，按歐拉公式得出公式得出 y*n+1=yn+1+ n+1 的計(jì)算過程不再有新的誤差，的計(jì)算過程不再有新的誤差，則擾動(dòng)值滿足則擾動(dòng)值滿足.)1(1nnh 可見擾動(dòng)值滿足原來的差

37、分方程可見擾動(dòng)值滿足原來的差分方程(4.9). 這樣，如果差這樣，如果差分方程的解是不增長(zhǎng)的，即有分方程的解是不增長(zhǎng)的，即有.1nnyy 則它就是穩(wěn)定的則它就是穩(wěn)定的. 這一論斷對(duì)于下面將要研究的其它這一論斷對(duì)于下面將要研究的其它方法同樣適用方法同樣適用.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)顯然，為要保證差分方程顯然，為要保證差分方程(4.9)的解是不增長(zhǎng)的，的解是不增長(zhǎng)的，只要選取只要選取h充分小，使充分小，使)10. 4(. 11 h在在 =h 的復(fù)平面上，這是以的復(fù)平面上，這是以(- -1,0)為圓心，為圓心，1為為半徑的單位圓內(nèi)部半徑的單位圓內(nèi)部. 稱為歐拉法的稱為歐拉法的

38、絕對(duì)穩(wěn)定域絕對(duì)穩(wěn)定域，一般，一般情形可由下面定義情形可由下面定義. 定義定義6 單步法單步法(4.1)用于解模型方程用于解模型方程y = y，若得，若得到的解到的解yn+1=E(h )yn，滿足，滿足|E(h )|1，則稱方法，則稱方法(4.1)是是絕對(duì)穩(wěn)定絕對(duì)穩(wěn)定的的. 在在 =h 的平面上的平面上, 使使|E(h )|1的變量的變量圍成的區(qū)域，稱為圍成的區(qū)域，稱為絕對(duì)穩(wěn)定區(qū)域絕對(duì)穩(wěn)定區(qū)域，它與實(shí)軸的交稱，它與實(shí)軸的交稱為為絕對(duì)穩(wěn)定區(qū)間絕對(duì)穩(wěn)定區(qū)間.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)對(duì)歐拉法對(duì)歐拉法E(h )=1+h ，其，其絕對(duì)穩(wěn)定域絕對(duì)穩(wěn)定域?yàn)闉閨1+h |1,絕對(duì)穩(wěn)定區(qū)間

39、絕對(duì)穩(wěn)定區(qū)間為為- -2h 0,在例在例5中中 =- -100,- -2- -100h0,即即0h2/100=0.02為為穩(wěn)定區(qū)間穩(wěn)定區(qū)間，在例，在例4中取中取h=0.025，故它是不穩(wěn)定的，當(dāng)取故它是不穩(wěn)定的，當(dāng)取h=0.005時(shí)它是穩(wěn)定的時(shí)它是穩(wěn)定的. 對(duì)二階對(duì)二階R- -K方法方法，解模型方程，解模型方程(4.1)可得到可得到,2)(121nnyhhy 故故.2)(1)(2 hhhE 絕對(duì)穩(wěn)定域絕對(duì)穩(wěn)定域由由|E(h )|1得到，于是可得得到，于是可得絕對(duì)穩(wěn)定區(qū)間絕對(duì)穩(wěn)定區(qū)間為為- -2h 0，即，即0h2/ .上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 類似可得三階及四階類似可得

40、三階及四階R- -K方法方法的的E(h )分別為分別為.!3)(!2)(1)(32 hhhhE .!4)(!3)(!2)(1)(432 hhhhhE 由由|E(h )|1可得到相應(yīng)的可得到相應(yīng)的絕對(duì)穩(wěn)定域絕對(duì)穩(wěn)定域. 當(dāng)當(dāng) 為實(shí)數(shù)時(shí)，為實(shí)數(shù)時(shí)，則得則得絕對(duì)穩(wěn)定區(qū)間絕對(duì)穩(wěn)定區(qū)間，它們分別為，它們分別為(圖形復(fù)雜圖形復(fù)雜p296) 三階顯式三階顯式R- -K方法方法： 四階顯式四階顯式R- -K方法方法：./51. 20, 051. 2 hh即即./78. 20, 078. 2 hh即即從以上討論可知顯式從以上討論可知顯式R- -K方法方法的的絕對(duì)穩(wěn)定域絕對(duì)穩(wěn)定域均為均為有限域，都對(duì)步長(zhǎng)有限域，都

41、對(duì)步長(zhǎng)h有限制有限制. 如果如果h不在所給的不在所給的絕對(duì)穩(wěn)絕對(duì)穩(wěn)定區(qū)間定區(qū)間內(nèi)，方法就不穩(wěn)定內(nèi)，方法就不穩(wěn)定.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 例例5 分別取分別取h=0.1及及h=0.2，用經(jīng)典的四階，用經(jīng)典的四階R- -K方方法法(3.13)計(jì)算計(jì)算初值問題初值問題 . 1)0()10(,20yxyy 解解 本例本例 =- -20, ,h 分別為分別為- -2及及- -4，前者在，前者在絕對(duì)絕對(duì)穩(wěn)定區(qū)間穩(wěn)定區(qū)間內(nèi)，后者則不在，用四階內(nèi)，后者則不在，用四階R- -K方法方法計(jì)算其計(jì)算其誤差見下表誤差見下表xn0.20.40.60.81.0h=0.1h=0.20.9310-

42、-14.980.1210- -125.00.1410- -2125.00.1510- -3625.00.1710- -43125.0從以上結(jié)果看到，如果步長(zhǎng)從以上結(jié)果看到，如果步長(zhǎng)h不滿足不滿足絕對(duì)穩(wěn)定條絕對(duì)穩(wěn)定條件件，誤差增長(zhǎng)很快，誤差增長(zhǎng)很快. .上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)對(duì)隱式單步法對(duì)隱式單步法, ,可以同樣討論方法的可以同樣討論方法的絕對(duì)穩(wěn)定性絕對(duì)穩(wěn)定性, ,例如對(duì)例如對(duì)后退歐拉法后退歐拉法，用它解模型方程可得，用它解模型方程可得,111nnyhy 故故.11)( hhE 由由|E(h )|1，這是以，這是以(1,0)為圓心，為圓心，1為半徑的單位圓外部為半徑的單

43、位圓外部. 故方法的故方法的絕對(duì)穩(wěn)定區(qū)絕對(duì)穩(wěn)定區(qū)間間為為- - h 0. 當(dāng)當(dāng) 0時(shí)，則時(shí)，則0h，即對(duì)任何步長(zhǎng)均，即對(duì)任何步長(zhǎng)均為穩(wěn)定的為穩(wěn)定的.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 對(duì)對(duì)隱式梯形法隱式梯形法，它用于解模型方程，它用于解模型方程(4.8)得得,21211nnyhhy 故故.2121)( hhhE 對(duì)對(duì)Re( )0有有|E(h )|1，故，故絕對(duì)穩(wěn)定域絕對(duì)穩(wěn)定域?yàn)闉?=h 的左半的左半平面，平面，絕對(duì)穩(wěn)定區(qū)間絕對(duì)穩(wěn)定區(qū)間為為- - h 0，即，即0h 時(shí)時(shí)隱式梯隱式梯形法形法均是穩(wěn)定的均是穩(wěn)定的.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)分析得到分析得到隱式歐

44、拉法隱式歐拉法與與梯形方法梯形方法的的絕對(duì)穩(wěn)定域絕對(duì)穩(wěn)定域均均為為h | Re(h )0，在具體計(jì)算中步長(zhǎng)，在具體計(jì)算中步長(zhǎng)h的選取只需考的選取只需考慮計(jì)算精度及迭代收斂性要求而不必考慮穩(wěn)定性，具慮計(jì)算精度及迭代收斂性要求而不必考慮穩(wěn)定性，具有這種特點(diǎn)的方法需特別重視有這種特點(diǎn)的方法需特別重視, 由此給出下面的定義由此給出下面的定義. 定義定義7 如果數(shù)值方法的絕對(duì)穩(wěn)定域包含了集合如果數(shù)值方法的絕對(duì)穩(wěn)定域包含了集合h | Re(h )0，那么稱此方法是，那么稱此方法是A- -穩(wěn)定穩(wěn)定的的. .由定義知由定義知A- -穩(wěn)定穩(wěn)定方法對(duì)步長(zhǎng)方法對(duì)步長(zhǎng)h沒有限制沒有限制. .上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下

45、頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)9.5 線性多步法線性多步法 在逐步推進(jìn)的求解過程中在逐步推進(jìn)的求解過程中, 計(jì)算計(jì)算yn+1之前事實(shí)上之前事實(shí)上已經(jīng)求出了一系列的近似值已經(jīng)求出了一系列的近似值y0,y1,yn, 如果充分利如果充分利用前面多步的信息來預(yù)測(cè)用前面多步的信息來預(yù)測(cè)yn+1, 則可以期望會(huì)獲得較則可以期望會(huì)獲得較高的精度高的精度. 這就是構(gòu)造所得這就是構(gòu)造所得線性多步法線性多步法的基本思想的基本思想. 構(gòu)造構(gòu)造多步法多步法的主要途徑基于數(shù)值積分方法和基的主要途徑基于數(shù)值積分方法和基于泰勒展開方法，前者可直接由方程于泰勒展開方法，前者可直接由方程(1.1)兩端積分兩端積分后利用插值求積公式得到

46、后利用插值求積公式得到. 本節(jié)主要介紹基于泰勒本節(jié)主要介紹基于泰勒展開的構(gòu)造方法展開的構(gòu)造方法.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)9.5.1 線性多步法的一般公式線性多步法的一般公式 如果計(jì)算如果計(jì)算yn+k時(shí)，除用時(shí)，除用yn+k- -1的值，還要用到的值，還要用到y(tǒng)n+i (i=0,1,k- -2)的值，則稱此方法為的值，則稱此方法為線性多步法線性多步法. 一般一般的的線性多步法公式線性多步法公式可表示為可表示為) 1 . 5(,010 kiinikiiniknfhyy 其中其中yn+1為為y(xn+1)的近似，的近似，fn+i=f(xn+i, yn+i), 這里這里xn+i

47、=xn+ih, i, i為常數(shù)為常數(shù), 0及及 0不全為零不全為零, 則稱則稱(5.1)為為線性線性k步法步法, 計(jì)算時(shí)需先給出前面計(jì)算時(shí)需先給出前面k個(gè)近似值個(gè)近似值y0,y1,yk- -1, 再由再由(5.1)逐次求出逐次求出 yk,yk+1,.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè) 如果如果 k=0，則，則(5.1)稱為稱為顯式顯式k步法步法，這時(shí)，這時(shí)yn+k可可直接由直接由(5.1)算出；如果算出；如果 k0, 則則(5.1)稱為稱為隱式隱式k步法步法,求解時(shí)與梯形法求解時(shí)與梯形法(2.7)相同相同, 要用迭代法方可算出要用迭代法方可算出yn+k. (5.1)中系數(shù)中系數(shù)

48、i及及 i可根據(jù)方法的局部截?cái)嗾`差及階確可根據(jù)方法的局部截?cái)嗾`差及階確定，其定義為定，其定義為 定義定義8 設(shè)設(shè)y(x)是初值問題是初值問題(1.1), (1.2)的準(zhǔn)確解，的準(zhǔn)確解，線性多步法線性多步法(5.1)在在xn+k上上局部截?cái)嗾`差局部截?cái)嗾`差為為);(hxyLTnkn ) 2 . 5(. )()()(010 kiinikiiniknxyhxyxy 若若Tn+k=O(hp+1)，則稱方法，則稱方法(5.1)是是p階的階的，如果，如果p 1, 則則稱方法稱方法(5.1)與微分方程與微分方程(1.1)是是相容的相容的.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)由定義由定義8，對(duì)，對(duì)

49、Tn+k在在xn處泰勒展開，由于處泰勒展開，由于,)(! 3)()(! 2)()()()(32 nnnnnxyihxyihxyihxyihxy.)(! 2)()()()(2 nnnnxyihxyihxyihxy代入代入(5.2)得得) 3 . 5(,)()()()()(2210 npppnnnknxyhcxyhcxyhcxycT上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)其中其中., 3 , 2)2()!1(1) 1(2(!11110121 qkqkkqckqqkqqqq ),(11100 kc ),() 1(2101211kkkkc ) 4 . 5( 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)

50、下頁(yè)下頁(yè)下頁(yè)若在公式若在公式(5.1)中選擇系數(shù)中選擇系數(shù) i 及及 i ，使它滿足，使它滿足. 0, 0110 ppcccc由定義可知此時(shí)所構(gòu)造的多步法是由定義可知此時(shí)所構(gòu)造的多步法是p階的，且階的，且) 5 . 5().()(2) 1(11 pnpppknhOxyhcT稱右端第一項(xiàng)為稱右端第一項(xiàng)為局部截?cái)嗾`差主項(xiàng)局部截?cái)嗾`差主項(xiàng), cp+1稱為稱為誤差常數(shù)誤差常數(shù).根據(jù)相容性定義根據(jù)相容性定義, p 1，即，即c0=c1=0，由，由(5.4)得得) 6 . 5(., 1011110 kikiikiik 故線性多步法故線性多步法(5.1)與微分方程與微分方程(1.1)相容的充分必要條相容的充

51、分必要條件是件是(5.6)成立成立.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)顯然，當(dāng)顯然，當(dāng)k=1時(shí)，若時(shí)，若 1=0，則由，則由(5.6)可求得可求得 0=1， 0=1.此時(shí)公式此時(shí)公式(5.1)為為,1nnnhfyy 即為歐拉法即為歐拉法. 從從(5.4)可求得可求得c2=1/20，故方法為，故方法為1階精階精度，且局部截?cái)嗾`差為度，且局部截?cái)嗾`差為).()(21321hOxyhTnn 這和第這和第2節(jié)給出的定義及結(jié)果是一致的節(jié)給出的定義及結(jié)果是一致的.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)對(duì)對(duì)k=1，若，若 10，此時(shí)方法為隱式公式，為了確，此時(shí)方法為隱式公式，為了確定

52、系數(shù)定系數(shù) 0, 0, 1, 可由可由c0=c1=c2=0解得解得 0=1, 0= 1=1/2. 于是得到公式于是得到公式).(211 nnnnffhyy即為梯形公式即為梯形公式.由由(5.4)可求得可求得c2=- -1/12，故，故p=2，所以梯形法是二，所以梯形法是二階方法，其局部截?cái)嗾`差主項(xiàng)是階方法，其局部截?cái)嗾`差主項(xiàng)是).(123nxyh 這與這與9.2節(jié)中討論是一致的節(jié)中討論是一致的.對(duì)對(duì)k 2的多步法公式都可利用的多步法公式都可利用(5.4)確定系數(shù)確定系數(shù) i, i,并由并由(5.5)式給出局部截?cái)嗾`差，下面只就若干常用的式給出局部截?cái)嗾`差，下面只就若干常用的多步法導(dǎo)出具體公式多

53、步法導(dǎo)出具體公式.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)9.5.2 阿當(dāng)姆斯顯式與隱式公式阿當(dāng)姆斯顯式與隱式公式p299自學(xué)自學(xué).9.5.3 米爾尼方法與辛普森方法米爾尼方法與辛普森方法p301自學(xué)自學(xué).9.5.4 漢明方法漢明方法p302自學(xué)自學(xué).9.5.5 預(yù)測(cè)預(yù)測(cè)-校正方法校正方法p303自學(xué)自學(xué).9.5.6 構(gòu)造多步法公式的注記和例構(gòu)造多步法公式的注記和例p305自學(xué)自學(xué).上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)9.6 線性多步法的收斂性與穩(wěn)定性線性多步法的收斂性與穩(wěn)定性p306自看自看.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)9.7 一階方程組與剛性方程組一

54、階方程組與剛性方程組9.7.1 一階方程組一階方程組 前面我們研究了單個(gè)方程前面我們研究了單個(gè)方程y = =f 的數(shù)值解，只要的數(shù)值解，只要把把y 和和f 理解為向量，那么，所提供的各種計(jì)算公式理解為向量，那么，所提供的各種計(jì)算公式即可應(yīng)用到一階方程組的情形即可應(yīng)用到一階方程組的情形. 考察考察一階方程組一階方程組), 2 , 1(),(21NiyyyxfyNii 的初值問題，初始條件給為的初值問題，初始條件給為)., 2 , 1(,)(00Niyxyii 上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)若采用向量的記號(hào)，記若采用向量的記號(hào)，記(向量向量),),(21TNyyyy ,),(00

55、2010TNyyyy .),(21TNffff 則上述方程組的初值問題可表示為則上述方程組的初值問題可表示為) 1 . 7(.)(),(00 yxyyxfy上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)求解這一初值問題的求解這一初值問題的四階龍格四階龍格- -庫(kù)塔公式庫(kù)塔公式為為(向量向量),22(643211kkkkhyynn ).,(),2,2(),2,2(),(3423121hkyhxfkkhyhxfkkhyhxfkyxfknnnnnnnn式中式中(向量向量)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)或表示為或表示為(分量分量), 2 , 1()22(643211,NiKKKKhyyiiiiinni 其中其中(分量分量) ).,(),2,2,2,2(),2,2,2,2(),(323213142222121312121112211NNnnnniiNNnnnniiNNnnnniiNnnnniihKyhKyhKyhxfKKhyKhyKhyhxfKKhyKhyKhyhxfKyyyxfK這里這里yin是第是第i個(gè)因變量個(gè)因變量yi(x)在節(jié)點(diǎn)在節(jié)點(diǎn)xn=x0+nh的近似值的近似值.上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)上頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)下頁(yè)為了幫助理解這一公式的計(jì)算過程，