
a=0; b=3; N=100; h=(b-a)/N;

t = linspace(a,b,N+1);

u(1) = 1;  % y(a) = y(0) = 1

for i=1:N,
    u(i+1) = u(i) + h * ode1( t(i), u(i) );
end







