Neither Newton's method nor the secant method can guarantee that an existing solution will be found (see Exercise 73: Understand why Newton's method can fail and Exercise 74: See if the secant method fails). The bisection method, however, does that. However, if there are several solutions present, it finds only one of them, just as Newton's method and the secant method. The bisection method is slower than the other two methods, so reliability comes with a cost of speed.
To solve \( x^2 - 9 = 0 \), \( x \in \left[0, 1000\right] \), with the bisection method, we reason as follows. The first key idea is that if \( f(x) = x^2 - 9 \) is continuous on the interval and the function values for the interval endpoints (\( x_L = 0 \), \( x_R =1000 \)) have opposite signs, \( f(x) \) must cross the \( x \) axis at least once on the interval. That is, we know there is at least one solution.
The second key idea comes from dividing the interval in two equal parts, one to the left and one to the right of the midpoint \( x_M = 500 \). By evaluating the sign of \( f(x_M) \), we will immediately know whether a solution must exist to the left or right of \( x_M \). This is so, since if \( f(x_M) \ge 0 \), we know that \( f(x) \) has to cross the \( x \) axis between \( x_L \) and \( x_M \) at least once (using the same argument as for the original interval). Likewise, if instead \( f(x_M) \le 0 \), we know that \( f(x) \) has to cross the \( x \) axis between \( x_M \) and \( x_R \) at least once.
In any case, we may proceed with half the interval only. The exception is if \( f(x_M) \approx 0 \), in which case a solution is found. Such interval halving can be continued until a solution is found. A "solution" in this case, is when \( |f(x_M)| \) is sufficiently close to zero, more precisely (as before): \( |f(x_M)| < \epsilon \), where \( \epsilon \) is a small number specified by the user.
The sketched strategy seems reasonable, so let us write a reusable function that can solve a general algebraic equation \( f(x)=0 \) (bisection_method.py):
def bisection(f, x_L, x_R, eps, return_x_list=False): f_L = f(x_L) if f_L*f(x_R) > 0: print "Error! Function does not have opposite \ signs at interval endpoints!" sys.exit(1) x_M = float(x_L + x_R)/2.0 f_M = f(x_M) iteration_counter = 1 if return_x_list: x_list =  while abs(f_M) > eps: if f_L*f_M > 0: # i.e. same sign x_L = x_M f_L = f_M else: x_R = x_M x_M = float(x_L + x_R)/2 f_M = f(x_M) iteration_counter += 1 if return_x_list: x_list.append(x_M) if return_x_list: return x_list, iteration_counter else: return x_M, iteration_counter def f(x): return x**2 - 9 a = 0; b = 1000 solution, no_iterations = bisection(f, a, b, eps=1.0e-6) print "Number of function calls: %d" % (1 + 2*no_iterations) print "A solution is: %f" % (solution)
Note that we first check if \( f \) changes sign in \( [a,b] \), because that is a requirement for the algorithm to work. The algorithm also relies on a continuous \( f(x) \) function, but this is very challenging for a computer code to check.
We get the following printout to the screen when
bisection_method.py is run:
Number of function calls: 61 A solution is: 3.000000
We notice that the number of function calls is much higher than with the previous methods.
As with the two previous methods, the function
placed in the file
nonlinear_solvers.py for easy import and use.