常微分方程的數(shù)值解法

第五章常微分方程的數(shù)值解法主要內(nèi)容：1、引言2、歐拉方法3、龍格－庫塔方法4、單步法的收斂性和穩(wěn)定性5、線性多步法6、一階方程組與高階方程1.第一節(jié)引言●在常微分方程課程里面討論的是一些典型方程求解解析解的根本方法?！袢欢谏a(chǎn)實踐和科學研究中遇到的微分方程往往比較復雜，在很多情況下，不能給出解的解析表達式；有時候即時能用解析表達式來表示，又因為計算量太大而不實用，有時候一些是已經(jīng)有了求解的根本方法的典型方程，但實際使用時也是有困難的。●以上情況說明用求解解析解的根本方法來求微分方程的解往往是不適宜的，甚至很難辦到?！駥嶋H問題中，對于求解微分方程，一般只要求得到解的假設干個點上的近似值或者解的便于計算的近似表達式?！癖菊卵芯课⒎址匠痰臄?shù)值解法，而且著重討論微分方程中最簡單的一類問題——一階方程的初值問題。2.第一節(jié)引言1、一階方程的初值問題假定上式在區(qū)間[a,b]上存在唯一且足夠光滑的解y(x)?！袼^數(shù)值解法就是尋求解y(x)在一系列離散點，也稱為節(jié)點處的值：要計算出解函數(shù)y(x)在一系列節(jié)點a=x0<x1<…<xn=b處的近似值3.第一節(jié)引言●節(jié)點間距，即步長為：通常采用等距節(jié)點，即hi=h(常數(shù))●等間距節(jié)點●在這些節(jié)點上采用離散化方法〔通常用數(shù)值積分、微分、泰勒展開等〕將上述初值問題化成關于離散變量的相應問題。把這個相應問題的解yn作為y(xn)的近似值。這樣求得的yn就是上述初值問題在節(jié)點xn上的數(shù)值解。一般說來，不同的離散化導致不同的方法。4.第二節(jié)歐拉方法一、歐拉法Euler1、向前差商近似導數(shù)記為5.第二節(jié)歐拉方法2、舉例例1

用歐拉法求初值問題當h=0.02時在區(qū)間[0,0.10]上的數(shù)值解。解：根據(jù)歐拉公式可以得到：此外，可以得到方程的真解：6.第二節(jié)歐拉方法求解過程如下：nxnyny(xn)

n=y(xn)-

