% Function applying the simple Euler method to a given equation
%                   y' = fun(x,y) 
%
% The user must specify:
%  -  the function  fun; 
%  -  the integration span  xspan (a 2-component vector);
%  -  the integration step size  h;  and
%  -  the initial value,  y_init,  for y.

function [X,Y]=hw3_AdamsBashforth3(fun,xspan,h,y_0,y_1,y_2);

Y(1:3)=[y_0 y_1 y_2];
X=xspan(1):h:xspan(2);

f_nminus1 = feval(fun,X(2),Y(2));
f_nminus2 = feval(fun,X(1),Y(1));

for n=3:length(X)-1
    f_n = feval(fun,X(n),Y(n));
    Y(n+1) = Y(n) + (h/12)*( 23*f_n - 16*f_nminus1 + 5*f_nminus2 );
    f_nminus2 = f_nminus1;
    f_nminus1 = f_n;
end