The Bisection Method With Python Code

in #mathematics7 years ago

Hi there. This page is about the bisection method with Python code for finding approximate roots or x-intercepts to polynomial equations.


Image Source

Sections


  • References
  • What Is The Bisection Method?
  • The Intermediate Value Theorem From Calculus
  • The Algorithm With Pseudocode
  • Python Code Example

References


The original posts are from my website here and here. Here are the references.

What Is The Bisection Method?


The Bisection method is a numerical method which finds approximate solutions to polynomial equations with the use of midpoints. Numerical methods provide approaches to certain mathematical problems when finding the exact numeric answers is not possible. Instead of finding exact answers, we find approximate solutions with a certain loss of precision. Some approximate answers are close to the real answer by 0.01 and some are as close as one out of a million.

The Intermediate Value Theorem


To help us understand the bisection method, I think it is important to go over the Intermediate Value Theorem (IVT) from calculus.

Suppose there is a continuous function f(x) in a closed interval [a, b] with f(a) and $f(b). Then there is some numbercin the interval [a, b] such that we have a ≤ c ≤ b and min(f(a), f(b)) ≤ f(c) ≤ max(f(a), f(b)). Thisf(c)would be a horizontal linek` or even 0.

Note that I have put in the min and max parts above to deal with both increasing and decreasing functions. If f(x) is an increasing function in the interval [a, b] then I have f(a) ≤ f(c) ≤ f(b)$. With the other case I would have f(a) ≥ f(c) ≥ f(b) if f(x) is a decreasing function in the interval [a, b].


Image Source: IVT

In the bisection method, we try to look for the number c such that f(c) = 0.

The Algorithm With Pseudocode


Here is a rough idea of what the bisection method function would look like. At each iteration, we want to shrink the interval [a, b].

# Pseudocode For Bisection Method

if f(a)*f(b) > 0

   print("No root found") # Both f(a) and f(b) are the same sign

else

   while (b - a) / 2 > tolerance 

       c = (b + a) /2 # c is like a midpoint

       if f(c) == 0

          return(midpt)  # The midpt is the root such that f(midpt) = 0

       else if f(a)*f(c) < 0

          b = c  # Shrink interval from right.

       else

          a = c

return(c)

In f(a)f(c) < 0, either is f(a) is negative or f(c) is negative. If f(a) is positive then the line has to cross the x-axis to get to f(c) which is negative. In this case, the new b would be c. In the increasing case, f(a) would be negative and f(c) would be positive. The choice here would be having the new b as c.

With the f(a)f(c) > 0 case, you would choose the new a as c.

This process is repeated until (b - a)/2 is small enough (or the interval is very small).

Some Picture Examples



Image Source

Suppose c is the midpoint in the [a,b] interval and the function f(x) is increasing. In this case, f(c) is negative, f(a) is negative and f(b) is positive. You would not choose c as the new b since both f(c) and f(a) are negative and there would be no root in the interval [a, c]. Instead, you choose c as the new a to obtain the new subinterval as [c, b]. You trim the left side of the [a. b] interval.


Bisection Method: Image Source

The midpoint x_s would be the new b giving the interval [a, x_s] because f(a) is positive and f(x_s) is negative.

Python Code Example


Here is one example done in Python. I find an approximate solution to three decimal places to the polynomial . Instead of c, I use the variable midpoint.

The user has the choice of choosing the [a, b] interval. If the real solution is not in this interval, no solution would be found from the bisection method. Be careful.

A tolerance or margin of error can be chosen. This is represented by tol.

# Bisection Method In Python

# Goal: Finding Roots to Polynomials
# Inputs: Function, a and b for interval [a, b], a guess in [a, b]
# Output: Roots/Zeroes To Polynomials

# Reference: Tim Sauer - Numerical Analysis Second Edition
# http://code.activestate.com/recipes/578417-bisection-method-in-python/

# https://stackoverflow.com/questions/6289646/python-function-as-a-function-argument

def f(x):
    return(x**2 - 11)

def bisection_method(a, b, tol):
    if f(a)*f(b) > 0:
        #end function, no root.
        print("No root found.")
    else:
        while (b - a)/2.0 > tol:
            midpoint = (a + b)/2.0
            if f(midpoint) == 0:
                return(midpoint) #The midpoint is the x-intercept/root.
            elif f(a)*f(midpoint) < 0: # Increasing but below 0 case
                b = midpoint
            else:
                a = midpoint
        
        return(midpoint)

answer = bisection_method(-1, 5, 0.0001)

print("Answer:", round(answer, 3))

Answer: 3.317

To three decimal places, the approximate solution to is 3.317.

Sort:  

This is awesome with programming.
Thanks

You can also implement this method in MATLAB and in R.