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

1、第六章常微分方程的數(shù)值解法§6.0 引言§6.1 算法構(gòu)造的主要途徑 §6.2 Runge-Kutta Method算法§6.3 線性多步法§6.4 線性多步法的一般形式§6.5 算法的穩(wěn)定性、收斂性§6.0 引 言1主要考慮如下的一階常微分方程初值問(wèn)題的求解:微分方程的解就是求一個(gè)函數(shù)y=y(x)，使得該函數(shù)滿足微分方程并且符合初值條件。2. 例如微分方程:xy'-2y=4x ；初始條件: y(1)=-3。于是可得一階常微分方程的初始問(wèn)題。顯然函數(shù)y(x)=x2-4x滿足以上條件,因而是該初始問(wèn)題的微分方程的解。3

2、. 但是,只有一些特殊類型的微分方程問(wèn)題能夠得到用解析表達(dá)式表示的函數(shù)解，而大量的微分方程問(wèn)題很難得到其解析解，有的甚至無(wú)法用解析表達(dá)式來(lái)表示。因此，只能依賴于數(shù)值方法去獲得微分方程的數(shù)值解。4 微分方程的數(shù)值解：設(shè)微分方程問(wèn)題的解y(x)的存在區(qū)間是a,b，初始點(diǎn)x0=a，將a,b進(jìn)行劃分得一系列節(jié)點(diǎn)x0 , x1 ,.,xn，其中a= x0< x1<< xn =b。y(x)的解析表達(dá)式不容易得到或根本無(wú)法得到，我們用數(shù)值方法求得y(x)在每個(gè)節(jié)點(diǎn)xk的近似值y(xk)，即yy(xk)，這樣y0 , y1 ,.,yn稱為微分方程的數(shù)值解。如圖所示:ab x0 x1 x2 .

3、 xn-1 xn§6.1 算法構(gòu)造的主要途徑1 歐拉公式1.1 構(gòu)造的思想：微分方程初值問(wèn)題: 利用差商代替一階導(dǎo)數(shù)，即，則。于是，可求出y(x1)的近似值y1，同樣地，可利用x1處的微分方程可得：一般地，利用在xn處的微分方程可得： 此式稱為歐拉公式。1.2 幾何意義：對(duì)于微分方程y'=2(x+1),其通解是y=(x+1)2+c,是一個(gè)曲線族，當(dāng)給定初值條件y(0)=2，其特解為y=(x+1)2+1。如圖所示：由y(0)=2，過(guò)該曲線上一點(diǎn)(0,2)作曲線的切線，其斜率，切線為：,因此可計(jì)算出y1,如此，可根據(jù)：，故歐拉法又稱歐拉折線法。1.3 算例：例：解：h=0.2 ,

4、 xi=1+ih y1= y0+hf(x0，y0)=-1+0.2× y2= y1+hf(x1，y1)=-1+0.2× y3= y2+hf(x2，y2)= -0.9333+0.2×y4= y3+hf(x3，y3)=0.8+0.2×y5= y4+hf(x4，y4)=0.6+0.2×y6= y5+hf(x5，y5)=0.3333+0.2×精確解為：y=x2-2xxky(xk)ykek1.2-0.96-10.041.4-0.84-0.93330.09331.6-0.64-0.80.161.8-0.36-0.60.242.00-0.33330.

