clear;
%程序C8-3  求解常微分方程組的四階R-K法
y(:,1)=[1,0]';t(1)=0;
h=0.5;n=1;
while t<=4
    K1=h*f_um(y(:,n));
    K2=h*f_um(y(:,n)+0.5*K1);
    K3=h*f_um(y(:,n)+0.5*K2);
    K4=h*f_um(y(:,n)+K3);
    y(:,n+1)=y(:,n)+(K1+2*K2+2*K3+K4)/6;
    t(n+1)=n*h;
    n=n+1;
end
plot(t,y(1,:),'-',t,y(2,:),':')
y