function [T,W]=taylor2(Fun,dFun,t0,tlast,w0,N);
%     function [T,W]=taylor2(Fun,dFun,t0,tlast,w0,N);
%%%%% Taylor's method of order 2  for initial value problems %%%
%%
%% Solve y'= f(t,y),   t0<= t <=tlast,  y(t0)=w0
%%
%% It is assumed that f(t,y), f'(t,y), ... are defined in suitable m-files
% INPUT:
% Fun   - String containing name of function m-file (enclose in '  ')
% dFun  - String containing the name of the m-file (enclosed in ' ') for
%         the function f'(t,w)
% t0    - Initial value of t.
% tlast - Last value of t.
% w0    - Initial value.
% N     - The number of steps
%
% OUTPUT:
% T  - time points (column-vector).
% W  - Solution column-vector for corresponding T-values.
%
% The result can be displayed by: plot(T,W)

%% Step 1
h=(tlast-t0)/N;t=t0;w=w0;
T=t;W=w;

%% Step 2
for i=1:N;
    w = w + h*feval(Fun, t, w) +(h^2/2)*feval(dFun, t, w); 
    t = t0 + i*h;
    T=[T;t]; W = [W;w];
end; 
