BISECTION is a fast and robust root-finding method that handles n-dim arrays.
[x,fVal,ExitFlag] = BISECTION(f,LB,UB,target,options) finds x +/- TolX (LB < x < UB) such that f(x) = target +/- TolFun.
x = BISECTION(f,LB,UB) finds the root(s) of function f on the interval [LB, UB], i.e. finds x such that f(x) = 0 where LB <= x <= UB. f will never be evaluated outside of the interval specified by LB and UB. f should have only one root and f(UB) and f(LB) must bound it. Elements of x are NaN for instances where a solution could not be found.
x = BISECTION(f,LB,UB,target) finds x such that f(x) = target.
x = BISECTION(f,LB,UB,target,TolX) will terminate the search when the search interval is smaller than TolX (TolX must be positive).
x = BISECTION(f,LB,UB,target,options) solves with the default parameters replaced by values in the structure OPTIONS, an argument created with the OPTIMSET function. Used options are TolX and TolFun. Note that OPTIMSET will not allow arrays for tolerances, so set the fields of the options structure manually for non-scalar TolX or TolFun.
[x,fVal] = BISECTION(f,...) returns the value of f evaluated at x.
[x,fVal,ExitFlag] = BISECTION(...) returns an ExitFlag that describes the exit condition of BISECTION. Possible values of elements of ExitFlag and the corresponding exit conditions are
1 Search interval smaller than TolX.
2 Function value within TolFun of target.
3 Search interval smaller than TolX AND function value within TolFun of target.
-1 No solution found.
Any or all of f(scalar), f(array), LB, UB, target, TolX, or TolFun may be scalar or n-dim arrays. All non-scalar arrays must be the same size. All outputs will be this size.
Default values are target = 0, TolX = 1e-6, and TolFun = 0.
There is no iteration limit. This is because BISECTION (with a TolX that won't introduce numerical issues) is guaranteed to converge if f is a continuous function on the interval [UB, LB] and f(x)-target changes sign on the interval.
The bisection method is very robust root-finding method. The absolute error is halved at each step so the method converges linearly. However, Brent's method (such as implemented in FZERO) can converge superlinearly and is as robust. FZERO also has more features and input checking, so use BISECTION in cases where either the optimization toolbox is unavailable or if FZERO would have to be implemented in a loop to solve multiple cases, in which case BISECTION will be much faster because of vectorization.
Define LB, UB, target, TolX, and TolFun for each specific application using great care for the following reasons:
- There is no iteration limit, so given an unsolvable task, such as TolX = TolFun = 0, BISECTION remains in an unending loop.
- Spacing between very large floating point numbers is likely to be greater than TolX.
- There is no initial check to make sure that f(x) - target changes sign between LB and UB.
- Very large or very small numbers can introduce numerical issues.
Example 1: Find cube root of array 'target' without using NTHROOT and compare speed to using FZERO (typical result is FZERO taking ~70 times as long as BISECTION).
options = optimset('TolX', 1e-9);
target = [(-100:.1:100)' (-1000:1:1000)'];
xfz = zeros(size(target));
for ii = 1:numel(target)
xfz(ii) = fzero(@(x) x.^3-target(ii), [-20 20], options);
fzero_time = toc
xbis = bisection(@(x) x.^3, -20, 20, target, options);
bisection_time = toc
fprintf('FZERO took %0.0f times longer than BISECTION.\n', fzero_time/bisection_time)
Example 2: Find roots by varying the function coefficients.
[A, B] = meshgrid(linspace(1,2,6), linspace(4,12,10));
f = @(x) A.*x.^0.2 + B.*x.^0.87 - 15;
xstar = bisection(f,0,5);