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

第五章 常微分方程的數(shù)值解法,主要內(nèi)容： 1、引言 2、歐拉方法 3、龍格庫塔方法 4、單步法的收斂性和穩(wěn)定性 5、線性多步法 6、一階方程組與高階方程,1,第一節(jié) 引言,在常微分方程課程里面討論的是一些典型方程求解解析解的基本方法。 然而在生產(chǎn)實(shí)踐和科學(xué)研究中遇到的微分方程往往比較復(fù)雜，在很多情況下，不能給出解的解析表達(dá)式；有時(shí)候即時(shí)能用解析表達(dá)式來表示，又因?yàn)橛?jì)算量太大而不實(shí)用，有時(shí)候一些是已經(jīng)有了求解的基本方法的典型方程，但實(shí)際使用時(shí)也是有困難的。 以上情況說明用求解解析解的基本方法來求微分方程的解往往是不適宜的，甚至很難辦到。 實(shí)際問題中，對(duì)于求解微分方程，一般只要求得到解的若干個(gè)點(diǎn)上的近似值或者解的便于計(jì)算的近似表達(dá)式。 本章研究微分方程的數(shù)值解法，而且著重討論微分方程中最簡單的一類問題一階方程的初值問題。,2,第一節(jié) 引言,1、一階方程的初值問題 假定上式在區(qū)間a, b上存在唯一且足夠 光滑的解y(x)。 所謂數(shù)值解法就是尋求解y(x)在一系列離散點(diǎn)，也稱為節(jié)點(diǎn)處的值： 要計(jì)算出解函數(shù) y(x) 在一系列節(jié)點(diǎn) a = x0 x1 xn= b 處的近似值,3,第一節(jié) 引言,節(jié)點(diǎn)間距，即步長為： 通常采用等距節(jié)點(diǎn)，即hi = h (常數(shù)) 等間距節(jié)點(diǎn) 在這些節(jié)點(diǎn)上采用離散化方法（通常用數(shù)值積分、微分、泰勒展開等）將上述初值問題化成關(guān)于離散變量的相應(yīng)問題。把這個(gè)相應(yīng)問題的解yn作為y(xn)的近似值。這樣求得的yn就是上述初值問題在節(jié)點(diǎn)xn上的數(shù)值解。一般說來，不同的離散化導(dǎo)致不同的方法。,4,第二節(jié) 歐拉方法,一、歐拉法Euler 1、向前差商近似導(dǎo)數(shù),5,第二節(jié) 歐拉方法,2、舉例 例1 用歐拉法求初值問題 當(dāng)h = 0.02時(shí)在區(qū)間0, 0.10上的數(shù)值解。 解：根據(jù)歐拉公式可以得到： 此外，可以得到方程的真解：,6,第二節(jié) 歐拉方法,求解過程如下：,7,第二節(jié) 歐拉方法,3、歐拉方法的幾何意義,根據(jù)已知條件：曲線y(x)上的點(diǎn)（x0,y0）及該點(diǎn)處曲線的導(dǎo)數(shù)f(x0,y0)，則可以得到過該點(diǎn)的直線：,該直線與xx1的交點(diǎn)P1，則P1的縱坐標(biāo)y1為：,就用y1作為y(x1)的近似值,逐次進(jìn)行后可以得到一條折線P0P1Pn，該折線看作是初值問題的積分曲線的近似，因此歐拉方法也稱為歐拉折線法,8,第二節(jié) 歐拉方法,從上述幾何意義上得知，由Euler法所得的折線明顯偏離了積分曲線，可見此方法非常粗糙即誤差太大。 4、歐拉法的局部截?cái)嗾`差 （1）截?cái)嗾`差定義 在假設(shè) yi = y(xi)，即第 i 步計(jì)算是精確的前提下，考慮的截?cái)嗾`差 Ri+1 = y(xi+1) yi+1 ，稱為局部截?cái)嗾`差 如圖所示：APi+1即為歐拉方法在xi+1點(diǎn)的截?cái)嗾`差 （2）如果某種方法的局部截?cái)嗾`差是 則稱該方法具有p階精度,9,第二節(jié) 歐拉方法,（3）則截?cái)嗾`差的大??？ 寫出y(xn+1)的泰勒展開式： 由歐拉方法可以得到： 則上面兩個(gè)公式相減得到：,具有1階精度,10,第二節(jié) 歐拉方法,二、改進(jìn)的歐拉法 一階方程的初值問題與如下積分方程是等價(jià)的： 當(dāng)x = x1時(shí) 可以借助于數(shù)值積分，求y(x1)的值 1、用矩形公式,11,第二節(jié) 歐拉方法,可以推導(dǎo)出： 用矩形法計(jì)算右端的積分與用歐拉法計(jì)出的結(jié)果完全相同 2、用梯形公式 則可以推導(dǎo)出：,12,第二節(jié) 歐拉方法,梯形公式的截?cái)嗾`差：,梯形公式具有二階精度，比歐拉方法有了進(jìn)步,13,第二節(jié) 歐拉方法,和歐拉公式相比較，梯形公式在計(jì)算yi+1時(shí)候也只用到前一步的值yi，但是若yi已知，將yi帶入公式求解時(shí)候，一般不能直接得到y(tǒng)i+1，而需要通過其他方法（比如迭代法）求解，所以梯形公式被稱為隱式公式。 3、改進(jìn)的歐拉方法 梯形公式是隱式的，一般用迭代法求解，計(jì)算量較大。實(shí)際中常將歐拉公式和梯形公式聯(lián)合使用，先用歐拉公式得出一個(gè)y(xi+1)的近似值 稱為預(yù)估值，然后對(duì)預(yù)估值使用梯形公式對(duì)它進(jìn)行精確化，得到較為精確的近似值yi+1，稱之為校正值，計(jì)算公式為： 這樣的預(yù)估校正系統(tǒng)稱為改進(jìn)的歐拉方法。,14,第二節(jié) 歐拉方法,為了便于編寫程序，常將上面的公式改寫為如下式：,15,第二節(jié) 歐拉方法,4、舉例 P90，例題5-1 在區(qū)間0, 1.5上，取h = 0.1，求解。 解：（1）用歐拉法計(jì)算公式如下： （2）用改進(jìn)歐拉法計(jì)算公式如下：,16,本題的精確解為 ， 可用來檢驗(yàn)近似解 的精確程度。 計(jì)算結(jié)果如表：,17,第二節(jié) 歐拉方法,P108 習(xí)題5-1,18,第二節(jié) 歐拉方法,5、歐拉兩步公式 中心差商,19,第三節(jié) 龍格庫塔方法,一、龍格庫塔法的基本思想 1、平均斜率 考察差商： 根據(jù)微分中值定理：在閉區(qū)間a, b上連續(xù)，開區(qū)間(a, b)上可導(dǎo)，則至少存在(a, b)上的一點(diǎn)，使得下式成立： 根據(jù)上面公式可以得到：,稱為區(qū)間xi,xi+1上的平均斜率,20,第三節(jié) 龍格庫塔方法,因此只要對(duì) K* 提供一種算法，就可以求得數(shù)值解，根據(jù)該觀點(diǎn)對(duì)歐拉法及改進(jìn)的歐拉進(jìn)行分析。 2、基于平均斜率對(duì)歐拉法和改進(jìn)的歐拉法進(jìn)行分析 （1）在歐拉公式中，是取了一個(gè)點(diǎn) xi上的斜率值f(xi,yi) 作為平均斜率K*的近似值的，已經(jīng)知道其精度較低。 （2）對(duì)于改進(jìn)的歐拉公式：,可以看出它用兩個(gè)點(diǎn)xi 和 xi+1上的斜率K1和K2 的算術(shù)平均值作為平均斜率 K* 的近似值的，而xi+1處的斜率是由已知信息預(yù)測得到的。,21,第三節(jié) 龍格庫塔方法,根據(jù)上面的分析得到：如果在區(qū)間 xi, xi+1上多測幾個(gè)點(diǎn)的斜率值，然后取其加權(quán)平均作為平均斜率的近似值，有可能構(gòu)造出具有更高精度的計(jì)算公式，此即為龍格庫塔算法的基本思想。 二、二階龍格庫塔方法 1、推廣改進(jìn)的歐拉方法，考察區(qū)間xi, xi+1上的一點(diǎn)： 用xi 和xp 兩個(gè)點(diǎn)上的斜率值K1和K2 的加權(quán)平均作為平均斜率 K*的近似值： 即取： 與改進(jìn)的歐拉法類似，有： 如何得到xp的斜率？,22,第三節(jié) 龍格庫塔方法,如何得到xp的斜率？ 根據(jù)改進(jìn)的歐拉法，可以利用歐拉法預(yù)測 的值： 則可以得到點(diǎn)xi+p斜率K2： 則可以得到算法的具體表達(dá)式：,23,第三節(jié) 龍格庫塔方法,2、選擇參數(shù)使得算法具有2階精度 計(jì)算上面公式的局部截?cái)嗾`差： 根據(jù)泰勒公式有： 另外：,24,25,第三節(jié) 龍格庫塔方法,考慮到： 則有：,二元泰勒展開,兩式相減可得到：,26,第三節(jié) 龍格庫塔方法,兩式相減可得到： 則要使得上式滿足二階精度，即Ri+1=O(h3)，只需要： 共有3個(gè)參數(shù)，但只需要滿足2個(gè)條件，因此滿足該式的參數(shù)不止一組，而是一簇，所有滿足條件的公式通稱為二階龍格庫塔公式。,27,第三節(jié) 龍格庫塔方法,該公式的特別情況： 有： 此時(shí)，二階龍格庫塔公式就是改進(jìn)的歐拉公式。,28,第三節(jié) 龍格庫塔方法,如果p=1/2，則有： 此時(shí)二階龍格庫塔公式變成： 為變形的歐拉公式。,29,第三節(jié) 龍格庫塔方法,3、高階龍格庫塔公式 為了進(jìn)一步提高精度，可以考慮在區(qū)間xi, xi+1上，除了xi，xp外，再增加幾個(gè)點(diǎn)。如再增加一個(gè)點(diǎn)：xi+m 利用該三點(diǎn)處的斜率的加權(quán)平均值作為K* 近似值： 此時(shí)計(jì)算公式為： 式中K1和K2的取法與前相同，如何獲得K3？,30,第三節(jié) 龍格庫塔方法,為了得到xi+m處的斜率K3，需要確定該點(diǎn)的函數(shù)值y(xi+m)，可以用二階龍格庫塔公式估計(jì)y(xi+m)的近似值yi+m： 然后通過計(jì)算函數(shù)值，就可以得到K3： 由此得到的 3階龍格庫塔 公式為：,31,第三節(jié) 龍格庫塔方法,通過選擇未知參數(shù)使得公式具有3階精度，采用前面相同的方法，可以得到參數(shù)需滿足如下條件：,32,第三節(jié) 龍格庫塔方法,常用的三階龍格庫塔公式為：,其中,33,第三節(jié) 龍格庫塔方法,經(jīng)典的四階龍格庫塔公式為：,34,第三節(jié) 龍格庫塔方法,35,