5、33332.20.4400.44 可以看出誤差隨著計(jì)算在積累。1.4 Euler法的特點(diǎn)和誤差迭代格式：特點(diǎn)：(1)單步方法；(2)顯式格式；(3)局部截?cái)嗾`差為。局部截?cái)嗾`差：當(dāng)時(shí)，由按照歐拉方法計(jì)算來(lái)的的誤差稱為局部截?cái)嗾`差。即，是局部截?cái)嗾`差。如：歐拉法得：因此，局部截?cái)嗾`差是。2改進(jìn)Euler法2.1方法構(gòu)造在微分方程初值問(wèn)題,對(duì)其從到進(jìn)行定積分得:將右端的定積分用梯形公式來(lái)進(jìn)行近似計(jì)算。用和來(lái)分別代替和得計(jì)算格式:這就是改進(jìn)歐拉方法。2.2 顯式格式和隱式格式在歐拉式中每一步計(jì)算已知,直接用格式可以計(jì)算出,此類格式稱為顯式格式。而在改進(jìn)歐拉方法中在每一步計(jì)算中是未知,待求的,未知量在

6、中這是一個(gè)方程,如是非線形或超越函數(shù),此方程是無(wú)法直接解出來(lái)(要依靠迭代法才能解出)。這類格式稱為隱式格式。2.3 算例例: 用改進(jìn)歐拉方法求解。解：解得:注意：由于,是線形函數(shù)可以從隱式格式中解出。問(wèn)題的精確解是歐拉方法誤差3預(yù)測(cè)校正方法由于改進(jìn)Euler法是隱式格式，無(wú)法從格式中直接求出必須要解方程。下面用預(yù)測(cè)校正方法來(lái)求隱式格式中的。 預(yù)測(cè)值: 校正值: 此式相當(dāng)于對(duì)隱式格式求時(shí)采用迭代的方法,用歐拉格式得到的作為初始值迭代公式迭代一次而已,此公式代入后得:如改寫成平均的形式為:§6.2 龍格-庫(kù)塔法Runge-Kutta Method1 龍格-庫(kù)塔法的思想1.1 考慮微分方程

7、的初值問(wèn)題:根據(jù)微分中值定理有：，其中0<<1。于是即我們稱為y(x)在區(qū)間xk, xk+1上的平均斜率，記作K，其中，是存在但是未知的。因此，如何對(duì)平均斜率K進(jìn)行近似計(jì)算，相應(yīng)地就得到一種近似公式，或稱為微分方程的一種計(jì)算格式。1.2 例如：用f(xk,yk)作為平均斜率K的近似值就得到歐拉格式：。用作為平均斜率K的近似值就得到比歐拉格式高一階精度的格式，即，改進(jìn)歐拉格式的預(yù)測(cè)-校正方法:1.3 啟發(fā)(Motivative)能否在二維平面中xxk, xk+1, yyk, yk+1上多找一些f(x,y)點(diǎn)，在這些點(diǎn)上作函數(shù)值的平均數(shù)并以此作為平均斜率K的近似值。由于自由度的增加，使

8、得的p能夠提高，從而達(dá)到提高精度的目標(biāo)，這就是龍格-庫(kù)塔法的基本思想。1.4 R-K的一般形式 在二維平面區(qū)域上取N個(gè)點(diǎn)，得到公式:，其中Kj是二元函數(shù)f(x,y)在這N個(gè)點(diǎn)上的值。其中是待定常數(shù)。2 二階龍格-庫(kù)塔法2.1計(jì)算格式: 應(yīng)用Taylor公式求參數(shù)：。將f(x,y)在(xk , yk)上展開(kāi)又因 以上兩式代入中，又由得比較兩邊得四個(gè)未知量,三個(gè)方程,有無(wú)窮多組解，但局部截?cái)嗾`差均為,是二階精度的。例如:(1) 取,則,就是改進(jìn)歐拉公式的預(yù)測(cè)-校正方法(或二階龍格-庫(kù)塔法)。(2) 取，則，則,得二階中點(diǎn)法：3 三階龍格庫(kù)塔方法3.1 計(jì)算格式: 參數(shù)有,共有八個(gè)參數(shù)。將f和均在點(diǎn)

9、Tayor展開(kāi),有利用,按展開(kāi)式代入(高于的不計(jì))得:以上展開(kāi)式代入中比較得:(1)(2)(3)(4)(5)(6)(7)(8)上述八個(gè)方程有六個(gè)獨(dú)立方程。八個(gè)未知量,六個(gè)方程有無(wú)窮多組解,但其截?cái)嗾`差均為,具有三階精度。例如:，得:此方法稱為三階Kutta方法。例如:,得:此方法稱為三階Heun法。4 四階龍格庫(kù)塔方法4.1 顯式四階龍格-庫(kù)塔方法的一般形式是:利用與三階龍格-庫(kù)塔類似的處理方法，經(jīng)過(guò)推導(dǎo)，可得古典四階龍格-庫(kù)塔算法:其特點(diǎn)是：截?cái)嗾`差為，階數(shù)為四階精度，為顯式四階龍格庫(kù)塔方法。計(jì)算均要計(jì)算四個(gè)f的值。一般一階常微分方程初值問(wèn)題的解均用四階龍格-庫(kù)塔方法進(jìn)行計(jì)算，其精度滿足實(shí)際

10、問(wèn)題的精度要求。數(shù)值例子：例：初值問(wèn)題用四階古典RungeKutta方法，。0.21.18322931.18321600.41.34168031.34164080.61.48328381.48323970.81.61251721.61245151.01.73214631.7320508例： §6.3 線性多步法1 單步方法和多步方法單步方法：歐拉方法、改進(jìn)歐拉方法、龍格-庫(kù)塔方法均是單步方法，即在每一步計(jì)算時(shí)，只要前面一個(gè)值已知的條件下就可以計(jì)算出。單步方法特點(diǎn)：(1)可以自成系統(tǒng)進(jìn)行直接計(jì)算，因?yàn)槌跏紬l件只有一個(gè)已知，由可以計(jì)算，不必借助于其它方法，這種單步方法是自開(kāi)始的。(2)如

11、果格式簡(jiǎn)單如歐拉方法，則只有一階精度，如果提高精度，則計(jì)算很復(fù)雜，如RungeKutta方法。(3)公式的構(gòu)造推導(dǎo)也很復(fù)雜。多步方法：利用前面已知計(jì)算出來(lái)的，由已經(jīng)計(jì)算好的個(gè)值來(lái)計(jì)算，這樣可以提高算法的精度，該方法稱為多步方法，利用k個(gè)值計(jì)算，稱為k步方法。多步方法的特點(diǎn)：(1)因初始條件只有一個(gè)，運(yùn)用多步方法要借助高階的單步方法來(lái)開(kāi)始。例如，已知用單步的四階Runge-Kutta方法計(jì)算，再計(jì)算，再由計(jì)算，用單步方法有后運(yùn)用四階的四步方法，由計(jì)算；由計(jì)算；由計(jì)算；一直下去， 可以用多步方法，并且始終達(dá)到四階精度。(2)多步方法比較簡(jiǎn)單，只要在這四個(gè)點(diǎn)的函數(shù)值的線性組合，而且每步中后三個(gè)函數(shù)值

