clear all;
format long;

%
% A sample matlab code to calculate f'(x) at x=x0=0 if f(x) = sin(x)
%

% INPUT
x0 = 0;

DevExact = cos(x0); % Exact value of f'(0)

eps = 1.d-8; % Tolerance to stop the loop
step = 1; % Inital stepsize to calculate the derivative

Dev0 = 0;
Dev1 = (sin(x0+step) - sin(x0))/(x0 + step - x0); % Definition of derivative

n = 0; % Initialize the counter
figure(1);clf; % visualize the convergence
axes('FontSize',16);
hold on;
while abs(Dev1-Dev0) >= eps 
    step = step/2;
    Dev0 = Dev1;
    Dev1 = (sin(x0+step) - sin(x0))/(x0 + step - x0);
    n = n + 1;
    plot(n, Dev1, 'r*','LineWidth', 2, 'MarkerSize', 16);
    size(n) = step;
    value(n) = Dev1 - DevExact;
end
xlabel('Counter')
ylabel('Approximated derivative')

disp('The value of derivative of sin(x) at x0=0');
disp(Dev1);

% Check the convergence rate of the algorithm
figure(2);clf;
axes('FontSize',16);
loglog(size, value, 'r-o','LineWidth', 2, 'MarkerSize', 16);
hold on;
loglog(size, -size.^2, 'k-o','LineWidth', 2, 'MarkerSize', 16);
xlabel('h')
legend('Derivative difference','h^2')
% Second-order convergence