Bisection method in matlab

Question

function r=bisection(f,a,b,tol,nmax)
% function r=bisection(f,a,b,tol,nmax)
% inputs: f: function handle or string
% a,b: the interval where there is a root
% tol: error tolerance
% nmax: max number of iterations
% output: r: a root
c=(a+b)/2;
nit=1;
if f(a)*f(b)>0
    r=NaN;
    fprintf("The bisection method failed \n")
else
    while(abs(f(c))>=tol && nit<nmax)
        if (f(c)*f(a))<0
            c=(a+c)/2;
        elseif (f(c)*f(b))<0
            c=(a+b)/2;
        elseif f(c)==0
            break;
        end
        nit=nit+1;
    end
    r=c;
end

Above are my code for the Bisection method. I am confused about why that code don't work well. The result of f(c) is repeated every three times when running this. Does anyone can tell me why this code won't work?

Answer 1

In your solution, you forgot to consider that you need to reset one of the 2 extremes a and b of the interval to c at each iteration.

function r=bisection(f,a,b,tol,nmax)
% function r=bisection(f,a,b,tol,nmax)
% inputs: f: function handle or string
% a,b: the interval where there is a root
% tol: error tolerance
% nmax: max number of iterations
% output: r: a root
c=(a+b)/2;
nit=1;
if f(a)*f(b)>0
    r=NaN;
    fprintf("The bisection method failed \n")
else
    while(abs(f(c))>=tol && nit<nmax)
        if (f(c)*f(a))<0
            b=c;                % new line
            c=(a+c)/2;            
        elseif (f(c)*f(b))<0
            a=c;                % new line
            c=(c+b)/2;
        elseif f(c)==0
            break;
        end
        nit=nit+1;
    end
    r=c;
end

Answer 2

I think you need to update the boundaries for next round of bisection (within your while loop) like below

function r = bisection(f,a,b,tol,nmax)
c=mean([a,b]);
nit=1;
if f(a)*f(b)>0
    r=NaN;
    fprintf("The bisection method failed \n")
else
    while(abs(f(c))>=tol && nit<nmax)
        if (f(c)*f(a))<0
            b=c;                      
        elseif (f(c)*f(b))<0
            a=c;
        elseif f(c)==0
            break;
        end
        c=mean([a,b]);
        nit=nit+1;
    end
    r=c;
end
end