clear all
close all

% Set the sequence of differences to be used:
h = 10.^[-1 -2 -3 -4 -5 -6 -7 -8 -9 -10];

% Set the different evaluation points:
x0vec = [0 1 10 50 100];
nx0 = length(x0vec);

figure()
% Compute the derivatives for all x0s:
for i = 1:nx0
    % Get the next point:
    x0 = x0vec(i);
    % Compute the analytic derivative:
    df0 = 1/(1+x0^2);
    % Compute all the finite-difference approximations:
    df = (atan(x0 + h) - atan(x0))./h;
    % Plot the relative error:
    loglog(h,abs(df-df0)./abs(df0))
    hold all
end
% Annotate the plot:
xlabel('h')
ylabel('Relative Error')
legend('x0 = 0','x0 = 1','x0 = 10','x0 = 50','x0 = 100','Location','southeast')