11

Is there a bisection method I can find online, specifically for python?

For example, given these equations how can I solve them using the bisection method?

x^3 = 9  
3 * x^3 + x^2 = x + 5  
cos^2x + 6 = x
Hooked
  • 84,485
  • 43
  • 192
  • 261
bbnn
  • 3,505
  • 10
  • 50
  • 68
  • I wish my numerical methods course had used Python. :/ This is really instructive to implement yourself; just read Wikipedia's description for the algorithm. – Derrick Turk Dec 01 '10 at 17:12
  • Best to use something that's already been in use by many people than to try to write it yourself. [75%-90% of binary search implementations are incorrect.](https://en.wikipedia.org/wiki/Binary_search_algorithm#Implementation_issues) – endolith Apr 16 '13 at 21:21

2 Answers2

16

Using scipy.optimize.bisect:

import scipy.optimize as optimize
import numpy as np

def func(x):
    return np.cos(x)**2 + 6 - x

# 0<=cos(x)**2<=1, so the root has to be between x=6 and x=7
print(optimize.bisect(func, 6, 7))
# 6.77609231632

optimize.bisect calls _zeros._bisect, which is implemented in C.

unutbu
  • 842,883
  • 184
  • 1,785
  • 1,677
  • how to get the a and b value? and also the number of loop? – bbnn Dec 01 '10 at 17:00
  • like in the example http://math.fullerton.edu/mathews/n2003/BisectionMod.html it is trying to solve this equation (x^3+4x^2-10=0) and it uses function Bisection[1,2,30] how to get the number 1 2 and 30? is it from the equation? – bbnn Dec 01 '10 at 17:04
  • @bn: To use `bisect`, you must supply `a` and `b` such that `func(a)` and `func(b)` have opposite signs, thus guaranteeing that there is a root in `[a,b]` since `func` is required to be continuous. You could try to guess the values for `a` and `b`, use a bit of analysis, or if you want to do it programmatically, you could devise some method of generating candidate `a` and `b` until you find two that have opposite signs. That's all beyond the scope of the simple `bisect` method however. – unutbu Dec 01 '10 at 17:06
  • @bn: Regarding the number of iterations, `scipy.optimize.bisect` has a `maxiter` argument. For example, `so.bisect(func,6,8,maxiter=500)`. See http://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.bisect.html#scipy.optimize.bisect. – unutbu Dec 01 '10 at 17:08
  • @unutbu, thanks! for (equation x^3 = 9 ) i tried value of 0 for a and 1000 for b (2.10571289062) . and it gives me different answer with 0 for a and 5 for b (2.0802307129) according to wolfram, this 2.0802307129 is the answer http://www.wolframalpha.com/input/?i=x^3+%3D+9++ – bbnn Dec 01 '10 at 17:16
  • @bn: For a odd-order polynomial like `f(x) = x^3 + 4x^2 - 10`, the `x^3` term is going to dominate when `x` is large and negative, and when `x` is large and positive. When `x` is large and negative, `f(x)` is going to be negative, and similarly, when `x` is large and positive, `f(x)` is going to be positive. So it is easy to guess `a` and `b` in this case: Just choose a large negative value (like -100) for `a` and a large positive value (like 100) for `b`. As long as the number you choose is larger (in magnitude) than the coefficients, it should work. – unutbu Dec 01 '10 at 17:17
  • @bn: The `bisect` method quits when the root is known to lie within `xtol` of the returned value. You can try to make the answer more accurate by making the `xtol` parameter smaller. For example, `bisect(func,0,5,xtol=1e-15)`. – unutbu Dec 01 '10 at 17:24
  • @bn: Correction: the number you choose should be *significantly* larger than the coefficients. Being just slightly larger than the largest coefficient might not cut it. – unutbu Dec 01 '10 at 17:26
  • @bn: If you're looking for the `a` and `b` values just plot the function over the range you are interested in. This is especially important if the function has a lot of zeros. – dtlussier Dec 01 '10 at 18:37
  • Please remove the link to the **execrable** implementation in Python. Thanks. – Olivier Verdier Jun 05 '12 at 12:08
0

This could help you!

import numpy as np

def fn(x):
    # This the equation to find the root
    return (x**3 - x - 1) #x**2 - x - 1

def find_root_interval():
    for x in range(0, 1000):
        if fn(x) < 0:
            lower_interval = x
            if fn(x+1) > 0:
                higher_interval = x + 1
                return lower_interval, higher_interval
    return False

def bisection():
    a,b = find_root_interval()
    print("Interval: [{},{}]".format(a,b))
    # Create a 1000 equally spaced values between interval
    mid = 0
    while True:
        prev_mid = mid
        mid = (a+b)/2
        print("Mid value: "+str(mid))
        # 0.0005 is set as the error range
        if abs(mid-prev_mid) < 0.0005:
            return mid
        elif fn(mid) > 0:
            b = mid
        else:
            a = mid

root = bisection()
print(root)
Shamanth Muroor
  • 3
  • 3