%
% Luke Olson
% 20060206
%

function [k,x,err,X] = fixedpointplot(g,x0,tol,maxit)
% Input 
%       o g = iteration matrix
%       o x0 = initial guess
%       o tol = tolerance
%       o maxit = maximum number of iterations
% Output 
%       o k = number of iterations
%       o x = approx fixed-point
%       o err = error in the approximation
%       o X = iteration history {x_n}

X(1)= x0;

figure(1);clf;
figure(2);clf;
minx = .5*x0;
maxx = 1.5*x0;

for k=2:maxit
   X(k)=feval(g,X(k-1));
   err=abs(X(k)-X(k-1));
   relerr=err/(abs(X(k))+eps);
   x=X(k);
   %
   %
   % plot
   %---------------------------------------------------------------------
   figure(1);hold on;
   maxx = max(maxx,x);
   minx = min(minx,x);
   xx=linspace(.9*minx,1.1*maxx,50);
   gx=feval(g,xx);
   plot(xx,xx,'k-','LineWidth',1);
   plot(xx,gx,'b-','LineWidth',1);
   plot(xx,0*xx,'k:');
   if(k==2)
    plot([X(k-1) X(k)],[X(k) X(k)],'ro-','LineWidth',2,'MarkerSize',5,'MarkerEdgeColor','r','MarkerFaceColor','k');
   else
    plot([X(k-1) X(k-1)], [X(k-1) X(k)],'ro-','LineWidth',2,'MarkerSize',5,'MarkerEdgeColor','r','MarkerFaceColor','k');
    plot([X(k-1) X(k)],[X(k) X(k)],'ro-','LineWidth',2,'MarkerSize',5,'MarkerEdgeColor','r','MarkerFaceColor','k');
     delete(h1);
     delete(h2);
   end
   minX=min(X);
   maxX=max(X);
   d=abs(minX-maxX);
   b1=-0.9*d+minX;
   b2=1.1*d+maxX;
   c1 = min(b1,feval('gfunc',b1)) ;
   c2 = max(b2,feval('gfunc',b2));
   axis([b1 b2 c1 c2]);
   %axis([min(xx) max(xx) min(-0.11,min(gx)) max(0.11,max(gx))])
   h1 = text(X(k-1),c1+(c2-c1)/10,['\uparrow x_{' num2str(k-2) '}'],...
    'HorizontalAlignment','left',...
    'VerticalAlignment','top',...
    'FontSize',16);

   h2 = text(X(k),c1+(c2-c1)/5,['\uparrow x_{' num2str(k-1) '}'],...
    'HorizontalAlignment','left',...
    'VerticalAlignment','top',...
    'FontSize',16);
   title(['Error(' num2str(k-1) ') = ' num2str(err)])
   xlabel('x');
   %---------------------------------------------------------------------

   %---------------------------------------------------------------------
   % plot zoom
   p1=gca;
   p2=figure(2);clf;
   tmp=copyobj(p1,p2);
   minX=min(X(k-1:k));
   maxX=max(X(k-1:k));
   d=abs(minX-maxX);
   b1=-0.9*d+minX;
   b2=1.1*d+maxX;
   axis([b1 b2 min(b1,feval('gfunc',b1)) max(b2,feval('gfunc',b2))])
   %---------------------------------------------------------------------
   %
   pause;
   if (err<tol) | (relerr<tol),break;end
end

if k == maxit
   disp('maximum number of iterations exceeded')
end
X=X';