% This is another simple implementation of Algorithm2.1

% Read in the relevant parameters

fprintf(1,'Enter tol, a, b, nmax (on separate lines)\n');

tol = input(' ');
a = input(' ');
b = input(' ');
nmax  = input(' ');

% Specify the function f(x)

 fprintf(1,'Input the function f(x) in terms of x\n');
 s = input(' ','s');
 f = inline(s,'x');

p=b;
false = 1;
true = 0;
finish = false;
iteration = 0;

while (finish == false & iteration <= nmax );

   iteration = iteration + 1;
   pold=p;
%   p=(a+b)/2;
    p=a + (b-a)/2;

   fa=f(a);
   fb=f(b);
   fp=f(p);
   if (abs(p-pold)/abs(p) < tol) | (fp==0.0) ;

      % finished !

      p_final = p;
      finish = true;

   else;

      % not converged ! calculate new a and b

      if (fa*fp > 0);

         %  root is between p and b   
         a=p;

      else;

         %  root is between a and p
         b=p;

      end
    end

end

if (iteration < nmax);

  fprintf(1,'result after iteration number %3d ',iteration);
  fprintf(1,'is %15.8e \n',p_final);

else

  fprintf(1,'no convergence after iteration number %3d ',nmax);

end
