资源描述:
《数学实验题目3 四阶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
3、-1],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
4、('x.^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)'))T
5、estRK4('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(
6、'ode7',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(
7、fun,xEnd,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]=siz
8、e(ini