function p=bisection_table(f,a,b,tol,Niter)
% Like bisection.m, but printing a convergence table, up to 100 iterates.
% Ex: 
% f = @(x) x.^2 - 5;
% [p,iter] = bisection_table(f,0,5,1e-6)

if nargin == 4
   Niter = 100; 
end

fprintf(' n         a             b             p           f(p)    \n');
fprintf('-----------------------------------------------------------\n');
for n=1:Niter
  p=(a+b)/2;
  if p-a<tol, break; end
  fprintf('%2d  %12.8f  %12.8f  %12.8f  %12.8f\n',n,a,b,p,f(p));
  if f(a)*f(p)>0
    a=p;
  else
    b=p;
  end
end
