常微分方程初值問題數(shù)值解法.ppt

第9章 常微分方程初值問題數(shù)值解法,9.1 引言 9.2 簡單的數(shù)值方法與基本概念 9.3 龍格-庫塔方法 9.4 單步法的收斂性與穩(wěn)定性 9.5 線性多步法 9.6 方程組和高階方程,9.1 引 言,科學(xué)技術(shù)中常常需要求解常微分方程的定解問題. 這類問題最簡單的形式，是本章將著重考察的一階方程的初值問題,我們知道，只有f(x, y)適當(dāng)光滑譬如關(guān)于y滿足利普希茨(Lipschitz)條件,理論上就可以保證初值問題的解yf(x)存在并且唯一.,雖然求解常微分方程有各種各樣的解析方法，但解析方法只能用來求解一些特殊類型的方程，實際問題中歸結(jié)出來的微分方程主要靠數(shù)值解法.,所謂數(shù)值解法, 就是尋求解y(x)在一系列離散節(jié)點,上的近似值 y1,y2,yn,yn+1,. 相鄰兩個節(jié)點的間距hn=xn+1-xn稱為步長. 今后如不特別說明，總是假定 hi=h(i=1,2,)為定數(shù), 這時節(jié)點為xn=x0+nh(i=0,1,2,) (等距節(jié)點).,初值問題的數(shù)值解法有個基本特點，他們都采取“步進式”，即求解過程順著節(jié)點排列的次序一步一步地向前推進. 描述這類算法，只要給出用已知信息yn,yn-1,yn-2,計算yn+1的遞推公式.,首先，要對微分方程離散化，建立求解數(shù)值解的遞推公式. 一類是計算yn+1時只用到前一點的值yn，稱為單步法. 另一類是用到y(tǒng)n+1前面 k 點的值yn,yn-1, yn-k+1，稱為k步法. 其次，要研究公式的局部截斷誤差和階，數(shù)值解yn與精確解y(xn)的誤差估計及收斂性，還有遞推公式的計算穩(wěn)定性等問題.,9.2 簡單的數(shù)值方法與基本概念,9.2.1 歐拉法與后退歐拉法,我們知道，在xy平面上，微分方程(1.1)式的解y=f(x)稱作它的積分曲線，積分曲線上一點(x, y)的切線斜率等于函數(shù)f(x, y)的值. 如果按f(x, y)在xy平面上建立一個方向場，那么，積分曲線上每一點的切線方向均與方向場在該點的方向相一致.,基于上述幾何解釋，我們從初始點P0(x0, y0)出發(fā),先依方向場在該點的方向推進到x=x1上一點P1，然后再從P1點依方向場在該點的方向推進到 x=x2 上一點P2 , 循環(huán)前進做出一條折線P0 P1 P2.,一般地，設(shè)已做出該折線的頂點Pn,過Pn(xn, yn)依方向場的方向再推進到Pn+1(xn+1, yn+1)，顯然兩個頂點Pn,Pn+1的坐標有關(guān)系,這就是著名的(顯式)歐拉(Euler)公式. 若初值y0已知，則依公式(2.1)可逐次逐步算出各點數(shù)值解.,即,例1 用歐拉公式求解初值問題,解 取步長h=0.1，歐拉公式的具體形式為,其中xn=nh=0.1n (n=0,1,10), 已知y0 =1, 由此式可得,依次計算下去，部分計算結(jié)果見下表.,與準確解 相比，可看出歐拉公式的計算結(jié)果精度很差.,歐拉公式具有明顯的幾何意義, 就是用折線近似代替方程的解曲線，因而常稱公式(2.1)為歐拉折線法.,還可以通過幾何直觀來考察歐拉方法的精度.假設(shè)yn=y(xn),即頂點Pn落在積分曲線y=y(x)上，那么，,按歐拉方法做出的折線PnPn+1便是y=y(x)過點Pn的切線.從圖形上看,這樣定出的頂點Pn+1顯著地偏離了原來的積分曲線，可見歐拉方法是相當(dāng)粗糙的.,為了分析計算公式的精度，通常可用泰勒展開將y(xn+1)在xn處展開，則有,在yn=y(xn)的前提下，f(xn,yn )=f(xn,y(xn)=y(xn).于是可得歐拉法(2.1)的公式誤差為,稱為此方法的局部截斷誤差.,如果對方程(1.1)從xn到xn+1積分，得,右端積分用左矩形公式hf(xn,y(xn)近似，再以yn代替y(xn)，yn+1代替y(xn+1)也得到歐拉公式(2.1)，局部截斷誤差也是(2.3).,稱為(隱式)后退的歐拉公式.,如果右端積分用右矩形公式hf(xn+1,y(xn+1)近似，則得到另一個公式,后退的歐拉公式與歐拉公式有著本質(zhì)的區(qū)別, 后者是關(guān)于yn+1的一個直接計算公式，這類公式稱作是顯式的；前者公式的右端含有未知的yn+1，它實際上是關(guān)于yn+1的一個函數(shù)方程,這類方程稱作是隱式的.,顯式與隱式兩類方法各有特點，考慮到數(shù)值穩(wěn)定性等其他因素，人們有時需要選用隱式方法，但使用顯式算法遠比隱式方便.,隱式方程通常用迭代法求解，而迭代過程的實質(zhì)是逐步顯式化.,設(shè)用歐拉公式,給出迭代初值 ，用它代入(2.5)式的右端，使之轉(zhuǎn)化為顯式，直接計算得,然后再用 代入(2.5)式，又有,如此反復(fù)進行，得,由于f(x, y)對y滿足Lipschitz條件(1.3). 由(2.6)減(2.5)得,由此可知，只要hL1，迭代法(2.6)就收斂到解.關(guān)于后退歐拉方法的公式誤差，從積分公式看到它與歐拉法是相似的.,9.2.2 梯形方法,為得到比歐拉法精度高的計算公式，在等式(2.4) 右端積分用梯形求積公式近似, 并用yn代替y(xn), yn+1代替y(xn+1)，則得,稱為矩形方法.,矩形方法是隱式單步法，用迭代法求解，同后退的歐拉方法一樣，仍用歐拉法提供迭代初值，則矩形迭代公式為,為了分析迭代過程的收斂性, 將(2.7)與(2.8)相減, 得,于是有,使得,則當(dāng)k時有 , 這說明迭代過程(2.8)是收斂的.,9.2.3 單步法的局部截斷誤差與階,初值問題(1.1),(1.2)的單步法可用一般形式表示為,其中多元函數(shù)與f(x, y )有關(guān)，當(dāng)含有yn+1時，方法是隱式的，若不含yn+1則為顯式方法，所以顯式單步法可表示為,(x, y, h)稱為增量函數(shù)，例如對歐拉法(2.1)有,它的局部截斷誤差已由(2.3)給出, 對一般顯式單步法則可如下定義.,定義1 設(shè)y(x)是初值問題(1.1),(1.2)的準確解, 稱,為顯式單步法(2.10)的局部截斷誤差.,Tn+1之所以稱為局部的，是假設(shè)在xn前各步?jīng)]有誤差.當(dāng)yn=y(xn)時，計算一步，則有,所以，局部截斷誤差可理解為用方法(2.10)計算一步的誤差，也即公式(2.10)中用準確解y(x)代替數(shù)值解產(chǎn)生的公式誤差. 根據(jù)定義, 顯然歐拉法的局部截斷誤差為,即為(2.3)的結(jié)果. 這里 稱為局部截斷誤差主項. 顯然Tn+1=O(h2). 一般情形的定義如下,定義2 設(shè)y (x)是初值問題的準確解，若存在最大整數(shù)p使顯式單步法(2.10)的局部截斷誤差滿足,則稱方法(2.10)具有p階精度.,若將(2.10)展開式寫成,則 稱為局部截斷誤差主項.,以上定義對隱式單步法(2.9)也是適用的.例如，對后退歐拉法(2.5)其局部截斷誤差為,這里p=1是1階方法，局部截斷誤差主項為,同樣對梯形法(2.7)有,所以梯形方法(2.7)是二階的. 其局部截斷誤差主項為,9.2.4 改進的歐拉公式,我們看到，梯形方法雖然提高了精度，但其算法復(fù)雜，在應(yīng)用迭代公式(2.9)進行實際計算時，每迭代一次，都要重新計算函數(shù)f(x, y )的值，而迭代又要反復(fù)進行若干次，計算量很大，而且往往難以預(yù)測. 為了控制計算量，通常只迭代一兩次就轉(zhuǎn)入下一步的計算，這就簡化了算法.,具體地說，我們先用歐拉公式求得一個初步的近似值 ，稱之為預(yù)測值，此預(yù)測值 的精度可能很差，再用梯形公式(2.7)將它校正一次，即按(2.8)式迭代一次，這個結(jié)果稱之為校正值.,這樣建立的預(yù)測校正系統(tǒng)通常稱為改進的歐拉公式:,或表為下列平均化形式,(2.13),預(yù)測,校正,例2 用改進的歐拉法解例1中的初值問題(2.2).,解 仍取步長h=0.1,改進的歐拉公式為,部分計算結(jié)果見下表,同例1中的歐拉法的計算結(jié)果比較，明顯改善了精度.,例 (兩種方法的精度比較),用歐拉公式和改進的歐拉公式解下述初值問題，并與其準確解y=e-x+x進行比較.,解 取步長h=0.1，xk=kh(k=0,1,6).用兩種方法進行計算對應(yīng)結(jié)果及絕對誤差見下表,9.3 龍格庫塔方法,對許多實際問題來說，歐拉公式與改進歐拉公式精度還不能滿足要求，為此從另一個角度來分析這兩個公式的特點，從而探索一條構(gòu)造高精度方法的途徑.,9.3.1 顯式龍格庫塔法的一般形式,上節(jié)給出了顯式單步法的表達式(2.10), 其局部截斷誤差為O(hp+1)，對歐拉法Tn+1=O(h2)，即方法為p=1階，若用改進歐拉法(2.13)，它可表為,此時增量函數(shù)為,它比歐拉法的(xn, yn, h)=f(xn, yn), 增加了計算一個右函數(shù)f 的值，可望 p=2.若要使得到的公式階數(shù)p更大, 就必須包含更多的f 值. 實際上從方程(1.1)等價的積分形式(2.4) ，即,若要使公式階數(shù)提高，就必須使右端積分的數(shù)值求積公式精度提高，它必然要增加求積節(jié)點，為此可將(3.3)的右端用求積公式表示為,一般說來，點數(shù)r 越多，精度越高，上式右端相當(dāng)于增量函數(shù)(xn, yn, h)，為得到便于計算的顯式方法，可類似于改進歐拉法(3.1),(3.2)，將公式表示為,其中,這里ci, i, ij均為常數(shù). (3.4)和(3.5)稱為r級顯式龍格-庫塔(Runge-Kutta)法, 簡稱R-K方法.,當(dāng)r=1, (xn, yn, h)=f(xn, yn)時，就是歐拉法，此時方法的階為p=1. 當(dāng)r=2時，改進歐拉法(3.1),(3.2)是其中一種，下面將證明其階p=2. 要使公式(3.4),(3.5)具有更高的階p，就要增加點數(shù)r. 下面我們只就r=2推導(dǎo)R-K方法. 并給出 r=3,4 時的常用公式，其推導(dǎo)方法與r=2時類似，只是計算較復(fù)雜.,9.3.2 二階顯式R-K方法,對r=2的R-K方法，由(3.4),(3.5)式可得如下計算公式,這里 c1, c2, 2, 21 均為待定常數(shù)，我們希望適當(dāng)選取這些系數(shù)，使公式階數(shù) p 盡量高. 根據(jù)局部截斷誤差定義，推導(dǎo)出(3.6)的局部截斷誤差為,其中,這里yn=y(xn), yn+1=y(xn+1). 為得到Tn+1的階p，要將上式各項在(xn, yn)處做泰勒展開，由于f(x, y )是二元函數(shù)，故要用二元泰勒展開，各項展開式為,將以上結(jié)果代入(3.7)，則有,要使公式(3.6)具有p=2階，必須使,即,(3.9)的解是不唯一的. 可令c2=a0，則得,這樣得到的公式稱為二階R-K方法.,則由此可以看出在改進的歐拉公式中相當(dāng)于取(xn,yn), (xn+1,yn+1)兩點處斜率的平均值，近似代替平均斜率，其精度比歐拉公式提高了.,如取a=1/2，則c1= c2=1/2, 2=21=1. 這就是改進的歐拉公式(3.1).,稱為中點公式(變形的歐拉公式)，相當(dāng)于數(shù)值積分的中矩形公式.也可以表示為,如取a=1，則c1=0, c2=1, 2=21=1/2. 得計算公式,對r=2的R-K公式(3.6)能否使局部誤差提高到O(h4)? 為此 需把K2多展開一項，從(3.8)的 看到展開式中的項 是不能通過選擇參數(shù)削掉的，實際上要使 h3 的項為零，需增加3個方程，要確定4個參數(shù)c1, c2, 2, 21，這是不可能的. 故r2的顯式R-K方法的階只能是p=2，而不能得到三階公式.,9.3.3 三階與四階顯式R-K方法,要得到三階顯式R-K方法，必須r=3. 此時計算(3.4), (3.5)的公式表示為,其中c1, c2, c3及2, 21, 3, 31, 32均為待定常數(shù)，公式(3.11)的局部截斷誤差為,只要K1, K2將按二元泰勒展開，使Tn+1O(h4)，可得待定參數(shù)滿足方程,這是8個未知數(shù)6個方程的方程組，解不是唯一的. 可以得到很多公式. 滿足條件(3.12)的公式(3.11)統(tǒng)稱為三階R-K公式. 下面只給出其中一個常見的公式.,此公式稱為庫塔三階方法.,繼續(xù)上述過程，經(jīng)過較復(fù)雜的數(shù)學(xué)演算，可以導(dǎo)出各種四階R-K公式，下列經(jīng)典公式是其中常用的一個：,四階R-K方法的每一步需要計算四次函數(shù)值f，可以證明其局部截斷誤差為O(h5). (例3見書p349),然而值得指出的是，龍格-庫塔方法的推導(dǎo)基于泰勒展開方法，因而它要求所求的解具有較好的光滑性質(zhì). 反之，如果解的光滑性差，那么，使用龍格-庫塔方法求得的數(shù)值解，其精度可能反而不如改進的歐拉方法. 實際計算時，我們應(yīng)當(dāng)針對問題的具體特點選擇合適的算法.,9.3.4 變步長的龍格-庫塔方法,單從每一步看，步長越小，截斷誤差就越小，但隨著步長的縮小，在一定求解范圍內(nèi)所要完成的步數(shù)就增加了. 步數(shù)的增加不但引起計算量的增大，而且可能導(dǎo)致舍入誤差的嚴重積累. 因此同積分的數(shù)值計算一樣，微分方程的數(shù)值解法也有個選擇步長的問題.,在選擇步長時，需要考慮兩個問題： 1. 怎樣衡量和檢驗計算結(jié)果的精度？ 2. 如何依據(jù)所獲得的精度處理步長？,我們考察四階R-K公式(3.13) ，從節(jié)點xn出發(fā)，先以h為步長求出一個近似值，記為 ，由于公式的局部截斷誤差為O(h5)，故,然后將步長折半，即取為步長 ,從xn跨兩步到xn+1，再求得一個近似值 ，每跨一步的局部截斷誤差是 ，因此有,比較(3.14)式和(3.15)式我們看到，步長折半后，誤差大約減少到1/16，即有,由此易得下列事后估計式,這樣，我們可以通過檢查步長，折半前后兩次計算結(jié)果的偏差,來判定所選的步長是否合適，具體地說，將區(qū)分以下兩種情況處理：,這種通過加倍或折半處理步長的方法稱為變步長方法.表面上看，為了選擇步長，每一步的計算量增加了，但總體考慮往往是合算的.,1.對于給定的精度,如果,我們反復(fù)將步長折半計算,直至為止,這時取最終得到的 作為結(jié)果；,2.如果為止，這時再將步長折半計算一次，就得到所要的結(jié)果.,9.4 單步法的收斂性與穩(wěn)定性,9.4.1 收斂性與相容性,數(shù)值解法的基本思想是，通過某種離散化手段將微分方程轉(zhuǎn)化為差分方程，如單步法,它在點xn處的解為yn，而初值問題在點xn處的精確解為y(xn)，記en=y(xn)-yn稱為整體截斷誤差. 收斂性就是討論當(dāng) x=xn 固定且 時en0的問題.,定義3 若一種數(shù)值方法對于固定的xn=x0+nh, 當(dāng)h0時有yny(xn)，其中y(x)是(1.1),(1.2)的準確解，則稱該方法是收斂的.,顯然數(shù)值方法收斂是指en=y(xn)-yn0，對單步法(4.1)有下述收斂性定理：,定理1 假設(shè)單步法(4.1)具有p階精度，且增量函數(shù)(x, y ,h)關(guān)于y滿足利普希次條件,又設(shè)初值y0是準確的, 即y0=f(x0), 則其整體截斷誤差,證明 設(shè)以yn+1表示取yn=y(xn)用公式(4.1)求得的結(jié)果，即,則y(xn)-yn+1為局部截斷誤差，由于所給方法具有p階精度，按定義2，存在定數(shù)C ，使,又由式(4.4)與(4.1)，得,利用利普希次條件(4.2)，有,從而有,即對整體截斷誤差en=y(xn)-yn成立下列遞推關(guān)系式,據(jù)此不等式反復(fù)遞推，可得,由此可以斷定，如果初值是準確的，即e0=0，則(4.3)式成立. 定理證畢.,再注意到當(dāng)x=x0+nhT時,最終得下列估計式,依據(jù)這一定理，判斷單步法(4.1)的收斂性，歸結(jié)為驗證增量函數(shù)能否滿足利普希次條件(4.2).,對于歐拉方法，由于其增量函數(shù) 就是f(x, y), 故當(dāng)f(x, y)關(guān)于y滿足利普希次條件時它是收斂的.,再考察改進的歐拉方法，其增量函數(shù)已由(3.2)式給出，假定f(x, y)關(guān)于y滿足利普希次條件，即,這時有,即,設(shè)限定步長hh0(h0為定數(shù))，上式表明關(guān)于y的利普希次常數(shù)為,因此改進的歐拉方法也是收斂的.,類似地, 不難驗證其它龍格-庫塔方法的收斂性.,定理1表明p1時單步法收斂, 并且當(dāng)y(x)是初值問題(1.1),(1.2)的解, (4.1)具有p階精度時, 則有展開式,所以p1的充分必要條件是 ，而 ，于是可給出如下定義：,定義4 若單步法(4.1)的增量函數(shù)滿足,以上討論表明p階方法(4.1)當(dāng)p1時與(1.1), (1.2)相容，反之相容方法至少是1階的.,于是由定理1可知方法(4.1)收斂的充分必要條件是此方法是相容的.,則稱單步法(4.1)與初值問題(1.1),(1.2)相容.,9.4.2 絕對穩(wěn)定性與絕對穩(wěn)定域,前面關(guān)于收斂性的討論有個前提，必須假定數(shù)值方法本身的計算是準確的. 實際情形并不是這樣，差分方程的求解還會有計算誤差. 譬如由于數(shù)字舍入而引起的小擾動. 這類小擾動在傳播過程中會不會惡性增長，以至于“淹沒”了差分方程的“真解”呢？這就是差分方程的穩(wěn)定性問題. 在實際計算時，我們希望某一步產(chǎn)生的擾動值，在后面的計算中能夠被控制，甚至是逐步衰減的.,定義5 若一種數(shù)值方法在節(jié)點值yn上大小為的擾動，于以后各節(jié)點值ym(mn)上產(chǎn)生的偏差均不超過，則稱該方法是穩(wěn)定的.,下面以歐拉法為例考察計算穩(wěn)定性.,例4 用歐拉公式求解初值問題,解 用歐拉法解方程y=-100y 得,其準確解 是一個按指數(shù)曲線衰減很快的函數(shù).,若取步長h=0.025，則歐拉公式的具體形式為,計算結(jié)果見表, 明顯計算過程不穩(wěn)定, 但取h=0.005, yn+1=0.5yn, 則計算過程穩(wěn)定.,對后退的歐拉公式，取h=0.025時，則計算公式為yn+1=(1/3.5)yn .計算結(jié)果見表, 這時計算過程是穩(wěn)定的.,例題表明穩(wěn)定性不但與方法有關(guān)，也與步長h有關(guān)，當(dāng)然與方程中的f(x, y)有關(guān). 為了只考察數(shù)值方法本身，通常只檢驗數(shù)值方法用于解模型方程的穩(wěn)定性，模型方程為,其中為復(fù)數(shù)，這個方程分析較簡單，對一般方程可以通過局部線性化化為這種形式，例如在(x, y)的鄰域，可展開為,略去高階項，再做變換即可得到的形式 u=u. 對于m個方程的方程組, 可線性化為y=Ay, 這里A為mm雅可比矩陣(fi/yj)，若A有m個特征值1,2,m，其中i可能是復(fù)數(shù)，所以，為了使模型方程結(jié)果能推廣到方程組，方程(4.8)中為復(fù)數(shù). 為保證微分方程本身的穩(wěn)定性，還應(yīng)假定Re()0.,下面先研究歐拉方法的穩(wěn)定性. 模型方程y=y的歐拉公式為,設(shè)在節(jié)點yn上有一擾動值n，它的傳播使節(jié)點值yn+1產(chǎn)生大小為的擾動值n+1，假設(shè)用yn*=yn+n，按歐拉公式得出yn+1*=yn+1+n+1的計算過程不再有新的誤差，則擾動值滿足,可見擾動值滿足原來的差分方程(4.9). 這樣，如果差分方程的解是不增長的，即有,則它就是穩(wěn)定的. 這一論斷對于下面將要研究的其它方法同樣適用.,顯然，為要保證差分方程(4.9)的解是不增長的，只要選取h充分小，使,在=h的復(fù)平面上，這是以(-1,0)為圓心，1為半徑的單位圓. 稱為歐拉法的絕對穩(wěn)定域，一般情形可由下面定義.,定義6 單步法(4.1)用于解模型方程y=y，若得到的解yn+1=E(h)yn，滿足|E(h)|1，則稱方法(4.1)是絕對穩(wěn)定的. 在=h的平面上, 使|E(h)|1的變量圍成的區(qū)域，稱為絕對穩(wěn)定區(qū)域，它與實軸的交稱為絕對穩(wěn)定區(qū)間.,對歐拉法E(h)=1+h，其絕對穩(wěn)定域為|1+h|1,絕對穩(wěn)定區(qū)間為-20,在例5中=-100,-2-100h0,即0h2/100=0.02為穩(wěn)定區(qū)間，在例4中取h=0.025，故它是不穩(wěn)定的，當(dāng)取h=0.005時它是穩(wěn)定的.,對二階R-K方法，解模型方程(4.1)可得到,故,絕對穩(wěn)定域由|E(h)|1得到，于是可得絕對穩(wěn)定區(qū)間為-2h0，即0h2/.,類似可得三階及四階R-K方法的E(h)分別為,由|E(h)|1可得到相應(yīng)的絕對穩(wěn)定域. 當(dāng)為實數(shù)時，則得絕對穩(wěn)定區(qū)間，它們分別為,三階顯式R-K方法：,四階顯式R-K方法：,從以上討論可知顯式R-K方法的絕對穩(wěn)定域均為有限域，都對步長h有限制. 如果h不在所給的絕對穩(wěn)定區(qū)間內(nèi)，方法就不穩(wěn)定.,例4 分別取h=0.1及h=0.2，用經(jīng)典的四階R-K方法(3.1)計算初值問題,解 本例=-20,h分別為-2及-4，前者在絕對穩(wěn)定區(qū)間內(nèi)，后者則不在，用四階R-K方法計算其誤差見下表,從以上結(jié)果看到，如果步長h不滿足絕對穩(wěn)定條件，誤差增長很快.,對隱式單步法,可以同樣討論方法的絕對穩(wěn)定性,例如對后退歐拉法，用它解模型方程可得,故,由|E(h)|1，這是以(1,0)為圓心，1為半徑的單位圓外部. 故方法的絕對穩(wěn)定區(qū)間為-h0. 當(dāng)0時，則0h，即對任何步長均為穩(wěn)定的.,對隱式梯形法，它用于解模型方程(4.8)得,故,對Re()0有|E(h)|1，故絕對穩(wěn)定域為=h的左半平面，絕對穩(wěn)定區(qū)間為-h0，即0h時隱式梯形法均是穩(wěn)定的.,9.5 線性多步法,在逐步推進的求解過程中，計算yn+1之前事實上已經(jīng)求出了一系列的近似值y0,y1,yn，如果充分利用前面多步的信息來預(yù)測yn+1，則可以期望會獲得較高的精度. 這就是構(gòu)造所得線性多步法的基本思想.,構(gòu)造多步法的主要途徑基于數(shù)值積分方法和基于泰勒展開方法，前者可直接由方程(1.1)兩端積分后利用插值求積公式得到. 本節(jié)主要介紹基于泰勒展開的構(gòu)造方法.,9.5.1 線性多步法的一般公式,如果計算yn+k時，除用yn+k-1的值，還要用到y(tǒng)n+i (i=0,1,k-2)的值，則稱此方法為線性多步法. 一般的線性多步法公式可表示為,其中yn+1為y(xn+1)的近似，fn+i=f(xn+i, yn+i), 這里xn+i=xn+ih，i, i為常數(shù)， 0及0不全為零，則稱(5.1)為線性k步法，計算時需先給出前面k個近似值y0,y1,yk-1，再由(5.1)逐次求出yk,yk+1,.,如果k=0，則(5.1)稱為顯式k步法，這時yn+k可直接由(5.1)算出；如果k0, 則(5.1)稱為隱式k步法,求解時與梯形法(2.7)相同, 要用迭代法方可算出yn+k. (5.1)中系數(shù)i及i可根據(jù)方法的局部截斷誤差及階確定，其定義為,定義7 設(shè)y(x)是初值問題(1.1), (1.2)的準確解，線性多步法(5.1)在xn+k上局部截斷誤差為,若Tn+k=O(hp+1)，則稱方法(5.1)是p階的，p1則稱方法(5.1)與方程(1.1)是相容的.,由定義7，對Tn+k在xn處泰勒展開，由于,代入(5.2)得,其中,若在公式(5.1)中選擇系數(shù)i及i，使它滿足,由定義可知此時所構(gòu)造的多步法是p階的，且,稱右端第一項為局部截斷誤差主項, cp+1稱為誤差常數(shù).,根據(jù)相容性定義, p1，即c0=c1=0，由(5.4)得,故線性多步法(5.1)與微分方程(1.1)相容的充分必要條件是(5.6)成立.,顯然，當(dāng)k=1時，若1=0，則由(5.6)可求得,0=1，0=1.,此時公式(5.1)為,即為歐拉法. 從(5.4)可求得c2=1/20