function r = fixedpoint(fun, xinit, tol, nmax, verbose)

% fixedpoint  Uses fixed point iteration to find a root 
%             of the equation  x = fun(x)  
%             starting the iteration at xinit
%
% Input:  fun = the function that we iterate
%         xinit = the value from which we start iterating
%         tol = x tolerance
%         nmax = maximum number of iterations
%         verbose = 1 in case all iterations are needed, 0 otherwise
%
% Output: r = the root of fun(x) = x with the desired tolerance

if tol < 3*eps       %  The tolerance is too small
  warning(sprintf('the tolerance is too small\n'));
  return;
end

if verbose == 1         %  Printing the headings
  fprintf('\niters         x                 x-xold      log_10(|x-xold|)\n\n');
end


%  The main part of the program starts here

xold = xinit;

for i = 1:nmax
  x = feval(fun,xold);       %  Assigning value fun(x) to x

  if verbose == 1         %  Printing all iterations
    fprintf('%3d  %18.14f  %18.14f  %10.5f\n', ...
	    i, x, x-xold,log(abs(x-xold))/log(10.) );
  end

  if abs(x-xold) < tol
    r = x;    return;   %  We found the root
  end

  xold = x;

end

warning(sprintf('root not found within tol after %d iters\n',i));