Bisection and False position - Executable Engineering: Interactive Numerical Methods in Python

import numpy as np
import plotly.graph_objects as go
from scipy.optimize import bisect

try:
    import executable_engineering as exe
except ImportError:
    %pip install -q executable_engineering
    import executable_engineering as exe

from executable_engineering import RootFinderClosed

Order of convergence¶

Iterative search methods are characterized by their order of convergence, which measure how successive guesses approach the true root. Given the true root $x$ , we can check how successive guesses approach it by calculating the error:

x^{i+1}-x \propto [x^{i}-x]^k

(1)

where $k$ is the order.

Bisection methods¶

Bisection methods are essentially a binary search for the root. If $f(x)$ is continuous between bounds $a$ and $b>a$ and $f(a)$ and $f(b)$ are opposite signs, there must be a point $c$ where $f(c)=0$ .

The algorithm is:

Given brackets $a$ and $b$
Calculate the midpoint $c = (b-a)/2$ .
If $f(c) \approx 0$ or $a \approx b$ : exit.
If $f(c)>0$ : set $a = c$
If $f(c)<0$ : set $b = c$
repeat

Graphically this is:

The error of the Bisection Method generaly follows:

x^{i+1}-x = \frac{1}{2} [x^{i}-x]

(2)

and therefore has a linear order of convergence ( $k=1$ ).

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

print('Setting x tolerance: \n')
xtols = [1e-1, 1e-2, 1e-3, 1e-4, 1e-5]
for xtol in xtols:
  print(bisect(f, 1, 2, maxiter = 100, xtol = xtol))

print('\nSetting relative tolerance \n', bisect(f, 1, 2, maxiter = 100, rtol = 1e-12))

Method of False Position (Regula falsi)¶

Let’s use more information! In bisection we are only interested in $f(a,b)$ switching signs, but it stands to reason the larger $f(b)$ is compared to $f(a)$ , the further the root is from $b$ !

The method of false position uses this information to determine the next candidate solution:

c= b - f(b) \frac{b-a}{f(b)-f(a)}

(3)

The same algorithm as for bisection is then applied to replace either $a$ or $b$ with $c$ so that the root remains bracketted.

Note that the line joining the brackets approximates the tangent of the curve.

The method of False Position has an order of convergence of 1.618 (the Golden Ratio!):

x^{i+1}-x \propto [x^i-x]^{1.618}

(4)

This is considered superlinear and a good thing!

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

# Example values for a and b
a = 1
b = 2

# Calculate f(a) and f(b)
fa = f(a)
fb = f(b)

# Calculate the next candidate solution c using the method of false position
c = b - fb * (b - a) / (fb - fa)
fc = f(c)

# Generate x values for the plot
x = np.linspace(a - 0.5, b + 0.5, 100)

fig = go.Figure()
fig.add_trace(go.Scatter(x=x, y=f(x), mode='lines', name='f(x)'))
fig.add_trace(go.Scatter(x=[a, b], y=[fa, fb], mode='markers', marker=dict(size=10, color='blue'), name='Initial Bracket'))
fig.add_trace(go.Scatter(x=[c], y=[fc], mode='markers', marker=dict(size=10, color='red'), name='New Candidate (c)'))
fig.add_trace(go.Scatter(x=[a, b], y=[fa, fb], mode='lines', line=dict(dash='dash', color='gray'), name='Linear interpolation'))
fig.add_hline(y=0, line_dash='dot', line_color='black')
fig.update_layout(title='Method of False Position (Regula falsi) - One Iteration', xaxis_title='x', yaxis_title='f(x)')
fig.show()

Comparison of bracketting methods¶

We can compare bisection vs False position directly.

def f(x):
    return np.cos(x) - x

rf = RootFinderClosed(f, a=0, b=1, tol=1e-6, max_iter=20)
rf.show_toggle_closed()

Summary of bracketting methods¶

Bracketting methods are usually robust but slow to converge and generlization to N-D is not trivial.

Through incorporating the additional information of linear interpolation, the method of False Position achieves superlinear convergence.

But bracketting is complicated in N-D, so let’s try to remove it.