function xc = bisect(f,a,b,tol) ind = b-a; while ind > tol xx = (a+b)/2; if f(a)*f(xx) < 0 b = xx; else a = xx; end ind = b - a; end xc = xx;