function Taylor(n)
% Solves the ivp y'=y - 1/(1+x)^2 - 1/(1+x), y(0)=2 
% by the third order Taylor method
%
fprintf('\n')
%
for iit=1:9
%
nn=2^(iit-1)*n;
uu(iit)=1/nn;
%
h=1./nn;
%
y=2.;
u=0.;
%
for i=0:nn
t=i*h;
t1=(i+1)*h;
ex=exp(t)+(1/(1+t));
er(i+1)=abs(ex-y);
%y=y+h*f1(t,y)+(h^2/2)*fx1(t,y)+(h^3/6)*fxx1(t,y);
y=y+h*f1(t,y)+(h^2/2)*fx1(t,y)+(h^3/6)*fxx1(t,y)+(h^4/24)*fxxx1(t,y);
end
u1(iit)=max(er);
fprintf('h =%12.6f        maxerr = %15.6e\n',h,u1(iit))
end
fprintf('\n')
for iit=1:8
      p=log(u1(iit)/u1(iit+1))/log(2.);
      fprintf('h =%12.6f        rate %15.6e\n',uu(iit),p)
end
u1';
loglog(uu,u1,'b+')
xlabel('h');ylabel('Error');
function f=f1(x,y)
f=y-(1/(1+x)^2)-(1/(1+x));
function fx=fx1(x,y)
fx=y-(1/(1+x))+(2/(1+x)^3);
function fxx=fxx1(x,y)
fxx=y-(1/(1+x))-(6/(1+x)^4);
function fxxx=fxxx1(x,y)
fxxx=y-(1/(1+x))+(24/(1+x)^5);

