function [ a, b, it ] = bisection ( a, b, tol, f )

%*****************************************************************************80
%
%% bisection() carries out the bisection method to find the zero of a function.
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    04 December 2023
%
%  Author:
%
%    John Burkardt
%
%  Input:
%
%    real A, B: the endpoints of a change of sign interval.
%
%    real TOL: the interval size tolerance.  Once |B-A| < TOL,
%    the iteration will stop.
%
%    pointer value=F(x): a pointer to a function.  This can be
%    an anonymous function or the name of a MATLAB M-file function.
%
%  Output:
%
%    real A, B: the new endpoints that constitute an
%    change of sign interval no larger than TOL.
%
%    integer IT: the number of bisections.
%
  fa = f ( a );
  fb = f ( b );
  it = 0;

  while ( tol < abs ( b - a ) )
    c = ( a + b ) / 2.0;
    fc = f(c);
    it = it + 1;

    if ( 100 < it )
      error ( 'bisection(): Too many iterations!' );
    end

    if ( sign ( fc ) == sign ( fa ) )
      a = c;
      fa = fc;
    else
      b = c;
      fb = fc;
    end
  end

  return
end

