第九章常微分方程數(shù)值解

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

在自然界與工程界中,很多問題的數(shù)學(xué)表述都可以歸結(jié)為解常微分方程和偏微分方程問題.但是,只有極少數(shù)特殊的方程有解析解。對(duì)于絕大部分的微分方程求其解析解是很困難的。因此研究微分方程的數(shù)值解是很重要的。

本章討論常微分方程的數(shù)值解法歐拉法,后退的歐拉法,梯形公式;顯式的,隱式的方法,單步法多步法,局部截?cái)嗾`差;常用的四階R-K方法另外,還有收斂性分析,穩(wěn)定性分析,剛性方程,高性能算法等.對(duì)于一個(gè)常微分方程：通常會(huì)有無窮個(gè)解。例如：因此，我們要加入一個(gè)限定條件。通常會(huì)在端點(diǎn)出給出，如下面的初值問題：為了使解存在唯一，一般，要加限制條件在f上，要求f對(duì)y滿足Lipschitz條件：如何求解解析解法：（常微分方程理論）只能求解極少一類常微分方程；實(shí)際中給定的問題不一定是解析表達(dá)式，而是函數(shù)表，無法用解析解法。計(jì)算解函數(shù)y(x)在一系列節(jié)點(diǎn)a=x0<x1<…<xn=b

處的近似值節(jié)點(diǎn)間距為步長，通常采用等距節(jié)點(diǎn)，即取hi=h

(常數(shù))。數(shù)值解法:求解所有的常微分方程步進(jìn)式：根據(jù)已知的或已求出的節(jié)點(diǎn)上的函數(shù)值計(jì)算當(dāng)前節(jié)點(diǎn)上的函數(shù)值，一步一步向前推進(jìn)。因此只需建立由已知的或已求出的節(jié)點(diǎn)上的函數(shù)值求當(dāng)前節(jié)點(diǎn)函數(shù)值的遞推公式即可。例：我們對(duì)區(qū)間做等距分割：設(shè)解函數(shù)在節(jié)點(diǎn)的近似為由數(shù)值微分公式，我們有，則：向前差商公式可以看到，給出初值，就可以用上式求出所有的這種方法，稱為數(shù)值離散方法。在一系列離散點(diǎn)列上，求未知函數(shù)y在這些點(diǎn)上的值的近似。為了考察數(shù)值方法提供的數(shù)值解，是否有實(shí)用價(jià)值，需要知道如下幾個(gè)結(jié)論：①步長充分小時(shí)，所得到的數(shù)值解能否逼近問題得真解；即收斂性問題②誤差估計(jì)③產(chǎn)生得舍入誤差，在以后得各步計(jì)算中，是否會(huì)無限制擴(kuò)大；即穩(wěn)定性問題單步法：在計(jì)算yi+1

時(shí)只利用yi多步法：在計(jì)算yi+1

時(shí)不僅利用yi,還要利用yi?1,yi?2,…,k步法：在計(jì)算yi+1

時(shí)要用到y(tǒng)i,yi?1,…,yi?k+1顯式計(jì)算公式可寫成:yk+1=yk+hΦf(xk，yk;h）隱式格式:yk+1=yk+hΦf（xk，yk,yk+1;h）它每步求解yk+1需要解一個(gè)隱式方程1.Euler公式做等距分割，利用數(shù)值微分代替導(dǎo)數(shù)項(xiàng)，向前差商公式于是得到Euler公式稱為局部截?cái)嗾`差。顯然，這個(gè)誤差在逐步計(jì)算過程中會(huì)傳播，積累。因此還要估計(jì)這種積累定義在假設(shè)yi=y(xi)，即第

i

步計(jì)算是精確的前提下，考慮的截?cái)嗾`差Ri=y(xi+1)

yi+1稱為局部截?cái)嗾`差/*localtruncationerror*/。定義若某算法的局部截?cái)嗾`差為O(hp+1)，則稱該算法有p

階精度。Euler公式一階顯式的2.后退的Euler公式是隱式的，要迭代求解后退的Euler公式一階隱式的預(yù)測-校正系統(tǒng):用顯式算法預(yù)測,

再用隱式算法校正PQR兩式相加平均得到梯形公式預(yù)測-校正系統(tǒng):用顯式算法預(yù)測,

再用隱式算法校正2階隱式的改進(jìn)的歐拉格式3.梯形公式2階隱式的單步法幾何解釋xnxn+1ABPn+1=(A+B)/2尤拉法后退尤拉法梯形法

基于數(shù)值積分的構(gòu)造法將在上積分，得到只要近似地算出右邊的積分，則可通過近似y(xn+1)

。而選用不同近似式Ik，可得到不同的計(jì)算公式。Anotherpointofview所以，有格式為：局部截?cái)嗾`差梯形公式隱式的對(duì)微分方程做積分，則：例如公式局部截?cái)嗾`差精度顯隱穩(wěn)定性步數(shù)尤拉顯式公式1階顯差單步尤拉隱式公式1階隱好單步梯形公式2階隱差單步中點(diǎn)法2階顯好二步summary4.Runge－Kutta法由Taylor展開記為所以，可以構(gòu)造格式這種格式使用到了各階偏導(dǎo)數(shù)，使用不便。從另一個(gè)角度看，取(x,y)及其附近的點(diǎn)做線性組合，表示F，問題就好辦了。當(dāng)然，要求此時(shí)的展開精度相同。這種方法稱為Runge－Kutta法。在(x,y)處展開，比較以2階為例，設(shè)有：改進(jìn)的Euler公式常用的四階Runge－Kutta方法四階Runge－Kutta方法的每一步需要四次計(jì)算函數(shù)值f,可以證明其截?cái)嗾`差為O(h5)§2Runge-KuttaMethod注：

龍格-庫塔法的主要運(yùn)算在于計(jì)算Ki

的值，即計(jì)算f

的值。Butcher于1965年給出了計(jì)算量與可達(dá)到的最高精度階數(shù)的關(guān)系：753可達(dá)到的最高精度642每步須算Ki的個(gè)數(shù)

由于龍格-庫塔法的導(dǎo)出基于泰勒展開，故精度主要受解函數(shù)的光滑性影響。對(duì)于光滑性不太好的解，最好采用低階算法而將步長h

取小。方程組和高階方程的數(shù)值解法寫成向量的形式：各種方法都可以直接運(yùn)用過來。Euler公式以兩個(gè)方程的方程組為例Runge-Kutta公式1、2、確定方法，然后求解（0.202760.0881157）（0.2130070.0934037）（0.2237630.0988499）（0.2350520.104437）（0.2469020.110146）4階Runge-Kutta法，h=1例高階方程則有：令例：考察初值問題在區(qū)間[0,0.5]上的解。分別用歐拉法、后退的歐拉法和改進(jìn)的歐拉法計(jì)算數(shù)值解。0.00.10.20.30.40.5精確解改進(jìn)的歐拉法

后退的歐拉法歐拉顯式

節(jié)點(diǎn)xi

1.0000

2.00004.0000

8.00001.6000101

3.2000101

1.00002.5000101

6.25001021.56251023.90631039.