function res = illinois_method(func,xim1,xi,tol)
%
% Based from the paper by Dowell and Jerrett (1971): 
% 
% http://link.springer.com/article/10.1007/BF01934364#page-1
% 
% Written by:
% --
% John L. Weatherwax                2009-04-21
%
% email: wax@alum.mit.edu
%
% Please send comments and especially bug reports to the
% above email address.
%
% Problem on EPage 62; Solutions on EPage 161
%
%-----

res=xi; % initial best guess at the roots location (in case of quick exit)
it=0; 
fim1 = feval(func,xim1);
fi = feval(func,xi);
assert( fim1*fi <= 0, 'inputs xim1 and xi dont bracket a root' );
denom = ( fi - fim1 )/(xi - xim1);
d = fi/denom;
while abs(d)>tol
  xip1 = xi-d;
  fip1 = feval(func,xip1);
  
  % Update (xim1,fim1): 
  if fi * fip1 < 0 % the new interval brackets our root
    xim1 = xi;
    fim1 = fi;
  else % we no longer bracket our root
    fim1 = 0.5*fim1;
  end
  
  % Update (xi,fi):
  xi = xip1;
  fi = fip1;
  
  assert( fi * fim1 < 0, 'both points should bracket root' );

  it=it+1;
  denom = ( fi - fim1 )/(xi - xim1);
  d = fi/denom;
end;
res = xi;