% ----------------------------
% euler2.m
% Euler method 
% ----------------------------
clear all;close all;
% initial condition
x0  = 0;
f   = 0 ;        %   f(x)=0
p   = 1;         %  f'(x)=1
g   = 0;         % f''(x)=0
gp  =  -g-5*x0 ; % f'''(x) = -f''(x)-5x
%
h = 0.01;
x = x0:h:1;
res =[x(1) f p g (-g-5*x(1))];
N=size(x,2);
for i = 2:N
    gh = g + gp*h;       % g(x+h)
    ph = p + g*h;        % p(x+h)
    fh = f + p*h;        % f(x+h)
    gph = (-g-5*x(i));    % g'(x)
    res = [res ; [x(i) fh ph gh gph]];
    % move to next point
    gp = gph;
    g  = gh;
    p  = ph;
    f  = fh;
end
func = plot(res(:,1),res(:,2),'bo');hold on;
der1 = plot(res(:,1),res(:,3),'r-');
der2 = plot(res(:,1),res(:,4),'g-');
der3 = plot(res(:,1),res(:,5),'m-');
legend('f','df/dx','d^2f/dx^2','d^3f/dx^3');
disp('result; f(x) and higher derivatives');
res