rk4.m

常微分方程的龙格库塔法系列程序

function rk4(f,u0,a,b,h)
%Runge_Kutta法,s=4
%调用格式:rk4(function_handle,u0,a,b,h)
%test: du/dt=u/2;u0=1;h=0.1;a=0;b=1;
%调用:先定义微分方程f=@(t,u) u/2;rk4(f,1,0,1,0.1)

%2007-03-15
tic;
len=1+(b-a)/h;
u=zeros(1,len);
t=a:h:b;
u(1)=u0;
for i=2:len,
    K1=f(t(i-1),u(i-1));
    K2=f(t(i-1),u(i-1)+h*K1/2);
    K3=f(t(i-1),u(i-1)+h*K2/2);
    K4=f(t(i-1),u(i-1)+h*K3);
    u(i)=u(i-1)+h*(K1+2*K2+2*K3+K4)/6;
end
plot(t,u,'b+');
hold on;
line(t,u,'Color','b');
fun=@(t) exp(t/2);
fplot(fun,[0,1],'r');
error=sqrt(sum(abs(fun(t)-u).^2))
toc