【正文】 W Y 阜師院數科院 第八章 常微分方程數值解法 81 第八章 常微分方程數值解法 W Y 阜師院數科院 第八章 常微分方程數值解法 82 第八章目錄 167。 1 歐拉（ Euler）方法 Euler法及其簡單改進 改進的 Euler法 167。 2 龍格庫塔（ Rungekutta）方法 龍格 庫塔方法的基本思想 二階龍格 庫塔公式 高階 RK公式 變步長 RK法 167。 3 線性多步法 167。 4 一階方程組與高階方程初值問題 167。 5 收斂性與穩(wěn)定性 W Y 阜師院數科院 第八章 常微分方程數值解法 83 第八章 序 許多科學技術問題，例如天文學中的星體運動，空間 技術中的物體飛行，自動控制中的系統(tǒng)分析，力學中的振 動，工程問題中的電路分析等，都可歸結為常微分方程的 初值問題。 所謂初值問題，是函數及其必要的導數在積分的起始 點為已知的一類問題，一般形式為： ),()( )1()()( ?????? nnn yyyxfyxfy ?或我們先介紹 簡單的一階問題： )18()( ))(,()(0?????????yaybxaxyxfxyyyLyxfyxfL i p sc h i t zyxf??? ),(),(:))((),( 條件李卜希茲滿足只要W Y 阜師院數科院 第八章 常微分方程數值解法 84 第八章 序 由 常微分方程 的理論可知：上述問題的解唯一存在。 常微分方程 求解求什么？應求一滿足初值問題（ 8— 1） 的解函數 y = y(x) ，如對下列微分方程： 222022224444122221221)0(xyxyxycycxyx d xydyxydxdyyxyyx???????????????????????????????《 高等數學 》 中， 微分方程求解，如對一階微分方程： y? =f(x,y)是求解解函數 y = y(x) ，使?jié)M足上述方程。但能夠 求出準確的解析函數 y(x)的微分方程是很少的， 《 高數 》 中研究微分方程的求解，是 分門別類討論 ，對不同類型的 微分方程，求解方法不一樣，因此，要求解微分方程， 首 先必須認清類型 。 W Y 阜師院數科院 第八章 常微分方程數值解法 85 微分方程 數值解 而 常微分方程 初值問題的數值解法，是要尋求解函數 y(x) 在一系列點 y(xi ) （離散點） : 上 y(xi )的 近似值 yi（ i=1,2,…,n ），并且還可由這些（ xi,yi) （ i=1,2,…,n ）構造插值函數作為近似函數。上述離散點相 鄰兩點間的距離 hi=xi1xi 稱為步長，若 hi 都相等為一定數 h, 則稱為定步長，否則為變步長。 ?? ???? nxxx 10 由于在實際問題和科學研究中遇到的微分方程往往很 復雜， 絕大多數很難，甚至不可能求出解析 函數 y(x)，因 此只能考慮求其數值解。 本章重點討論如下 一階微分方程： 在此基礎上介紹一階微分方程組與 高階微分方程的數值解法。 ????????0)( ))(,()(yaybxaxyxfxyW Y 阜師院數科院 第八章 常微分方程數值解法 86 167。 1 歐拉（ Euler）法 以 Euler法及其改進方法為例 ,說明 常微分方程 初值問題數值解法的一般概 念， Euler法很簡單，準確度也不高， 介紹此方法的目的，是由于對它的分析 討論能夠比較清楚地顯示出方法的一些 特點，而這些特點及基本方法反映了其 它方法的特點。 Euler法用于求 解一階微分方 程初值問題： ??? ? ???? 0)( ))(,()( yay bxaxyxfxyW Y 阜師院數科院 第八章 常微分方程數值解法 87 Euler法及其簡單改進 Euler公 式為： 由 x0出發(fā) ?x1,x2,…, xN，而利用此式可算出對應的 y1,y2,…, yN，式（ 82）稱為 差分方程 （序列 {yn}滿足的方 程）。 下面是 Euler公式的推導 ： 一 、從幾何意義出發(fā)： y?=f (x,y)的解函數 y=y (x) 在 xoy平 面上是一族解曲線， 而初值問題則是其中一條積分曲線， 假定 y = y(x)的曲線如圖 81從給定的點 P0(x0,y0) 出發(fā)，以P0為切點，作切線，切線斜率為曲線 y(x)的切線斜率 y? =f (x0,y0)，因此可 得切線：（點斜式） )(,( 0000 xxyxfyy ??? ）P1 P2 y(x) P0 x2 x1 x0 NabhNnnhxxn????? ),1,0(0 ?對等距節(jié)點)28(,2,1 )( ),(001 ????????? Nnyxyyxhfyy nnnn ?W Y 阜師院數科院 第八章 常微分方程數值解法 88 Euler公式的推導 (續(xù) 1） )],[(),())(,(),())(,(),(),(),())(,(11211111211122222111111111000010001上在近似曲線而用可得：相交于與為斜率，可得切線：為切點，以再以yxPPyxhfyxxyxfyyyxPxxxxyxfyyyxfyxPyxhfyxxyxfyy??????????????),( ))(,(:))(,(:)(),(111 nnnnnnnnnnnnnnnnnnyxhfyxxyxfyyxxxxyxfyyxyyyx????????????? 相交得與可作直線一般地假定已求出????,, . . . ), . . . ,2,1)(()(,1010NiNyyyNixyxyyxxx的近似值的初值問題的與所對應的微分方程出在重復上述過程，即可求??幾何意義 ：用折線近似解曲線 y = y(x)，折線不會偏離太遠 ，因為每項以 f (x, y)（斜率）修正。 切線與 x = x1交于 P1(x1,y1)，在 [x0,x1]上以切線 近似曲線， 10 PP：時，可求，當近似以 1111 )( yxxxyy ?W Y 阜師院數科院 第八章 常微分方程數值解法 89 Euler公式 的推導（續(xù) 2） 二 、利用 Taylor級數：將 y(x)在 xn處展開： )()(2),(),( ),(),(211局部的截斷誤差稱為nnnnnnnnnnnnxyhhxEyxhfyyhxEyxhfyy?????????????112),(2)(,)()(!2))(,()()(!2)()()(???????????? ? ? ? ?????????nnnnnnnnyxfynnnnyxyyxyhxyhxyxhfxyxyhxyhxyhxy，并以的線性部分作為近似式取??W Y 阜師院數科院 第八章 常微分方程數值解法 810 Euler公式的推導（續(xù) 3） 3)(8),(),()(),(),()(1),()()(1)(2))()((1)()(2))()((1)()(2)(1)()(2)(1)(1111111,1111011010???????????????????????????????????????????????????????????????????????????nnnnnnnnnnnnnnnnnnnnnnnnnnyxhfyyyxhfyxyyxfhxEyyhyxfhxEyyhyhxyxyhxfyhxyxyhxffhyyhxffhyyhxf????公式（ 83）稱為后退 Euler公式 )(2),(2nn yhhxE ???? )(2),(2nn yhhxE ????? 所謂局部載斷誤差是指以 yn代替 y(xn)而得到 y(xn+1)的近 似值 yn+1的誤差（只估計這一步的誤差）。 三、利用數值微分公式：利用兩點公式并且 Euler公式的局部截斷誤差為 : 后退 Euler公式的 局部截斷誤差為 : W Y 阜師院數科院 第八章 常微分方程數值解法 811 Euler公式的推導（續(xù) 4） )(3)()(3),(2))(,()(6))()((21)()()(6)(21)(:33112112021nnnnnnnnnnnnnyhhxEyhyxhfyyxyxfyhxyxyhxyfhyyhxf???????????????????????????????????，去掉誤差項，也稱為中點法公式還可利用三點公式4)(8 ),(211 nnnn yxhfyy ??? ??W Y 阜師院數科院 第八章 常微分方程數值解法 812 Euler公式的推導（續(xù) 5） 四 、 利用數值積分公式：在 [xn,xn+1]上對 y?(x)=f (x,y(x)) 積分 ?????????????????1111))(,()()( ))(,()()())(,()(11nnnnnnnnxxnnxxnnxxxxdxxyxfxyxydxxyxfxyxydxxyxfdxxy所以 對右端積分項采用不同的數值積分公式，便可得到各種 不同的求解 dE初值問題的計算公式。 如，以矩形面積代替曲邊梯形面積 ? ? 1 ))(,(nnxx dxxyxf1） 以左矩形面積代替曲邊梯形面積如圖 82，亦即以 ))(,(),( nnn xyxfyxf ?公式為則有E u l e ryxhfyyxfxxyynnnnnnnnn),(),()( 11????? ??y f(x, y) xn x xn+1 圖 82 W Y 阜師院數科院 第八章 常微分方程數值解法 813 y f(x, y) xn x xn+1 圖 83 y f(x, y) xn x xn+1 圖 84 3）以梯形公式（面積）代替曲邊形如圖 84則有 式（ 85）稱為求 dE初值問題的梯形公式，梯形公式看作 是以（ xn,yn） (xn+1,yn+1)構造的插值多項式代替被積函數 得到的，而 Euler公式 則是以左端點