% This demo computes the solution on the interval [1,2] to
%  y' = (2y-2)/t
%  y(1) = 2
% using the Euler's method.
% We know the exact solution to this problem is y = t^2+1

% In the first test, we set N=5, that is, h=(b-a)/N = 0.2
N = 5;
% plot the exact solution
figure(1);
clf;
ezplot('t^2+1',[1,2]);
hold on;
% Use the Euler's method to compute an approximate solution
h = 1/N;
t = 1;
w = 2;
plot(t,w,'*');
for i=1:N
    w = w + h*(2*w-2)/t;
    t = 1 + i*h;
    plot(t,w,'*');
end
% Finaly, compute the error of y - w at t=2
error = (2^2+1) - w;
fprintf('When h = %2.2f, the error at t=2 is %5.4f\n',h, error);
hold off

% In the second test, we set N=10, that is, h=(b-a)/N = 0.1
N = 10;
% plot the exact solution
figure(2);
clf;
ezplot('t^2+1',[1,2]);
hold on;
% Use the Euler's method to compute an approximate solution
h = 1/N;
t = 1;
w = 2;
plot(t,w,'*');
for i=1:N
    w = w + h*(2*w-2)/t;
    t = 1 + i*h;
    plot(t,w,'*');
end
% Finaly, compute the error of y - w at t=2
error = (2^2+1) - w;
fprintf('When h = %2.2f, the error at t=2 is %5.4f\n',h, error);
hold off;

% In the second test, we set N=20, that is, h=(b-a)/N = 0.05
N = 20;
% plot the exact solution
figure(3);
clf;
ezplot('t^2+1',[1,2]);
hold on;
% Use the Euler's method to compute an approximate solution
h = 1/N;
t = 1;
w = 2;
plot(t,w,'*');
for i=1:N
    w = w + h*(2*w-2)/t;
    t = 1 + i*h;
    plot(t,w,'*');
end
% Finaly, compute the error of y - w at t=2
error = (2^2+1) - w;
fprintf('When h = %2.2f, the error at t=2 is %5.4f\n',h, error);
hold off;
When h = 0.20, the error at t=2 is 0.3333
When h = 0.10, the error at t=2 is 0.1818
When h = 0.05, the error at t=2 is 0.0952