12、下一步還可使用。2顯式Adams方法：2.1 構(gòu)造的思想考慮微分方程初值問(wèn)題將微分方程在上積分，2.2 顯式Adams方法若已知來(lái)計(jì)算，簡(jiǎn)記，用的拉格朗日插值多項(xiàng)式代替f， 截?cái)嗖糠?，用等距步長(zhǎng)，上面積分很簡(jiǎn)單，得到的方法就是顯式四階Adams方法。以上可以看到該方法的局部截?cái)嗾`差是因而是四階精度的。例如：解：取，首先用四階RungeKutta方法來(lái)起步，計(jì)算出，,下面不必用RungeKutta方法，而開(kāi)始用四階Adams方法。(1)、求(2)、求只要補(bǔ)算(3)、求只要補(bǔ)算現(xiàn)列表看用Adams方法求出的誤差，精解為0.81.61142311.61245151.0210-31.01.729840

13、31.73205082.2110-31.21.84066161.84390893.2410-32.3 隱式Adams方法用作為插值結(jié)點(diǎn)，由于也是插值結(jié)點(diǎn)，必帶來(lái)從而導(dǎo)致是隱式格式。用插值多項(xiàng)式來(lái)代替積分中的得：截掉得近似公式：得：，從而得四階隱式Adams方法。因，而是未知的，故這是隱式格式。隱式格式的解法用預(yù)測(cè)-校正法：用顯式格式作為預(yù)測(cè)值，再用隱式格式來(lái)校正。預(yù)測(cè)值：校正值：§6.4 線性多步法的一般形式微分方程初值問(wèn)題: 1 k步線性法的一般結(jié)構(gòu)是常數(shù)線性多步法 , k步法不同時(shí)為零 k步法 ,隱式結(jié)構(gòu);否則為顯式結(jié)構(gòu)。2 k步p階線性法的推導(dǎo)由: 定義算子：由Taylor展開(kāi)

14、右端有代入上式：其中，因此，若有次可微，令,則因而算子對(duì)應(yīng)的線性多步法為p階k步方法。3 算例例：考慮2步法記于是：故一般的二步法為：同時(shí)可以算出：當(dāng)時(shí)，上述公式為三階的。當(dāng)時(shí)，上述公式為四階的。 §6.5 收斂性與穩(wěn)定性1 定義收斂性： 若某算法對(duì)于任意固定的 x = xn= x0 + n h，當(dāng) h®0 ( 同時(shí) n ® ¥) 時(shí),有 yn® y( xn)，則稱該算法是收斂的。 1.2 算例例：對(duì)于初值問(wèn)題 考察歐拉顯式格式的收斂性。解：該問(wèn)題的精確解為 歐拉公式為即：對(duì)任意固定的 x = xn = n h ，有由 2 穩(wěn)定性 例：考察初值

15、問(wèn)題 在區(qū)間0, 0.5上的解。分別用歐拉顯式、隱式格式和改進(jìn)的歐拉格式計(jì)算數(shù)值解。節(jié)點(diǎn) xi歐拉顯式歐拉隱式改進(jìn)歐拉法精確解0.00.10.20.30.40.51.0000-2.00004.0000-8.0000 1.6000´101-3.2000´1011.00002.5000´10-16.2500´10-21.5625´10-23.9063´10-39.7656´10-41.00002.50006.25001.5626´1013.9063´1019.7656´1011.00004.9787

16、´10-22.4788´10-31.2341´10-46.1442´10-63.0590´10-7 誤差包括截?cái)嗾`差(算法理論誤差)和 舍入誤差兩個(gè)部分。后者由計(jì)算機(jī)字長(zhǎng)等決定，屬于穩(wěn)定性問(wèn)題。若某算法在計(jì)算過(guò)程中任一步產(chǎn)生的誤差在以后的計(jì)算中都逐步衰減，則稱該算法是絕對(duì)穩(wěn)定的 討論數(shù)值方法的穩(wěn)定性，通常只考慮試驗(yàn)方程 當(dāng)步長(zhǎng)取為 h 時(shí)，將某算法應(yīng)用于上式，并假設(shè)只在初值產(chǎn)生誤差 ，則若此誤差以后逐步衰減，就稱該算法相對(duì)于絕對(duì)穩(wěn)定，的全體構(gòu)成絕對(duì)穩(wěn)定區(qū)域。有時(shí)，我們稱算法A 比算法B 穩(wěn)定，就是指 A 的絕對(duì)穩(wěn)定區(qū)域比 B 的絕對(duì)穩(wěn)定區(qū)域大。2.1 算例：考察顯式歐拉法算法的穩(wěn)定性解：由 得：于是：因此，要保證初始誤差e0 以后逐步衰減，必須滿足：0-1-2ReIm例：考察隱式歐拉法 可見(jiàn)