function newton()
% An implementation of Newton's method for solving f(x) = 0

  xn = -1.8;     % This is the initial condition
  maxit = 75;    % Allow a maximum of maxit iterations
  tol = 1.0e-12; % Stop when |x_{n+1} - x_n| < tol
%
  err = 1.0;     % Initial error
  iter = 0;      % Iteration counter
  fail = 0;      % set to one if f' is "small"
%
% Write header
%
  fprintf(1, ' n            xn                     f(xn)                  f''(xn)                  error\n');
  fprintf(1, '--+-----------------------+-----------------------+-----------------------+-----------------------\n');
%
% The Newton iterations
%
  while iter < maxit && err > tol
%
%   Compute f and f' at the current guess
%
    [fn, dn] = f(xn);
%
%   Print them along with the current error
%
    fprintf(1, '%2d %23.16e %23.16e %23.16e %23.16e\n', iter, xn, fn, dn, err);
%
%   Check if f' vanishes at the current guess
%
    if abs(dn) < 1.0e-12
       fail = 1;
       break;
    end
%
%   Otherwise compute the next iterate
%
    xnew = xn - fn / dn;
    err = abs(xnew - xn);
%
    iter = iter + 1;
    xn = xnew;
  end
%
  if fail == 1
     fprintf(1, 'A near zero derivative found at step %d\n', iter);
  elseif iter >= maxit
     fprintf(1, 'Newton''s method failed to converge in %d\n steps', maxit);
  else
    fprintf(1, '%2d %23.16e %23.16e %23.16e %23.16e\n', iter, xn, f(xn), d(xn), err);
  end
%
end

function [fv, dv] = f(x)
   fv = x^3 -   x^2 -   x - 1;
   dv =       3*x^2 - 2*x - 1;
end
