% Euler and 2nd Order Taylor Example
% to solve dy/dt = f(t,y), y(a)=y0, on [a,b]

a=0;
b=1;
N=10;
h=(b-a)/N;

% Slope function

%f = inline('t+y','t','y');
%sol = inline('2*exp(t)-t-1','t');

f = inline('t^2+y','t','y');
sol = inline('-2-2*t-t^2+3*exp(t)','t');

%f = inline('t^3+y','t','y');
%sol = inline('-6-6*t-3*t^2-t^3+7*exp(t)','t');

% Initial Condition
t(1)=0;
y(1)=1;

% Euler's Method
for i=2:N+1
    y(i)=y(i-1)+h*f(t(i-1),y(i-1));
    t(i)=t(i-1)+h;
    %disp(['(' num2str(t(i)) ',' num2str(y(i)) ')'])
end
Err1 = abs(y(N+1)-sol(1))

% Taylor's Method, Order 2
y2(1)=1;
h2 = h^2/2;
for i=2:N+1
    %y2(i)=y2(i-1)+h*f(t(i-1),y2(i-1))+h2*(1+t(i-1)+y2(i-1));
    y2(i)=y2(i-1)+h*f(t(i-1),y2(i-1))+h2*(2*t(i-1)+t(i-1)^2+y2(i-1));
    %y2(i)=y2(i-1)+h*f(t(i-1),y2(i-1))+h2*(3*t(i-1)^2+t(i-1)^3+y2(i-1));
    t(i)=t(i-1)+h;
end
Err2 = abs(y2(N+1)-sol(1))

% Taylor's Method, Order 3
y3(1)=1;
h3 = h^2*(1/2+h/6);
for i=2:N+1
    y3(i)=y3(i-1)+h*f(t(i-1),y3(i-1))+h3*(1+t(i-1)+y3(i-1));
    t(i)=t(i-1)+h;
    %disp(['(' num2str(t(i)) ',' num2str(y3(i)) ')'])
    %disp([num2str(y3(i))])
end
Err3 = abs(y3(N+1)-sol(1))

% Taylor's Method, Order 4
y4(1)=1;
h4 = h^2*(1/2+h/6+h^2/24);
for i=2:N+1
    y4(i)=y4(i-1)+h*f(t(i-1),y4(i-1))+h4*(1+t(i-1)+y4(i-1));
    t(i)=t(i-1)+h;
end
Err4 = abs(y4(N+1)-sol(1))

% Midpoint Method
y5(1)=1;
for i=2:N+1
    k1=h/2*f(t(i-1),y5(i-1));
    k2=h*f(t(i-1)+h/2,y5(i-1)+k1);
    y5(i)=y5(i-1)+k2;
    t(i)=t(i-1)+h;
end
Err5 = abs(y5(N+1)-sol(1))


for i=1:N+1
    disp([num2str(y2(i)) ' & ' num2str(y5(i)) ' & ' num2str(sol(t(i))) '\\'])
end

for i=1:N+1
    disp([num2str(sol(t(i))) ' & ' num2str(y2(i)) ' & ' num2str(sol(t(i))-y2(i)) ' & ' num2str(y5(i)) ' & ' num2str(sol(t(i))-y5(i)) '\\' ]) %    ])
end