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

第七章常微分方程初值問題數(shù)值解法

數(shù)值分析2023/8/101NumericalAnalysis第七章常微分方程初值問題數(shù)值解法

數(shù)值分析2023/8/本章內(nèi)容歐拉法歐拉公式兩步歐拉公式梯形法改進歐拉法龍格-庫塔法基本思路二階、三階龍格-庫塔法經(jīng)典龍格-庫塔法隱式龍格-庫塔法線性多步法亞當斯法亞當斯預(yù)報-校正公式誤差修正法收斂性與穩(wěn)定性微分方程組和高階微分方程2023/8/102NumericalAnalysis本章內(nèi)容歐拉法線性多步法2023/8/12Numerical簡單的數(shù)值方法與基本概念科學(xué)技術(shù)中常常需要求解常微分方程的定解問題.這類問題最簡單的形式，是本章將著重考察的一階方程的初值問題我們知道，只有f(x,y)適當光滑—譬如關(guān)于y滿足利普希茨(Lipschitz)條件理論上就可以保證初值問題的解y＝f(x)存在并且唯一.2023/8/103NumericalAnalysis簡單的數(shù)值方法與基本概念科學(xué)技術(shù)中常常需要雖然求解常微分方程有各種各樣的解析方法，但解析方法只能用來求解一些特殊類型的方程，實際問題中歸結(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é)點).2023/8/104NumericalAnalysis雖然求解常微分方程有各種各樣的解析方法，但解初值問題的數(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)定性等問題.2023/8/105NumericalAnalysis初值問題的數(shù)值解法有個基本特點，他們都采取歐拉法2023/8/106NumericalAnalysis歐拉法2023/8/16NumericalAnalysis2023/8/107NumericalAnalysis2023/8/17NumericalAnalysis2023/8/108NumericalAnalysis2023/8/18NumericalAnalysis2023/8/109NumericalAnalysis2023/8/19NumericalAnalysis2歐拉法數(shù)學(xué)推導(dǎo)

：數(shù)值微分（用差商代替導(dǎo)數(shù)）設(shè)等距，步長令x=xn,x+h=xn+1,

y(xn)≈yn

，y(xn+1)≈yn+1,初值問題離散化為(歐拉公式)

2023/8/1010NumericalAnalysis2歐拉法數(shù)學(xué)推導(dǎo)：數(shù)值微分（用差商代替導(dǎo)數(shù)）設(shè)等距，2023/8/1011NumericalAnalysis2023/8/111NumericalAnalysis

例7-1

用歐拉公式求解初值問題

解取步長h=0.1，歐拉公式的具體形式為其中xn=nh=0.1n(n=0,1,

,10),已知y0=1,由此式可得2023/8/1012NumericalAnalysis例7-1用歐拉公式求解初值問題解取步長h依次計算下去，部分計算結(jié)果見下表.與準確解相比，可看出歐拉公式的計算結(jié)果精度很差.

xn

歐拉公式數(shù)值解yn準確解y(xn)

誤差0.20.40.60.81.01.1918181.3582131.5089661.6497831.7847701.1832161.3416411.4832401.6124521.7320510.0086020.0165720.0257260.0373310.0527192023/8/1013NumericalAnalysis依次計算下去，部分計算結(jié)果見下表.與準確解

局部截斷誤差和階：數(shù)值公式的精度定義局部截斷誤差：假設(shè)第n步是準確的，即y(xn)=yn，將y(xn+1)-yn+1定義為數(shù)值方法的局部截斷誤差。由于實際上yn不是準確值，因此它的誤差會傳播下去。實際計算時，每一步都可能產(chǎn)生舍入誤差。定義若局部截斷誤差為O(hp+1),p為正整數(shù)，則稱數(shù)值公式是p階公式。

2023/8/1014NumericalAnalysis局部截斷誤差和階：2023/8/114Numeric歐拉公式的截斷誤差是O(h2)，公式是1階的。二階泰勒公式

兩式相減，由設(shè)yn=y(xn)，有

