Introduction

The bisection method is simple, robust, and straight-forward: take an interval [a, b] such that f(a) and f(b) have opposite signs, find the midpoint of [a, b], and then decide whether the root lies on [a, (a + b)/2] or [(a + b)/2, b]. Repeat until the interval is sufficiently small.

Background

Useful background for this topic includes:

3. Iteration
8. Bracketing

References

Bradie, Section 2.1, The Bisection Method, p.58.
Mathews, Section 2.2, Bisection Method of Bolzano, p.53.
Weisstein, http://mathworld.wolfram.com/Bisection.html.

Theory

Intermediate Value Theorem (IVT)

Given a continuous real-valued function f(x) defined on an interval [a, b], then if y is a point between the values of f(a) and f(b), then there exists a point r such that y = f(r).

As an example, consider the function f(x) = sin(x) defined on [1, 6]. The function is continuous on this interval, and the point 0.5 lies between the values of sin(1) ≈ 0.841 and sin(6) ≈ -0.279. Thus, there is at least one point r (there may be more) on the interval [1, 6] such that sin(r) = 0.5. In this case, r ≈ 2.61799. This is shown in Figure 1.

Figure 1. The IVT applied to sin(x) on [1, 6] with y = 0.5.

Using the IVT to Bound a Root

Suppose we have a function f(x) and an interval [a, b] such that either the case that f(a) > 0 and f(b) < 0 or the case that f(a) < 0 and f(b) > 0, that is, f(a) and f(b) have opposite signs. Then, the value 0 lies between f(a) and f(b), and therefore, there must exist a point r on [a, b] such that f(r) = 0.

We may refine our approximation to the root by dividing the interval into two: find the midpoint c = (a + b)/2. In any real world problem, it is very unlikely that f(c) = 0, however if we are that lucky, then we have found a root. More likely, if f(a) and f(c) have opposite signs, then a root must lie on the interval [a, c]. The only other possibility is that f(c) and f(b) have opposite signs, and therefore the root must lie on the interval [c, b].

We may repeat this process numerous times, each time halving the size of the interval.

Example

An example of bisecting is shown in Figure 2. With each step, the midpoint is shown in blue and the portion of the function which does not contain the root is shaded in grey. As the iteration continues, the interval on which the root lies gets smaller and smaller. The first two bisection points are 3 and 4.

Figure 2. The bisection method applied to sin(x) starting with the interval [1, 5].

HOWTO

Problem

Given a function of one variable, f(x), find a value r (called a root) such that f(r) = 0.

Assumptions

We will assume that the function f(x) is continuous.

Tools

We will use sampling, bracketing, and iteration.

Initial Requirements

We have an initial bound [a, b] on the root, that is, f(a) and (b) have opposite signs.

Iteration Process

Given the interval [a, b], define c = (a + b)/2. Then

if f(c) = 0 (unlikely in practice), then halt, as we have found a root,
if f(c) and f(a) have opposite signs, then a root must lie on [a, c], so assign b = c,
else f(c) and f(b) must have opposite signs, and thus a root must lie on [c, b], so assign a = c.

Halting Conditions

There are three conditions which may cause the iteration process to halt:

As indicated, if f(c) = 0.
We halt if both of the following conditions are met:
- The width of the interval (after the assignment) is sufficiently small, that is b - a < ε_step, and
- The function evaluated at one of the end point |f(a)| or |f(b)| < ε_abs.
If we have iterated some maximum number of times, say N, and have not met Condition 1, we halt and indicate that a solution was not found.

If we halt due to Condition 1, we state that c is our approximation to the root. If we halt according to Condition 2, we choose either a or b, depending on whether |f(a)| < |f(b)| or |f(a)| > |f(b)|, respectively.

If we halt due to Condition 3, then we indicate that a solution may not exist (the function may be discontinuous).

Error Analysis

