function [root,iter] = bisect(f, a, b, tol)
%
% Call by: 
%          [xstar,e]=bisect(@f,1,1.5,1e-6)
% or
%          [xstar,e]=bisect(@(x) exp(x)-2*x-1,1,1.5,1e-6)
%
fa = f(a); fb = f(b);
if fa * fb/abs(fb) > 0
    fprintf('\n   Wrong initial interval (f(a)*f(b) > 0).\n')
    return;
end    
iter = 0;
while abs((b - a)/2) > tol
    mid = a + 0.5 * (b - a); 
    fmid = f(mid);
    if fmid == 0
        a = mid; b = mid;
    elseif fa * (fmid/abs(fmid)) < 0
        b = mid; fb = fmid;
    else
        a = mid; fa = fmid;
    end;
    iter = iter + 1;
end;
root = (a + b) / 2;