歐拉公式的局部截斷誤差和階2023/8/1015NumericalAnalysis歐拉公式的截斷誤差是O(h2)，公式是1隱式（后退）歐拉公式2023/8/1016NumericalAnalysis隱式（后退）歐拉公式2023/8/116Numerical兩步歐拉公式2023/8/1017NumericalAnalysis兩步歐拉公式2023/8/117NumericalAnal2023/8/1018NumericalAnalysis2023/8/118NumericalAnalysis梯形法對微分方程y′=f(x,y)兩邊求xn到xn+1的定積分，有利用梯形公式計算積分，有

將y(xn)

、y(xn+1)分別用yn、yn+1

代替，構(gòu)造數(shù)值公式2023/8/1019NumericalAnalysis梯形法對微分方程y′=f(x,y)兩邊求xn到xn+1改進的歐拉法(預(yù)報-校正公式)

歐拉法，顯式，計算量小，精度低。梯形公式是隱式公式,計算量大，精度高。

實際計算時，將二者綜合之，先用歐拉公式計算出yn+1作為初始值,初始值精度不高，取作預(yù)報值，代入梯形公式，得到校正值yn+1。寫成預(yù)報-校正公式

2023/8/1020NumericalAnalysis改進的歐拉法(預(yù)報-校正公式)2023/8/120Nume預(yù)報-校正公式又常常寫成一步嵌套顯式形式或?qū)懗善骄问筋A(yù)報-校正公式的局部截斷誤差y(xn+1)-yn+1=O(h3)2023/8/1021NumericalAnalysis預(yù)報-校正公式或?qū)懗善骄问筋A(yù)報-校正公式的局部截斷誤差y2023/8/1022NumericalAnalysis2023/8/122NumericalAnalysis取步長h=0.1，有n=02023/8/1023NumericalAnalysis取步長h=0.1，有2023/8/123Numerical龍格-庫塔（Runge-Kutta)法

基本思想

二階龍格-庫塔法經(jīng)典龍格-庫塔法2023/8/1024NumericalAnalysis龍格-庫塔（Runge-Kutta)法2023/8/124N龍格-庫塔（Runge-Kutta)法2023/8/1025NumericalAnalysis龍格-庫塔（Runge-Kutta)法2023/8/125N2023/8/1026NumericalAnalysis2023/8/126NumericalAnalysis2023/8/1027NumericalAnalysis2023/8/127NumericalAnalysis2023/8/1028NumericalAnalysis2023/8/128NumericalAnalysis2023/8/1029NumericalAnalysis2023/8/129NumericalAnalysis經(jīng)典（四階）龍格-庫塔法仿照上述的討論，可導(dǎo)出四階經(jīng)典龍格－庫塔公式經(jīng)典龍格－庫塔公式每步要四次計算函數(shù)值，具有四階精度，局部截斷誤差是O(h5).

2023/8/1030NumericalAnalysis經(jīng)典（四階）龍格-庫塔法經(jīng)典龍格－庫塔公式每步要四次計算函數(shù)2023/8/1031NumericalAnalysis2023/8/131NumericalAnalysis然而值得指出的是，龍格-庫塔方法的推導(dǎo)基于泰勒展開方法，因而它要求所求的解具有較好的光滑性質(zhì).反之，如果解的光滑性差，那么，使用龍格-庫塔方法求得的數(shù)值解，其精度可能反而不如改進的歐拉方法.實際計算時，我們應(yīng)當針對問題的具體特點選擇合適的算法.2023/8/1032NumericalAnalysis然而值得指出的是，龍格-庫塔方法的推導(dǎo)基于泰2023/8/1033NumericalAnalysis2023/8/133NumericalAnalysis線性多步法(n=0,1,···)其中,xn+i=x0+(n+i)h,fn+i=f(xn+i,yn+i)局部載斷誤差A(yù)damas顯格式:yn+2=yn+1+h(3fn+1-fn)/2yn+3=yn+2+h(23fn+2-16fn+1+5fn)/122023/8/1034NumericalAnalysis