% Script to illustrate how forward difference approximation of the
% first-order derivative works.
y = @(t) atan(t);   % a (rather arbitrary) function to illustrate the concept
yprime_exact = @(t) 1/(1+t^2);   % its exact derivative (known from calculus)
t0 = 2;   % point at which we evaluate the derivative
h = 0.1;  % initial stepsize
fwddiff = @(t,h) (y(t+h)-y(t))/h;   % forward difference approximation
% Repeat forward difference approximation many times for smaller and
% smaller stepsize h
disp(sprintf('         h         approx. deriv.       error        rate'));
for i=1:20
   err = abs(yprime_exact(t0)-fwddiff(t0,h));  % absolute error of forward difference approximation
   if i>1
      rate = log(old_err/err)/log(2);    % exponent of h in the truncation error
      disp(sprintf('%15.12f   %15.12f   %e   %f', h, fwddiff(t0,h), err, rate));
   else
      disp(sprintf('%15.12f   %15.12f   %e', h, fwddiff(t0,h), err));
   end
   h = h/2;  % halve stepsize for next iteration
   old_err = err;
end
