MATH566

MATH566 Lesson 14: Root refinement

Overview

Motivation:

Root refinement: apply function approximation
- linear approximants: secant, Newton
- quadratic approximants: Halley, Brent

Root refinement

Bisection has linear order of convergence. Are there faster methods?
Observations:
- Bisection only uses function values at interval endpoints
- Function behavior over the interval is not used
Idea: approximate the function over an interval
$f (x) = 0 \leftrightarrow f (x) ≅ g (x) = 0$
Choices:
- $g$ first-degree polynomial interpolant
- $g$ second-degree polynomial interpolant
- $g$ first-degree Taylor polynomial (series truncation)
- $g$ second-degree Taylor polynomial (series truncation)

Secant method

Instead of $f (x) = 0$ solve $g (x) = 0$ with $f ≅ g$ , $g$ linear approximant
Goal construct approximation sequence ${x_{n}}_{n \in ℕ}$ , $x_{n} \to z$ , $f (z) = 0$
Linear approximant = polynomial of degree $m = 1$ , requires 2 data points
$𝒟 = {(x_{i}, y_{i} = f (x_{i})), i = n - 1, n}$ $p (t) = \sum_{i = 0}^{1} y_{n - 1 + i} ℓ_{n - 1 + i} (t) = y_{n - 1} \frac{t - x_{n}}{x_{n - 1} - x_{n}} + y_{n} \frac{t - x_{n - 1}}{x_{n} - x_{n - 1}}$
Find next term in approximation sequence ${x_{n}}_{n \in ℕ}$ by setting $p (x_{n + 1}) = 0 \Rightarrow$
$y_{n - 1} (x_{n + 1} - x_{n}) - y_{n} (x_{n + 1} - x_{n - 1}) = 0 \Rightarrow$ $x_{n + 1} = \frac{x_{n - 1} y_{n} - x_{n} y_{n - 1}}{y_{n} - y_{n - 1}} = x_{n} - \frac{y_{n}}{y_{n} - y_{n - 1}} (x_{n} - x_{n - 1})$
The above is known as the secant method, and requires 1 function evaluation per iteration, same as bisection.

Secant method analysis

Error analysis: define $e_{n} = x_{n} - z$
Replace in secant iteration formula $x_{n + 1} = x_{n} - \frac{y_{n}}{y_{n} - y_{n - 1}} (x_{n} - x_{n - 1})$

