function [xvect,xdif,fx,nit]=fixpoint(x0,nmax,toll,fun,phi)err=toll+1; nit=0; xvect=x0; x=x0; fx=eval(fun); xdif=[];while (nit < nmax & err > toll),    nit=nit+1; x=xvect(nit); xn=eval(phi); err=abs(xn-x);     xdif=[xdif; err]; x=xn; xvect=[xvect;x]; fx=[fx;eval(fun)];end;