4階runge-kutta法求解一階常微分方程

《4階runge-kutta法求解一階常微分方程》由會員上傳分享，免費(fèi)在線閱讀，更多相關(guān)內(nèi)容在行業(yè)資料-天天文庫。

1、《MATLAB語言及應(yīng)用》大作業(yè)姓名：學(xué)號：學(xué)院：班級：題目編號：2013年10月134階Runge-Kutta法求解一階常微分方程。一、Runge-Kutta法的數(shù)學(xué)理論龍格-庫塔(Runge-Kutta)方法是一種在工程上應(yīng)用廣泛的高精度單步算法。由于此算法精度高，采取措施對誤差進(jìn)行抑制，所以其實現(xiàn)原理也較復(fù)雜。該算法是構(gòu)建在數(shù)學(xué)支持的基礎(chǔ)之上的。龍格庫塔方法的理論基礎(chǔ)來源于泰勒公式和使用斜率近似表達(dá)微分，它在積分區(qū)間多預(yù)計算出幾個點(diǎn)的斜率，然后進(jìn)行加權(quán)平均，用做下一點(diǎn)的依據(jù)，從而構(gòu)造出了精度更高的數(shù)值積分計算方法。如果預(yù)先求兩個點(diǎn)的斜率就是二階龍格庫塔法，如果預(yù)先取四個

2、點(diǎn)就是四階龍格庫塔法。一階常微分方程可以寫作：y'=f(x,y),使用差分概念。(Yn+1-Yn)/h=f(Xn,Yn)推出（近似等于，極限為Yn'）Yn+1=Yn+h*f(Xn,Yn)另外根據(jù)微分中值定理，存在0

3、/2)*K2);K4=f(Xn+h,Yn+h*K3);Yn+1=Yn+h*(K1+2K2+2K3+K4)*(1/6)一、Runge-Kutta的算法和流程圖在龍格-庫塔法中，四階龍格-庫塔法的局部截斷誤差約為0(h5),被廣泛應(yīng)用于解微分方程的初值問題。其算法公式為:流程圖：(1)、四階龍格-庫塔方法流程圖:（2）、實例求解流程圖：一、Runge-Kutta的Matlab實現(xiàn)function[x,y]=runge_kutta1(ufunc,y0,h,a,b)%參數(shù)表順序依次是微分方程組的函數(shù)名稱，初始值向量，步長，時間起點(diǎn)，時間終點(diǎn)（參數(shù)形式參考了ode45函數(shù)）n=floo

4、r((b-a)/h);%求步數(shù)x(1)=a;%時間起點(diǎn)y(:,1)=y0;%賦初值，可以是向量，但是要注意維數(shù)forii=1:nx(ii+1)=x(ii)+h;k1=ufunc(x(ii),y(:,ii));k2=ufunc(x(ii)+h/2,y(:,ii)+h*k1/2);k3=ufunc(x(ii)+h/2,y(:,ii)+h*k2/2);k4=ufunc(x(ii)+h,y(:,ii)+h*k3);y(:,ii+1)=y(:,ii)+h*(k1+2*k2+2*k3+k4)/6;%按照龍格庫塔方法進(jìn)行數(shù)值求解end一、Runge-Kutta的算例實現(xiàn)例:求解常微分方程d

5、(y)/d(x)=-2*y+2*x*x+2*x，0≤x≤0.5，y(0)=1.編輯m文件Untitled3：fun=inline('-2*y+2*x*x+2*x');[x,y]=ode45(fun,[0,0.5],1)>>Untitled3x=00.01250.02500.03750.05000.06250.07500.08750.10000.11250.12500.13750.15000.16250.17500.18750.20000.21250.22500.23750.25000.26250.27500.28750.30000.31250.32500.33750.3500

6、0.36250.37500.38750.40000.41250.42500.43750.45000.46250.47500.48750.5000y=1.00000.97550.95190.92910.90730.88640.86630.84710.82870.81120.79440.77850.76330.74890.73530.72240.71030.69890.68830.67830.66900.66050.65260.64540.63880.63290.62770.62310.61910.61570.61300.61090.60930.60840.60800.60830

7、.60910.61040.61240.61480.6179