歐拉法、梯形法和龍格-庫塔解微分方程

1、- - 歐拉法、梯形法和龍格 - 庫塔一、解方程： ? x( ) (0）=1 二、算出方程的解析解為：y 2 - ?-4?2三、實驗原理 : 1. 歐拉法原理 : 將區(qū)間a,b 分成 n 段，那么方程在第xi點有y （xi）= f(xi,y（x） ）, 再用向前差商近似代替導數(shù)則為：，在這里,h 是步長，即相鄰兩個結(jié)點間的距離。因此可以根據(jù)xi 點和 i點的數(shù)值計算出+來: ,i=0,1 ，,n 這就是歐拉格式，若初值yi+ 1是已知的，則可依據(jù)上式逐步算出數(shù)值解1,y.yn 2. 梯形法原理：將向前歐拉公式中的導數(shù) （i,y )改為微元兩端 f(xi,yi）和 f （i+ ,yi+ )的平均

2、 , 即梯形公式。3.龍格庫塔方法的基本思想: - - 在區(qū)間 xn ， n+ 內(nèi)多取幾個點 ,將他們的斜率加權(quán)平均 , 作為導數(shù)的近似。令初值問題表述如下。則, 對于該問題的 r4 由如下方程給出：其中這樣, 下一個值 (n+1)由現(xiàn)在的值 (）加上時間間隔 (h) 和一個估算的斜率的乘積決定。該斜率是以下斜率的加權(quán)平均：k是時間段開始時的斜率；?2 是時間段中點的斜率 ,通過歐拉法采用斜率 1 來決定 y 在點 tn + h/2 的值；? k3也是中點的斜率 , 但是這次采用斜率k2 決定 y 值; ?是時- - 間段終點的斜率，其y 值用 k3 決定。當四個斜率取平均時 , 中點的斜率有

3、更大的權(quán)值：r4 法是四階方法，也就是說每步的誤差是h階，而總積累誤差為h階。四、歐拉法、梯形法和龍格庫塔的實現(xiàn)代碼：h=0.; =0:h：1; y1zer s（size(x); y1(1)= ；=zro (ze(x ）) ；y(1) 1; y3=zeo(size(x)）; y3( ）=; or i =2:len th(x ）y1( ) = y1( 1） + h*x(i1-1)*(2-y( 1) ） ； %歐拉法1= 8*x(i1-1)*(2y2（11) ；% 梯形法- - m2 =8* (i )(2-y2(i-1) h*m1); y2(）=y2(i -1)+h*( 1m )/2; % 梯形法

4、公式k1=8*x （11) (2- 2(i1-1 ）; 龍格庫塔 k2=8*(x（1-1)+h/2 ）( (y (i1 ) *k /2 ）); 3=*( （1-1)+h/2)*(2（2(i -1)+h k/2 ）; k4=8(x(i1-1)+h）（ -( (i1- )+h*k ) ；y3(i ) y3（i1-1) (1+2*k2+2*k3+k4)*h/ ； % 龍格- 庫塔公式=-e p( 4*(x. ) ；解析解ot(x,y ,x ，y,x ， , , ) % 解析解與數(shù)值解圖像eged(y1 , y2 , y3 ,y4 ) plot(x,4-y1,x,y4 ,x ，y4-y3) 解析解與數(shù)

5、值解誤差圖像leg n(y4 1, y4-y2，y4y3 ) 五、圖像 : 1. 解析解 y4 與各個數(shù)值解 y1,y ，y3 的圖像：- - (1）當步長 h=0.1 時：()當步長 =0.05 時: - - (3) 當步長 h=. 1 時: - - 結(jié)論: 三個圖中，方程的解析解都是y4。 y1表示歐拉法 , y2 表示梯形法 ,y3 表示龍格庫塔法。從上面三個圖像可以看出， y與 y4 的擬合度比其他兩個都好, 即龍格庫塔法得到的數(shù)值解與真實解的擬合度比歐拉法和梯形法都高。同時, 從上圖可以看出，步長 h=0.1 比 h=.5 和 h=1 所得到的的解更接近于真實解，這也說明了, 步長越小，區(qū)間內(nèi)的點越多, 估計解也就越接近與真實解。2. 解析解 y4 與各個數(shù)值解 1,y ,y3 的誤差圖像 : (）當步長 h=.1 時- - （2）當步長 h0.05 - - (3) 當步長 =0. 1 時- - 結(jié)論: 同樣的 , 上面三個誤差圖可以看出:y -y3 表示龍格- 庫塔法與真實解的誤差,y 2 表示梯形法與真實解的誤差,y 1 表示歐拉法與真實解的誤差。從