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

第五章常微分方程的數(shù)值解法許多工程實(shí)際問(wèn)題的數(shù)學(xué)模型可以用常微分方程來(lái)描述，但是，除了常系數(shù)線性微分方程和一些特殊的微分方程可以用解析方法求解以外，在生產(chǎn)實(shí)踐和科學(xué)研究中遇到的絕大多數(shù)微分方程往往比較復(fù)雜，不能給出解的解析表達(dá)式；有時(shí)候即使能用解析表達(dá)式來(lái)表示，又因?yàn)橛?jì)算量太大而不實(shí)用，有些是已經(jīng)有了求解的基本方法的典型方程，但實(shí)際使用時(shí)也是有困難的。很多實(shí)際問(wèn)題只需要求某些點(diǎn)上的函數(shù)值，不必求出函數(shù)表達(dá)式。

因此，我們有必要研究常微分方程的數(shù)值近似解法。本章主要研究一階常微分方程初值問(wèn)題的幾個(gè)常用解法。例：某傘兵與降落傘質(zhì)量為90kg，當(dāng)傘張開(kāi)時(shí)以12m/s速度垂直下降。假設(shè)空氣阻力與傘的下降速度成正比，在速度為6m/s時(shí)，測(cè)得的空氣阻力為353N。試求傘兵開(kāi)傘后，第1秒末、第2秒末、直到第5秒末時(shí)各時(shí)刻的速度。解：設(shè)阻力為，則由力學(xué)知識(shí)可得，本章主要內(nèi)容一階常微分方程初值問(wèn)題歐拉方法、梯形格式、改進(jìn)的歐拉方法龍格－庫(kù)塔（Runge-Kutta）方法第一節(jié)引言一階方程的初值問(wèn)題：由常微分方程理論，在滿(mǎn)足一定條件時(shí)存在唯一解函數(shù)。要計(jì)算出解函數(shù)y(x)在一系列節(jié)點(diǎn)a=x0<x1<…<xn=b處的近似值。常微分方程的數(shù)值解節(jié)點(diǎn)間距為步長(zhǎng)，通常采用等距節(jié)點(diǎn)，即取hi=h

(常數(shù))。在這些節(jié)點(diǎn)上采用離散化方法，（通常用數(shù)值積分、微分、泰勒展開(kāi)等）將上述初值問(wèn)題化成關(guān)于離散變量的相應(yīng)問(wèn)題。把這個(gè)相應(yīng)問(wèn)題的解

作為的近似值。這樣求得的就是上述初值問(wèn)題在節(jié)點(diǎn)上的數(shù)值解。一般說(shuō)來(lái)，不同的離散化導(dǎo)致不同的方法。第二節(jié)歐拉方法1、歐拉方法的導(dǎo)出對(duì)于上面的初值問(wèn)題，若將函數(shù)在處的導(dǎo)數(shù)用差商來(lái)表示，則有，…解：根據(jù)歐拉公式可以得到：此外，方程的真解為：2、舉例例1用歐拉法求初值問(wèn)題取步長(zhǎng)h=0.2。i012345xiyi

01.00000.20001.20000.40001.37330.60001.53150.80001.68111.00001.8269求解結(jié)果如下：y(xi)

1.00001.18321.34161.48321.61251.7321

i=y(xi)-

yi

00.01680.03170.04830.06860.0949根據(jù)已知條件：曲線y(x)上的點(diǎn)（x0,y0）及該點(diǎn)處曲線的導(dǎo)數(shù)f(x0,y0)，則可以得到過(guò)該點(diǎn)的直線：該直線與x＝x1的交點(diǎn)P1，則P1的縱坐標(biāo)y1為：就用y1作為y(x1)的近似值…逐次進(jìn)行后可以得到一條折線P0P1…Pn，該折線看作是初值問(wèn)題的積分曲線的近似，因此歐拉方法也稱(chēng)為歐拉折線法。3、歐拉方法的幾何意義歐拉法雖然形式簡(jiǎn)單，計(jì)算方便。但從上述幾何意義上得知，由Euler法所得的折線明顯偏離了積分曲線，可見(jiàn)此方法非常粗糙,精度也低，即誤差太大。特別當(dāng)?shù)那€曲率較大時(shí)，歐拉法的效果更差。二、梯形格式對(duì)微分方程兩邊積分得，歐拉方法左矩形公式后退的歐拉方法右矩形公式梯形方法梯形公式和歐拉公式相比較，梯形公式精度有所提高；但梯形公式是隱式的，在計(jì)算yi+1時(shí)候也只用到前一步的值yi，若yi已知，將yi代入公式求解時(shí)，一般不能直接得到y(tǒng)i+1，而需要通過(guò)其他方法（比如迭代法）求解，計(jì)算量較大。實(shí)際中，將歐拉公式和梯形公式聯(lián)合使用，從而得到改進(jìn)的歐拉方法。三、改進(jìn)的歐拉方法將歐拉公式和梯形公式聯(lián)合使用，先用歐拉公式得出一個(gè)y(xi+1)的粗糙近似值，稱(chēng)為預(yù)估值，然后對(duì)預(yù)估值使用梯形公式對(duì)它進(jìn)行精確化，得到較為精確的近似值yi+1，稱(chēng)之為校正值，計(jì)算公式為：這樣的預(yù)估校正系統(tǒng)稱(chēng)為改進(jìn)的歐拉方法為了便于編寫(xiě)程序，常將上面的公式改寫(xiě)為如下式：

改進(jìn)的歐拉方法與梯形方法具有同樣的精度級(jí)別，但改進(jìn)的歐拉方法為顯式格式，使用更方便。解：根據(jù)改進(jìn)的歐拉公式得：此外，方程的真解為：舉例例1用改進(jìn)的歐拉法求初值問(wèn)題取步長(zhǎng)h=0.2。i012345

01.00000.20001.18670.40001.34830.60001.49370.80001.62791.00001.7542xiyi求解結(jié)果如下：y(xi)

1.00001.18321.34161.48321.61251.7321

i=y(xi)-

yi

00.00350.00670.01050.01540.0222例：某傘兵與降落傘質(zhì)量為90kg，當(dāng)傘張開(kāi)時(shí)以12m/s速度垂直下降。假設(shè)空氣阻力與傘的下降速度成正比，在速度為6m/s時(shí)，測(cè)得的空氣阻力為353N。試求傘兵開(kāi)傘后，第1秒末、第2秒末、直到第5秒末時(shí)各時(shí)刻的速度。012.00001.000013.96002.000014.63953.000014.87504.000014.95675.000014.9850