function [t,a] = euler(h, npts, t0, y0)
    t(1) = t0;
    a(1) = y0;
    for i = 2: npts+1;
        a(i) = a(i - 1) + h*f(t(i-1), a(i - 1));
        t(i) = t(i-1) + h;
    end;