function [x,n] = fixed_point(func,x0,atol)
%
% This function implements the fixed point iteration
%      x_k+1=func(x_k)
% written by Jessica Bosch
%
% Assuming ga = func(a), gb = func(b), ga >= a, and gb <= b,
% there is a fixed point in [a,b] such that func(fixed point) = fixed point.
% x0 is the initial guess.
% This function returns in x a value such that 
%      | x_k - x_k-1 | < atol,
% in n the number of iterations required, and in iterates the sequence of
% iterates.
%
    
    
  % initialization
    x_old=x0+1;                     % x_old will be x_k-1
    x_new=x0;                       % x_new will be x_k
    n=0;                            % iteration counter

    
  % fixed point iteration  
    while (abs(x_new-x_old)>=atol)
        
        x_old=x_new;                % update previous iterate
        x_new = feval(func,x_new);  % compute next iterate
        n=n+1;                      % update iteration counter

        % break if it takes too long
        if (n>=1e5)
           break;
        end
        
    end

    x=x_new;                        % approximate fixed point

end

