% sample function to find the zero of a function func(x) using Bisection
% [xout,fx,iters] = bisect1(xin,tol,debug)
% Inputs:
%    xin(1:2) - lower and upper limits of interval to search for a zero
%                func must have different signs at xin(1),xin(2)
%    tol - stop bisecting interval when its width is <= tol
%    debug - print results of each iteration if debug>0
% Outputs:
%    xout(1:2) - output interval "guaranteed" to contain a zero
%    fx(1:2) - value of func at xout(1),xout(2)
%    iters - number of iterations to converge
function [xout,fx,iters] = bisect1(xin,tol,debug)
xlower = xin(1);
xupper = xin(2);
flower = func(xlower); 
fupper = func(xupper);
width = xupper-xlower;
iters=0;
if (debug>0),
   disp('iteration            interval width')
   disp('xlower               xupper')
   disp('f(xlower)            f(xupper) = ')
   format long e
   [[i,width];[xlower,xupper];[flower,fupper]]
   disp('pause - hit return to continue'), pause
end
while (width > tol),
   xmid = (xlower + xupper)/2;
   fmid = func(xmid) ;
   if (fmid*flower <= 0),
      xupper = xmid;
      fupper = fmid;
   else
      xlower = xmid;
      flower = fmid;
   end
   width = xupper - xlower;
   iters=iters+1;
   if (debug>0),
      disp('iteration            interval width')
      disp('xlower               xupper')
      disp('f(xlower)            f(xupper) = ')
      format long e
      [[iiters,width];[xlower,xupper];[flower,fupper]]
      disp('pause - hit return to continue'), pause
   end
end
xout=[xlower,xupper];
fx=[flower,fupper];