14 KiB
# Tutorial 4: Root-Finding Methods
In this tutorial, we study classical root-finding algorithms for nonlinear equations. We will:
- Define the root-finding problem mathematically
- Derive several algorithms (bisection, fixed-point, Newton, secant)
- Discuss convergence conditions and error behavior
- Compare methods with worked examples using the
numethodspackage.
1. Problem Setup and Notation
We seek to solve a nonlinear scalar equation: \[ f(x) = 0, \quad f: \mathbb{R} \to \mathbb{R}. \] with a continuously differentiable function.
Root, residual, and error
- A root \(r\) satisfies \(f(r)=0\).
- Absolute error: \((e_k = |x_k - x^\star|)\).
- Residual: \((r_k = |f(x_k)|)\). Note that small residual does not always imply small error.
Multiplicity
A root \(r\) has multiplicity \(m\) if \[ f(r) = f'(r) = \dots = f^{(m-1)}(r) = 0, \quad f^{(m)}(r) \neq 0. \] If \(x^\star\) satisfies \(f(x^\star)=0\) and \(f'(x^\star)\ne 0\), we say the root is simple (multiplicity 1).
If \(f'(x^\star)=\cdots=f^{(m-1)}(x^\star)=0\) and \(f^{(m)}(x^\star)\ne 0\), we say the root has multiplicity (m).
- Simple roots (\(m=1\)): most methods converge rapidly.
- Multiple roots (\(m>1\)): convergence often slows.
2. Bisection Method
Assumption (Intermediate Value Theorem): If f is continuous on ([a,b]) and (f(a),f(b) < 0), then there exists \(x^\star\) in (a,b) with \(f(x^\star)=0\).
- Assumes \(f\) is continuous on \([a,b]\) with \(f(a)f(b)<0\).
- Repeatedly bisect interval and select subinterval containing the root.
Iteration: \[ c_k = \frac{a_k+b_k}{2}, \quad f(c_k). \]
Error bound: interval length halves each step: \[ |c_k-r| \le \frac{b-a}{2^k}. \]
- Convergence: linear, guaranteed.
from numethods import Bisection, FixedPoint, NewtonRoot, Secant, print_trace
import math
# Example function: f(x) = x^2 - 2
f = lambda x: x**2 - 2
df = lambda x: 2*x
# Bisection
bisect = Bisection(f, 0, 2, tol=1e-8)
root_b = bisect.solve()
print('Bisection root:', root_b)
steps = bisect.trace()
print("\nBisection Method Trace (x^2 - 2):")
print_trace(steps)Bisection root: 1.4142135605216026
Bisection Method Trace (x^2 - 2):
iter | a | b | c | f(a) | f(b) | f(c) | interval
--------------------------------------------------------------------------------------------------------
0 | 0 | 2 | 1 | -2 | 2 | -1 | 2
1 | 1 | 2 | 1.5 | -1 | 2 | 0.25 | 1
2 | 1 | 1.5 | 1.25 | -1 | 0.25 | -0.4375 | 0.5
3 | 1.25 | 1.5 | 1.375 | -0.4375 | 0.25 | -0.109375 | 0.25
4 | 1.375 | 1.5 | 1.4375 | -0.109375 | 0.25 | 0.0664062 | 0.125
5 | 1.375 | 1.4375 | 1.40625 | -0.109375 | 0.0664062 | -0.0224609 | 0.0625
6 | 1.40625 | 1.4375 | 1.42188 | -0.0224609 | 0.0664062 | 0.0217285 | 0.03125
7 | 1.40625 | 1.42188 | 1.41406 | -0.0224609 | 0.0217285 | -0.000427246 | 0.015625
8 | 1.41406 | 1.42188 | 1.41797 | -0.000427246 | 0.0217285 | 0.0106354 | 0.0078125
9 | 1.41406 | 1.41797 | 1.41602 | -0.000427246 | 0.0106354 | 0.00510025 | 0.00390625
10 | 1.41406 | 1.41602 | 1.41504 | -0.000427246 | 0.00510025 | 0.00233555 | 0.00195312
11 | 1.41406 | 1.41504 | 1.41455 | -0.000427246 | 0.00233555 | 0.000953913 | 0.000976562
12 | 1.41406 | 1.41455 | 1.41431 | -0.000427246 | 0.000953913 | 0.000263274 | 0.000488281
13 | 1.41406 | 1.41431 | 1.41418 | -0.000427246 | 0.000263274 | -8.20011e-05 | 0.000244141
14 | 1.41418 | 1.41431 | 1.41425 | -8.20011e-05 | 0.000263274 | 9.06326e-05 | 0.00012207
15 | 1.41418 | 1.41425 | 1.41422 | -8.20011e-05 | 9.06326e-05 | 4.31482e-06 | 6.10352e-05
16 | 1.41418 | 1.41422 | 1.4142 | -8.20011e-05 | 4.31482e-06 | -3.88434e-05 | 3.05176e-05
17 | 1.4142 | 1.41422 | 1.41421 | -3.88434e-05 | 4.31482e-06 | -1.72643e-05 | 1.52588e-05
18 | 1.41421 | 1.41422 | 1.41421 | -1.72643e-05 | 4.31482e-06 | -6.47477e-06 | 7.62939e-06
19 | 1.41421 | 1.41422 | 1.41421 | -6.47477e-06 | 4.31482e-06 | -1.07998e-06 | 3.8147e-06
20 | 1.41421 | 1.41422 | 1.41421 | -1.07998e-06 | 4.31482e-06 | 1.61742e-06 | 1.90735e-06
21 | 1.41421 | 1.41421 | 1.41421 | -1.07998e-06 | 1.61742e-06 | 2.68718e-07 | 9.53674e-07
22 | 1.41421 | 1.41421 | 1.41421 | -1.07998e-06 | 2.68718e-07 | -4.05632e-07 | 4.76837e-07
23 | 1.41421 | 1.41421 | 1.41421 | -4.05632e-07 | 2.68718e-07 | -6.84571e-08 | 2.38419e-07
24 | 1.41421 | 1.41421 | 1.41421 | -6.84571e-08 | 2.68718e-07 | 1.0013e-07 | 1.19209e-07
25 | 1.41421 | 1.41421 | 1.41421 | -6.84571e-08 | 1.0013e-07 | 1.58366e-08 | 5.96046e-08
26 | 1.41421 | 1.41421 | 1.41421 | -6.84571e-08 | 1.58366e-08 | -2.63102e-08 | 2.98023e-08
27 | 1.41421 | 1.41421 | 1.41421 | -2.63102e-08 | 1.58366e-08 | -5.23681e-09 | 1.49012e-08
3. Fixed-Point Iteration
- Rewrite equation as \(x=g(x)\).
- Iterate: \[ x_{k+1} = g(x_k). \]
Convergence theorem (Banach fixed-point): If g is continuously differentiable near \((x^\star)\) and \[ |g'(x_\star)| < 1, \] then for initial guesses \(x_0\) sufficiently close to \(x^\star\), the iterates converge linearly to \(x^\star\) with asymptotic rate \(|g'(x^\star)|\).
Choice of g. Different rearrangements yield different g's with different convergence properties. A poor choice (with \((|g'|\ge 1))\) can diverge.
- Rate: linear with factor \(|g'(r)|\).
# Fixed point: g(x)=sqrt(2) ~ rewriting
g = lambda x: (2/x) # not always convergent, but demonstrates
try:
fp = FixedPoint(g, 1.0, tol=1e-8)
root_fp = fp.solve()
print('Fixed-point root:', root_fp)
except Exception as e:
print('Fixed-point failed:', e)
g = lambda x: 0.5 * (x + 2 / x)
steps = FixedPoint(g, 1.0).trace()
print("\nFixed-Point Iteration Trace (x^2 - 2):")
print_trace(steps)Fixed-point failed: Fixed-point iteration did not converge
Fixed-Point Iteration Trace (x^2 - 2):
iter | x | x_new | error
----------------------------------------------------
0 | 1 | 1.5 | 0.5
1 | 1.5 | 1.41667 | 0.0833333
2 | 1.41667 | 1.41422 | 0.00245098
3 | 1.41422 | 1.41421 | 2.1239e-06
4 | 1.41421 | 1.41421 | 1.59495e-12
4. Newton’s Method
From Taylor expansion: \[ f(x) \approx f(x_k) + f'(x_k)(x-x_k). \] Set \(f(x)=0\) to solve for next iterate: \[ x_{k+1} = x_k - \frac{f(x_k)}{f'(x_k)}. \]
- Quadratic convergence for simple roots.
- For multiple roots: drops to linear.
- Requires derivative, sensitive to initial guess.
# Newton
newton = NewtonRoot(f, df, 1.0, tol=1e-12)
root_n = newton.solve()
print('Newton root:', root_n)
steps = newton.trace()
print("Newton Method Trace (x^2 - 2):")
print_trace(steps)Newton root: 1.4142135623730951
Newton Method Trace (x^2 - 2):
iter | x | f(x) | df(x) | x_new | error
------------------------------------------------------------------------------
0 | 1 | -1 | 2 | 1.5 | 0.5
1 | 1.5 | 0.25 | 3 | 1.41667 | 0.0833333
2 | 1.41667 | 0.00694444 | 2.83333 | 1.41422 | 0.00245098
3 | 1.41422 | 6.0073e-06 | 2.82843 | 1.41421 | 2.1239e-06
4 | 1.41421 | 4.51061e-12 | 2.82843 | 1.41421 | 1.59472e-12
5. Secant Method
Avoids derivative by approximating slope with finite difference: \[ x_{k+1} = x_k - f(x_k) \frac{x_k - x_{k-1}}{f(x_k)-f(x_{k-1})}. \]
Convergence order: \(\approx 1.618\) (superlinear).
More efficient than Newton if derivative expensive.
# Secant
sec = Secant(f, 0, 2, tol=1e-12)
root_s = sec.solve()
print('Secant root:', root_s)
steps = sec.trace()
print("\nSecant Method Trace (x^2 - 2):")
print_trace(steps)Secant root: 1.414213562373095
Secant Method Trace (x^2 - 2):
iter | x0 | x1 | x2 | f(x0) | f(x1) | error
-------------------------------------------------------------------------------------------
0 | 0 | 2 | 1 | -2 | 2 | 1
1 | 2 | 1 | 1.33333 | 2 | -1 | 0.333333
2 | 1 | 1.33333 | 1.42857 | -1 | -0.222222 | 0.0952381
3 | 1.33333 | 1.42857 | 1.41379 | -0.222222 | 0.0408163 | 0.0147783
4 | 1.42857 | 1.41379 | 1.41421 | 0.0408163 | -0.00118906 | 0.000418335
5 | 1.41379 | 1.41421 | 1.41421 | -0.00118906 | -6.00729e-06 | 2.12421e-06
6 | 1.41421 | 1.41421 | 1.41421 | -6.00729e-06 | 8.93146e-10 | 3.15775e-10
7 | 1.41421 | 1.41421 | 1.41421 | 8.93146e-10 | -8.88178e-16 | 2.22045e-16
6. Stopping Criteria
We stop iteration when:
- \(|f(x_k)| \le \varepsilon\) (residual small), or
- \(|x_{k+1}-x_k| \le \varepsilon(1+|x_{k+1}|)\) (relative step small).
7. Comparison of Methods
| Method | Requires derivative | Convergence rate | Guarantee? |
|---|---|---|---|
| Bisection | No | Linear | Yes (if sign change) |
| Fixed-Point | No | Linear | Not always |
| Newton | Yes | Quadratic | Locally (good guess) |
| Secant | No | ~1.618 (superlinear) | Locally (good guess) |
8. Exercises
- Apply all four methods to \(f(x)=\cos x - x\).
- Try Newton’s method on \(f(x)=(x-1)^2\) and compare convergence rate with \(f(x)=x^2-2\).
- Modify Secant to stop if denominator too small.