function [x,H]=bisez(f,a0,b0,sfa0,tol)
%BISEZ bisection method 
%X=BISEZ(F,A0,B0,SFA0,TOL) gives by the bisection method an approximation
%X of the zero of the continuous strictly monotone function F located in the interval
%[A0,B0] with sign(F(A0))*sign(F(B0))<0 and sign(F(A0))=SFA.
%The approximation X is within a tolerance TOL from the zero.
%
%[X,H]=BISEZ(F,A0,B0,SFA0,TOL) also gives a three columns matrix H, whose
%rows are the ends and the middle point of the location intervals
%generated by the method.


a=a0;
b=b0;
c=(a+b)/2;
sfa=sfa0;
n=ceil(log2((b0-a0)/tol))-1;
H=zeros(n+1,3);
H(1,:)=[a b c];
for i=1:n
    sfc=sign(f(c));
    if sfc==0
        x=c;
        return
    end
    if sfc*sfa<0
        b=c;
    else
        a=c;
        sfa=sfc;
    end
    c=(a+b)/2;
    H(i+1,:)=[a b c];
end
x=c;