非線形方程式を数値的に解く
非線形方程式を計算機で解く方法について解説する。 齊藤宣一『数値解析』の第1章の読書ノートである。
1.1 ニュートン法
非線形方程式の解を求める方法として、ニュートン法がある。ニュートン法は、二分探索で解の範囲を狭めていく二分法よりも収束が速いことが知られている。
問
\(f \colon \mathbb{R} \to\mathbb{R}, C^2\)級とする。 方程式\(f(x)=0\)の解\(a\in\mathbb{R}\)を数値的に求めたい。
答
反復列\(\{x_k\}_{k=0}^\infty\)であって、\(x_k \to a(k \to \infty)\)となるものを構成する。
1. ニュートン法
\[x_{k+1} = x_k - \frac{f(x_k)}{f'(x_k)}\]2. 簡易ニュートン法
\[x_{k+1} = x_k - \frac{f(x_k)}{f'(x_0)}\]3. セカント法
\[x_{k+2} = x_{k+1} - \frac{x_{k+1}- x_k}{f(x_{k+1})-f(x_k)}f(x_{k+1})\]4. 緩和反復法
\[x_{k+1} = x_k - \beta f(x_k)\] \[(k = 0, 1, \ldots)\]Rem.
上の4つの方法のうち、ニュートン法が最も高速なので、ニュートン法にしぼって解説する。
1.2. 縮小写像の定理
ニュートン法が収束することを証明する。閉区間の完備性や、縮小写像の原理などを使う。このあたりの詳しいことは、内田伏一『集合と写像』を参照。
Th. （ニュートン法の収束）
\(f: \mathbb{R} \to \mathbb{R}, C^2 class\)で、解を\(a \in \mathbb{R}, f(a) = 0, f'(a) \neq 0\)とする。
このとき、初期値\(x_0 \in \mathbb{R}\)を解\(a\)の十分近くにとれば、
\[x_k \to a(k \to \infty)\]となる。
proof.
次の項を与える写像を
\[\phi: \mathbb{R} \setminus Z \to \mathbb{R}, x\mapsto x-\frac{f(x)}{f'(x)}\]とおく。ただし\(Z=\{x\colon f'(x)=0\}\)とする。このとき、
\[\varphi(a) = a, \varphi'(a) = 0\]である。
\(\varphi\)の連続性から、
\[\exists \delta > 0 |x-a| \leq \delta \implies |\varphi'(x) - \varphi'(a)| = |\varphi'(x)| \leq 1/2 =: \lambda\]となる。必要なら\(\delta\)を小さく取り直して、
\[J := [a-\delta, a+\delta] \subset \mathbb{R}\setminus Z\]とおく。
このとき、
\[x\in J \implies\varphi(x) \in J\]が成り立つ。（\(\varphi(a) = a\)と\(\varphi\)に対する平均値の定理を使う）
したがって、
\[\varphi: J \to J\]となる。初期値を
\[x_0\in J\]となるように選べば、
\[\{x_k\} \subset J\]となる。
縮小写像の原理から、
\[\{x_k\}\]は\(J\)上の値 \(b\in J\)に収束する。\(a=b\)を示す。
\(\varphi\)の連続性から
\[b = \lim x_{k+1} = \lim \varphi(x_k) = \varphi(\lim x_k) = \varphi(b)\]また、平均値の定理から
\[|a-b| = |\varphi(\xi)||a-b|\]いま、 \(\| \varphi (\xi )\| \leq \lambda = 1/2\) なので、
\[a = b\]したがって、反復列\(\{x_k\}\)は（初期値\(x_0\)を\(J\)内にとれば）解\(a\)に収束する。（証明終わり）
Rem.
\(\varphi(x) = \frac{f(x)f''(x)}{f'(x)^2}\) が\(0\)に近い範囲が大きいほど、ニュートン法は収束しやすい。
\(f''(x)\)が0に近いのは\(y=f(x)\)のグラフが直線に近いとき。グラフが直線ならば、ニュートン法は1回の反復で回に到達する。
\(f'(x)\)の絶対値が大きいのは接線の傾きが大きいとき。接線の傾きが大きいほど\(\|x_{k+1}-x_k\|\)が小さくなり、収束しやすい。
1.3 収束の速さ
Th 1.13. （ニュートン法の収束の速さ）
ニュートン法で、
\[\lim_{k\to\infty} \frac{x_{k+1}-a}{(x_k-a)^2} = \frac{1}{2}\frac{f''(a)}{f'(a)}\]が成り立つ。
proof.
テイラーの定理から
\[f(a) = f(x) + f'(x)(a-x) + \frac{1}{2}f''(\xi)(a-x)^2\]であり、いま\(f(a)=0\)である。ニュートン法の漸化式の定義と、上の式から、
\[x_{k+1} - a = x_k - \frac{f(x_k)}{f'(x_k)} - a = \frac{f''(\xi)}{2f'(x_k)}(a-x_k)^2\]となる。両辺を\((x_k-a)^2\)で割って\(k\to\infty\)の極限をとればよい。（証明終わり）
Rem.
\(f\)が\(C^3\)級ならば、
\[\frac{1}{2}\frac{f''(a)}{f'(a)} = \frac{1}{2}\varphi''(a)\]である。
1.4 停止条件
誤差\(e_k := x_k-a\)が\(\| e_k \| < \epsilon\)となったところで反復をやめたい。\(\mu = \frac{1}{2}\frac{f''(a)}{f'(a)}\)とおく。さっき証明した定理と、\(e_{k-1}^2 \ll e_{k-1}\)から、
\[x_k - x_{k-1} = e_k - e_{k-1} \simeq \mu e_{k-1}^2 + e_{k-1} \simeq e_{k-1}\]であるから、
\(\mu > R\)となるような\(R\)を一つ取り、
\[|x_k - x_{k-1}| < \left|\frac{\epsilon}{R}\right|^{1/2}\]となったところで反復をやめればよい。
Rem.
この式によって、計算した数値の精度の目安を与えることができる。
1.5 連立非線型方程式
連立非線型方程式の場合は、1変数のニュートン法での微分係数の逆数のかわりに、ヤコビ行列の逆行列をかければよい。
1.6 代数方程式に対する連立法
本記事の執筆時間の都合で省略。（連立法、アーバスの方法、etc.）