%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Program 1.2 Bisection Method
%Computes approximate solution of f(x)=0
%Input: a,b such that f(a)*f(b)<0, and number of steps n
%Output: Approximate solution xc
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function bisect(a, b, n)


eb=(b-a)/2^(n+1);  %error bound
if f(a)*f(b) >= 0
  disp('f(a)f(b)<0 not satisfied!')
end
fa=f(a);
fb=f(b);
for k=1:n
    a1(k)=a;b1(k)=b;
  c=(a+b)/2;c1(k)=c;
  fc=f(c);
  if fc == 0       %c is a solution, done
    %k
    break
  end
  if fc*fa<0       %a and c make the new interval
    b=c;fb=fc;
  else             %c and b make the new interval
    a=c;fa=fc;
  end
end
xc=(a+b)/2         %new midpoint is best estimate
f(xc,1);
fprintf(1,'The solution is %16.10f +- %16.10f \n',xc,eb)