$e_{n + 1} = e_{n} - \frac{y_{n}}{y_{n} - y_{n - 1}} (e_{n} - e_{n - 1})$ (1)
Assume $f \in C^{2}$ , and use Taylor series expansions around root $z$ ( $f (z) = 0$ )
$y_{n} = f^{'} (z) e_{n} + \frac{1}{2} f^{''} (z) e_{n}^{2} + \dots, y_{n - 1} = f^{'} (z) e_{n - 1} + \frac{1}{2} f^{''} (z) e_{n - 1}^{2} + \dots$
Replace in (1) (use notation $f^{'} = f^{'} (z)$ , $f^{''} = f^{''} (z)$ , $c = f^{''} / (2 f^{'})$ )
$\frac{y_{n}}{y_{n} - y_{n - 1}} ≅ \frac{f^{'} (z) e_{n} + \frac{1}{2} f^{''} (z) e_{n}^{2}}{f^{'} (z) (e_{n} - e_{n - 1}) + \frac{1}{2} f^{''} (z) {(e_{n} - e_{n - 1})}^{2}} \Rightarrow$ $e_{n + 1} = e_{n} [1 - \frac{f^{'} + \frac{1}{2} f^{''} e_{n}}{f^{'} + \frac{1}{2} (e_{n} + e_{n - 1}) f^{''}}] = e_{n} [1 - \frac{1 + c e_{n}}{1 + (e_{n} + e_{n - 1}) c}] \Rightarrow$

Superlinear convergence of secant method

Recall that for small $δ$
$\frac{1}{1 + δ} = 1 - δ + δ^{2} - \dots$
Apply to with $δ = e_{n - 1} / (1 + c e_{n})$
$e_{n + 1} = e_{n} [1 - \frac{1 + c e_{n}}{1 + (e_{n} + e_{n - 1}) c}] = e_{n} [1 - \frac{1}{1 + c δ}] ≅ e_{n} c δ = \frac{c e_{n - 1} e_{n}}{1 + c e_{n}}$
Assume $e_{n}$ is small such that $1 + c e_{n} ≅ 1$ to obtain $e_{n + 1} = c e_{n - 1} e_{n}$
At what order does ${e_{n}}$ converge to zero? For large enough $n$
$| e_{n} | ≅ r {| e_{n - 1} |}^{p}, | e_{n + 1} | ≅ r {| e_{n} |}^{p} = > r^{p + 1} {| e_{n - 1} |}^{p^{2}} = | c | r {| e_{n - 1} |}^{p + 1}$
Since $e_{n} \to 0$ , obtain $p^{2} - p - 1 = 0$ with solution $p = ϕ = (1 + \sqrt{5}) / 2 ≅ 1.62$
Secant converges at rate $p$ faster than linear, slower than quadratic

Newton method

Instead of $f (x) = 0$ solve $g (x) = 0$ with $f ≅ g$ , $g$ linear approximant
Goal construct approximation sequence ${x_{n}}_{n \in ℕ}$ , $x_{n} \to z$ , $f (z) = 0$
Linear approximant = Taylor polynomial of degree $m = 1$
$p (t) = f (x_{n}) + f^{'} (x_{n}) (t - x_{n})$
Find next term in approximation sequence ${x_{n}}_{n \in ℕ}$ by setting $p (x_{n + 1}) = 0 \Rightarrow$
$f (x_{n}) + f^{'} (x_{n}) (x_{n + 1} - x_{n}) = 0 \Rightarrow$ $x_{n + 1} = x_{n} - \frac{f (x_{n})}{f^{'} (x_{n})}$
The above is known as the Newton method, and requires 1 function evaluation and 1 derivative evaluation per iteration, more expensive than secant or bisection.

Newton method analysis

Error analysis: define $e_{n} = x_{n} - z$
Replace in Newton iteration formula $x_{n + 1} = x_{n} - \frac{f (x_{n})}{f^{'} (x_{n})}$

$e_{n + 1} = e_{n} - \frac{f (x_{n})}{f^{'} (x_{n})}$ (2)
Assume $f \in C^{2}$ , and use Taylor series expansions around root $z$ ( $f (z) = 0$ )
$f (x_{n}) = f^{'} e_{n} + \frac{1}{2} f^{''} e_{n}^{2} + \dots, f^{'} (x_{n}) = f^{'} + f^{''} (z) e_{n} + \frac{1}{2} f^{'''} (z) e_{n}^{2} + \dots$
Replace in (2) to obtain
$e_{n + 1} = e_{n} - \frac{f^{'} e_{n} + \frac{1}{2} f^{''} e_{n}^{2} + \dots +}{f^{'} + f^{''} e_{n} + \dots} ≅ \frac{1}{2} \frac{f^{''}}{f^{'}} e_{n}^{2}$
Newton's method is second order (as long as $f^{'} \neq 0$ )
If $f^{'} = f^{'} (z) = 0$ (double roots) convergence drops to first order

Muller's method

Instead of $f (x) = 0$ solve $g (x) = 0$ with $f ≅ g$ , $g$ quadratic approximant
Goal construct approximation sequence ${x_{n}}_{n \in ℕ}$ , $x_{n} \to z$ , $f (z) = 0$
Quadratic approximant = polynomial of degree $m = 2$ , Newton form
$p (t) = y_{n - 2} + [y_{n - 1}, y_{n - 2}] (t - x_{n - 2}) + [y_{n}, y_{n - 1}, y_{n - 2}] (t - x_{n - 2}) (t - x_{n - 1})$
Find next term in approximation sequence ${x_{n}}_{n \in ℕ}$ by setting $p (x_{n + 1}) = 0 \Rightarrow$
$x_{n + 1} = x_{n} - \frac{2 y_{n}}{w \pm \sqrt{w^{2} - 4 y_{n} [y_{n}, y_{n - 1}, y_{n - 2}]}},$
The $w$ quantity involves divided differences
$w = [y_{n}, y_{n - 1}] + [y_{n}, y_{n - 2}] - [y_{n - 1}, y_{n - 2}]$
This is Muller's method, and requires 1 function evaluation per iteration and has order of convergence 1.84 (superlinear, better than secant 1.62)

Halley's method

Instead of $f (x) = 0$ solve $g (x) = 0$ with $f ≅ g$ , $g$ quadratic approximant
Quadratic approximant = Taylor polynomial of degree $m = 2$
$p (t) = f (x_{n}) + f^{'} (x_{n}) (t - x_{n}) + \frac{1}{2} f^{''} (x_{n}) {(t - x_{n})}^{2}$
Find next term in approximation sequence ${x_{n}}_{n \in ℕ}$ by setting $p (x_{n + 1}) = 0 \Rightarrow$
$x_{n + 1} = x_{n} - \frac{2 f (x_{n}) f^{'} (x_{n})}{2 {[f^{'} (x_{n})]}^{2} - f (x_{n}) f^{''} (x_{n})}$
The above is known as the Halley's method, and requires 1 function evaluation and 2 derivative evaluations per iteration, more expensive than Newton, but has cubic order of convergence.