四階龍格庫(kù)塔實(shí)驗(yàn)報(bào)告

1、三、四階 Runge-Kutta 法求解常微分方程一、龍格庫(kù)塔法的思想根據(jù)第九章的知識(shí)可知道，Euler 方法的局部截?cái)嗾`差是 ，而當(dāng)用 2 ( ) O hEuler 方法估計(jì)出1 , ( ) (1) n n n n y y hf x y ? ? ? ??????????????再用梯形公式1 1 1 [ ( , ) ( , )] (2) 2n n n n n nh y y f x y f x y ? ? ? ? ? ? ???????

2、?進(jìn)行校正，即采用改進(jìn) Euler 方法得出數(shù)值解的截?cái)嗾`差為 。 3 ( ) O h由 Lagrange 微分中值定理'1 1 ( ) ( ) ( )( ) ( ) ( , ( )) (3) n n n n n y x y x y x x y x hf y ? ? ? ? ? ? ? ? ? ? ?????記 ，得到 * ( , ( )) k hf y ? ? ?*1 ( ) ( ) (4) n n y x y x k ? ?

3、 ? ??????????這樣只要給出一種計(jì)算 的算法，就能得到相應(yīng)的計(jì)算公式。 * k用這種觀點(diǎn)的來(lái)分析 Euler 方法和改進(jìn) Euler 方法，Euler 方法的迭代公式可改寫(xiě)為1 11 ( , )n nn ny y kk hf x y? ? ??改進(jìn) Euler 方法的預(yù)報(bào)-校正公式可改寫(xiě)為1 1 21 2 11 ( ) 2( , ), ( , )n nn n n ny y k kk hf x y k hf x h y k? ?

4、 ? ?? ? ? ?Euler 方法實(shí)際上是用一個(gè)點(diǎn)處的值 近似 ，而改進(jìn) Euler 方法是用兩個(gè)點(diǎn)處 1 k * k的值 ，和 ，做算術(shù)平均值近似 自然改進(jìn) Euler 方法要優(yōu)于 Euler 方法。 1 k 2 k * k因此，可以想到假如在 內(nèi)多預(yù)報(bào)幾個(gè)點(diǎn)值 ，并用他們的加權(quán)平均值作 1 [ , ] n n x x ? i k為 的近似值，則有可能構(gòu)造出具有更高精度的計(jì)算公式，這就是 Runge-Kutta * k法的基本思想

5、。二、四階龍格庫(kù)塔法3.2 程序解釋及使用該算法可以對(duì)一階微分方程，一階微分方程組進(jìn)行有效的求解。ydot_fun 為一階微分方程的函數(shù)，x0 為初始點(diǎn)，y0 為初始向量，h 為步長(zhǎng)，N 為區(qū)間的 等分?jǐn)?shù)，x 為 Xn 構(gòu)成的向量，y 為 Yn 構(gòu)成的矩陣。程序調(diào)用方法：1，先編寫(xiě)要求解的一階微分方程或方程組的函數(shù)文件 ydot_fun.m 文件，將該文件和 Runge 文件放到同一個(gè)目錄下。 2. 調(diào)用求解程序，[x,y]=Runge

6、(@dot_fun,x0,y0,h,N)，運(yùn)行后即可得出結(jié)果?；蛘哂脙?nèi)部函數(shù)調(diào)用： 輸入：ydot_fun=(x,y)[][x,y]= Runge(ydot_fun,x0,y0,h,N)3.3 實(shí)例求解課本 304 頁(yè) 9.2 題目：用標(biāo)準(zhǔn) 4 級(jí) 4 階 R-K 法求解， ，''' '' '' ''2 3(0) 1, (0) 3, (0) 2y y y y xy

7、y y? ? ? ? ? ? ? ? ? ? ? ? ? ?取步長(zhǎng) h=0.1，計(jì)算 的近似值，并與解析解 作比較。 (1) y ( ) 2 1 x y x xe x ? ? ?解：首先將三階方程改寫(xiě)成微分方程組的形式：令 得如下微分方程組' ' '' '1 2 1 3 2 , , y y y y y y y y ? ? ? ? ?' 1 2' 2 3' ' 

8、9;3 2 1 1' ''1 2 32 3(0) 1, (0) 3, (0) 2y yy yy y y y xy y y? ?? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?在 ydot_fun.m 文件中編寫(xiě)待求解微分方程組，調(diào)用計(jì)算程序，保留 5 位小數(shù) 得：表 3-1 三階微分方程求解結(jié)果x 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1 y -1.00

9、000 -0.68948 -0.35572 0.00496 0.39673 0.82436 1.29327 1.80962 2.38042 3.01363 3.71827 2 y 3.00000 3.21569 3.46568 3.75481 4.08855 4.47308 4.91538 5.42337 6.00596 6.67323 7.43655 3 y 2.00000 2.32086 2.68708 3.10467 3.5803