function ShowDD()
% Rounding vs Calculus error
% Use divided difference to estimate the derivative of sin(x) at x = 1...
N = 100;
deltaVals = logspace(-1,-16,N);
errorVals = zeros(1,N);
for k=1:N
h = deltaVals(k);
Dsin = (sin(1+h)-sin(1))/h;
errorVals(k) = abs(Dsin - cos(1));
end
loglog(deltaVals,errorVals)
xlabel('h','fontsize',14)
ylabel('error','fontsize',14)
shg

function [d,err] = Derivative(f,a,delta,M2)
% f is a handle that references a function f(x) whose derivative
% at x = a is sought. delta is the absolute error associated with
% an f-evaluation and M2 is an estimate of the second derivative
% magnitude near a. d is an approximation to f'(a) and err is an estimate
% of its absolute error.
%
% Usage:
% [d,err] = Derivative(@f,a)
% [d,err] = Derivative(@f,a,delta)
% [d,err] = Derivative(@f,a,delta,M2)

if nargin <= 3
% No derivative bound supplied, so assume the
% second derivative bound is 1.
M2 = 1;
end
if nargin == 2
% No function evaluation error supplied, so
% set delta to eps.
delta = eps;
end
% Compute optimum h and divided difference
hopt = 2*sqrt(delta/M2);
d = (f(a+hopt) - f(a))/hopt;
err = 2*sqrt(delta*M2);