Given that we an initial bound on the problem [a, b], then the maximum error of using either a or b as our approximation ish = b − a. Because we halve the width of the interval with each iteration, the error is reduced by a factor of 2, and thus, the error after n iterations will be h/2ⁿ.

Therefore, thus, if ε_step is fixed, then we may immediately know how many steps are required, after which we are assured that the absolute error is less than ε_step. The inequality

h / 2ⁿ

may be solved for an integer value of n by finding:

For example, suppose that our initial interval is [0.7, 1.5]. If we have an ε_step value of 1e-5, then we require a minimum of ⌈log₂( 0.8/1e-5 )⌉ = 17 steps.

Examples

Example 1

Consider finding the root of f(x) = x² - 3. Let ε_step = 0.01, ε_abs = 0.01 and start with the interval [1, 2].

Table 1. Bisection method applied to f(x) = x² - 3.

a	b	f(a)	f(b)	c = (a + b)/2	f(c)	Update	new b − a
1.0	2.0	-2.0	1.0	1.5	-0.75	a = c	0.5
1.5	2.0	-0.75	1.0	1.75	0.062	b = c	0.25
1.5	1.75	-0.75	0.0625	1.625	-0.359	a = c	0.125
1.625	1.75	-0.3594	0.0625	1.6875	-0.1523	a = c	0.0625
1.6875	1.75	-0.1523	0.0625	1.7188	-0.0457	a = c	0.0313
1.7188	1.75	-0.0457	0.0625	1.7344	0.0081	b = c	0.0156
1.71988/td>	1.7344	-0.0457	0.0081	1.7266	-0.0189	a = c	0.0078

Thus, with the seventh iteration, we note that the final interval, [1.7266, 1.7344], has a width less than 0.01 and |f(1.7344)| < 0.01, and therefore we chose b = 1.7344 to be our approximation of the root.

Example 2

Consider finding the root of f(x) = e^-x(3.2 sin(x) - 0.5 cos(x)) on the interval [3, 4], this time with ε_step = 0.001, ε_abs = 0.001.

Table 1. Bisection method applied to f(x) = e^-x(3.2 sin(x) - 0.5 cos(x)).

a	b	f(a)	f(b)	c = (a + b)/2	f(c)	Update	new b − a
3.0	4.0	0.047127	-0.038372	3.5	-0.019757	b = c	0.5
3.0	3.5	0.047127	-0.019757	3.25	0.0058479	a = c	0.25
3.25	3.5	0.0058479	-0.019757	3.375	-0.0086808	b = c	0.125
3.25	3.375	0.0058479	-0.0086808	3.3125	-0.0018773	b = c	0.0625
3.25	3.3125	0.0058479	-0.0018773	3.2812	0.0018739	a = c	0.0313
3.2812	3.3125	0.0018739	-0.0018773	3.2968	-0.000024791	b = c	0.0156
3.2812	3.2968	0.0018739	-0.000024791	3.289	0.00091736	a = c	0.0078
3.289	3.2968	0.00091736	-0.000024791	3.2929	0.00044352	a = c	0.0039
3.2929	3.2968	0.00044352	-0.000024791	3.2948	0.00021466	a = c	0.002
3.2948	3.2968	0.00021466	-0.000024791	3.2958	0.000094077	a = c	0.001
3.2958	3.2968	0.000094077	-0.000024791	3.2963	0.000034799	a = c	0.0005

Thus, after the 11th iteration, we note that the final interval, [3.2958, 3.2968] has a width less than 0.001 and |f(3.2968)| < 0.001 and therefore we chose b = 3.2968 to be our approximation of the root.

Example 3

Apply the bisection method to f(x) = sin(x) starting with [1, 99], ε_step = ε_abs = 0.00001, and comment.

After 24 iterations, we have the interval [40.84070158, 40.84070742] and sin(40.84070158) ≈ 0.0000028967. Note however that sin(x) has 31 roots on the interval [1, 99], however the bisection method neither suggests that more roots exist nor gives any suggestion as to where they may be.

Questions

Question 1

