function runge_kutta_stability

h=2.;
t=0:h:100;

[l,m]=size(t);
x(1)=1;

for i=2:m
    F1=-h*x(i-1);
    F2=-h*(x(i-1)+1/2*F1);
    F3=-h*(x(i-1)+1/2*F2);
    F4=-h*(x(i-1)+F3);
    x(i)=x(i-1)+1/6*(F1+2*F2+2*F3+F4);
end               

M=max(abs(x));
if (M>1)
    M=10;
    Mi=-10;
else
    Mi=-1;
    M=1;
end

y=exp(-t);
figure(1)
plot(t,y,'r',t,x,'b')
title('runge-kutta method')
axis([0 100 Mi M])