% Example 1.1 
% Table 1.1 and Figure 1.2 of  http://www.amath.washington.edu/~rjl/fdmbook

% Order of the derivative
k = 1;

% Point at which the derivative will be approximated
xb = 1;

% Exact value
uptrue = fp(xb);

% Indices of grid points
a = -2;
b = 1;
j = a:b;

% Vector containing Values of step size "h"
 hvals = [1e-1 5e-2 1e-2 5e-3 1e-3];
Epu = [];  Emu = [];  E0u = [];  E3u = [];

disp('')
disp('Approximations to u^(k)(xb) using different FD approximations')
% table headings:
disp(' ')
disp('       h              Epu            Emu             E0u             E3u')

for i=1:length(hvals)
   h = hvals(i);
   % approximations to u'(xb):
   Dpu = (f(xb+h) - f(xb))/h;
   Dmu = (f(xb) - f(xb-h))/h;
   D0u = (f(xb+h) - f(xb-h))/(2*h);
   xpts = xb + h*j';
   D3u = fdcoeffF(k,xb,xpts) * f(xpts) ;
   % errors:
   Epu(i) = Dpu - uptrue;
   Emu(i) = Dmu - uptrue;
   E0u(i) = D0u - uptrue;
   E3u(i) = D3u - uptrue;

   % print line of table:
   disp(sprintf('%13.4e   %13.4e   %13.4e   %13.4e  %13.4e',...
                 h,Epu(i),Emu(i),E0u(i),E3u(i)))
   end

% plot errors:
clf
loglog(hvals,abs(Epu),'o-')
axis([1e-4 .2 1e-12 1])
%Another possible scale for the axes.
%axis([1e-8 .2 1e-13 1])
xlabel('log(h)');
ylabel('log(|E(h)|)');
hold on
loglog(hvals,abs(E0u),'o-')
loglog(hvals,abs(E3u),'o-')
hold off