yn001.00001.0000010.020.98200.98250.000520.040.96500.96600.000530.060.94890.95030.001440.080.93360.93540.001850.100.91920.9230.00217.第二節(jié)歐拉方法3、歐拉方法的幾何意義根據(jù)條件：曲線y(x)上的點〔x0,y0〕及該點處曲線的導數(shù)f(x0,y0)，那么可以得到過該點的直線：該直線與x＝x1的交點P1，那么P1的縱坐標y1為：就用y1作為y(x1)的近似值…逐次進行后可以得到一條折線P0P1…Pn，該折線看作是初值問題的積分曲線的近似，因此歐拉方法也稱為歐拉折線法8.第二節(jié)歐拉方法從上述幾何意義上得知，由Euler法所得的折線明顯偏離了積分曲線，可見此方法非常粗糙即誤差太大。4、歐拉法的局部截斷誤差〔1〕截斷誤差定義在假設yi=y(xi)，即第i步計算是精確的前提下，考慮的截斷誤差Ri+1=y(xi+1)yi+1，稱為局部截斷誤差如下圖：APi+1即為歐拉方法在xi+1點的截斷誤差〔2〕如果某種方法的局部截斷誤差是那么稱該方法具有p階精度A9.第二節(jié)歐拉方法〔3〕那么截斷誤差的大?。繉懗鰕(xn+1)的泰勒展開式：由歐拉方法可以得到：那么上面兩個公式相減得到：具有1階精度10.第二節(jié)歐拉方法二、改進的歐拉法一階方程的初值問題與如下積分方程是等價的：當x=x1時可以借助于數(shù)值積分，求y(x1)的值1、用矩形公式11.第二節(jié)歐拉方法可以推導出：用矩形法計算右端的積分與用歐拉法計出的結(jié)果完全相同2、用梯形公式那么可以推導出：12.第二節(jié)歐拉方法梯形公式的截斷誤差：梯形公式具有二階精度，比歐拉方法有了進步13.第二節(jié)歐拉方法和歐拉公式相比較，梯形公式在計算yi+1時候也只用到前一步的值yi，但是假設yi，將yi帶入公式求解時候，一般不能直接得到y(tǒng)i+1，而需要通過其他方法〔比方迭代法〕求解，所以梯形公式被稱為隱式公式。3、改進的歐拉方法梯形公式是隱式的，一般用迭代法求解，計算量較大。實際中常將歐拉公式和梯形公式聯(lián)合使用，先用歐拉公式得出一個y(xi+1)的近似值稱為預估值，然后對預估值使用梯形公式對它進行精確化，得到較為精確的近似值yi+1，稱之為校正值，計算公式為：這樣的預估校正系統(tǒng)稱為改進的歐拉方法。14.第二節(jié)歐拉方法為了便于編寫程序，常將上面的公式改寫為如下式：15.第二節(jié)歐拉方法4、舉例P90，例題5-1在區(qū)間[0,1.5]上，取h=0.1，求解。解：〔1〕用歐拉法計算公式如下：〔2〕用改進歐拉法計算公式如下：16.此題的精確解為，可用來檢驗近似解的精確程度。計算結(jié)果如表：xn歐拉法yn迭代一次改進歐拉法yn準確解01110.11.11.0959091.0954450.21.1918181.1840961.1832160.31.2774381.2602011.2649110.41.3582131.3433601.3416410.51.4351331.4161021.4142140.61.5089661.4829561.4832400.71.5803381.5525151.5491930.81.6497831.6164761.6124520.91.7177791.6781681.6733201.01.7847701.7378691.7320511.11.851181.7958221.7888541.21.9174641.8522421.8439091.31.9840461.9073231.8973671.42.0514041.9612531.9493591.52.1200522.0142072.00000017.第二節(jié)歐拉方法P108習題5-118.第二節(jié)歐拉方法5、歐拉兩步公式中心差商19.第三節(jié)龍格－庫塔方法一、龍格庫塔法的根本思想1、平均斜率考察差商：根據(jù)微分中值定理：在閉區(qū)間[a,b]上連續(xù)，開區(qū)間(a,b)上可導，那么至少存在(a,b)上的一點，使得下式成立：根據(jù)上面公式可以得到：稱為區(qū)間[xi,xi+1]上的平均斜率20.第三節(jié)龍格－庫塔方法因此只要對K*提供一種算法，就可以求得數(shù)值解，根據(jù)該觀點對歐拉法及改進的歐拉進行分析。2、基于平均斜率對歐拉法和改進的歐拉法進行分析〔1〕在歐拉公式中，是取了一個點xi上的斜率值f(xi,yi)作為平均斜率K*的近似值的，已經(jīng)知道其精度較低?！?〕對于改進的歐拉公式：可以看出它用兩個點xi和xi+1上的斜率K1和K2的算術平均值作為平均斜率K*的近似值的，而xi+1處的斜率是由信息預測得到的。21.第三節(jié)龍格－庫塔方法根據(jù)上面的分析得到：如果在區(qū)間[xi,xi+1]上多測幾個點的斜率值，然后取其加權平均作為平均斜率的近似值，有可能構造出具有更高精度的計算公式，此即為龍格庫塔算法的根本思想。二、二階龍格庫塔方法1、推廣改進的歐拉方法，考察區(qū)間[xi,xi+1]上的一點：用xi和xp兩個點上的斜率值K1和K2的加權平均作為平均斜率K*的近似值：即?。号c改進的歐拉法類似，有：如何得到xp的斜率？22.第三節(jié)龍格－庫塔方法如何得到xp的斜率？根據(jù)改進的歐拉法，可以利用歐拉法預測的值：那么可以得到點xi+p斜率K2：那么可以得到算法的具體表達式：23.第三節(jié)龍格－庫塔方法2、選擇參數(shù)使得算法具有2階精度計算上面公式的局部截斷誤差：根據(jù)泰勒公式有：另外：24.25.第三節(jié)龍格－庫塔方法考慮到：那么有：二元泰勒展開兩式相減可得到：26.第三節(jié)龍格－庫塔方法兩式相減可得到：那么要使得上式滿足二階精度，即Ri+1=O(h3)，只需要：共有3個參數(shù)，但只需要滿足2個條件，因此滿足該式的參數(shù)不止一組，而是一簇，所有滿足條件的公式通稱為二階龍格庫塔公式。27.第三節(jié)龍格－庫塔方法該公式的特別情況：有：此時，二階龍格庫塔公式就是改進的歐拉公式。28.第三節(jié)龍格－庫塔方法如果p=1/2，那么有：此時二階龍格庫塔公式變成：為變形的歐拉公式。29.第三節(jié)龍格－庫塔方法3、高階龍格庫塔公式為了進一步提高精度，可以考慮在區(qū)間[xi,xi+1]上，除了xi，xp外，再增加幾個點。如再增加一個點：xi+m利用該三點處的斜率的加權平均值作為K*近似值：此時計算公式為：式中K1和K2的取法與前相同，如何獲得K3？30.第三節(jié)龍格－庫塔方法為了得到xi+m處的斜率K3，需要確定該點的函數(shù)值y(xi+m)，可以用二階龍格庫塔公式估計