%
% Compute errors on successively refined grids
% for a variety of finite difference approximations
% to the derivative.
%
% Date:   January 19, 2011
% Author: John Stockie
%         Department of Mathematics
%         Simon Fraser University

function fderrors( )

x0      = 0.0;     % point to evaluate derivative
nlevels = 25;      % levels of refinement
dx      = 1.0 * 0.5.^[0:nlevels];
dferror = 0 * dx;
errplus   = dferror; 
errminus  = dferror;  
errcenter = dferror;

for i = 1 : nlevels+1,
  dfplus  = (func(x0+dx(i)) - func(x0)) / dx(i);
  dfminus = (func(x0) - func(x0-dx(i))) / dx(i);
  dfcenter= (func(x0+dx(i)) - func(x0-dx(i))) / (2*dx(i));
  dfminus2= (func(x0-2*dx(i)) - 4*func(x0-dx(i)) + 3*func(x0)) / (2*dx(i));
  
  dfexact = dfunc(x0);
  errplus(i)   = abs(dfplus   - dfexact);
  errminus(i)  = abs(dfminus  - dfexact);
  errcenter(i) = abs(dfcenter - dfexact);
  errminus2(i) = abs(dfminus2 - dfexact);
end
min(dx)

figure(1)
loglog(dx, errplus,   'b-',  ...
       dx, errminus,  'r--', ...
       dx, errcenter, 'g-.', ...
       dx, errminus2, 'm-')
axis tight
grid on
set(gca, 'XTick', 10.^[-10:2:0],'XMinorGrid','off')
set(gca, 'YTick', 10.^[-20:2:0],'YMinorGrid','off')
legend( 'Forward', 'Backward', 'Centered', 'Backward-2', 'Location', 'SouthEast' )

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function f = func(x)
f = sin(x) .* exp(-x);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function df = dfunc(x)
df = (cos(x) - sin(x)) .* exp(-x);
