Geek Answers Handbook

How to make a Bisection method in Python

I want to create a Python program that will run the bisection method to determine the root:

f(x) = -26 + 85x - 91x2 +44x3 -8x4 + x5

The Bisection method is a numerical method for estimating the roots of the polynomial f (x).

Is there any pseudo code, algorithms or libraries available that I could use to tell me the answer?

+4
source share
4 answers

Here is the code showing the basic technique:

 >>> def samesign(a, b): return a * b > 0 >>> def bisect(func, low, high): 'Find root of continuous function where f(low) and f(high) have opposite signs' assert not samesign(func(low), func(high)) for i in range(54): midpoint = (low + high) / 2.0 if samesign(func(low), func(midpoint)): low = midpoint else: high = midpoint return midpoint >>> def f(x): return -26 + 85*x - 91*x**2 +44*x**3 -8*x**4 + x**5 >>> x = bisect(f, 0, 1) >>> print x, f(x) 0.557025516287 3.74700270811e-16
+6
source

You can see the solution in the previous question here , which uses scipy.optimize.bisect . Or, if your goal is learning, the bisection pseudo-code in a Wikipedia entry is a good guide for doing your own implementation in Python, as suggested by a commenter on an earlier question.

+4
source

My implementation is more general and yet simpler than other solutions: (and a public domain)

 def solve(func, x = 0.0, step = 1e3, prec = 1e-10): """Find a root of func(x) using the bisection method. The function may be rising or falling, or a boolean expression, as long as the end points have differing signs or boolean values. Examples: solve(lambda x: x**3 > 1000) to calculate the cubic root of 1000. solve(math.sin, x=6, step=1) to solve sin(x)=0 with x=[6,7). """ test = lambda x: func(x) > 0 # Convert into bool function begin, end = test(x), test(x + step) assert begin is not end # func(x) and func(x+step) must be on opposite sides while abs(step) > prec: step *= 0.5 if test(x + step) is not end: x += step return x
+1
source

With the admission:

 # there is only one root def fn(x): return x**3 + 5*x - 9 # define bisection method def bisection( eq, segment, app = 0.3 ): a, b = segment['a'], segment['b'] Fa, Fb = eq(a), eq(b) if Fa * Fb > 0: raise Exception('No change of sign - bisection not possible') while( b - a > app ): x = ( a + b ) / 2.0 f = eq(x) if f * Fa > 0: a = x else: b = x return x #test it print bisection(fn,{'a':0,'b':5}, 0.00003) # => 1.32974624634

Live: http://repl.it/k6q

0
source

More articles:


All Articles