Approximate the root of f(x) = x³ - 3 with the bisection method starting with the interval [1, 2] and use ε_step = 0.1 and ε_abs= 0.1 .

Answer: 1.4375

Question 2

Approximate the root of f(x) = x² - 10 with the bisection method starting with the interval [3, 4] and use ε_step = 0.1 and ε_abs= 0.1 .

Answer: 3.15625 (you need a few extra steps for ε_abs)

Applications to Engineering

To be completed.

Matlab

The bisection method in Matlab is quite straight-forward. Assume a file f.m with contents

function y = f(x)
   y = x.^3 - 2;

exists. Then:

>> format long
>> eps_abs = 1e-5;
>> eps_step = 1e-5;
>> a = 0.0;
>> b = 2.0;
>> while (b - a >= eps_step || ( abs( f(a) ) >= eps_abs && abs( f(b) )  >= eps_abs ) )
    c = (a + b)/2;
    if ( f(c) == 0 )
       break;
    elseif ( f(a)*f(c) < 0 )
       b = c;
    else
       a = c;
    end
end
>> [a b]
ans = 1.259918212890625  1.259925842285156
>> abs(f(a))
ans =    0.0000135103601622
>> abs(f(b))
ans =    0.0000228224229404

Thus, we would choose 1.259918212890625 as our approximation to the cube-root of 2, which has an actual value (to 16 digits) of 1.259921049894873.

An implementation of a function would be

function [ r ] = bisection( f, a, b, N, eps_step, eps_abs )
    % Check that that neither end-point is a root
    % and if f(a) and f(b) have the same sign, throw an exception.

    if ( f(a) == 0 )
 r = a;
 return;
    elseif ( f(b) == 0 )
 r = b;
 return;
    elseif ( f(a) * f(b) > 0 )
        error( 'f(a) and f(b) do not have opposite signs' );
    end

    % We will iterate N times and if a root was not
    % found after N iterations, an exception will be thrown.

    for k = 1:N
        % Find the mid-point
        c = (a + b)/2;

        % Check if we found a root or whether or not
        % we should continue with:
        %          [a, c] if f(a) and f(c) have opposite signs, or
        %          [c, b] if f(c) and f(b) have opposite signs.

        if ( f(c) == 0 )
            r = c;
            return;
        elseif ( f(c)*f(a) < 0 )
            b = c;
        else
            a = c;
        end

        % If |b - a| < eps_step, check whether or not
        %       |f(a)| < |f(b)| and |f(a)| < eps_abs and return 'a', or
        %       |f(b)| < eps_abs and return 'b'.

        if ( b - a < eps_step )
            if ( abs( f(a) ) < abs( f(b) ) && abs( f(a) ) < eps_abs )
                r = a;
                return;
            elseif ( abs( f(b) ) < eps_abs )
                r = b;
                return;
            end
        end
    end

    error( 'the method did not converge' );
end

Maple

The bisection method in Maple is quite straight-forward:

> epsilon[abs] := 1e-5:
> epsilon[step] := 1e-5:
> f := x -> x^3 - 2:
> a := 0.0:
> b := 2.0:
> while b - a >= epsilon[step] or ( abs( f(a) ) >= epsilon[abs] and abs( f(b) ) >= epsilon[abs] ) do
    c := (a + b)/2;

    if f(c) = 0 then
       break;
    elif f(a)*f(c) < 0 then
       b := c;
    else
       a := c;
    end if;
end do:
> [a, b];
                               [1.259918214, 1.259925843]
> abs( f(a) );
                                      0.000013505
> abs( f(b) );
                                      0.000022826

Thus, we would choose 1.259918214 as our approximation to the cube-root of 2, which has an actual value (to 10 digits) of 1.259921050.

==========================================

SOURCE: https://ece.uwaterloo.ca/~dwharder/NumericalAnalysis/10RootFinding/bisection/

Engineer Scraps

Sunday, October 23, 2011

Numerical Analysis (The Bisection Method)