求解四阶龙格库塔方程问题附方程和程序代码

gigicool333

用公式y(k+1)=y(k)+h/6(k1+2k2+2k3+k4)计算[a,b]上的初值问题y'=f(t,y),y(a)=y0的近似解
function R=rk4(f,a,b,ya,M)
h=(b-a)/M;
T=zeros(1,M+1);
Y=zeros(1,M+1);
T=a:h:b;
Y(1)=ya;
for j=1:M
   k1=h*feval(f,T(j),Y(j));
   k2=h*feval(f,T(j)+h/2,Y(j)+k1/2);
   k3=h*feval(f,T(j)+h/2,Y(j)+k2/2);
   k4=h*feval(f,T(j)+h,Y(j)+k3);
   Y(j+1)=Y(j)+(k1+2*k2+2*k3+k4)/6;
end
R=[T' Y'];
书上这样写的程序为啥我运行不了提示错误；
正确的话应该怎样写；
小弟渴望好心人的帮助!

haisujiang

我不太懂,不过感觉用C写也行吧!

eight

原帖由 gigicool333 于 2007-5-12 00:56 发表
用公式y(k+1)=y(k)+h/6(k1+2k2+2k3+k4)计算[a,b]上的初值问题y'=f(t,y),y(a)=y0的近似解
function R=rk4(f,a,b,ya,M)
h=(b-a)/M;
T=zeros(1,M+1);
Y=zeros(1,M+1);
T=a:h:b;
Y(1)=ya;
for j=1:M
   k1= ...

请先看看置顶贴：聚宝盆，然后把你的问题叙述清楚，特别是出错信息