四阶Runge-Kutta方法

《四阶Runge-Kutta方法》由会员上传分享，免费在线阅读，更多相关内容在教育资源-天天文库。

1、实验题目3四阶Runge-Kutta方法摘要一阶常微分方程初值问题（6.1）的数值解法是近似计算中很重要的部分。常微分方程初值问题的数值解法是求方程（6.1）的解在点列上的近似值，这里是到的步长，一般略去下标记为。常微分方程初值问题的数值解法一般分为两大类：（1）单步法：这类方法在计算时，只用到、和，即前一步的值。因此，在有了初值以后就可以逐步往下计算。典型方法如龙格–库塔方法。（2）多步法：这类方法在计算时，除用到、和以外，还要用，即前面步的值。典型方法如Adams方法。经典的方法是一个四阶的方法，它的计算公式是：（6.2）方法的优点是：单步法、精度高，计算过程便

2、于改变步长，缺点是计算量较大，每前进一步需要计算四次函数值。前言利用四阶龙格-库塔方法求解微分方程的初值问题程序设计流程开始输入x0,y0,h,Nx1=x0+hk1=f(x0,y0),k2=f(x0+h/2,y0+hk1/2)k3=f(x0+h/2,y0+hk2/2),k4=f(x1,y0+hk3)y1=y0+h(k1+2k2+2k3+k4)/6n=1输出x1,y1n=N?结束n=n+1x0=x1,y0=y1否是龙格-库塔法流程图问题1（1）TestRK4('ode1',1,[0-1],5,inline('-x-1'))TestRK4('ode1',1,[0-1],

3、10,inline('-x-1'))TestRK4('ode1',1,[0-1],20,inline('-x-1'))（2）TestRK4('ode2',1,[01],5,inline('1./(x+1)'))TestRK4('ode2',1,[01],10,inline('1./(x+1)'))TestRK4('ode2',1,[01],20,inline('1./(x+1)'))问题2（1）TestRK4('ode3',3,[10],5,inline('x.^2.*(exp(x)-x)'))TestRK4('ode3',3,[10],10,inline('x.^

4、2.*(exp(x)-x)'))TestRK4('ode3',3,[10],20,inline('x.^2.*(exp(x)-x)'))（1）TestRK4('ode4',3,[1-2],5,inline('2*x./(1-2*x)'))TestRK4('ode4',3,[1-2],10,inline('2*x./(1-2*x)'))TestRK4('ode4',3,[1-2],20,inline('2*x./(1-2*x)'))问题3（1）TestRK4('ode5',1,[01/3],5,inline('x.^2+1/3*exp(-20*x)'))TestRK4

5、('ode5',1,[01/3],10,inline('x.^2+1/3*exp(-20*x)'))TestRK4('ode5',1,[01/3],20,inline('x.^2+1/3*exp(-20*x)'))（1）TestRK4('ode6',1,[01],5,inline('exp(-20*x)+sin(x)'))TestRK4('ode6',1,[01],10,inline('exp(-20*x)+sin(x)'))TestRK4('ode6',1,[01],20,inline('exp(-20*x)+sin(x)'))（1）TestRK4('ode7',

6、1,[00],5,inline('exp(x).*sin(x)'))TestRK4('ode7',1,[00],10,inline('exp(x).*sin(x)'))TestRK4('ode7',1,[00],20,inline('exp(x).*sin(x)'))实验所用函数function[x,y]=RK4ODE(fun,xEnd,ini,h)%RK4ODE用四阶Runge-Kutta法解初值问题dy/dx=f(x,y),y(x0)=y0，在x处y的值%%Synopsis:[x,y]=RK4ODE(fun,xEnd)%[x,y]=RK4ODE(fun,xEnd

7、,ini)%[x,y]=RK4ODE(fun,xEnd,ini,h)%%Input:fun=(string)初值问题的函数%xEnd=使用Euler法的截止点%ini=(optional)初始条件[x0y0]，默认为[00]%h=(optional)步长，默认为0.05%%Output:y=初值问题在x处y的近似值ifnargin<3ini=[00];%若未给初始条件，将初始条件设为[00]endifnargin<4h=0.05;%若未给出步长，将步长设为0.05endini=ini(:);%将初始条件转为列向量，便于判断是否正确[m,n]=size(ini);