ch常微分方程初值問題的數(shù)值解法PPT課件

1、第六章 常微分方程初值問題的數(shù)值解法6.1 歐拉方法 6.2 計算公式的誤差分析 6.3 龍格庫塔方法 6.4 向一階方程組與高階方程的推廣第1頁/共20頁問題：) 15(,)(),(000bxxyxyyxfy求解數(shù)值求解方法：的近似值的函數(shù)值然后逐個求解出節(jié)點(diǎn)上生成節(jié)點(diǎn)作等距分割對區(qū)間)(,),(),(,), 1 , 0(,21000NixyxyxyNxbhniihxxbx第2頁/共20頁 歐拉公式與改進(jìn)歐拉公式歐拉公式與改進(jìn)歐拉公式)(,()()(:)(,()()()(11iiiiiiiiixyxfhxyxyxyxfxyhxyxy得由)26()1, 1 ,0()(),(001Nixyyyx

2、fhyyiiii由此可建立計算公式6.1 歐拉方法歐拉方法這稱為歐拉公式歐拉公式例例6.1 以 h=0.1為步長，用歐拉法求常微分方程初值問題第3頁/共20頁)(,()()(:)(,()()()(,1111111iiiiiiiiixyxfhxyxyxyxfxyhxyxy得由類似地)36()1, 1 ,0()(),(00111Nixyyyxfhyyiiii由此可建立另一公式后退歐拉公式是一個隱式公式，通常采用迭代法求解。這稱為后退歐拉公式后退歐拉公式第4頁/共20頁 梯形公式與改進(jìn)歐拉公式梯形公式與改進(jìn)歐拉公式歐拉公式與后退歐拉公式也可采用積分近似的方法推出)46()(,()()()(,()()

3、()(),(11111iiiiiixxiixxxxiidxxyxfxyxydxxyxfdxxyxyxyyxfy化為積分方程問題從而把微分方程問題轉(zhuǎn)得對微分方程積分從而可導(dǎo)出歐拉公式得由)(,()()()(,()(,()(,(111iiiiiixxiixxxyxfhxyxyxyxfhdxxyxfdxxyxfiiii第5頁/共20頁式從而可導(dǎo)出后退歐拉公得由)(,()()()(,()(,()(,(111111111iiiiiixxiixxxyxfhxyxyxyxfhdxxyxfdxxyxfiiii)56()(),(),(2)(,()(,(2)()()(,()(,(2)(,(00111111111x

4、yyyxfyxfhyyxyxfxyxfhxyxyxyxfxyxfhdxxyxfiiiiiiiiiiiiiiiixxii從而導(dǎo)出梯形公式得由梯形公式梯形公式也是隱式單步法公式第6頁/共20頁11)(1)1(1)(1)(,|,iikikikiyxyyyyk的近似值值作為函數(shù)取時使得迭代計算到某個用梯形公式計算時，通常取歐拉公式的解作為迭代初值進(jìn)行迭代計算，即采用下式), 1 ,0(),(),(2),()(11)1(1)0(1kyxfyxfhyyyxfhyykiiiiikiiiii) 1, 1 , 0(),(),(2),(,11111)1(1Niyxfyxfhyyyxfhyyyyiiiiiiiiii

5、ii則得到顯式計算公式代替取直接簡單地這稱為改進(jìn)歐拉公式改進(jìn)歐拉公式第7頁/共20頁例例 仍取步長h ，采用改進(jìn)歐拉法重新計算例 6.1 的常微分方程初值問題。這時改進(jìn)歐拉公式為計算結(jié)果見表6-2（書125頁）解解第8頁/共20頁6.2 計算公式的誤差分析計算公式的誤差分析 定義定義 若 yi+1 是 yi=y(xi) 從計算得到的近似解，則稱y(xi+1) yi+1為所用公式的局部截斷誤差局部截斷誤差，簡稱為截斷誤差截斷誤差。 定理定理6.1 若單步法 yi+1 = yi+h (xi , yi , h) 的局部截斷誤差為 O (h p+1) ，且增量函數(shù) (x , y , h) 關(guān)于 y 滿

6、足李普希茲條件，即存在常數(shù) L0，使對 成立不等式y(tǒng)y,| ),(),(|yyLhyxhyx則其整體截斷誤差 y(xi) yi=O(hp) 第9頁/共20頁截斷誤差的估計(基本假設(shè)： yi = y( xi ) )設(shè) y(x)C 3 x0 , b , 則 )76()()(2)()()()(321 hOxyhxyhxyhxyxyiiiii（1）對歐拉公式，有)86()()()(2)()()()(,()(232111 hOhOxyhyxyxyhxyxyxfhxyyiiiiiiiii故因此，歐拉公式的局部截斷誤差為 O (h2)（2）對后退歐拉公式，有)96()()()(2)()(,()(232111

7、11 hOhOxyhyxyxyxfhxyyiiiiiii故因此，后退歐拉公式的局部截斷誤差為 O (h2)第10頁/共20頁（3）對梯形公式，注意到其公式可改寫為),(),(21111iiiiiiiyxfhyyxhfyy故由式（6-9）和（6-9）得)106()()()(2)()(221),()(),()(21)(33232111111 hOhOxyhhOxyhyxfhyxyyxfhyxyyxyiiiiiiiiiiii因此，梯形公式的局部截斷誤差為 O ( h3 )第11頁/共20頁（4）對改進(jìn)歐拉公式，有),(),(2)()(1111iiiiiiiiyxfyxfhyxyhxyy而由 ，故有)

8、,(),(),(yxfyyxfyyxfyyx 得)()(2)()()()()()()(,()()(,(),(32122)(,(11hOxyhxyhxyyhOxyhxyhOfyhfhxyxfxyhxyhxfyxfiiiiiixyxyxiiiiiiiii 所以與式（6-7）比較得 y(xi+1) yi+1 = O ( h3 ) 因此，改進(jìn)歐拉公式的局部截斷誤差為 O ( h3 ) 第12頁/共20頁 定義定義 若一種求解常微分方程初值問題的數(shù)值計算方法的局部截斷誤差為 O ( hp+1 ) ，則稱該方法為 p階精度階精度，或稱該方法為 p階方法階方法。 由此定義知，歐拉方法與后退歐拉方法為一階精度

9、，梯形法與改進(jìn)歐拉方法為二階精度。第13頁/共20頁6.3 龍格龍格-庫塔方法庫塔方法由中值定理，有),(,)(,()()()()(111iiiiiixxyfhyxxxyxy 因此，以上介紹的各種單步法本質(zhì)上都是對平均斜率 f( , y( ) 進(jìn)行近似，龍格-庫塔據(jù)之提出了適當(dāng)選取若干點(diǎn)上的斜率值作近似以構(gòu)造高精度計算公式的方法，其基本思想是基于泰勒展式的待定系數(shù)法。第14頁/共20頁 二階二階R-KR-K公式公式問題問題：建立二階精度的計算格式形為在 y(xi) = yi 的假設(shè)下，有)()()()(,()(,()()(,(2)(,(2)(,(1121hOfybfahxyhOfbhKfahx

10、yxfbhKxyahxfKxyxyxfKiiiixyxyxixyxyxiiiiiii)()()()(3)(,(222211hOfybfahxyhxyyiixyxyxiii故解解第15頁/共20頁)()()()(3)(,(222211hOfybfahxyhxyyiixyxyxiii而)(2121)()()()()2/()()()(3)(,(2321hOfyfhxyhxyhOxyhxyhxyxyiixyxyxiiiiii 根據(jù)格式為二階精度，即 y(xi+1) yi+1 = O(h3) 比較兩式系數(shù)得)136(2/12/112221ba 系數(shù)滿足(6-13)的形為(6-12)計算格式統(tǒng)稱為二階R-

11、K公式。當(dāng)令1=1/2時，解得 2=1/2 ，a=b=1，即為改進(jìn)歐拉公式。若令 1=0，解得 2=1，a=b=1/2，則得另一計算公式)146()2/, 2/(),(12121hKyhxfKyxfKhKyyiiiiii變形歐拉公式變形歐拉公式第16頁/共20頁 四階四階 R-K R-K 公式公式每一步需計算的 f 值的個數(shù)1234567n8精度階1234456n-2 1965年，Butcher研究發(fā)現(xiàn)顯式R-K公式的精度與需要組合的斜率值的個數(shù)具有如下關(guān)系 可見，超過四階精度的R-K公式效率并不高，實際計算通常選用如下四階格式)156(),()2/, 2/()2/, 2/(),()22(63

12、42312143211hKyhxfKhKyhxfKhKyhxfKyxfKKKKKhyyiiiiiiiiii經(jīng)典經(jīng)典R-KR-K公式公式第17頁/共20頁這時經(jīng)典R-K公式為 例例 取步長h ，采用經(jīng)典R-K法計算例 6.1 的常微分方程初值問題。)()()2()2()2()2()22(63)(42)2/(31)2/(2143211hKyehxKKhyehxKKhyehxKyexKKKKKhyyihxiihxiihxiixiiiiiii 取 h=0.2 計算得到表6-4(書133頁）。 與例和例比較可見，用經(jīng)典R-K法計算得到的解比用歐拉法和改進(jìn)歐拉法所得到的解精確得多。解解第18頁/共20頁 步長的自動選擇步長的自動選擇 對于 p 階精度的計算格式，當(dāng)取步長為 h 時，記 為從 y(xi) 計算得到的 y (