Root Finding
Definitions
Formulas
Bisection Method
Properties
- Guaranteed convergence: Always finds a root if one exists in the bracket
- Linear convergence: Error halved each iteration (order 1)
- Error bound: |x_n − x*| ≤ (b − a) / 2^n
- Derivative-free: Only requires function evaluations
- Drawback: Slow; cannot find roots without a bracket; cannot find multiple roots simultaneously
Newton-Raphson Method
Properties
- Quadratic convergence: Error squared each iteration near a simple root (order 2)
- Requires derivative: f'(x) must be known and computable
- Requires good initial guess: Can diverge if x_0 is far from the root
- Division by zero: Fails when f'(x_n) = 0
- Multiple roots: Convergence degrades to linear for roots with multiplicity > 1
Secant Method
Properties
- Superlinear convergence: Order φ ≈ 1.618 (golden ratio)
- No derivative needed: Uses finite difference approximation
- Two initial points: Requires x_0 and x_1
- Not guaranteed to converge: Can diverge or cycle if initial guesses are poor
- No bracket required: Unlike bisection, does not need a sign change
Fixed-Point Iteration
Choosing g(x)
For f(x) = x³ − x − 2 = 0, possible reformulations:
- g₁(x) = x³ − 2 -> g₁'(x) = 3x², |g₁'(1.5)| = 6.75 > 1 (diverges)
- g₂(x) = (x + 2)^(1/3) -> g₂'(x) = 1/[3(x+2)^(2/3)], |g₂'(1.5)| = 0.21 < 1 (converges)
- g₃(x) = x − f(x)/f'(x) -> This is Newton's method
Examples
Theorems
Python Implementation
import numpy as np
from typing import Callable, Tuple, List, Optional
class RootFinder:
"""Comprehensive root-finding toolkit."""
@staticmethod
def bisection(f: Callable, a: float, b: float,
tol: float = 1e-10, max_iter: int = 100) -> Tuple[float, List]:
"""Bisection method with full history."""
if f(a) * f(b) > 0:
raise ValueError("No sign change in [a, b]")
if f(a) * f(b) == 0:
return a if f(a) == 0 else b, []
history = []
for i in range(max_iter):
c = (a + b) / 2
fc = f(c)
history.append({'iter': i, 'x': c, 'f(x)': fc,
'interval': (a, b), 'width': b - a})
if abs(fc) < tol or (b - a) / 2 < tol:
return c, history
if f(a) * fc < 0:
b = c
else:
a = c
return (a + b) / 2, history
@staticmethod
def newton(f: Callable, df: Callable, x0: float,
tol: float = 1e-10, max_iter: int = 100) -> Tuple[float, List]:
"""Newton-Raphson method."""
history = []
x = x0
for i in range(max_iter):
fx = f(x)
dfx = df(x)
history.append({'iter': i, 'x': x, 'f(x)': fx, 'df(x)': dfx})
if abs(fx) < tol:
return x, history
if abs(dfx) < 1e-15:
raise RuntimeError(f"Derivative near zero at x={x}")
x = x - fx / dfx
return x, history
@staticmethod
def secant(f: Callable, x0: float, x1: float,
tol: float = 1e-10, max_iter: int = 100) -> Tuple[float, List]:
"""Secant method."""
history = []
f0, f1 = f(x0), f(x1)
for i in range(max_iter):
history.append({'iter': i, 'x': x1, 'f(x)': f1})
if abs(f1) < tol:
return x1, history
if abs(f1 - f0) < 1e-15:
raise RuntimeError("Denominator too small")
x_new = x1 - f1 * (x1 - x0) / (f1 - f0)
x0, f0 = x1, f1
x1 = x_new
f1 = f(x1)
return x1, history
@staticmethod
def fixed_point(g: Callable, x0: float,
tol: float = 1e-10, max_iter: int = 100) -> Tuple[float, List]:
"""Fixed-point iteration."""
history = []
x = x0
for i in range(max_iter):
x_new = g(x)
history.append({'iter': i, 'x': x, 'x_new': x_new,
'diff': abs(x_new - x)})
if abs(x_new - x) < tol:
return x_new, history
x = x_new
return x, history
@staticmethod
def brentq(f: Callable, a: float, b: float,
tol: float = 1e-10, max_iter: int = 100) -> Tuple[float, List]:
"""Brent's method — combines bisection with inverse quadratic interpolation."""
if f(a) * f(b) > 0:
raise ValueError("No sign change")
fa, fb = f(a), f(b)
history = []
for i in range(max_iter):
if abs(fb) < tol:
return b, history
# Try inverse quadratic interpolation
c, fc = a, fa
d = b - a
e = d
if abs(fa) > abs(fb):
a, b = b, a
fa, fb = fb, fa
m = 0.5 * (a - b)
if abs(m) < tol or fb == 0:
return b, history
if abs(e) >= tol and abs(fa) > abs(fb):
s = fb / fa
if a == c:
# Secant
p = 2 * m * s
q = 1 - s
else:
# Inverse quadratic
q = fa / fc
r = fb / fc
p = s * (2 * m * q * (q - r) - (b - a) * (r - 1))
q = (q - 1) * (r - 1) * (s - 1)
if p > 0:
q = -q
else:
p = -p
if 2 * p < min(3 * m * q - tol * abs(q),
abs(e * q)):
e = d
d = p / q
else:
d = m
e = m
else:
d = m
e = m
a, fa = b, fb
if abs(d) > tol:
b = b + d
else:
b = b - tol if m > 0 else b + tol
fb = f(b)
history.append({'iter': i, 'x': b, 'f(x)': fb})
return b, history
def demonstrate_methods():
"""Compare all methods on several test functions."""
print("=" * 60)
print("ROOT FINDING METHOD COMPARISON")
print("=" * 60)
tests = [
("x^3 - x - 2 = 0",
lambda x: x**3 - x - 2,
lambda x: 3*x**2 - 1,
(1, 2), 1.5),
("cos(x) - x = 0",
lambda x: np.cos(x) - x,
lambda x: -np.sin(x) - 1,
(0, 1), 0.5),
("e^x - 3 = 0",
lambda x: np.exp(x) - 3,
lambda x: np.exp(x),
(0, 2), 1.0),
]
finder = RootFinder()
for name, f, df, bracket, x0 in tests:
print(f"\n--- {name} ---")
root, _ = finder.bisection(f, bracket[0], bracket[1])
print(f" Bisection: {root:.12f}")
root, _ = finder.newton(f, df, x0)
print(f" Newton: {root:.12f}")
root, _ = finder.secant(f, bracket[0], bracket[1])
print(f" Secant: {root:.12f}")
demonstrate_methods()
Applications in AI/ML
Optimization Connection
| Root Finding | Optimization Equivalent |
|---|---|
| Find x where f(x) = 0 | Find x where ∇f(x) = 0 |
| Newton: x − f/f' | Newton: x − H⁻¹∇f |
| Secant approximation | BFGS Hessian approximation |
| Bisection | Golden section search |
Specific Applications
- Logistic regression: Solving for weights requires finding roots of the gradient
- Newton's method in neural networks: Second-order optimizers use Hessian inverse
- Brent's method in scipy.optimize: Used for line search in L-BFGS
- Fixed-point iteration: Expectation-Maximization (EM) algorithm is a fixed-point iteration
- Eigenvalue problems: Finding eigenvalues reduces to finding roots of det(A − λI) = 0
Common Mistakes
| Mistake | Problem | Solution |
|---|---|---|
| No bracket for bisection | May miss the root entirely | Verify f(a)·f(b) < 0 before starting |
| Starting Newton far from root | Method may diverge | Use bisection to get close, then switch to Newton |
| Ignoring f'(x) = 0 in Newton | Division by zero crashes | Check |
| Using fixed-point iteration blindly | May diverge if | g'(x*) |
| Not checking convergence | May return early or iterate forever | Use both |
| Assuming unique root | May find wrong root | Use bracketing methods to isolate roots first |
| Too tight tolerance | More iterations than needed | Set tolerance relative to problem scale |
| Ignoring floating point issues | Cancellation in secant method | Use robust methods like Brent's when possible |
Interview Questions
Q1: When would you use bisection vs Newton's method? A: Use bisection when you need guaranteed convergence and have a bracket, or when derivatives are unavailable. Use Newton's when you have a good initial guess and the derivative, for faster convergence. In practice, use a hybrid: bisection to bracket, Newton to refine.
Q2: Why does Newton's method fail for multiple roots? A: At a root of multiplicity m, f'(x*) = 0, so the tangent line is nearly horizontal. The iteration becomes x_{n+1} = x_n − m·f(x_n)/f'(x_n), converging linearly instead of quadratically. The modification requires knowing m, which is often unknown.
Q3: Explain the secant method's convergence order. A: The secant method has order p = (1+√5)/2 ≈ 1.618. This is proven by analyzing the error recurrence e_{n+1} ≈ C·e_n·e_{n-1}, which leads to the characteristic equation p² = p + 1, giving the golden ratio. It is faster than linear (bisection) but slower than quadratic (Newton).
Q4: How does Brent's method combine bisection and interpolation? A: Brent's method uses inverse quadratic interpolation when three points are available and the interpolation step stays within the current bracket. If the interpolation step would leave the bracket or fails to reduce the error sufficiently, it falls back to bisection. This gives the robustness of bisection with the speed of higher-order methods.
Q5: Give an example where fixed-point iteration diverges. A: For x² − 2 = 0, rewriting as g(x) = x² + 2 gives |g'(√2)| = 2√2 ≈ 2.83 > 1, so iteration diverges. Rewriting as g(x) = 2/x gives |g'(√2)| = 2/2 = 1, which is borderline. Only g(x) = (x + 2/x)/2 (Newton's method) has |g'(√2)| = 0, converging fastest.
Q6: How is root finding used in line search for optimization? A: Line search finds α that minimizes f(x + αd). This requires finding where d/dα f(x + αd) = 0, a one-dimensional root-finding problem. Brent's method is commonly used because it combines guaranteed convergence with speed. Wolfe conditions provide stopping criteria that ensure sufficient decrease.
Practice Problems
Quick Reference
Cross-References
- 086-Floating-Point — Rounding errors affect convergence criteria; machine epsilon determines achievable tolerance
- 088-Integration — Quadrature methods use function evaluations subject to floating point error
- 090-Error Analysis — Condition numbers determine sensitivity of root locations to perturbations in f
- Optimization — Root finding of ∇f = 0 is the foundation of optimization; Newton's method generalizes to Newton's method for optimization
- Linear Algebra — Eigenvalue computation requires finding roots of the characteristic polynomial
- Differential Equations — Implicit ODE solvers require solving nonlinear equations at each timestep