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