電卓で 2 \sqrt2 2 を求めるとき、中では「だいたいの答えを出し、それを少し直して、また直す」という繰り返しが行われています。式で解けない方程式でも、同じ発想で答えに近づけます。この記事では、そのやり方がいくつかあること、どれが速いか、そして本当に答えに近づいていると言える理由を、順に確かめます。読み終えると、状況に応じてどの方法を選べばよいかと、その根拠を説明できるようになります。
方程式 f ( x ) = 0 f(x) = 0 f ( x ) = 0 の解を求める反復法は、誤差 e k = x k − x ∗ e_k = x_k - x^* e k = x k − x ∗ がどんな速さで 0 0 0 に近づくか で比べます。速さは「次数」(誤差が次の回に何乗になるか)と「漸近誤差定数」(そのときの比例係数)の 2 つの数で表します。
二分法 は連続性だけで動き、k k k 回で誤差が ( b − a ) / 2 k + 1 (b-a)/2^{k+1} ( b − a ) / 2 k + 1 以下になります。遅いが必ず収束します。
不動点反復 x k + 1 = g ( x k ) x_{k+1} = g(x_k) x k + 1 = g ( x k ) は、∣ g ′ ( x ∗ ) ∣ < 1 |g'(x^*)| \lt 1 ∣ g ′ ( x ∗ ) ∣ < 1 なら 1 次収束し、g g g の p − 1 p-1 p − 1 階までの導関数が x ∗ x^* x ∗ で消えれば p p p 次収束します。証明は平均値の定理と Taylor の定理です。
Newton 法 は単根の近くで 2 次収束し、e k + 1 / e k 2 → f ′ ′ ( x ∗ ) / ( 2 f ′ ( x ∗ ) ) e_{k+1}/e_k^2 \to f''(x^*)/(2f'(x^*)) e k + 1 / e k 2 → f ′′ ( x ∗ ) / ( 2 f ′ ( x ∗ )) です。重根では 1 次収束に落ち、比は ( m − 1 ) / m (m-1)/m ( m − 1 ) / m です。凸関数 では初期値を解の右に取れば単調に収束します。
割線法 は導関数を使わずに、次数が黄金比 ( 1 + 5 ) / 2 ≈ 1.618 (1+\sqrt5)/2 \approx 1.618 ( 1 + 5 ) /2 ≈ 1.618 以上になります。
多変数の Newton 法も 2 次収束し、Kantorovich 型の定理 は初期値での計算だけから解の存在と一意性を保証します。証明は縮小写像 の原理に帰着します。
この記事で使う事前知識
平均値の定理・Taylor の定理 : 関数の変化量や 1 次近似の誤差を、途中のある点での導関数の値で表す定理。→ 平均値の定理とテイラーの定理
収束 : 点列が限りなくある点に近づくこと。→ 収束
縮小写像 : 2 点の距離を一定の割合 L < 1 L \lt 1 L < 1 以下に縮める写像。→ 縮小写像
凸関数 : グラフの 2 点を結ぶ線分が、グラフより上にある関数。→ 凸関数
x 2 = 2 x^2 = 2 x 2 = 2 の正の解、x = cos x x = \cos x x = cos x の解、ある固有値問題の特性方程式の根。解を式で書けない方程式は、近似解の列 x 0 , x 1 , x 2 , … x_0, x_1, x_2, \ldots x 0 , x 1 , x 2 , … を作って解に近づけるしかありません。ここで問うべきことは 3 つあります。
列は解に収束するか 。どんな初期値から出発すればよいか。
収束するならどれだけ速いか 。正しい桁数が 1 回の反復で 1 桁ずつ増えるのか、2 倍ずつ増えるのか。
計算をいつ止めてよいか 。止めたときの誤差をどう見積もるか。
二分法は 1 に強く、Newton 法は 2 に強い、という対照がこの章の軸です。二分法は区間の端で符号が違えば必ず収束しますが、1 回あたり 1 ビットしか進みません。Newton 法は解の近くでは 1 回ごとに正しい桁数がおよそ 2 倍になりますが、初期値が悪いと発散したり周期軌道に入ったりします。実用の算法(Brent 法など。二分法のような確実な方法と、速い方法を切り替えて使う算法)は両者を組み合わせています。
3 の問いには、不動点反復の事後誤差評価 と、Kantorovich 型の定理が答えます。前者は「直前の 2 つの反復値の差」から誤差の上界を出し、後者は初期値での値と導関数だけから、解が存在する球の半径を出します。後者は計算機援用証明(計算機の出力を数学の証明として使う方法)で使う区間 Newton 法(誤差の範囲を区間で追いながら行う Newton 法)の考え方の土台です(浮動小数と区間演算 )。
以下、実数列 ( x k ) (x_k) ( x k ) が x ∗ x^* x ∗ に収束するとき e k = x k − x ∗ e_k = x_k - x^* e k = x k − x ∗ と書きます。
次の定義のねらいは、「誤差が 0 に近づく速さ」を 1 つの数で言えるようにすることです。次の回の誤差が今の誤差の何乗に比例するかを見て、その指数を次数と呼びます。
定義 2.1 (収束の次数(Q 次数)と漸近誤差定数 ) ↗
x k → x ∗ x_k \to x^* x k → x ∗ で、すべての k k k で e k ≠ 0 e_k \ne 0 e k = 0 とする。p ≥ 1 p \ge 1 p ≥ 1 と 0 < C < ∞ 0 \lt C \lt \infty 0 < C < ∞ について
lim k → ∞ ∣ e k + 1 ∣ ∣ e k ∣ p = C \lim_{k\to\infty} \frac{|e_{k+1}|}{|e_k|^p} = C k → ∞ lim ∣ e k ∣ p ∣ e k + 1 ∣ = C が成り立ち、さらに p = 1 p = 1 p = 1 のときは C < 1 C \lt 1 C < 1 であるとき、列は p p p 次収束 する(Q 次数が p p p である)といい、C C C を漸近誤差定数 という。p = 1 , 2 , 3 p = 1, 2, 3 p = 1 , 2 , 3 をそれぞれ 1 次(線形)・2 次・3 次収束という。また ∣ e k + 1 ∣ / ∣ e k ∣ → 0 |e_{k+1}|/|e_k| \to 0 ∣ e k + 1 ∣/∣ e k ∣ → 0 のとき超 1 次収束 するという。
p = 2 p = 2 p = 2 で C C C がほどほどの大きさなら、∣ e k + 1 ∣ ≈ C ∣ e k ∣ 2 |e_{k+1}| \approx C|e_k|^2 ∣ e k + 1 ∣ ≈ C ∣ e k ∣ 2 から、正しい桁数(− log 10 ∣ e k ∣ -\log_{10}|e_k| − log 10 ∣ e k ∣ )が 1 回ごとにほぼ 2 倍になります。1 次収束では、桁数が 1 回ごとに − log 10 C -\log_{10} C − log 10 C ずつ増えるだけです。
Q 次数の定義は極限の存在を要求するので、誤差が規則正しく減らない列には使えません。そこで上界だけで速さを測る定義も用意します。2 種類あるのは、毎回きれいな比が取れる場合(Q 次数)と、上から押さえることしかできない場合(R 次数)の両方に使えるようにするためです。
定義 2.2 (R 次数 ) ↗
x k → x ∗ x_k \to x^* x k → x ∗ とする。p > 1 p \gt 1 p > 1 について、定数 A > 0 A \gt 0 A > 0 と 0 < q < 1 0 \lt q \lt 1 0 < q < 1 があってすべての k k k で
∣ e k ∣ ≤ A q p k |e_k| \le A\, q^{p^k} ∣ e k ∣ ≤ A q p k が成り立つとき、列は R 次数 p p p 以上 で収束するという。p = 1 p = 1 p = 1 では ∣ e k ∣ ≤ A q k |e_k| \le A q^k ∣ e k ∣ ≤ A q k を条件とする。
q p k q^{p^k} q p k は p k p^k p k が指数に乗るので、p = 2 p = 2 p = 2 なら q , q 2 , q 4 , q 8 , … q, q^2, q^4, q^8, \ldots q , q 2 , q 4 , q 8 , … と桁数が倍々に増える上界です。以下の定理では、多くの場合に両方の意味での次数を示します。
定理 3.1 (二分法の収束 ) ✓
a < b a \lt b a < b 、f : [ a , b ] → R f : [a,b] \to \mathbb{R} f : [ a , b ] → R は連続で f ( a ) f ( b ) < 0 f(a) f(b) \lt 0 f ( a ) f ( b ) < 0 とする。a 0 = a a_0 = a a 0 = a 、b 0 = b b_0 = b b 0 = b とし、k = 0 , 1 , 2 , … k = 0, 1, 2, \ldots k = 0 , 1 , 2 , … について c k = ( a k + b k ) / 2 c_k = (a_k + b_k)/2 c k = ( a k + b k ) /2 とおいて、f ( a k ) f ( c k ) < 0 f(a_k) f(c_k) \lt 0 f ( a k ) f ( c k ) < 0 なら ( a k + 1 , b k + 1 ) = ( a k , c k ) (a_{k+1}, b_{k+1}) = (a_k, c_k) ( a k + 1 , b k + 1 ) = ( a k , c k ) 、f ( c k ) f ( b k ) < 0 f(c_k) f(b_k) \lt 0 f ( c k ) f ( b k ) < 0 なら ( a k + 1 , b k + 1 ) = ( c k , b k ) (a_{k+1}, b_{k+1}) = (c_k, b_k) ( a k + 1 , b k + 1 ) = ( c k , b k ) と定める(f ( c k ) = 0 f(c_k) = 0 f ( c k ) = 0 なら c k c_k c k を解として止める)。止まらないとき、
各 [ a k , b k ] [a_k, b_k] [ a k , b k ] は f f f の零点を含み、b k − a k = ( b − a ) / 2 k b_k - a_k = (b-a)/2^k b k − a k = ( b − a ) / 2 k である。
( a k ) (a_k) ( a k ) と ( b k ) (b_k) ( b k ) は共通の極限 x ∗ x^* x ∗ に収束し、f ( x ∗ ) = 0 f(x^*) = 0 f ( x ∗ ) = 0 である。
すべての k k k で ∣ c k − x ∗ ∣ ≤ ( b − a ) / 2 k + 1 |c_k - x^*| \le (b-a)/2^{k+1} ∣ c k − x ∗ ∣ ≤ ( b − a ) / 2 k + 1 が成り立つ。
証明(定理 3.1) 1. 帰納法で f ( a k ) f ( b k ) < 0 f(a_k) f(b_k) \lt 0 f ( a k ) f ( b k ) < 0 と b k − a k = ( b − a ) / 2 k b_k - a_k = (b-a)/2^k b k − a k = ( b − a ) / 2 k を示します。k = 0 k = 0 k = 0 は仮定です。k k k で成り立つとき、f ( c k ) ≠ 0 f(c_k) \ne 0 f ( c k ) = 0 なので f ( c k ) f(c_k) f ( c k ) は f ( a k ) f(a_k) f ( a k ) と f ( b k ) f(b_k) f ( b k ) のどちらか一方と異符号であり、選んだ半分の区間の端点で f f f は異符号です。幅はちょうど半分になります。端点で異符号なら 0 0 0 は f ( a k ) f(a_k) f ( a k ) と f ( b k ) f(b_k) f ( b k ) の間にあるので、中間値の定理(定理 5.7[極限と連続性] )から [ a k , b k ] [a_k, b_k] [ a k , b k ] に零点があります。
2. 作り方から a k ≤ a k + 1 a_k \le a_{k+1} a k ≤ a k + 1 、b k + 1 ≤ b k b_{k+1} \le b_k b k + 1 ≤ b k 、a k < b k ≤ b a_k \lt b_k \le b a k < b k ≤ b です。( a k ) (a_k) ( a k ) は単調増加で上に有界なので、単調有界定理(定理 2.1[級数と収束判定] )から収束します。極限を x ∗ x^* x ∗ とします。b k = a k + ( b − a ) / 2 k b_k = a_k + (b-a)/2^k b k = a k + ( b − a ) / 2 k なので b k → x ∗ b_k \to x^* b k → x ∗ でもあります。a k a_k a k と b k b_k b k は [ a , b ] [a,b] [ a , b ] の点で、f f f は x ∗ x^* x ∗ で連続なので、数列による特徴づけ(定理 4.5[極限と連続性] 。a k = x ∗ a_k = x^* a k = x ∗ となる項では f ( a k ) = f ( x ∗ ) f(a_k) = f(x^*) f ( a k ) = f ( x ∗ ) は自明)から f ( a k ) → f ( x ∗ ) f(a_k) \to f(x^*) f ( a k ) → f ( x ∗ ) 、f ( b k ) → f ( x ∗ ) f(b_k) \to f(x^*) f ( b k ) → f ( x ∗ ) です。積の極限をとると、f ( a k ) f ( b k ) < 0 f(a_k) f(b_k) \lt 0 f ( a k ) f ( b k ) < 0 から f ( x ∗ ) 2 ≤ 0 f(x^*)^2 \le 0 f ( x ∗ ) 2 ≤ 0 、すなわち f ( x ∗ ) = 0 f(x^*) = 0 f ( x ∗ ) = 0 です。
3. a k a_k a k は増加して x ∗ x^* x ∗ に、b k b_k b k は減少して x ∗ x^* x ∗ に収束するので、a k ≤ x ∗ ≤ b k a_k \le x^* \le b_k a k ≤ x ∗ ≤ b k です。c k c_k c k はこの区間の中点なので、∣ c k − x ∗ ∣ ≤ ( b k − a k ) / 2 = ( b − a ) / 2 k + 1 |c_k - x^*| \le (b_k - a_k)/2 = (b-a)/2^{k+1} ∣ c k − x ∗ ∣ ≤ ( b k − a k ) /2 = ( b − a ) / 2 k + 1 です。
∎
誤差の上界が f f f によらず決まることが二分法の強みです。許容誤差 τ \tau τ を達成する回数は k ≥ log 2 ( ( b − a ) / τ ) − 1 k \ge \log_2((b-a)/\tau) - 1 k ≥ log 2 (( b − a ) / τ ) − 1 で、[ 0 , 1 ] [0,1] [ 0 , 1 ] から 10 − 12 10^{-12} 1 0 − 12 まで縮めるには約 40 回かかります。1 回で 1 ビットずつ、10 進では約 3.3 回で 1 桁ずつ進みます。Q 次数の言葉では、∣ c k − x ∗ ∣ |c_k - x^*| ∣ c k − x ∗ ∣ そのものは規則正しく減るとは限りませんが、上界 ( b − a ) / 2 k + 1 (b-a)/2^{k+1} ( b − a ) / 2 k + 1 が比 1 / 2 1/2 1/2 で減るので、R 次数 1 で収束します。
例 3.2 (二分法と Brent 法の反復回数 )
f ( x ) = cos x − x f(x) = \cos x - x f ( x ) = cos x − x を [ 0 , 1 ] [0,1] [ 0 , 1 ] で考えます。f ( 0 ) = 1 > 0 f(0) = 1 \gt 0 f ( 0 ) = 1 > 0 、f ( 1 ) = cos 1 − 1 < 0 f(1) = \cos 1 - 1 \lt 0 f ( 1 ) = cos 1 − 1 < 0 です。二分法と、SciPy の brentq(二分法・割線法・逆 2 次補間を組み合わせた Brent 法)で、区間幅の許容値ごとに反復回数を比べます。
# Python 3.12.10 / SciPy 1.18.1 / mpmath 1.3.0
from scipy.optimize import brentq
from mpmath import mp, findroot, cos as mcos
def bisect ( f , a , b , tol ) :
""" f(a)f(b) < 0。区間幅が tol 以下になるまで二分する """
f = lambda x : math. cos ( x ) - x # 根は x = cos x の解
ref = float ( findroot ( lambda x : mcos ( x ) - x , 0.7 ))
print ( " tol 二分法: 回数 誤差 brentq: 回数 誤差 " )
for tol in [ 1e-4 , 1e-8 , 1e-12 ] :
xb, nb = bisect ( f , 0.0 , 1.0 , tol )
xr, res = brentq ( f , 0.0 , 1.0 , xtol = tol , full_output = True )
print ( f " {tol :.0e } {nb :3d } { abs ( xb - ref ) :.2e } "
f " {res.iterations :2d } { abs ( xr - ref ) :.2e } " )
print ( " 理論の回数 ceil(log2(1/tol)): " ,
[ math. ceil ( math. log2 ( 1 / t )) for t in [ 1e-4 , 1e-8 , 1e-12 ]] )
tol 二分法: 回数 誤差 brentq: 回数 誤差
1e-04 14 2.01e-05 5 5.93e-06
1e-08 27 9.04e-10 6 1.35e-09
1e-12 40 5.17e-14 7 4.77e-15
理論の回数 ceil(log2(1/tol)): [14, 27, 40]
検算. (a) 参照値は mpmath の findroot で 30 桁まで求めた x = cos x x = \cos x x = cos x の解 0.739085 … 0.739085\ldots 0.739085 … です。(b) 二分法の誤差はどれも許容値の半分以下で、定理 3.1 の 3 の評価(最後の区間の中点なので幅の半分以下)の内側にあります。Brent 法は別の算法で、同じ解に許容値以下の誤差で達しています。(c) 二分法の回数は幅 2 − n ≤ τ 2^{-n} \le \tau 2 − n ≤ τ となる最小の n = ⌈ log 2 ( 1 / τ ) ⌉ n = \lceil \log_2(1/\tau) \rceil n = ⌈ log 2 ( 1/ τ )⌉ に一致し、許容値を 10 − 4 10^{-4} 1 0 − 4 倍するごとに約 13.3 回ずつ増えます。Brent 法は 1 回ずつしか増えず、超 1 次収束していることが読み取れます。
方程式 f ( x ) = 0 f(x) = 0 f ( x ) = 0 を x = g ( x ) x = g(x) x = g ( x ) の形に書き直し、x k + 1 = g ( x k ) x_{k+1} = g(x_k) x k + 1 = g ( x k ) と反復するのが不動点反復 です。g ( x ) = x − f ( x ) g(x) = x - f(x) g ( x ) = x − f ( x ) 、g ( x ) = x − f ( x ) / f ′ ( x ) g(x) = x - f(x)/f'(x) g ( x ) = x − f ( x ) / f ′ ( x ) など書き直し方はいくらでもあり、どれを選ぶかで収束の速さが変わります。
定理 4.1 (不動点反復の大域収束と誤差評価 ) ✓
a < b a \lt b a < b とし、g : [ a , b ] → [ a , b ] g : [a,b] \to [a,b] g : [ a , b ] → [ a , b ] が定数 0 ≤ L < 1 0 \le L \lt 1 0 ≤ L < 1 について
∣ g ( x ) − g ( y ) ∣ ≤ L ∣ x − y ∣ ( x , y ∈ [ a , b ] ) |g(x) - g(y)| \le L |x - y| \qquad (x, y \in [a,b]) ∣ g ( x ) − g ( y ) ∣ ≤ L ∣ x − y ∣ ( x , y ∈ [ a , b ]) を満たすとする。このとき g g g は [ a , b ] [a,b] [ a , b ] にただ 1 つの不動点 x ∗ x^* x ∗ をもち、任意の x 0 ∈ [ a , b ] x_0 \in [a,b] x 0 ∈ [ a , b ] から x k + 1 = g ( x k ) x_{k+1} = g(x_k) x k + 1 = g ( x k ) で作った列は x ∗ x^* x ∗ に収束して、k ≥ 1 k \ge 1 k ≥ 1 について次が成り立つ。
(事前評価)∣ x k − x ∗ ∣ ≤ L k 1 − L ∣ x 1 − x 0 ∣ |x_k - x^*| \le \dfrac{L^k}{1-L}|x_1 - x_0| ∣ x k − x ∗ ∣ ≤ 1 − L L k ∣ x 1 − x 0 ∣ 。
(事後評価)∣ x k − x ∗ ∣ ≤ L 1 − L ∣ x k − x k − 1 ∣ |x_k - x^*| \le \dfrac{L}{1-L}|x_k - x_{k-1}| ∣ x k − x ∗ ∣ ≤ 1 − L L ∣ x k − x k − 1 ∣ 。
∣ x k − x ∗ ∣ ≤ L k ∣ x 0 − x ∗ ∣ |x_k - x^*| \le L^k |x_0 - x^*| ∣ x k − x ∗ ∣ ≤ L k ∣ x 0 − x ∗ ∣ 。
g g g が [ a , b ] [a,b] [ a , b ] で微分可能で ∣ g ′ ( x ) ∣ ≤ L |g'(x)| \le L ∣ g ′ ( x ) ∣ ≤ L なら、上の Lipschitz 条件が成り立つ。
証明(定理 4.1) [ a , b ] [a,b] [ a , b ] は R \mathbb{R} R の空でない閉集合で、g g g はそれ自身への縮小率 L L L の縮小写像です。縮小写像の原理(定理 3.2[逆関数定理・陰関数定理と条件付き極値] 。n = 1 n = 1 n = 1 の場合)から、不動点がただ 1 つ存在し、列は x ∗ x^* x ∗ に収束して 1 が成り立ちます(一般の完備距離空間での形は 定理 6.2[距離空間] )。
3 は ∣ x j + 1 − x ∗ ∣ = ∣ g ( x j ) − g ( x ∗ ) ∣ ≤ L ∣ x j − x ∗ ∣ |x_{j+1} - x^*| = |g(x_j) - g(x^*)| \le L|x_j - x^*| ∣ x j + 1 − x ∗ ∣ = ∣ g ( x j ) − g ( x ∗ ) ∣ ≤ L ∣ x j − x ∗ ∣ を j = 0 , … , k − 1 j = 0, \ldots, k-1 j = 0 , … , k − 1 について繰り返せば得られます。
2 は三角不等式と 3 の 1 段分から
∣ x k − x ∗ ∣ ≤ ∣ x k − x k + 1 ∣ + ∣ x k + 1 − x ∗ ∣ = ∣ g ( x k − 1 ) − g ( x k ) ∣ + ∣ g ( x k ) − g ( x ∗ ) ∣ ≤ L ∣ x k − 1 − x k ∣ + L ∣ x k − x ∗ ∣ |x_k - x^*| \le |x_k - x_{k+1}| + |x_{k+1} - x^*| = |g(x_{k-1}) - g(x_k)| + |g(x_k) - g(x^*)| \le L|x_{k-1} - x_k| + L|x_k - x^*| ∣ x k − x ∗ ∣ ≤ ∣ x k − x k + 1 ∣ + ∣ x k + 1 − x ∗ ∣ = ∣ g ( x k − 1 ) − g ( x k ) ∣ + ∣ g ( x k ) − g ( x ∗ ) ∣ ≤ L ∣ x k − 1 − x k ∣ + L ∣ x k − x ∗ ∣ となり、∣ x k − x ∗ ∣ |x_k - x^*| ∣ x k − x ∗ ∣ について解けば得られます。
最後の主張は、x ≠ y x \ne y x = y について平均値の定理(定理 3.3[平均値の定理とテイラーの定理] )から g ( x ) − g ( y ) = g ′ ( c ) ( x − y ) g(x) - g(y) = g'(c)(x - y) g ( x ) − g ( y ) = g ′ ( c ) ( x − y ) となる c c c が x x x と y y y の間にあることから従います。
∎
事前評価は反復を始める前に必要な回数を見積もるのに、事後評価は計算しながら止めどきを決めるのに使います。事後評価の右辺は計算済みの量 ∣ x k − x k − 1 ∣ |x_k - x_{k-1}| ∣ x k − x k − 1 ∣ だけで書けているので、「差が τ ( 1 − L ) / L \tau(1-L)/L τ ( 1 − L ) / L 以下になったら止める」とすれば誤差は τ \tau τ 以下に保証されます。ただし L L L が 1 1 1 に近いと係数 L / ( 1 − L ) L/(1-L) L / ( 1 − L ) が大きくなり、差が小さくても誤差は小さいとは限りません。
例 4.2 (x = cos x の不動点反復 )
g ( x ) = cos x g(x) = \cos x g ( x ) = cos x は [ 0 , 1 ] [0,1] [ 0 , 1 ] を [ cos 1 , 1 ] ⊂ [ 0 , 1 ] [\cos 1, 1] \subset [0,1] [ cos 1 , 1 ] ⊂ [ 0 , 1 ] に写し、∣ g ′ ( x ) ∣ = ∣ sin x ∣ ≤ sin 1 ≈ 0.8415 = L |g'(x)| = |\sin x| \le \sin 1 \approx 0.8415 = L ∣ g ′ ( x ) ∣ = ∣ sin x ∣ ≤ sin 1 ≈ 0.8415 = L です。定理 4.1 から不動点 x ∗ ≈ 0.739085 x^* \approx 0.739085 x ∗ ≈ 0.739085 がただ 1 つあり、x 0 = 0 x_0 = 0 x 0 = 0 から反復すると x 1 = 1 x_1 = 1 x 1 = 1 です。事前評価で誤差 10 − 10 10^{-10} 1 0 − 10 を保証するには L k / ( 1 − L ) ≤ 10 − 10 L^k/(1-L) \le 10^{-10} L k / ( 1 − L ) ≤ 1 0 − 10 、すなわち k ≥ 145 k \ge 145 k ≥ 145 回が必要と見積もられます。実際に倍精度で反復すると、56 回目で誤差が 10 − 10 10^{-10} 1 0 − 10 を下回ります。差が出るのは、x ∗ x^* x ∗ の近くでの縮小率 ∣ g ′ ( x ∗ ) ∣ = sin x ∗ ≈ 0.674 |g'(x^*)| = \sin x^* \approx 0.674 ∣ g ′ ( x ∗ ) ∣ = sin x ∗ ≈ 0.674 が区間全体での最悪値 0.8415 0.8415 0.8415 より小さいためです。次の定理はこの「解の近くでの縮小率」を正確に述べます。
定理 4.3 (不動点反復の局所収束と次数 ) ✓
I I I を開区間、g : I → R g : I \to \mathbb{R} g : I → R を C 1 C^1 C 1 級、x ∗ ∈ I x^* \in I x ∗ ∈ I を g ( x ∗ ) = x ∗ g(x^*) = x^* g ( x ∗ ) = x ∗ を満たす点とする。
∣ g ′ ( x ∗ ) ∣ < 1 |g'(x^*)| \lt 1 ∣ g ′ ( x ∗ ) ∣ < 1 ならば δ > 0 \delta \gt 0 δ > 0 があって、x 0 ∈ [ x ∗ − δ , x ∗ + δ ] x_0 \in [x^* - \delta, x^* + \delta] x 0 ∈ [ x ∗ − δ , x ∗ + δ ] から作った x k + 1 = g ( x k ) x_{k+1} = g(x_k) x k + 1 = g ( x k ) はこの区間にとどまり、x ∗ x^* x ∗ に収束する。さらにすべての k k k で e k ≠ 0 e_k \ne 0 e k = 0 なら e k + 1 / e k → g ′ ( x ∗ ) e_{k+1}/e_k \to g'(x^*) e k + 1 / e k → g ′ ( x ∗ ) である。とくに g ′ ( x ∗ ) ≠ 0 g'(x^*) \ne 0 g ′ ( x ∗ ) = 0 なら 1 次収束で、漸近誤差定数は ∣ g ′ ( x ∗ ) ∣ |g'(x^*)| ∣ g ′ ( x ∗ ) ∣ である。
p ≥ 2 p \ge 2 p ≥ 2 とし、g g g が C p C^p C p 級で g ′ ( x ∗ ) = ⋯ = g ( p − 1 ) ( x ∗ ) = 0 g'(x^*) = \cdots = g^{(p-1)}(x^*) = 0 g ′ ( x ∗ ) = ⋯ = g ( p − 1 ) ( x ∗ ) = 0 ならば、1 の δ \delta δ について、すべての k k k で e k ≠ 0 e_k \ne 0 e k = 0 のとき
e k + 1 e k p → g ( p ) ( x ∗ ) p ! \frac{e_{k+1}}{e_k^p} \to \frac{g^{(p)}(x^*)}{p!} e k p e k + 1 → p ! g ( p ) ( x ∗ )
である。とくに g ( p ) ( x ∗ ) ≠ 0 g^{(p)}(x^*) \ne 0 g ( p ) ( x ∗ ) = 0 なら p p p 次収束で、漸近誤差定数は ∣ g ( p ) ( x ∗ ) ∣ / p ! |g^{(p)}(x^*)|/p! ∣ g ( p ) ( x ∗ ) ∣/ p ! である。
証明(定理 4.3) 1. ∣ g ′ ( x ∗ ) ∣ < L < 1 |g'(x^*)| \lt L \lt 1 ∣ g ′ ( x ∗ ) ∣ < L < 1 となる L L L を取ります。g ′ g' g ′ は x ∗ x^* x ∗ で連続なので、δ > 0 \delta \gt 0 δ > 0 を小さく取れば J = [ x ∗ − δ , x ∗ + δ ] ⊂ I J = [x^* - \delta, x^* + \delta] \subset I J = [ x ∗ − δ , x ∗ + δ ] ⊂ I かつ J J J 上で ∣ g ′ ∣ ≤ L |g'| \le L ∣ g ′ ∣ ≤ L にできます。x ∈ J x \in J x ∈ J 、x ≠ x ∗ x \ne x^* x = x ∗ なら平均値の定理(定理 3.3[平均値の定理とテイラーの定理] )から x x x と x ∗ x^* x ∗ の間の η \eta η があって
g ( x ) − x ∗ = g ( x ) − g ( x ∗ ) = g ′ ( η ) ( x − x ∗ ) g(x) - x^* = g(x) - g(x^*) = g'(\eta)(x - x^*) g ( x ) − x ∗ = g ( x ) − g ( x ∗ ) = g ′ ( η ) ( x − x ∗ ) です。η ∈ J \eta \in J η ∈ J なので ∣ g ( x ) − x ∗ ∣ ≤ L ∣ x − x ∗ ∣ ≤ δ |g(x) - x^*| \le L|x - x^*| \le \delta ∣ g ( x ) − x ∗ ∣ ≤ L ∣ x − x ∗ ∣ ≤ δ 、すなわち g ( J ) ⊂ J g(J) \subset J g ( J ) ⊂ J です。よって列は J J J にとどまり、∣ e k ∣ ≤ L k ∣ e 0 ∣ → 0 |e_k| \le L^k |e_0| \to 0 ∣ e k ∣ ≤ L k ∣ e 0 ∣ → 0 です。e k ≠ 0 e_k \ne 0 e k = 0 のとき上の式で x = x k x = x_k x = x k とすると e k + 1 / e k = g ′ ( η k ) e_{k+1}/e_k = g'(\eta_k) e k + 1 / e k = g ′ ( η k ) で、η k \eta_k η k は x k x_k x k と x ∗ x^* x ∗ の間にあるので η k → x ∗ \eta_k \to x^* η k → x ∗ です。g ′ g' g ′ の連続性から g ′ ( η k ) → g ′ ( x ∗ ) g'(\eta_k) \to g'(x^*) g ′ ( η k ) → g ′ ( x ∗ ) です。
2. g ′ ( x ∗ ) = 0 g'(x^*) = 0 g ′ ( x ∗ ) = 0 なので 1 が適用でき、列は J J J にとどまって x ∗ x^* x ∗ に収束します。e k ≠ 0 e_k \ne 0 e k = 0 のとき、Taylor の定理(定理 5.3[平均値の定理とテイラーの定理] 。a = x ∗ a = x^* a = x ∗ 、n = p − 1 n = p - 1 n = p − 1 )を x = x k x = x_k x = x k に使うと、x k x_k x k と x ∗ x^* x ∗ の間の ξ k \xi_k ξ k があって
g ( x k ) = g ( x ∗ ) + ∑ j = 1 p − 1 g ( j ) ( x ∗ ) j ! e k j + g ( p ) ( ξ k ) p ! e k p = x ∗ + g ( p ) ( ξ k ) p ! e k p g(x_k) = g(x^*) + \sum_{j=1}^{p-1} \frac{g^{(j)}(x^*)}{j!} e_k^j + \frac{g^{(p)}(\xi_k)}{p!} e_k^p = x^* + \frac{g^{(p)}(\xi_k)}{p!} e_k^p g ( x k ) = g ( x ∗ ) + j = 1 ∑ p − 1 j ! g ( j ) ( x ∗ ) e k j + p ! g ( p ) ( ξ k ) e k p = x ∗ + p ! g ( p ) ( ξ k ) e k p です。よって e k + 1 / e k p = g ( p ) ( ξ k ) / p ! e_{k+1}/e_k^p = g^{(p)}(\xi_k)/p! e k + 1 / e k p = g ( p ) ( ξ k ) / p ! で、ξ k → x ∗ \xi_k \to x^* ξ k → x ∗ と g ( p ) g^{(p)} g ( p ) の連続性から主張を得ます。
∎
この定理は「反復関数の導関数が解で消えるほど速い」ことを述べています。Newton 法はちょうど g ′ ( x ∗ ) = 0 g'(x^*) = 0 g ′ ( x ∗ ) = 0 となるように g g g を設計した反復と見ることができ、次の節の 2 次収束はその具体化です。
同じ方程式 x = cos x x = \cos x x = cos x を 2 通りに書き直した反復で、∣ g ′ ( x ∗ ) ∣ |g'(x^*)| ∣ g ′ ( x ∗ ) ∣ が 1 より小さいか大きいかの違いを見ます。左は g = cos g = \cos g = cos で、折れ線が不動点に巻き込むように近づきます。右は逆関数 g = arccos g = \arccos g = arccos で、同じ不動点から離れていきます。
x = cos x の不動点 x* ≈ 0.7391(星印)に向かう反復。左は |g'(x*)| ≈ 0.674 で収束し、右は |g'(x*)| ≈ 1.484 で、初期値 0.70 が x* のすぐ近くでも離れていく。
x k x_k x k での f f f の接線 y = f ( x k ) + f ′ ( x k ) ( x − x k ) y = f(x_k) + f'(x_k)(x - x_k) y = f ( x k ) + f ′ ( x k ) ( x − x k ) と x x x 軸の交点を次の近似とするのが Newton 法
x k + 1 = x k − f ( x k ) f ′ ( x k ) x_{k+1} = x_k - \frac{f(x_k)}{f'(x_k)} x k + 1 = x k − f ′ ( x k ) f ( x k )
です。f f f を x k x_k x k のまわりで 1 次式に置き換えて解く、と言っても同じです。置き換えの誤差は Taylor の定理の 2 次の剰余なので、誤差が 2 乗で効くことが予想できます。
x₀ = 2 x₁ = 1.5 x₂ √2 y = x² − 2 凸関数 f(x) = x² − 2 での Newton 法。x₀ = 2 での接線が x 軸と交わる点が x₁ = 1.5、x₁ での接線から x₂ ≈ 1.4167 が得られ、解 √2 ≈ 1.4142 に右から近づく。
次の定理は文字が多いですが、言っていることは 1 つです。解の近くの範囲 J J J で、m m m は ∣ f ′ ∣ |f'| ∣ f ′ ∣ の最小値、M M M は ∣ f ′ ′ ∣ |f''| ∣ f ′′ ∣ の最大値、K = M / ( 2 m ) K = M/(2m) K = M / ( 2 m ) はそこから決まる 1 つの定数です。K δ < 1 K\delta \lt 1 K δ < 1 は「範囲の幅が十分狭い」という条件で、これが成り立てば 1 回の反復で誤差が必ず小さくなり、2 乗の形で減っていきます。
定理 5.1 (Newton 法の局所 2 次収束 ) ✓
I I I を開区間、f : I → R f : I \to \mathbb{R} f : I → R を C 2 C^2 C 2 級とし、x ∗ ∈ I x^* \in I x ∗ ∈ I で f ( x ∗ ) = 0 f(x^*) = 0 f ( x ∗ ) = 0 、f ′ ( x ∗ ) ≠ 0 f'(x^*) \ne 0 f ′ ( x ∗ ) = 0 とする。
δ > 0 \delta \gt 0 δ > 0 が J = [ x ∗ − δ , x ∗ + δ ] ⊂ I J = [x^* - \delta, x^* + \delta] \subset I J = [ x ∗ − δ , x ∗ + δ ] ⊂ I と m = min J ∣ f ′ ∣ > 0 m = \min_J |f'| \gt 0 m = min J ∣ f ′ ∣ > 0 を満たすとし、M = max J ∣ f ′ ′ ∣ M = \max_J |f''| M = max J ∣ f ′′ ∣ 、K = M / ( 2 m ) K = M/(2m) K = M / ( 2 m ) とおく。K δ < 1 K\delta \lt 1 K δ < 1 ならば、任意の x 0 ∈ J x_0 \in J x 0 ∈ J について Newton 法の列は定義されて J J J にとどまり、
∣ e k + 1 ∣ ≤ K e k 2 , K ∣ e k ∣ ≤ ( K ∣ e 0 ∣ ) 2 k |e_{k+1}| \le K e_k^2, \qquad K|e_k| \le (K|e_0|)^{2^k} ∣ e k + 1 ∣ ≤ K e k 2 , K ∣ e k ∣ ≤ ( K ∣ e 0 ∣ ) 2 k
が成り立つ。とくに x k → x ∗ x_k \to x^* x k → x ∗ で、R 次数 2 以上で収束する。
1 の条件を満たす δ > 0 \delta \gt 0 δ > 0 は存在する。
1 の状況で、すべての k k k で e k ≠ 0 e_k \ne 0 e k = 0 ならば
e k + 1 e k 2 → f ′ ′ ( x ∗ ) 2 f ′ ( x ∗ ) \frac{e_{k+1}}{e_k^2} \to \frac{f''(x^*)}{2f'(x^*)} e k 2 e k + 1 → 2 f ′ ( x ∗ ) f ′′ ( x ∗ )
である。とくに f ′ ′ ( x ∗ ) ≠ 0 f''(x^*) \ne 0 f ′′ ( x ∗ ) = 0 なら 2 次収束で、漸近誤差定数は ∣ f ′ ′ ( x ∗ ) ∣ / ( 2 ∣ f ′ ( x ∗ ) ∣ ) |f''(x^*)|/(2|f'(x^*)|) ∣ f ′′ ( x ∗ ) ∣/ ( 2∣ f ′ ( x ∗ ) ∣ ) である。
証明(定理 5.1) m m m と M M M は、連続関数 ∣ f ′ ∣ |f'| ∣ f ′ ∣ 、∣ f ′ ′ ∣ |f''| ∣ f ′′ ∣ の閉区間 J J J 上の最小値・最大値として存在します(定理 2.3[平均値の定理とテイラーの定理] )。
1. x k ∈ J x_k \in J x k ∈ J と仮定します。∣ f ′ ( x k ) ∣ ≥ m > 0 |f'(x_k)| \ge m \gt 0 ∣ f ′ ( x k ) ∣ ≥ m > 0 なので x k + 1 x_{k+1} x k + 1 は定義されます。x k = x ∗ x_k = x^* x k = x ∗ なら f ( x k ) = 0 f(x_k) = 0 f ( x k ) = 0 で x k + 1 = x ∗ x_{k+1} = x^* x k + 1 = x ∗ となり、以後の主張は自明です。x k ≠ x ∗ x_k \ne x^* x k = x ∗ なら、Taylor の定理(定理 5.3[平均値の定理とテイラーの定理] 。a = x k a = x_k a = x k 、n = 1 n = 1 n = 1 )を x = x ∗ x = x^* x = x ∗ に使うと、x k x_k x k と x ∗ x^* x ∗ の間の ξ k \xi_k ξ k があって
0 = f ( x ∗ ) = f ( x k ) + f ′ ( x k ) ( x ∗ − x k ) + f ′ ′ ( ξ k ) 2 ( x ∗ − x k ) 2 0 = f(x^*) = f(x_k) + f'(x_k)(x^* - x_k) + \frac{f''(\xi_k)}{2}(x^* - x_k)^2 0 = f ( x ∗ ) = f ( x k ) + f ′ ( x k ) ( x ∗ − x k ) + 2 f ′′ ( ξ k ) ( x ∗ − x k ) 2 です。両辺を f ′ ( x k ) f'(x_k) f ′ ( x k ) で割って整理すると
e k + 1 = x k − f ( x k ) f ′ ( x k ) − x ∗ = f ′ ′ ( ξ k ) 2 f ′ ( x k ) e k 2 e_{k+1} = x_k - \frac{f(x_k)}{f'(x_k)} - x^* = \frac{f''(\xi_k)}{2f'(x_k)}\, e_k^2 e k + 1 = x k − f ′ ( x k ) f ( x k ) − x ∗ = 2 f ′ ( x k ) f ′′ ( ξ k ) e k 2 を得ます。ξ k ∈ J \xi_k \in J ξ k ∈ J なので ∣ e k + 1 ∣ ≤ K e k 2 ≤ ( K δ ) ∣ e k ∣ < ∣ e k ∣ ≤ δ |e_{k+1}| \le K e_k^2 \le (K\delta)|e_k| \lt |e_k| \le \delta ∣ e k + 1 ∣ ≤ K e k 2 ≤ ( K δ ) ∣ e k ∣ < ∣ e k ∣ ≤ δ となり、x k + 1 ∈ J x_{k+1} \in J x k + 1 ∈ J です。帰納法で列は J J J にとどまり、∣ e k + 1 ∣ ≤ K e k 2 |e_{k+1}| \le K e_k^2 ∣ e k + 1 ∣ ≤ K e k 2 がすべての k k k で成り立ちます。両辺に K K K を掛けると K ∣ e k + 1 ∣ ≤ ( K ∣ e k ∣ ) 2 K|e_{k+1}| \le (K|e_k|)^2 K ∣ e k + 1 ∣ ≤ ( K ∣ e k ∣ ) 2 で、これを繰り返すと K ∣ e k ∣ ≤ ( K ∣ e 0 ∣ ) 2 k K|e_k| \le (K|e_0|)^{2^k} K ∣ e k ∣ ≤ ( K ∣ e 0 ∣ ) 2 k です。K ∣ e 0 ∣ ≤ K δ < 1 K|e_0| \le K\delta \lt 1 K ∣ e 0 ∣ ≤ K δ < 1 なので右辺は 0 0 0 に収束します(K = 0 K = 0 K = 0 なら e 1 = 0 e_1 = 0 e 1 = 0 です)。K > 0 K \gt 0 K > 0 なら A = 1 / K A = 1/K A = 1/ K 、q = K δ q = K\delta q = K δ として R 次数 2 以上の定義(定義 2.2 )を満たします。
2. f ′ f' f ′ は x ∗ x^* x ∗ で連続で f ′ ( x ∗ ) ≠ 0 f'(x^*) \ne 0 f ′ ( x ∗ ) = 0 なので、δ 0 > 0 \delta_0 \gt 0 δ 0 > 0 を J 0 = [ x ∗ − δ 0 , x ∗ + δ 0 ] ⊂ I J_0 = [x^* - \delta_0, x^* + \delta_0] \subset I J 0 = [ x ∗ − δ 0 , x ∗ + δ 0 ] ⊂ I かつ J 0 J_0 J 0 上で ∣ f ′ ∣ ≥ ∣ f ′ ( x ∗ ) ∣ / 2 |f'| \ge |f'(x^*)|/2 ∣ f ′ ∣ ≥ ∣ f ′ ( x ∗ ) ∣/2 となるように取れます。M 0 = max J 0 ∣ f ′ ′ ∣ M_0 = \max_{J_0} |f''| M 0 = max J 0 ∣ f ′′ ∣ とします。0 < δ ≤ δ 0 0 \lt \delta \le \delta_0 0 < δ ≤ δ 0 なら m ≥ ∣ f ′ ( x ∗ ) ∣ / 2 m \ge |f'(x^*)|/2 m ≥ ∣ f ′ ( x ∗ ) ∣/2 、M ≤ M 0 M \le M_0 M ≤ M 0 なので K ≤ M 0 / ∣ f ′ ( x ∗ ) ∣ K \le M_0/|f'(x^*)| K ≤ M 0 /∣ f ′ ( x ∗ ) ∣ です。よって δ < ∣ f ′ ( x ∗ ) ∣ / M 0 \delta \lt |f'(x^*)|/M_0 δ < ∣ f ′ ( x ∗ ) ∣/ M 0 (M 0 = 0 M_0 = 0 M 0 = 0 なら任意)かつ δ ≤ δ 0 \delta \le \delta_0 δ ≤ δ 0 と取れば K δ < 1 K\delta \lt 1 K δ < 1 です。
3. 1 の式から e k + 1 / e k 2 = f ′ ′ ( ξ k ) / ( 2 f ′ ( x k ) ) e_{k+1}/e_k^2 = f''(\xi_k)/(2f'(x_k)) e k + 1 / e k 2 = f ′′ ( ξ k ) / ( 2 f ′ ( x k )) です。x k → x ∗ x_k \to x^* x k → x ∗ で、ξ k \xi_k ξ k は x k x_k x k と x ∗ x^* x ∗ の間にあるので ξ k → x ∗ \xi_k \to x^* ξ k → x ∗ です。f ′ f' f ′ 、f ′ ′ f'' f ′′ の連続性と f ′ ( x ∗ ) ≠ 0 f'(x^*) \ne 0 f ′ ( x ∗ ) = 0 から主張が従います。
∎
1 の評価は、解を知らなくても使える形にはなっていません(J J J が解を中心にしているため)。それでも「K ∣ e 0 ∣ < 1 K|e_0| \lt 1 K ∣ e 0 ∣ < 1 なら収束し、正しい桁数がほぼ倍々に増える」という見通しを与え、K K K が大きい(f ′ f' f ′ が小さい、または f ′ ′ f'' f ′′ が大きい)ほど収束域が狭いことを教えます。解を含む球を計算だけで保証する形は、後の節の Kantorovich 型の定理(定理 7.4 )で扱います。
例 5.2 (√2 の Newton 法 )
f ( x ) = x 2 − 2 f(x) = x^2 - 2 f ( x ) = x 2 − 2 では Newton 法は x k + 1 = ( x k + 2 / x k ) / 2 x_{k+1} = (x_k + 2/x_k)/2 x k + 1 = ( x k + 2/ x k ) /2 です。定理 5.1 の 3 が予言する極限は f ′ ′ ( x ∗ ) / ( 2 f ′ ( x ∗ ) ) = 2 / ( 4 2 ) = 1 / ( 2 2 ) ≈ 0.35355 f''(x^*)/(2f'(x^*)) = 2/(4\sqrt2) = 1/(2\sqrt2) \approx 0.35355 f ′′ ( x ∗ ) / ( 2 f ′ ( x ∗ )) = 2/ ( 4 2 ) = 1/ ( 2 2 ) ≈ 0.35355 です。倍精度では誤差がすぐ丸めの水準(10 − 16 10^{-16} 1 0 − 16 )で止まってしまうので、mpmath の 60 桁で反復して 2 次収束を最後まで追います。
# Python 3.12.10 / mpmath 1.3.0
from mpmath import mp, mpf, sqrt
mp.dps = 60 # 10 進 60 桁で反復する
r = sqrt ( mpf ( 2 )) # 参照値 √2
x, xd = mpf ( 2 ), 2.0 # x0 = 2(xd は倍精度で同じ反復)
print ( " k e_k(60桁) e_k/e_{k-1}^2 e_k(倍精度) " )
e, ed = abs ( x - r ), abs ( mpf ( xd ) - r )
q = f " { float ( e / ep ** 2 ) :.5f } " if ep else " - "
print ( f " {k :2d } { float ( e ) :.2e } {q} { float ( ed ) :.2e } " )
x = x - (x * x - 2 ) / ( 2 * x) # Newton 法
xd = xd - (xd * xd - 2 ) / ( 2 * xd)
print ( " 1/(2√2) = " , f " { float ( 1 / ( 2 * r) ) :.5f } " )
k e_k(60桁) e_k/e_{k-1}^2 e_k(倍精度)
1 8.58e-02 0.25000 8.58e-02
2 2.45e-03 0.33333 2.45e-03
3 2.12e-06 0.35294 2.12e-06
4 1.59e-12 0.35355 1.59e-12
5 8.99e-25 0.35355 9.67e-17
6 2.86e-49 0.35355 1.25e-16
検算. (a) この例では証明の式が閉じた形になり、e k + 1 = e k 2 / ( 2 x k ) e_{k+1} = e_k^2/(2x_k) e k + 1 = e k 2 / ( 2 x k ) です。表の比 0.25000 = 1 / ( 2 ⋅ 2 ) 0.25000 = 1/(2 \cdot 2) 0.25000 = 1/ ( 2 ⋅ 2 ) 、0.33333 = 1 / ( 2 ⋅ 1.5 ) 0.33333 = 1/(2 \cdot 1.5) 0.33333 = 1/ ( 2 ⋅ 1.5 ) はちょうど 1 / ( 2 x k − 1 ) 1/(2x_{k-1}) 1/ ( 2 x k − 1 ) に一致します。(b) 60 桁の計算と倍精度の計算は k = 4 k = 4 k = 4 まで表示の桁で一致し、倍精度は k = 5 k = 5 k = 5 以降、2 \sqrt2 2 に最も近い倍精度数の誤差 9.67 × 10 − 17 9.67 \times 10^{-17} 9.67 × 1 0 − 17 の近くで止まります(この先は丸め誤差の問題で、Newton 法の収束の問題ではありません)。(c) 誤差の指数は − 1 , − 3 , − 6 , − 12 , − 25 , − 49 -1, -3, -6, -12, -25, -49 − 1 , − 3 , − 6 , − 12 , − 25 , − 49 とほぼ倍々に増え、比 e k / e k − 1 2 e_k/e_{k-1}^2 e k / e k − 1 2 は理論値 1 / ( 2 2 ) 1/(2\sqrt2) 1/ ( 2 2 ) に近づきます。
定理 5.3 (凸関数での Newton 法の単調な大域収束 ) ✓
I I I を区間、f : I → R f : I \to \mathbb{R} f : I → R を微分可能で f ′ f' f ′ が広義単調増加(f f f が凸)な関数とし、x ∗ ∈ I x^* \in I x ∗ ∈ I で f ( x ∗ ) = 0 f(x^*) = 0 f ( x ∗ ) = 0 、f ′ ( x ∗ ) > 0 f'(x^*) \gt 0 f ′ ( x ∗ ) > 0 とする。x 0 ∈ I x_0 \in I x 0 ∈ I 、x 0 > x ∗ x_0 \gt x^* x 0 > x ∗ ならば、Newton 法の列は定義されて
x ∗ ≤ x k + 1 ≤ x k ( k ≥ 0 ) x^* \le x_{k+1} \le x_k \qquad (k \ge 0) x ∗ ≤ x k + 1 ≤ x k ( k ≥ 0 ) を満たし、x ∗ x^* x ∗ に収束する。
証明(定理 5.3) x k ∈ I x_k \in I x k ∈ I 、x k > x ∗ x_k \gt x^* x k > x ∗ と仮定します。f ′ f' f ′ は増加なので f ′ ( x k ) ≥ f ′ ( x ∗ ) > 0 f'(x_k) \ge f'(x^*) \gt 0 f ′ ( x k ) ≥ f ′ ( x ∗ ) > 0 で、x k + 1 x_{k+1} x k + 1 は定義されます。平均値の定理(定理 3.3[平均値の定理とテイラーの定理] )から f ( x k ) = f ( x k ) − f ( x ∗ ) = f ′ ( c ) ( x k − x ∗ ) f(x_k) = f(x_k) - f(x^*) = f'(c)(x_k - x^*) f ( x k ) = f ( x k ) − f ( x ∗ ) = f ′ ( c ) ( x k − x ∗ ) となる c ∈ ( x ∗ , x k ) c \in (x^*, x_k) c ∈ ( x ∗ , x k ) があり、0 < f ′ ( x ∗ ) ≤ f ′ ( c ) ≤ f ′ ( x k ) 0 \lt f'(x^*) \le f'(c) \le f'(x_k) 0 < f ′ ( x ∗ ) ≤ f ′ ( c ) ≤ f ′ ( x k ) です。よって
x k − x k + 1 = f ( x k ) f ′ ( x k ) = f ′ ( c ) f ′ ( x k ) ( x k − x ∗ ) ∈ ( 0 , x k − x ∗ ] x_k - x_{k+1} = \frac{f(x_k)}{f'(x_k)} = \frac{f'(c)}{f'(x_k)}(x_k - x^*) \in (0,\, x_k - x^*] x k − x k + 1 = f ′ ( x k ) f ( x k ) = f ′ ( x k ) f ′ ( c ) ( x k − x ∗ ) ∈ ( 0 , x k − x ∗ ] で、x ∗ ≤ x k + 1 < x k x^* \le x_{k+1} \lt x_k x ∗ ≤ x k + 1 < x k です。x k + 1 = x ∗ x_{k+1} = x^* x k + 1 = x ∗ なら以後 x j = x ∗ x_j = x^* x j = x ∗ で主張は自明です。そうでなければ x k + 1 ∈ ( x ∗ , x 0 ] ⊂ I x_{k+1} \in (x^*, x_0] \subset I x k + 1 ∈ ( x ∗ , x 0 ] ⊂ I で同じ議論を繰り返せます。
列 ( − x k ) (-x_k) ( − x k ) は単調増加で − x ∗ -x^* − x ∗ で上に有界なので、単調有界定理(定理 2.1[級数と収束判定] )から x k x_k x k はある ℓ ≥ x ∗ \ell \ge x^* ℓ ≥ x ∗ に収束します。f ′ ( x k ) f'(x_k) f ′ ( x k ) は [ f ′ ( x ∗ ) , f ′ ( x 0 ) ] [f'(x^*), f'(x_0)] [ f ′ ( x ∗ ) , f ′ ( x 0 )] に入るので、f ( x k ) = ( x k − x k + 1 ) f ′ ( x k ) f(x_k) = (x_k - x_{k+1}) f'(x_k) f ( x k ) = ( x k − x k + 1 ) f ′ ( x k ) から ∣ f ( x k ) ∣ ≤ f ′ ( x 0 ) ∣ x k − x k + 1 ∣ → 0 |f(x_k)| \le f'(x_0)|x_k - x_{k+1}| \to 0 ∣ f ( x k ) ∣ ≤ f ′ ( x 0 ) ∣ x k − x k + 1 ∣ → 0 です。f f f は連続なので(定理 4.5[極限と連続性] )f ( ℓ ) = 0 f(\ell) = 0 f ( ℓ ) = 0 です。一方 x > x ∗ x \gt x^* x > x ∗ なら上と同じく平均値の定理から f ( x ) = f ′ ( c ) ( x − x ∗ ) > 0 f(x) = f'(c)(x - x^*) \gt 0 f ( x ) = f ′ ( c ) ( x − x ∗ ) > 0 なので、ℓ = x ∗ \ell = x^* ℓ = x ∗ です。
∎
f ( x ) = x 2 − 2 f(x) = x^2 - 2 f ( x ) = x 2 − 2 は凸で、x 0 > 2 x_0 \gt \sqrt2 x 0 > 2 ならどこから始めても単調に減少して 2 \sqrt2 2 に収束します。0 < x 0 < 2 0 \lt x_0 \lt \sqrt2 0 < x 0 < 2 から始めても、1 回目で x 1 = ( x 0 + 2 / x 0 ) / 2 ≥ 2 x_1 = (x_0 + 2/x_0)/2 \ge \sqrt2 x 1 = ( x 0 + 2/ x 0 ) /2 ≥ 2 (相加相乗平均)となるので、2 回目以降はこの定理の状況に入ります。
f ′ ( x ∗ ) = 0 f'(x^*) = 0 f ′ ( x ∗ ) = 0 のとき、定理 5.1 は使えません。重根では Newton 法は収束はしますが、1 次収束に落ちます。
次の命題の h h h は、解 x ∗ x^* x ∗ を m m m 回掛けて取り出した残りの部分です。h ( x ∗ ) ≠ 0 h(x^*) \ne 0 h ( x ∗ ) = 0 なので、解の近くでは f f f は ( x − x ∗ ) m (x - x^*)^m ( x − x ∗ ) m とほぼ同じ形に振る舞います。たとえば ( x − 1 ) 2 (x-1)^2 ( x − 1 ) 2 なら m = 2 m = 2 m = 2 、h ≡ 1 h \equiv 1 h ≡ 1 です。
命題 5.4 (重根での Newton 法と修正 Newton 法 ) ✓
I I I を開区間、x ∗ ∈ I x^* \in I x ∗ ∈ I 、m ≥ 2 m \ge 2 m ≥ 2 を整数とし、f ( x ) = ( x − x ∗ ) m h ( x ) f(x) = (x - x^*)^m h(x) f ( x ) = ( x − x ∗ ) m h ( x ) (h : I → R h : I \to \mathbb{R} h : I → R は C 1 C^1 C 1 級、h ( x ∗ ) ≠ 0 h(x^*) \ne 0 h ( x ∗ ) = 0 )とする。
δ > 0 \delta \gt 0 δ > 0 があって、x 0 ∈ [ x ∗ − δ , x ∗ + δ ] x_0 \in [x^* - \delta, x^* + \delta] x 0 ∈ [ x ∗ − δ , x ∗ + δ ] から Newton 法の列は定義されて x ∗ x^* x ∗ に収束し、すべての k k k で e k ≠ 0 e_k \ne 0 e k = 0 なら e k + 1 / e k → ( m − 1 ) / m e_{k+1}/e_k \to (m-1)/m e k + 1 / e k → ( m − 1 ) / m である。すなわち 1 次収束で、漸近誤差定数は ( m − 1 ) / m (m-1)/m ( m − 1 ) / m である。
重複度 m m m を掛けた修正 Newton 法 x k + 1 = x k − m f ( x k ) / f ′ ( x k ) x_{k+1} = x_k - m f(x_k)/f'(x_k) x k + 1 = x k − m f ( x k ) / f ′ ( x k ) では、δ \delta δ を小さく取り直せば列は x ∗ x^* x ∗ に収束し、すべての k k k で e k ≠ 0 e_k \ne 0 e k = 0 なら
e k + 1 e k 2 → h ′ ( x ∗ ) m h ( x ∗ ) \frac{e_{k+1}}{e_k^2} \to \frac{h'(x^*)}{m\, h(x^*)} e k 2 e k + 1 → m h ( x ∗ ) h ′ ( x ∗ )
である。
証明(命題 5.4) 積の微分から f ′ ( x ) = ( x − x ∗ ) m − 1 D ( x ) f'(x) = (x - x^*)^{m-1} D(x) f ′ ( x ) = ( x − x ∗ ) m − 1 D ( x ) 、D ( x ) = m h ( x ) + ( x − x ∗ ) h ′ ( x ) D(x) = m h(x) + (x - x^*) h'(x) D ( x ) = mh ( x ) + ( x − x ∗ ) h ′ ( x ) です。D D D は連続で D ( x ∗ ) = m h ( x ∗ ) ≠ 0 D(x^*) = m h(x^*) \ne 0 D ( x ∗ ) = mh ( x ∗ ) = 0 なので、δ 0 > 0 \delta_0 \gt 0 δ 0 > 0 を J 0 = [ x ∗ − δ 0 , x ∗ + δ 0 ] ⊂ I J_0 = [x^* - \delta_0, x^* + \delta_0] \subset I J 0 = [ x ∗ − δ 0 , x ∗ + δ 0 ] ⊂ I 上で D D D も h h h も 0 0 0 にならないように取れます。このとき x ∈ J 0 x \in J_0 x ∈ J 0 、x ≠ x ∗ x \ne x^* x = x ∗ で f ′ ( x ) ≠ 0 f'(x) \ne 0 f ′ ( x ) = 0 であり、
f ( x ) f ′ ( x ) = ( x − x ∗ ) φ ( x ) , φ ( x ) = h ( x ) D ( x ) \frac{f(x)}{f'(x)} = (x - x^*)\,\varphi(x), \qquad \varphi(x) = \frac{h(x)}{D(x)} f ′ ( x ) f ( x ) = ( x − x ∗ ) φ ( x ) , φ ( x ) = D ( x ) h ( x ) と書けます。φ \varphi φ は J 0 J_0 J 0 上連続で φ ( x ∗ ) = 1 / m \varphi(x^*) = 1/m φ ( x ∗ ) = 1/ m です。
1. Newton 法では e k + 1 = e k − e k φ ( x k ) = e k ( 1 − φ ( x k ) ) e_{k+1} = e_k - e_k\varphi(x_k) = e_k(1 - \varphi(x_k)) e k + 1 = e k − e k φ ( x k ) = e k ( 1 − φ ( x k )) です(e k = 0 e_k = 0 e k = 0 なら止める)。1 − φ 1 - \varphi 1 − φ は連続で 1 − φ ( x ∗ ) = ( m − 1 ) / m < 1 1 - \varphi(x^*) = (m-1)/m \lt 1 1 − φ ( x ∗ ) = ( m − 1 ) / m < 1 なので、( m − 1 ) / m < L < 1 (m-1)/m \lt L \lt 1 ( m − 1 ) / m < L < 1 を取り、δ ≤ δ 0 \delta \le \delta_0 δ ≤ δ 0 を [ x ∗ − δ , x ∗ + δ ] [x^* - \delta, x^* + \delta] [ x ∗ − δ , x ∗ + δ ] 上で ∣ 1 − φ ∣ ≤ L |1 - \varphi| \le L ∣1 − φ ∣ ≤ L となるように取れます。すると ∣ e k + 1 ∣ ≤ L ∣ e k ∣ |e_{k+1}| \le L|e_k| ∣ e k + 1 ∣ ≤ L ∣ e k ∣ で列はこの区間にとどまり、e k → 0 e_k \to 0 e k → 0 です。さらに e k + 1 / e k = 1 − φ ( x k ) → 1 − 1 / m e_{k+1}/e_k = 1 - \varphi(x_k) \to 1 - 1/m e k + 1 / e k = 1 − φ ( x k ) → 1 − 1/ m です。
2. 修正 Newton 法では e k + 1 = e k ( 1 − m φ ( x k ) ) e_{k+1} = e_k(1 - m\varphi(x_k)) e k + 1 = e k ( 1 − m φ ( x k )) で、
1 − m φ ( x ) = D ( x ) − m h ( x ) D ( x ) = ( x − x ∗ ) h ′ ( x ) D ( x ) 1 - m\varphi(x) = \frac{D(x) - m h(x)}{D(x)} = \frac{(x - x^*)\, h'(x)}{D(x)} 1 − m φ ( x ) = D ( x ) D ( x ) − mh ( x ) = D ( x ) ( x − x ∗ ) h ′ ( x ) です。よって e k + 1 = e k 2 h ′ ( x k ) / D ( x k ) e_{k+1} = e_k^2\, h'(x_k)/D(x_k) e k + 1 = e k 2 h ′ ( x k ) / D ( x k ) です。ψ = h ′ / D \psi = h'/D ψ = h ′ / D は J 0 J_0 J 0 上連続なので ∣ ψ ∣ ≤ N |\psi| \le N ∣ ψ ∣ ≤ N と有界です(定理 2.3[平均値の定理とテイラーの定理] )。δ ≤ δ 0 \delta \le \delta_0 δ ≤ δ 0 を N δ < 1 N\delta \lt 1 N δ < 1 となるように取れば、∣ e k + 1 ∣ ≤ N δ ∣ e k ∣ |e_{k+1}| \le N\delta|e_k| ∣ e k + 1 ∣ ≤ N δ ∣ e k ∣ で列は区間にとどまり e k → 0 e_k \to 0 e k → 0 です。最後に e k + 1 / e k 2 = ψ ( x k ) → h ′ ( x ∗ ) / D ( x ∗ ) = h ′ ( x ∗ ) / ( m h ( x ∗ ) ) e_{k+1}/e_k^2 = \psi(x_k) \to h'(x^*)/D(x^*) = h'(x^*)/(m h(x^*)) e k + 1 / e k 2 = ψ ( x k ) → h ′ ( x ∗ ) / D ( x ∗ ) = h ′ ( x ∗ ) / ( mh ( x ∗ )) です。
∎
例 5.5 ((x − 1)² での遅い収束 )
f ( x ) = ( x − 1 ) 2 f(x) = (x-1)^2 f ( x ) = ( x − 1 ) 2 では Newton 法は x k + 1 = x k − ( x k − 1 ) / 2 x_{k+1} = x_k - (x_k - 1)/2 x k + 1 = x k − ( x k − 1 ) /2 、すなわち e k + 1 = e k / 2 e_{k+1} = e_k/2 e k + 1 = e k /2 です。誤差は 1 回で半分にしかならず、二分法と同じ速さです。命題 5.4 の比 ( m − 1 ) / m = 1 / 2 (m-1)/m = 1/2 ( m − 1 ) / m = 1/2 がちょうど実現しています。修正 Newton 法 x k + 1 = x k − 2 f ( x k ) / f ′ ( x k ) x_{k+1} = x_k - 2f(x_k)/f'(x_k) x k + 1 = x k − 2 f ( x k ) / f ′ ( x k ) では 1 回で x 1 = 1 x_1 = 1 x 1 = 1 に到達します(h ≡ 1 h \equiv 1 h ≡ 1 で h ′ = 0 h' = 0 h ′ = 0 )。実際の問題では重複度 m m m は未知のことが多く、e k + 1 / e k e_{k+1}/e_k e k + 1 / e k の観測値から m m m を推定するか、f / f ′ f/f' f / f ′ (単根しかもたない)に Newton 法を使う工夫をします。
重根は収束が遅いだけでなく、達成できる精度そのものが低い根でもあります。f f f の値の計算に含まれる誤差を f f f に加わる小さな摂動 ε q \varepsilon q εq とみなすと、多項式の単根は 1 次で − ε q ( x ∗ ) / f ′ ( x ∗ ) -\varepsilon q(x^*)/f'(x^*) − εq ( x ∗ ) / f ′ ( x ∗ ) だけ動きます(命題 3.1[条件数と安定性] )。∣ f ′ ( x ∗ ) ∣ |f'(x^*)| ∣ f ′ ( x ∗ ) ∣ が小さい根、つまり重根に近い根ほど、反復法をどう選んでも計算できる桁数は少なくなります。これは算法ではなく問題の性質(条件数)で、条件数と安定性 の章で扱います。
例 5.6 (周期軌道に入る Newton 法 )
f ( x ) = x 3 − 2 x + 2 f(x) = x^3 - 2x + 2 f ( x ) = x 3 − 2 x + 2 に x 0 = 0 x_0 = 0 x 0 = 0 から Newton 法を使います。f ( 0 ) = 2 f(0) = 2 f ( 0 ) = 2 、f ′ ( 0 ) = − 2 f'(0) = -2 f ′ ( 0 ) = − 2 なので x 1 = 1 x_1 = 1 x 1 = 1 、f ( 1 ) = 1 f(1) = 1 f ( 1 ) = 1 、f ′ ( 1 ) = 1 f'(1) = 1 f ′ ( 1 ) = 1 なので x 2 = 0 x_2 = 0 x 2 = 0 です。列は 0 , 1 , 0 , 1 , … 0, 1, 0, 1, \ldots 0 , 1 , 0 , 1 , … と周期 2 で回り続け、実数解 x ∗ ≈ − 1.769 x^* \approx -1.769 x ∗ ≈ − 1.769 には近づきません。f ′ f' f ′ は ± 2 / 3 \pm\sqrt{2/3} ± 2/3 で 0 0 0 になり、その間では接線が大きく外れた点に飛ばします。局所収束の定理の「解の十分近くから始めれば」という仮定は飾りではなく、実用では二分法で区間を絞ってから Newton 法に切り替える、減速(x k + 1 = x k − λ k f ( x k ) / f ′ ( x k ) x_{k+1} = x_k - \lambda_k f(x_k)/f'(x_k) x k + 1 = x k − λ k f ( x k ) / f ′ ( x k ) 、0 < λ k ≤ 1 0 \lt \lambda_k \le 1 0 < λ k ≤ 1 )を入れる、などの対策をとります。
Newton 法は毎回 f ′ ( x k ) f'(x_k) f ′ ( x k ) を必要とします。導関数の式がない、あるいは計算が高くつくときは、接線の傾きを直前の 2 点を結ぶ直線の傾きで置き換えます。これが割線法
x k + 1 = x k − f ( x k ) x k − x k − 1 f ( x k ) − f ( x k − 1 ) x_{k+1} = x_k - f(x_k)\,\frac{x_k - x_{k-1}}{f(x_k) - f(x_{k-1})} x k + 1 = x k − f ( x k ) f ( x k ) − f ( x k − 1 ) x k − x k − 1
です。1 回の反復で新しく計算する関数値は f ( x k ) f(x_k) f ( x k ) の 1 つだけです。解析には差分商 を使います。f [ a , b ] f[a,b] f [ a , b ] は 2 点 a , b a, b a , b を結ぶ線分の傾き、f [ a , b , c ] f[a,b,c] f [ a , b , c ] はその傾きの変化の割合で、導関数と 2 階導関数の「差分版」です。相異なる点 a , b , c a, b, c a , b , c について
f [ a , b ] = f ( b ) − f ( a ) b − a , f [ a , b , c ] = f [ a , b ] − f [ b , c ] a − c f[a,b] = \frac{f(b) - f(a)}{b - a}, \qquad f[a,b,c] = \frac{f[a,b] - f[b,c]}{a - c} f [ a , b ] = b − a f ( b ) − f ( a ) , f [ a , b , c ] = a − c f [ a , b ] − f [ b , c ]
と書きます。定義から f [ a , b ] = f [ b , a ] f[a,b] = f[b,a] f [ a , b ] = f [ b , a ] です。割線法は x k + 1 = x k − f ( x k ) / f [ x k − 1 , x k ] x_{k+1} = x_k - f(x_k)/f[x_{k-1}, x_k] x k + 1 = x k − f ( x k ) / f [ x k − 1 , x k ] と書けます。
次の補題のねらいは、割線法の誤差の式に出てくる差分商を、導関数の値に置き換えることです。1 つ目は平均値の定理の言い換えで、2 つ目は 3 点を通る 2 次式を作り、Rolle の定理を 2 回使って示します。
補題 6.1 (2 階の差分商の平均値定理 ) ✓
I I I を開区間、f : I → R f : I \to \mathbb{R} f : I → R を 2 回微分可能とし、a , b , c ∈ I a, b, c \in I a , b , c ∈ I を相異なる点とする。このとき a a a と b b b の間の η \eta η と、a , b , c a, b, c a , b , c の最小値と最大値の間の ξ \xi ξ があって
f [ a , b ] = f ′ ( η ) , f [ a , b , c ] = f ′ ′ ( ξ ) 2 f[a,b] = f'(\eta), \qquad f[a,b,c] = \frac{f''(\xi)}{2} f [ a , b ] = f ′ ( η ) , f [ a , b , c ] = 2 f ′′ ( ξ ) が成り立つ。
証明(補題 6.1) 1 つ目は平均値の定理(定理 3.3[平均値の定理とテイラーの定理] )そのものです。2 つ目のために 2 次式
p ( x ) = f ( a ) + f [ a , b ] ( x − a ) + f [ a , b , c ] ( x − a ) ( x − b ) p(x) = f(a) + f[a,b](x - a) + f[a,b,c](x - a)(x - b) p ( x ) = f ( a ) + f [ a , b ] ( x − a ) + f [ a , b , c ] ( x − a ) ( x − b ) を考えます。p ( a ) = f ( a ) p(a) = f(a) p ( a ) = f ( a ) 、p ( b ) = f ( a ) + f [ a , b ] ( b − a ) = f ( b ) p(b) = f(a) + f[a,b](b-a) = f(b) p ( b ) = f ( a ) + f [ a , b ] ( b − a ) = f ( b ) です。また定義から f [ a , b , c ] ( c − a ) = f [ b , c ] − f [ a , b ] f[a,b,c](c - a) = f[b,c] - f[a,b] f [ a , b , c ] ( c − a ) = f [ b , c ] − f [ a , b ] なので
f [ a , b , c ] ( c − a ) ( c − b ) = ( f [ b , c ] − f [ a , b ] ) ( c − b ) = f ( c ) − f ( b ) − f [ a , b ] ( c − b ) f[a,b,c](c-a)(c-b) = \bigl(f[b,c] - f[a,b]\bigr)(c - b) = f(c) - f(b) - f[a,b](c - b) f [ a , b , c ] ( c − a ) ( c − b ) = ( f [ b , c ] − f [ a , b ] ) ( c − b ) = f ( c ) − f ( b ) − f [ a , b ] ( c − b ) となり、p ( c ) = f ( a ) + f [ a , b ] ( c − a ) + f ( c ) − f ( b ) − f [ a , b ] ( c − b ) = f ( a ) + f [ a , b ] ( b − a ) + f ( c ) − f ( b ) = f ( c ) p(c) = f(a) + f[a,b](c - a) + f(c) - f(b) - f[a,b](c-b) = f(a) + f[a,b](b - a) + f(c) - f(b) = f(c) p ( c ) = f ( a ) + f [ a , b ] ( c − a ) + f ( c ) − f ( b ) − f [ a , b ] ( c − b ) = f ( a ) + f [ a , b ] ( b − a ) + f ( c ) − f ( b ) = f ( c ) です。
g = f − p g = f - p g = f − p は I I I 上 2 回微分可能で、3 点 a , b , c a, b, c a , b , c で 0 0 0 になります。3 点を小さい順に t 1 < t 2 < t 3 t_1 \lt t_2 \lt t_3 t 1 < t 2 < t 3 とすると、Rolle の定理(定理 3.1[平均値の定理とテイラーの定理] )を [ t 1 , t 2 ] [t_1, t_2] [ t 1 , t 2 ] と [ t 2 , t 3 ] [t_2, t_3] [ t 2 , t 3 ] に使って g ′ ( s 1 ) = g ′ ( s 2 ) = 0 g'(s_1) = g'(s_2) = 0 g ′ ( s 1 ) = g ′ ( s 2 ) = 0 となる t 1 < s 1 < t 2 < s 2 < t 3 t_1 \lt s_1 \lt t_2 \lt s_2 \lt t_3 t 1 < s 1 < t 2 < s 2 < t 3 を得ます。g ′ g' g ′ は微分可能なので連続で、もう一度 Rolle の定理を [ s 1 , s 2 ] [s_1, s_2] [ s 1 , s 2 ] に使うと g ′ ′ ( ξ ) = 0 g''(\xi) = 0 g ′′ ( ξ ) = 0 となる ξ ∈ ( s 1 , s 2 ) ⊂ ( t 1 , t 3 ) \xi \in (s_1, s_2) \subset (t_1, t_3) ξ ∈ ( s 1 , s 2 ) ⊂ ( t 1 , t 3 ) があります。p ′ ′ = 2 f [ a , b , c ] p'' = 2f[a,b,c] p ′′ = 2 f [ a , b , c ] なので f ′ ′ ( ξ ) = 2 f [ a , b , c ] f''(\xi) = 2f[a,b,c] f ′′ ( ξ ) = 2 f [ a , b , c ] です。
∎
一般の n n n 階の差分商についての同じ形の定理は、多項式補間 の章で補間多項式を使って示します。ここでは割線法に必要な 2 階の場合だけを使います。
定理 6.2 (割線法の局所収束 ) ✓
I I I を開区間、f : I → R f : I \to \mathbb{R} f : I → R を C 2 C^2 C 2 級とし、x ∗ ∈ I x^* \in I x ∗ ∈ I で f ( x ∗ ) = 0 f(x^*) = 0 f ( x ∗ ) = 0 、f ′ ( x ∗ ) ≠ 0 f'(x^*) \ne 0 f ′ ( x ∗ ) = 0 とする。δ > 0 \delta \gt 0 δ > 0 が J = [ x ∗ − δ , x ∗ + δ ] ⊂ I J = [x^* - \delta, x^* + \delta] \subset I J = [ x ∗ − δ , x ∗ + δ ] ⊂ I と m = min J ∣ f ′ ∣ > 0 m = \min_J |f'| \gt 0 m = min J ∣ f ′ ∣ > 0 を満たし、M = max J ∣ f ′ ′ ∣ M = \max_J |f''| M = max J ∣ f ′′ ∣ に対して K ≥ M / ( 2 m ) K \ge M/(2m) K ≥ M / ( 2 m ) かつ K δ < 1 K\delta \lt 1 K δ < 1 とする(このような δ \delta δ と K K K は存在する)。相異なる x 0 , x 1 ∈ J x_0, x_1 \in J x 0 , x 1 ∈ J から割線法を始め、f ( x k ) = 0 f(x_k) = 0 f ( x k ) = 0 となったら止める。このとき
列は定義されて J J J にとどまり、∣ e k + 1 ∣ ≤ K ∣ e k ∣ ∣ e k − 1 ∣ |e_{k+1}| \le K|e_k||e_{k-1}| ∣ e k + 1 ∣ ≤ K ∣ e k ∣∣ e k − 1 ∣ を満たして x ∗ x^* x ∗ に収束する。
止まらないならば、
e k + 1 e k e k − 1 → f ′ ′ ( x ∗ ) 2 f ′ ( x ∗ ) \frac{e_{k+1}}{e_k\, e_{k-1}} \to \frac{f''(x^*)}{2f'(x^*)} e k e k − 1 e k + 1 → 2 f ′ ( x ∗ ) f ′′ ( x ∗ )
である。
証明(定理 6.2) δ \delta δ と K K K の存在は、定理 5.1 の 2 と同じ議論で δ \delta δ を取り、K = M / ( 2 m ) K = M/(2m) K = M / ( 2 m ) とすれば得られます。
1. x k − 1 , x k ∈ J x_{k-1}, x_k \in J x k − 1 , x k ∈ J が相異なり、f ( x k − 1 ) ≠ 0 f(x_{k-1}) \ne 0 f ( x k − 1 ) = 0 、f ( x k ) ≠ 0 f(x_k) \ne 0 f ( x k ) = 0 と仮定します(k = 1 k = 1 k = 1 では仮定どおりで、止まっていなければ以下で帰納的に保たれます)。補題 6.1 から f [ x k − 1 , x k ] = f ′ ( η ) f[x_{k-1}, x_k] = f'(\eta) f [ x k − 1 , x k ] = f ′ ( η ) (η ∈ J \eta \in J η ∈ J )で、∣ f ′ ( η ) ∣ ≥ m > 0 |f'(\eta)| \ge m \gt 0 ∣ f ′ ( η ) ∣ ≥ m > 0 なので x k + 1 x_{k+1} x k + 1 は定義されます。f ( x k ) = f ( x k ) − f ( x ∗ ) = f [ x k , x ∗ ] e k f(x_k) = f(x_k) - f(x^*) = f[x_k, x^*]\, e_k f ( x k ) = f ( x k ) − f ( x ∗ ) = f [ x k , x ∗ ] e k なので
e k + 1 = e k − f [ x k , x ∗ ] f [ x k − 1 , x k ] e k = e k f [ x k − 1 , x k ] − f [ x k , x ∗ ] f [ x k − 1 , x k ] = e k e k − 1 f [ x k − 1 , x k , x ∗ ] f [ x k − 1 , x k ] e_{k+1} = e_k - \frac{f[x_k, x^*]}{f[x_{k-1}, x_k]}\, e_k = e_k\,\frac{f[x_{k-1}, x_k] - f[x_k, x^*]}{f[x_{k-1}, x_k]} = e_k\, e_{k-1}\,\frac{f[x_{k-1}, x_k, x^*]}{f[x_{k-1}, x_k]} e k + 1 = e k − f [ x k − 1 , x k ] f [ x k , x ∗ ] e k = e k f [ x k − 1 , x k ] f [ x k − 1 , x k ] − f [ x k , x ∗ ] = e k e k − 1 f [ x k − 1 , x k ] f [ x k − 1 , x k , x ∗ ] です。最後の等号は 2 階の差分商の定義そのものです(3 点 x k − 1 , x k , x ∗ x_{k-1}, x_k, x^* x k − 1 , x k , x ∗ は相異なります。f ( x k − 1 ) ≠ 0 f(x_{k-1}) \ne 0 f ( x k − 1 ) = 0 、f ( x k ) ≠ 0 f(x_k) \ne 0 f ( x k ) = 0 だからです)。補題 6.1 から分子は f ′ ′ ( ξ k ) / 2 f''(\xi_k)/2 f ′′ ( ξ k ) /2 (ξ k ∈ J \xi_k \in J ξ k ∈ J )なので、∣ e k + 1 ∣ ≤ M 2 m ∣ e k ∣ ∣ e k − 1 ∣ ≤ K ∣ e k ∣ ∣ e k − 1 ∣ |e_{k+1}| \le \frac{M}{2m}|e_k||e_{k-1}| \le K|e_k||e_{k-1}| ∣ e k + 1 ∣ ≤ 2 m M ∣ e k ∣∣ e k − 1 ∣ ≤ K ∣ e k ∣∣ e k − 1 ∣ です。∣ e k − 1 ∣ ≤ δ |e_{k-1}| \le \delta ∣ e k − 1 ∣ ≤ δ なので ∣ e k + 1 ∣ ≤ K δ ∣ e k ∣ < ∣ e k ∣ |e_{k+1}| \le K\delta|e_k| \lt |e_k| ∣ e k + 1 ∣ ≤ K δ ∣ e k ∣ < ∣ e k ∣ で、x k + 1 ∈ J x_{k+1} \in J x k + 1 ∈ J です。さらに x k + 1 = x k x_{k+1} = x_k x k + 1 = x k となるのは f ( x k ) = 0 f(x_k) = 0 f ( x k ) = 0 のときだけなので、x k ≠ x k + 1 x_k \ne x_{k+1} x k = x k + 1 です。帰納法ですべての段で成り立ち、∣ e k + 1 ∣ ≤ ( K δ ) k ∣ e 1 ∣ → 0 |e_{k+1}| \le (K\delta)^k |e_1| \to 0 ∣ e k + 1 ∣ ≤ ( K δ ) k ∣ e 1 ∣ → 0 です。
2. 1 の式から e k + 1 / ( e k e k − 1 ) = f ′ ′ ( ξ k ) / ( 2 f ′ ( η k ) ) e_{k+1}/(e_k e_{k-1}) = f''(\xi_k)/(2f'(\eta_k)) e k + 1 / ( e k e k − 1 ) = f ′′ ( ξ k ) / ( 2 f ′ ( η k )) です。η k \eta_k η k と ξ k \xi_k ξ k はいずれも x k − 1 , x k , x ∗ x_{k-1}, x_k, x^* x k − 1 , x k , x ∗ の最小値と最大値の間にあるので x ∗ x^* x ∗ に収束し、f ′ f' f ′ 、f ′ ′ f'' f ′′ の連続性から主張が従います。
∎
系 6.3 (割線法の R 次数は黄金比以上 ) ✓
定理 6.2 の状況で割線法が止まらないならば、φ = ( 1 + 5 ) / 2 \varphi = (1 + \sqrt5)/2 φ = ( 1 + 5 ) /2 について、列は R 次数 φ \varphi φ 以上で x ∗ x^* x ∗ に収束する。
証明の方針は、d k + 1 ≤ d k d k − 1 d_{k+1} \le d_k d_{k-1} d k + 1 ≤ d k d k − 1 という「前の 2 回の積」の形の不等式から、指数がどう増えるかを調べることです。指数が φ k \varphi^k φ k のとき、φ k + φ k − 1 = φ k + 1 \varphi^k + \varphi^{k-1} = \varphi^{k+1} φ k + φ k − 1 = φ k + 1 となって指数の足し算がぴったり合います。これが黄金比が出てくる理由で、q q q の取り方は最初の 2 項でこの形を満たすためのものです。
証明(系 6.3) K ~ = max ( K , 1 / ( 2 δ ) ) \tilde K = \max(K, 1/(2\delta)) K ~ = max ( K , 1/ ( 2 δ )) とおくと K ~ > 0 \tilde K \gt 0 K ~ > 0 、K ~ δ < 1 \tilde K\delta \lt 1 K ~ δ < 1 で、定理 6.2 の 1 の不等式は K K K を K ~ \tilde K K ~ に替えても成り立ちます。d k = K ~ ∣ e k ∣ d_k = \tilde K|e_k| d k = K ~ ∣ e k ∣ とおくと d k + 1 ≤ d k d k − 1 d_{k+1} \le d_k d_{k-1} d k + 1 ≤ d k d k − 1 、0 < d k ≤ K ~ δ < 1 0 \lt d_k \le \tilde K\delta \lt 1 0 < d k ≤ K ~ δ < 1 です。q = max ( d 0 , d 1 1 / φ ) q = \max(d_0, d_1^{1/\varphi}) q = max ( d 0 , d 1 1/ φ ) とすると 0 < q < 1 0 \lt q \lt 1 0 < q < 1 で、d 0 ≤ q = q φ 0 d_0 \le q = q^{\varphi^0} d 0 ≤ q = q φ 0 、d 1 ≤ q φ d_1 \le q^{\varphi} d 1 ≤ q φ です。d k − 1 ≤ q φ k − 1 d_{k-1} \le q^{\varphi^{k-1}} d k − 1 ≤ q φ k − 1 、d k ≤ q φ k d_k \le q^{\varphi^k} d k ≤ q φ k なら、φ 2 = φ + 1 \varphi^2 = \varphi + 1 φ 2 = φ + 1 を使って
d k + 1 ≤ q φ k + φ k − 1 = q φ k − 1 ( φ + 1 ) = q φ k + 1 d_{k+1} \le q^{\varphi^k + \varphi^{k-1}} = q^{\varphi^{k-1}(\varphi + 1)} = q^{\varphi^{k+1}} d k + 1 ≤ q φ k + φ k − 1 = q φ k − 1 ( φ + 1 ) = q φ k + 1 です。帰納法で ∣ e k ∣ ≤ K ~ − 1 q φ k |e_k| \le \tilde K^{-1} q^{\varphi^k} ∣ e k ∣ ≤ K ~ − 1 q φ k がすべての k k k で成り立ち、R 次数 φ \varphi φ 以上の定義(定義 2.2 )を満たします。
∎
1 回あたりの関数値の計算が Newton 法では f f f と f ′ f' f ′ の 2 回、割線法では f f f の 1 回です。f ′ f' f ′ の計算が f f f と同じ手間なら、割線法の 2 回分(次数 φ 2 ≈ 2.618 \varphi^2 \approx 2.618 φ 2 ≈ 2.618 )は Newton 法の 1 回分(次数 2)より速いことになります。
例 6.5 (割線法の次数を誤差から推定する )
f ( x ) = x 2 − 2 f(x) = x^2 - 2 f ( x ) = x 2 − 2 、x 0 = 2 x_0 = 2 x 0 = 2 、x 1 = 1.5 x_1 = 1.5 x 1 = 1.5 から割線法を mpmath の 300 桁で回し、誤差の比から次数を p k = log ( e k / e k − 1 ) / log ( e k − 1 / e k − 2 ) p_k = \log(e_k/e_{k-1})/\log(e_{k-1}/e_{k-2}) p k = log ( e k / e k − 1 ) / log ( e k − 1 / e k − 2 ) で推定します。
# Python 3.12.10 / mpmath 1.3.0
from mpmath import mp, mpf, sqrt, log
mp.dps = 300 # 誤差が 1e-250 程度まで追えるように
xs = [ mpf ( 2 ), mpf ( " 1.5 " ) ] # x0, x1
for k in range ( 10 ): # 割線法
xs. append ( b - f ( b ) * (b - a) / ( f ( b ) - f ( a ) ) )
e = [ abs ( x - r ) for x in xs ]
C = 2 / ( 2 * 2 * r) # f''(x*)/(2 f'(x*)) = 1/(2√2)
print ( " k e_k e_k/(e_{k-1}e_{k-2}) 推定次数 " )
for k in range ( 2 , len ( e )):
q = e[k] / (e[k - 1 ] * e[k - 2 ])
p = log ( e [ k ] / e [ k - 1 ]) / log ( e [ k - 1 ] / e [ k - 2 ])
s = f " { float ( e [ k ]) :.2e } " . ljust ( 10 )
print ( f " {k :2d } {s} { float ( q ) :.5f } { float ( p ) :.4f } " )
print ( f "1/(2√2) = { float ( C ) :.5f } , 黄金比 = { float ( ( 1 + sqrt ( 5 ) ) / 2 ) :.4f } " )
k e_k e_k/(e_{k-1}e_{k-2}) 推定次数
2 1.44e-02 0.28571 0.9305
3 4.21e-04 0.34146 1.9750
4 2.12e-06 0.35172 1.4980
5 3.16e-10 0.35350 1.6666
6 2.37e-16 0.35355 1.6000
7 2.65e-26 0.35355 1.6250
8 2.22e-42 0.35355 1.6154
9 2.08e-68 0.35355 1.6190
10 1.63e-110 0.35355 1.6176
11 1.20e-178 0.35355 1.6182
1/(2√2) = 0.35355, 黄金比 = 1.6180
検算. (a) 極限値 1 / ( 2 2 ) 1/(2\sqrt2) 1/ ( 2 2 ) は Newton 法の例と同じ f ′ ′ ( x ∗ ) / ( 2 f ′ ( x ∗ ) ) f''(x^*)/(2f'(x^*)) f ′′ ( x ∗ ) / ( 2 f ′ ( x ∗ )) で、表の第 3 列はそこに収束します(定理 6.2 の 2)。(b) x 2 − 2 x^2 - 2 x 2 − 2 では f [ a , b ] = a + b f[a,b] = a + b f [ a , b ] = a + b 、f [ a , b , c ] = 1 f[a,b,c] = 1 f [ a , b , c ] = 1 なので、証明の式は e k + 1 = e k e k − 1 / ( x k − 1 + x k ) e_{k+1} = e_k e_{k-1}/(x_{k-1} + x_k) e k + 1 = e k e k − 1 / ( x k − 1 + x k ) と閉じた形になり、第 3 列は 1 / ( x k − 2 + x k − 1 ) 1/(x_{k-2} + x_{k-1}) 1/ ( x k − 2 + x k − 1 ) に等しくなります(k = 2 k = 2 k = 2 で 1 / 3.5 = 0.28571 1/3.5 = 0.28571 1/3.5 = 0.28571 )。(c) 推定次数は 1.5 , 1.667 , 1.6 , 1.625 , 1.615 , … 1.5, 1.667, 1.6, 1.625, 1.615, \ldots 1.5 , 1.667 , 1.6 , 1.625 , 1.615 , … と黄金比 1.618 1.618 1.618 のまわりを振動しながら近づきます。これは誤差の指数が Fibonacci 数列に似た漸化式 log e k + 1 ≈ log e k + log e k − 1 \log e_{k+1} \approx \log e_k + \log e_{k-1} log e k + 1 ≈ log e k + log e k − 1 に従うためで、連続する Fibonacci 数の比が黄金比に近づく様子と同じです。
U ⊂ R n U \subset \mathbb{R}^n U ⊂ R n を開集合、F : U → R n F : U \to \mathbb{R}^n F : U → R n を全微分可能 な写像とし、ヤコビ行列を J ( x ) = J F ( x ) J(x) = J_F(x) J ( x ) = J F ( x ) と書きます。R n \mathbb{R}^n R n にはユークリッドノルム ∥ ⋅ ∥ \|\cdot\| ∥ ⋅ ∥ 、行列には作用素ノルム ∥ A ∥ = sup { ∥ A x ∥ : ∥ x ∥ ≤ 1 } \|A\| = \sup\{\|Ax\| : \|x\| \le 1\} ∥ A ∥ = sup { ∥ A x ∥ : ∥ x ∥ ≤ 1 } を使います(多変数関数の微分 と同じ記法です)。閉球を B ˉ ( a , r ) = { x : ∥ x − a ∥ ≤ r } \bar B(a, r) = \{x : \|x - a\| \le r\} B ˉ ( a , r ) = { x : ∥ x − a ∥ ≤ r } 、開球を B ( a , r ) B(a, r) B ( a , r ) と書きます。多変数の Newton 法 は
x k + 1 = x k − J ( x k ) − 1 F ( x k ) x_{k+1} = x_k - J(x_k)^{-1} F(x_k) x k + 1 = x k − J ( x k ) − 1 F ( x k )
です。実際の計算では逆行列を作らず、連立 1 次方程式 J ( x k ) d k = − F ( x k ) J(x_k)\, d_k = -F(x_k) J ( x k ) d k = − F ( x k ) を解いて x k + 1 = x k + d k x_{k+1} = x_k + d_k x k + 1 = x k + d k とします。証明の道具を 2 つ用意します。
補題 7.1 (逆行列の摂動 ) ✓
A , B A, B A , B を n n n 次の実正方行列とし、A A A は正則で ∥ A − 1 ( B − A ) ∥ ≤ q < 1 \|A^{-1}(B - A)\| \le q \lt 1 ∥ A − 1 ( B − A ) ∥ ≤ q < 1 とする。このとき B B B は正則で
∥ B − 1 ∥ ≤ ∥ A − 1 ∥ 1 − q \|B^{-1}\| \le \frac{\|A^{-1}\|}{1 - q} ∥ B − 1 ∥ ≤ 1 − q ∥ A − 1 ∥ である。
証明(補題 7.1) B x = 0 Bx = 0 B x = 0 とすると 0 = A − 1 B x = x + A − 1 ( B − A ) x 0 = A^{-1}Bx = x + A^{-1}(B - A)x 0 = A − 1 B x = x + A − 1 ( B − A ) x なので ∥ x ∥ = ∥ A − 1 ( B − A ) x ∥ ≤ q ∥ x ∥ \|x\| = \|A^{-1}(B - A)x\| \le q\|x\| ∥ x ∥ = ∥ A − 1 ( B − A ) x ∥ ≤ q ∥ x ∥ で、q < 1 q \lt 1 q < 1 から x = 0 x = 0 x = 0 です。斉次方程式が自明な解しかもたないので、正則性の判定(定理 7.1[行列式とその性質] の 4 ⇒ 2)から B B B は正則です。任意の y y y について x = B − 1 y x = B^{-1}y x = B − 1 y とおくと、上と同じ変形で x = A − 1 y − A − 1 ( B − A ) x x = A^{-1}y - A^{-1}(B - A)x x = A − 1 y − A − 1 ( B − A ) x なので
∥ x ∥ ≤ ∥ A − 1 ∥ ∥ y ∥ + q ∥ x ∥ , ∥ B − 1 y ∥ = ∥ x ∥ ≤ ∥ A − 1 ∥ 1 − q ∥ y ∥ \|x\| \le \|A^{-1}\|\,\|y\| + q\|x\|, \qquad \|B^{-1}y\| = \|x\| \le \frac{\|A^{-1}\|}{1 - q}\|y\| ∥ x ∥ ≤ ∥ A − 1 ∥ ∥ y ∥ + q ∥ x ∥ , ∥ B − 1 y ∥ = ∥ x ∥ ≤ 1 − q ∥ A − 1 ∥ ∥ y ∥ です。∥ y ∥ ≤ 1 \|y\| \le 1 ∥ y ∥ ≤ 1 について上限をとれば主張を得ます。
∎
これは Banach 空間でのノイマン級数(定理 6.1[有界作用素と共役作用素] )による評価の、有限次元での形です。有限次元では単射と全射が同値なので、級数を使わずに示せます。
次の補題は、多変数の場合に「1 次式で近似した誤差は距離の 2 乗程度」と言うためのものです。証明では、ベクトル値の量を v v v の向きの単位ベクトル w w w との内積で 1 つの実数にして、1 変数の単調性の議論に持ち込みます。
補題 7.2 (1 次近似の剰余の評価 ) ✓
U ⊂ R n U \subset \mathbb{R}^n U ⊂ R n を開集合、F : U → R n F : U \to \mathbb{R}^n F : U → R n を全微分可能、a , b ∈ U a, b \in U a , b ∈ U を結ぶ線分 [ a , b ] [a, b] [ a , b ] が U U U に含まれるとする。γ ≥ 0 \gamma \ge 0 γ ≥ 0 について、線分上のすべての点 z z z で ∥ J ( z ) − J ( a ) ∥ ≤ γ ∥ z − a ∥ \|J(z) - J(a)\| \le \gamma\|z - a\| ∥ J ( z ) − J ( a ) ∥ ≤ γ ∥ z − a ∥ ならば
∥ F ( b ) − F ( a ) − J ( a ) ( b − a ) ∥ ≤ γ 2 ∥ b − a ∥ 2 \|F(b) - F(a) - J(a)(b - a)\| \le \frac{\gamma}{2}\|b - a\|^2 ∥ F ( b ) − F ( a ) − J ( a ) ( b − a ) ∥ ≤ 2 γ ∥ b − a ∥ 2 である。
証明(補題 7.2) v = F ( b ) − F ( a ) − J ( a ) ( b − a ) v = F(b) - F(a) - J(a)(b - a) v = F ( b ) − F ( a ) − J ( a ) ( b − a ) とし、v ≠ 0 v \ne 0 v = 0 としてよいとします。w = v / ∥ v ∥ w = v/\|v\| w = v /∥ v ∥ とおきます。V = { s ∈ R : a + s ( b − a ) ∈ U } V = \{s \in \mathbb{R} : a + s(b - a) \in U\} V = { s ∈ R : a + s ( b − a ) ∈ U } は開集合 U U U の連続写像による逆像なので開集合で、[ 0 , 1 ] [0,1] [ 0 , 1 ] を含みます。V V V 上の実数値関数
ϕ ( s ) = w T ( F ( a + s ( b − a ) ) − s J ( a ) ( b − a ) ) \phi(s) = w^{\mathsf{T}}\bigl(F(a + s(b - a)) - s\,J(a)(b - a)\bigr) ϕ ( s ) = w T ( F ( a + s ( b − a )) − s J ( a ) ( b − a ) ) は、連鎖律(定理 9.4[多変数関数の微分と偏微分] )から微分可能で、z s = a + s ( b − a ) z_s = a + s(b-a) z s = a + s ( b − a ) とおくと ϕ ′ ( s ) = w T ( J ( z s ) − J ( a ) ) ( b − a ) \phi'(s) = w^{\mathsf{T}}\bigl(J(z_s) - J(a)\bigr)(b - a) ϕ ′ ( s ) = w T ( J ( z s ) − J ( a ) ) ( b − a ) です。s ∈ [ 0 , 1 ] s \in [0,1] s ∈ [ 0 , 1 ] では、コーシー・シュワルツの不等式(定理 4.2[内積空間とグラム・シュミット直交化] )と ∥ w ∥ = 1 \|w\| = 1 ∥ w ∥ = 1 、∥ z s − a ∥ = s ∥ b − a ∥ \|z_s - a\| = s\|b - a\| ∥ z s − a ∥ = s ∥ b − a ∥ から
∣ ϕ ′ ( s ) ∣ ≤ ∥ J ( z s ) − J ( a ) ∥ ∥ b − a ∥ ≤ γ s ∥ b − a ∥ 2 |\phi'(s)| \le \|J(z_s) - J(a)\|\,\|b - a\| \le \gamma s \|b - a\|^2 ∣ ϕ ′ ( s ) ∣ ≤ ∥ J ( z s ) − J ( a ) ∥ ∥ b − a ∥ ≤ γ s ∥ b − a ∥ 2 です。c = γ ∥ b − a ∥ 2 c = \gamma\|b - a\|^2 c = γ ∥ b − a ∥ 2 とし、ψ ( s ) = ϕ ( s ) − c s 2 / 2 \psi(s) = \phi(s) - c s^2/2 ψ ( s ) = ϕ ( s ) − c s 2 /2 とおくと [ 0 , 1 ] [0,1] [ 0 , 1 ] 上で ψ ′ ( s ) = ϕ ′ ( s ) − c s ≤ 0 \psi'(s) = \phi'(s) - cs \le 0 ψ ′ ( s ) = ϕ ′ ( s ) − cs ≤ 0 なので、導関数の符号と単調性(系 3.4[平均値の定理とテイラーの定理] 。− ψ -\psi − ψ に適用)から ψ ( 1 ) ≤ ψ ( 0 ) \psi(1) \le \psi(0) ψ ( 1 ) ≤ ψ ( 0 ) 、すなわち ϕ ( 1 ) − ϕ ( 0 ) ≤ c / 2 \phi(1) - \phi(0) \le c/2 ϕ ( 1 ) − ϕ ( 0 ) ≤ c /2 です。一方 ϕ ( 1 ) − ϕ ( 0 ) = w T v = ∥ v ∥ \phi(1) - \phi(0) = w^{\mathsf{T}} v = \|v\| ϕ ( 1 ) − ϕ ( 0 ) = w T v = ∥ v ∥ なので、∥ v ∥ ≤ c / 2 \|v\| \le c/2 ∥ v ∥ ≤ c /2 です。
∎
平均値不等式(定理 9.7[多変数関数の微分と偏微分] )をそのまま使うと係数は γ \gamma γ になり、1 / 2 1/2 1/2 が出ません。この 1 / 2 1/2 1/2 が、Kantorovich 型の定理の条件 h < 1 / 2 h \lt 1/2 h < 1/2 の 1 / 2 1/2 1/2 の出どころです。
定理 7.3 (多変数の Newton 法の局所 2 次収束 ) ✓
U ⊂ R n U \subset \mathbb{R}^n U ⊂ R n を開集合、F : U → R n F : U \to \mathbb{R}^n F : U → R n を全微分可能、x ∗ ∈ U x^* \in U x ∗ ∈ U で F ( x ∗ ) = 0 F(x^*) = 0 F ( x ∗ ) = 0 、J ( x ∗ ) J(x^*) J ( x ∗ ) は正則とし、β = ∥ J ( x ∗ ) − 1 ∥ \beta = \|J(x^*)^{-1}\| β = ∥ J ( x ∗ ) − 1 ∥ とおく。r > 0 r \gt 0 r > 0 と γ > 0 \gamma \gt 0 γ > 0 が B ˉ ( x ∗ , r ) ⊂ U \bar B(x^*, r) \subset U B ˉ ( x ∗ , r ) ⊂ U と
∥ J ( x ) − J ( y ) ∥ ≤ γ ∥ x − y ∥ ( x , y ∈ B ˉ ( x ∗ , r ) ) \|J(x) - J(y)\| \le \gamma\|x - y\| \qquad (x, y \in \bar B(x^*, r)) ∥ J ( x ) − J ( y ) ∥ ≤ γ ∥ x − y ∥ ( x , y ∈ B ˉ ( x ∗ , r )) を満たすとし、ρ = min ( r , 1 / ( 2 β γ ) ) \rho = \min(r, 1/(2\beta\gamma)) ρ = min ( r , 1/ ( 2 β γ )) とする。このとき ∥ x 0 − x ∗ ∥ ≤ ρ \|x_0 - x^*\| \le \rho ∥ x 0 − x ∗ ∥ ≤ ρ なら Newton 法の列は定義されて B ˉ ( x ∗ , ρ ) \bar B(x^*, \rho) B ˉ ( x ∗ , ρ ) にとどまり、e k = x k − x ∗ e_k = x_k - x^* e k = x k − x ∗ について
∥ e k + 1 ∥ ≤ β γ ∥ e k ∥ 2 , β γ ∥ e k ∥ ≤ ( β γ ∥ e 0 ∥ ) 2 k ≤ 2 − 2 k \|e_{k+1}\| \le \beta\gamma\|e_k\|^2, \qquad \beta\gamma\|e_k\| \le (\beta\gamma\|e_0\|)^{2^k} \le 2^{-2^k} ∥ e k + 1 ∥ ≤ β γ ∥ e k ∥ 2 , β γ ∥ e k ∥ ≤ ( β γ ∥ e 0 ∥ ) 2 k ≤ 2 − 2 k が成り立つ。とくに x k → x ∗ x_k \to x^* x k → x ∗ である。
証明(定理 7.3) ∥ e k ∥ ≤ ρ \|e_k\| \le \rho ∥ e k ∥ ≤ ρ と仮定します。∥ J ( x ∗ ) − 1 ( J ( x k ) − J ( x ∗ ) ) ∥ ≤ β γ ∥ e k ∥ ≤ 1 / 2 \|J(x^*)^{-1}(J(x_k) - J(x^*))\| \le \beta\gamma\|e_k\| \le 1/2 ∥ J ( x ∗ ) − 1 ( J ( x k ) − J ( x ∗ )) ∥ ≤ β γ ∥ e k ∥ ≤ 1/2 なので、補題 7.1 (A = J ( x ∗ ) A = J(x^*) A = J ( x ∗ ) 、B = J ( x k ) B = J(x_k) B = J ( x k ) 、q = 1 / 2 q = 1/2 q = 1/2 )から J ( x k ) J(x_k) J ( x k ) は正則で ∥ J ( x k ) − 1 ∥ ≤ 2 β \|J(x_k)^{-1}\| \le 2\beta ∥ J ( x k ) − 1 ∥ ≤ 2 β です。F ( x ∗ ) = 0 F(x^*) = 0 F ( x ∗ ) = 0 を使うと
e k + 1 = x k − x ∗ − J ( x k ) − 1 F ( x k ) = J ( x k ) − 1 [ F ( x ∗ ) − F ( x k ) − J ( x k ) ( x ∗ − x k ) ] e_{k+1} = x_k - x^* - J(x_k)^{-1}F(x_k) = J(x_k)^{-1}\bigl[F(x^*) - F(x_k) - J(x_k)(x^* - x_k)\bigr] e k + 1 = x k − x ∗ − J ( x k ) − 1 F ( x k ) = J ( x k ) − 1 [ F ( x ∗ ) − F ( x k ) − J ( x k ) ( x ∗ − x k ) ] です。閉球は凸なので線分 [ x k , x ∗ ] [x_k, x^*] [ x k , x ∗ ] は B ˉ ( x ∗ , r ) \bar B(x^*, r) B ˉ ( x ∗ , r ) に含まれ、その上で ∥ J ( z ) − J ( x k ) ∥ ≤ γ ∥ z − x k ∥ \|J(z) - J(x_k)\| \le \gamma\|z - x_k\| ∥ J ( z ) − J ( x k ) ∥ ≤ γ ∥ z − x k ∥ です。補題 7.2 (a = x k a = x_k a = x k 、b = x ∗ b = x^* b = x ∗ )から括弧の中のノルムは ( γ / 2 ) ∥ e k ∥ 2 (\gamma/2)\|e_k\|^2 ( γ /2 ) ∥ e k ∥ 2 以下なので、
∥ e k + 1 ∥ ≤ 2 β ⋅ γ 2 ∥ e k ∥ 2 = β γ ∥ e k ∥ 2 ≤ β γ ρ ∥ e k ∥ ≤ 1 2 ∥ e k ∥ \|e_{k+1}\| \le 2\beta \cdot \frac{\gamma}{2}\|e_k\|^2 = \beta\gamma\|e_k\|^2 \le \beta\gamma\rho\,\|e_k\| \le \frac12\|e_k\| ∥ e k + 1 ∥ ≤ 2 β ⋅ 2 γ ∥ e k ∥ 2 = β γ ∥ e k ∥ 2 ≤ β γ ρ ∥ e k ∥ ≤ 2 1 ∥ e k ∥ です。よって x k + 1 ∈ B ˉ ( x ∗ , ρ ) x_{k+1} \in \bar B(x^*, \rho) x k + 1 ∈ B ˉ ( x ∗ , ρ ) で、帰納法ですべての k k k で成り立ちます。β γ ∥ e k + 1 ∥ ≤ ( β γ ∥ e k ∥ ) 2 \beta\gamma\|e_{k+1}\| \le (\beta\gamma\|e_k\|)^2 β γ ∥ e k + 1 ∥ ≤ ( β γ ∥ e k ∥ ) 2 を繰り返せば 2 つ目の不等式を得ます。
∎
この定理も解 x ∗ x^* x ∗ を中心に述べているので、解が存在することは前提です。次の定理は、解を知らないまま、初期値 x 0 x_0 x 0 で計算できる 3 つの量 β , η , γ \beta, \eta, \gamma β , η , γ から、解の存在と一意性と、解を含む球の半径を結論します。
定理 7.4 (Kantorovich 型の存在と一意性(簡易 Newton 法) ) ✓
U ⊂ R n U \subset \mathbb{R}^n U ⊂ R n を開集合、F : U → R n F : U \to \mathbb{R}^n F : U → R n を全微分可能、x 0 ∈ U x_0 \in U x 0 ∈ U で J 0 = J ( x 0 ) J_0 = J(x_0) J 0 = J ( x 0 ) は正則とする。正の数 β , η , γ , R \beta, \eta, \gamma, R β , η , γ , R が
∥ J 0 − 1 ∥ ≤ β , ∥ J 0 − 1 F ( x 0 ) ∥ ≤ η , ∥ J ( x ) − J 0 ∥ ≤ γ ∥ x − x 0 ∥ ( x ∈ B ˉ ( x 0 , R ) ) \|J_0^{-1}\| \le \beta, \qquad \|J_0^{-1}F(x_0)\| \le \eta, \qquad \|J(x) - J_0\| \le \gamma\|x - x_0\| \quad (x \in \bar B(x_0, R)) ∥ J 0 − 1 ∥ ≤ β , ∥ J 0 − 1 F ( x 0 ) ∥ ≤ η , ∥ J ( x ) − J 0 ∥ ≤ γ ∥ x − x 0 ∥ ( x ∈ B ˉ ( x 0 , R )) と B ˉ ( x 0 , R ) ⊂ U \bar B(x_0, R) \subset U B ˉ ( x 0 , R ) ⊂ U を満たすとし、L = β γ L = \beta\gamma L = β γ 、h = L η h = L\eta h = L η とおく。h < 1 / 2 h \lt 1/2 h < 1/2 かつ
t − = 1 − 1 − 2 h L ≤ R t_- = \frac{1 - \sqrt{1 - 2h}}{L} \le R t − = L 1 − 1 − 2 h ≤ R ならば次が成り立つ。
F F F は B ˉ ( x 0 , t − ) \bar B(x_0, t_-) B ˉ ( x 0 , t − ) に零点 x ∗ x^* x ∗ をもち、B ˉ ( x 0 , t − ) \bar B(x_0, t_-) B ˉ ( x 0 , t − ) の中の零点はそれだけである。
簡易 Newton 法 x k + 1 = x k − J 0 − 1 F ( x k ) x_{k+1} = x_k - J_0^{-1}F(x_k) x k + 1 = x k − J 0 − 1 F ( x k ) の列は B ˉ ( x 0 , t − ) \bar B(x_0, t_-) B ˉ ( x 0 , t − ) にとどまって x ∗ x^* x ∗ に収束し、q = L t − = 1 − 1 − 2 h < 1 q = L t_- = 1 - \sqrt{1 - 2h} \lt 1 q = L t − = 1 − 1 − 2 h < 1 について ∥ x k − x ∗ ∥ ≤ q k η / ( 1 − q ) \|x_k - x^*\| \le q^k\eta/(1 - q) ∥ x k − x ∗ ∥ ≤ q k η / ( 1 − q ) である。
F F F は開球 B ( x 0 , min ( R , 1 / L ) ) B(x_0, \min(R, 1/L)) B ( x 0 , min ( R , 1/ L )) の中に x ∗ x^* x ∗ 以外の零点をもたない。
証明(定理 7.4) T ( x ) = x − J 0 − 1 F ( x ) T(x) = x - J_0^{-1}F(x) T ( x ) = x − J 0 − 1 F ( x ) とおきます。J 0 − 1 J_0^{-1} J 0 − 1 は正則なので、T ( x ) = x T(x) = x T ( x ) = x と F ( x ) = 0 F(x) = 0 F ( x ) = 0 は同値です。
縮小性. 0 < t ≤ R 0 \lt t \le R 0 < t ≤ R とします。T T T は全微分可能で J T ( z ) = I − J 0 − 1 J ( z ) = J 0 − 1 ( J 0 − J ( z ) ) J_T(z) = I - J_0^{-1}J(z) = J_0^{-1}(J_0 - J(z)) J T ( z ) = I − J 0 − 1 J ( z ) = J 0 − 1 ( J 0 − J ( z )) なので、z ∈ B ˉ ( x 0 , t ) z \in \bar B(x_0, t) z ∈ B ˉ ( x 0 , t ) で ∥ J T ( z ) ∥ ≤ β γ ∥ z − x 0 ∥ ≤ L t \|J_T(z)\| \le \beta\gamma\|z - x_0\| \le Lt ∥ J T ( z ) ∥ ≤ β γ ∥ z − x 0 ∥ ≤ L t です。x , y ∈ B ˉ ( x 0 , t ) x, y \in \bar B(x_0, t) x , y ∈ B ˉ ( x 0 , t ) を結ぶ線分は閉球に含まれるので、平均値不等式(定理 9.7[多変数関数の微分と偏微分] )から
∥ T ( x ) − T ( y ) ∥ ≤ L t ∥ x − y ∥ \|T(x) - T(y)\| \le Lt\,\|x - y\| ∥ T ( x ) − T ( y ) ∥ ≤ L t ∥ x − y ∥ です。
閉球を自分自身に写すこと. x ∈ B ˉ ( x 0 , t ) x \in \bar B(x_0, t) x ∈ B ˉ ( x 0 , t ) について
T ( x ) − x 0 = − J 0 − 1 [ F ( x ) − F ( x 0 ) − J 0 ( x − x 0 ) ] − J 0 − 1 F ( x 0 ) T(x) - x_0 = -J_0^{-1}\bigl[F(x) - F(x_0) - J_0(x - x_0)\bigr] - J_0^{-1}F(x_0) T ( x ) − x 0 = − J 0 − 1 [ F ( x ) − F ( x 0 ) − J 0 ( x − x 0 ) ] − J 0 − 1 F ( x 0 ) です。仮定から線分 [ x 0 , x ] [x_0, x] [ x 0 , x ] 上で ∥ J ( z ) − J 0 ∥ ≤ γ ∥ z − x 0 ∥ \|J(z) - J_0\| \le \gamma\|z - x_0\| ∥ J ( z ) − J 0 ∥ ≤ γ ∥ z − x 0 ∥ なので、補題 7.2 (a = x 0 a = x_0 a = x 0 、b = x b = x b = x )から
∥ T ( x ) − x 0 ∥ ≤ β ⋅ γ 2 ∥ x − x 0 ∥ 2 + η ≤ L 2 t 2 + η \|T(x) - x_0\| \le \beta\cdot\frac{\gamma}{2}\|x - x_0\|^2 + \eta \le \frac{L}{2}t^2 + \eta ∥ T ( x ) − x 0 ∥ ≤ β ⋅ 2 γ ∥ x − x 0 ∥ 2 + η ≤ 2 L t 2 + η です。t − t_- t − は 2 次方程式 L 2 t 2 − t + η = 0 \frac{L}{2}t^2 - t + \eta = 0 2 L t 2 − t + η = 0 の小さいほうの根(h < 1 / 2 h \lt 1/2 h < 1/2 なので実根)なので、t = t − t = t_- t = t − のとき右辺はちょうど t − t_- t − です。よって T T T は B ˉ ( x 0 , t − ) \bar B(x_0, t_-) B ˉ ( x 0 , t − ) を自分自身に写します。
1・2. B ˉ ( x 0 , t − ) \bar B(x_0, t_-) B ˉ ( x 0 , t − ) は R n \mathbb{R}^n R n の空でない閉集合で、T T T はそこで縮小率 q = L t − = 1 − 1 − 2 h q = Lt_- = 1 - \sqrt{1 - 2h} q = L t − = 1 − 1 − 2 h の縮小写像です。h < 1 / 2 h \lt 1/2 h < 1/2 から q < 1 q \lt 1 q < 1 です。縮小写像の原理(定理 3.2[逆関数定理・陰関数定理と条件付き極値] )から T T T は B ˉ ( x 0 , t − ) \bar B(x_0, t_-) B ˉ ( x 0 , t − ) にただ 1 つの不動点 x ∗ x^* x ∗ をもち、簡易 Newton 法の列(T T T の反復そのもの)はそこに収束して ∥ x k − x ∗ ∥ ≤ q k ∥ x 1 − x 0 ∥ / ( 1 − q ) \|x_k - x^*\| \le q^k\|x_1 - x_0\|/(1-q) ∥ x k − x ∗ ∥ ≤ q k ∥ x 1 − x 0 ∥/ ( 1 − q ) です。∥ x 1 − x 0 ∥ = ∥ J 0 − 1 F ( x 0 ) ∥ ≤ η \|x_1 - x_0\| = \|J_0^{-1}F(x_0)\| \le \eta ∥ x 1 − x 0 ∥ = ∥ J 0 − 1 F ( x 0 ) ∥ ≤ η から 2 の評価を得ます。不動点の一意性は、この閉球の中の零点の一意性です。
3. y y y を ∥ y − x 0 ∥ < min ( R , 1 / L ) \|y - x_0\| \lt \min(R, 1/L) ∥ y − x 0 ∥ < min ( R , 1/ L ) を満たす零点とし、t = max ( ∥ y − x 0 ∥ , t − ) t = \max(\|y - x_0\|, t_-) t = max ( ∥ y − x 0 ∥ , t − ) とおきます。t − ≤ R t_- \le R t − ≤ R と L t − = q < 1 Lt_- = q \lt 1 L t − = q < 1 から t ≤ R t \le R t ≤ R 、L t < 1 Lt \lt 1 L t < 1 です。x ∗ x^* x ∗ と y y y はともに B ˉ ( x 0 , t ) \bar B(x_0, t) B ˉ ( x 0 , t ) の中の T T T の不動点なので、縮小性から ∥ x ∗ − y ∥ = ∥ T ( x ∗ ) − T ( y ) ∥ ≤ L t ∥ x ∗ − y ∥ \|x^* - y\| = \|T(x^*) - T(y)\| \le Lt\|x^* - y\| ∥ x ∗ − y ∥ = ∥ T ( x ∗ ) − T ( y ) ∥ ≤ L t ∥ x ∗ − y ∥ となり、x ∗ = y x^* = y x ∗ = y です。
∎
条件 h = β γ η < 1 / 2 h = \beta\gamma\eta \lt 1/2 h = β γ η < 1/2 は、「初期値での Newton 補正の大きさ η \eta η 」が「ヤコビ行列の変化の速さ γ \gamma γ 」と「ヤコビ行列の逆の大きさ β \beta β 」に比べて十分小さい、という意味です。左辺の 3 つの量はすべて x 0 x_0 x 0 と F F F の情報から計算でき、とくに区間演算で上界を計算すれば、計算機の出力がそのまま「この球の中に解がただ 1 つある」ことの証明になります。
x₀ x*(零点) t₋ 破線の球の中の零点は x* だけ
Kantorovich 型の定理の状況(模式図)。初期値 x₀ を中心とする半径 t₋ の閉球の中に零点 x* があり、より大きい開球(半径 min(R, 1/L))の中では零点は x* だけである。
演習 8.1 易
f f f が [ 1 , 2 ] [1, 2] [ 1 , 2 ] で連続、f ( 1 ) f ( 2 ) < 0 f(1) f(2) \lt 0 f ( 1 ) f ( 2 ) < 0 とします。二分法の中点 c k c_k c k の誤差が 10 − 6 10^{-6} 1 0 − 6 以下であることが定理 3.1 から保証される最小の k k k を求めてください。
解答 誤差の上界は ( 2 − 1 ) / 2 k + 1 = 2 − ( k + 1 ) (2 - 1)/2^{k+1} = 2^{-(k+1)} ( 2 − 1 ) / 2 k + 1 = 2 − ( k + 1 ) です。2 − ( k + 1 ) ≤ 10 − 6 2^{-(k+1)} \le 10^{-6} 2 − ( k + 1 ) ≤ 1 0 − 6 は k + 1 ≥ 6 log 2 10 ≈ 19.93 k + 1 \ge 6\log_2 10 \approx 19.93 k + 1 ≥ 6 log 2 10 ≈ 19.93 、すなわち k ≥ 19 k \ge 19 k ≥ 19 と同値です。よって k = 19 k = 19 k = 19 で、c 0 , … , c 19 c_0, \ldots, c_{19} c 0 , … , c 19 の 20 回の関数値の計算で保証されます。
演習 8.2 標準
a > 0 a \gt 0 a > 0 とします。f ( x ) = 1 / x − a f(x) = 1/x - a f ( x ) = 1/ x − a に Newton 法を使うと x k + 1 = x k ( 2 − a x k ) x_{k+1} = x_k(2 - a x_k) x k + 1 = x k ( 2 − a x k ) となり、割り算なしで 1 / a 1/a 1/ a を計算できることを確かめてください。さらに r k = 1 − a x k r_k = 1 - a x_k r k = 1 − a x k とおいて r k + 1 = r k 2 r_{k+1} = r_k^2 r k + 1 = r k 2 を示し、列が 1 / a 1/a 1/ a に収束するための x 0 x_0 x 0 の必要十分条件と、e k + 1 / e k 2 e_{k+1}/e_k^2 e k + 1 / e k 2 の値を求めてください。
解答 f ′ ( x ) = − 1 / x 2 f'(x) = -1/x^2 f ′ ( x ) = − 1/ x 2 なので x k + 1 = x k + x k 2 ( 1 / x k − a ) = x k ( 2 − a x k ) x_{k+1} = x_k + x_k^2(1/x_k - a) = x_k(2 - a x_k) x k + 1 = x k + x k 2 ( 1/ x k − a ) = x k ( 2 − a x k ) です。r k + 1 = 1 − a x k ( 2 − a x k ) = 1 − 2 a x k + a 2 x k 2 = ( 1 − a x k ) 2 = r k 2 r_{k+1} = 1 - a x_k(2 - a x_k) = 1 - 2a x_k + a^2 x_k^2 = (1 - a x_k)^2 = r_k^2 r k + 1 = 1 − a x k ( 2 − a x k ) = 1 − 2 a x k + a 2 x k 2 = ( 1 − a x k ) 2 = r k 2 なので r k = r 0 2 k r_k = r_0^{2^k} r k = r 0 2 k です。x k → 1 / a x_k \to 1/a x k → 1/ a は r k → 0 r_k \to 0 r k → 0 と同値で、それは ∣ r 0 ∣ < 1 |r_0| \lt 1 ∣ r 0 ∣ < 1 、すなわち 0 < x 0 < 2 / a 0 \lt x_0 \lt 2/a 0 < x 0 < 2/ a と同値です(∣ r 0 ∣ = 1 |r_0| = 1 ∣ r 0 ∣ = 1 なら r k = 1 r_k = 1 r k = 1 で x k = 0 x_k = 0 x k = 0 にとどまり、∣ r 0 ∣ > 1 |r_0| \gt 1 ∣ r 0 ∣ > 1 なら ∣ r k ∣ → ∞ |r_k| \to \infty ∣ r k ∣ → ∞ です)。e k = x k − 1 / a = − r k / a e_k = x_k - 1/a = -r_k/a e k = x k − 1/ a = − r k / a なので e k + 1 = − r k 2 / a = − a e k 2 e_{k+1} = -r_k^2/a = -a e_k^2 e k + 1 = − r k 2 / a = − a e k 2 、すなわち e k + 1 / e k 2 = − a e_{k+1}/e_k^2 = -a e k + 1 / e k 2 = − a です。これは 定理 5.1 の極限 f ′ ′ ( x ∗ ) / ( 2 f ′ ( x ∗ ) ) = ( 2 a 3 ) / ( 2 ⋅ ( − a 2 ) ) = − a f''(x^*)/(2f'(x^*)) = (2a^3)/(2 \cdot (-a^2)) = -a f ′′ ( x ∗ ) / ( 2 f ′ ( x ∗ )) = ( 2 a 3 ) / ( 2 ⋅ ( − a 2 )) = − a と一致します。
演習 8.3 標準
a > 0 a \gt 0 a > 0 とし、g ( x ) = x ( x 2 + 3 a ) / ( 3 x 2 + a ) g(x) = x(x^2 + 3a)/(3x^2 + a) g ( x ) = x ( x 2 + 3 a ) / ( 3 x 2 + a ) とします。g ( a ) = a g(\sqrt a) = \sqrt a g ( a ) = a と
g ( x ) − a = ( x − a ) 3 3 x 2 + a g(x) - \sqrt a = \frac{(x - \sqrt a)^3}{3x^2 + a} g ( x ) − a = 3 x 2 + a ( x − a ) 3 を示し、x 0 x_0 x 0 が a \sqrt a a に十分近いとき不動点反復が 3 次収束することと、漸近誤差定数を求めてください。
解答 s = a s = \sqrt a s = a とおくと ( x − s ) 3 = x 3 − 3 s x 2 + 3 a x − a s (x - s)^3 = x^3 - 3sx^2 + 3ax - as ( x − s ) 3 = x 3 − 3 s x 2 + 3 a x − a s です。一方 g ( x ) − s = ( x 3 + 3 a x − s ( 3 x 2 + a ) ) / ( 3 x 2 + a ) = ( x 3 − 3 s x 2 + 3 a x − a s ) / ( 3 x 2 + a ) g(x) - s = \bigl(x^3 + 3ax - s(3x^2 + a)\bigr)/(3x^2 + a) = (x^3 - 3sx^2 + 3ax - as)/(3x^2 + a) g ( x ) − s = ( x 3 + 3 a x − s ( 3 x 2 + a ) ) / ( 3 x 2 + a ) = ( x 3 − 3 s x 2 + 3 a x − a s ) / ( 3 x 2 + a ) なので等式が成り立ち、とくに g ( s ) = s g(s) = s g ( s ) = s です。e k ≠ 0 e_k \ne 0 e k = 0 なら e k + 1 / e k 3 = 1 / ( 3 x k 2 + a ) e_{k+1}/e_k^3 = 1/(3x_k^2 + a) e k + 1 / e k 3 = 1/ ( 3 x k 2 + a ) です。g g g は x > 0 x \gt 0 x > 0 で C 1 C^1 C 1 級で、等式から g ′ ( s ) = 0 g'(s) = 0 g ′ ( s ) = 0 なので、定理 4.3 の 1 により s s s の近くから始めた列は s s s に収束します。よって e k + 1 / e k 3 → 1 / ( 3 a + a ) = 1 / ( 4 a ) e_{k+1}/e_k^3 \to 1/(3a + a) = 1/(4a) e k + 1 / e k 3 → 1/ ( 3 a + a ) = 1/ ( 4 a ) で、3 次収束、漸近誤差定数は 1 / ( 4 a ) 1/(4a) 1/ ( 4 a ) です(定理 4.3 の 2 で p = 3 p = 3 p = 3 としても同じ結論で、g ′ ′ ′ ( s ) / 3 ! = 1 / ( 4 a ) g'''(s)/3! = 1/(4a) g ′′′ ( s ) /3 ! = 1/ ( 4 a ) です)。
演習 8.4 難
F ( x , y ) = ( x 2 + y 2 − 1 , x − y ) F(x, y) = (x^2 + y^2 - 1,\ x - y) F ( x , y ) = ( x 2 + y 2 − 1 , x − y ) とします。(i) ( 1 , 0 ) (1, 0) ( 1 , 0 ) から多変数の Newton 法を 2 回行い、( 3 / 4 , 3 / 4 ) (3/4, 3/4) ( 3/4 , 3/4 ) を得ることを確かめてください。(ii) x 0 = ( 3 / 4 , 3 / 4 ) x_0 = (3/4, 3/4) x 0 = ( 3/4 , 3/4 ) で定理 7.4 の β , η , γ \beta, \eta, \gamma β , η , γ (R R R は任意)を求め、h h h と t − t_- t − を計算して、零点 ( 1 / 2 , 1 / 2 ) (1/\sqrt2, 1/\sqrt2) ( 1/ 2 , 1/ 2 ) が B ˉ ( x 0 , t − ) \bar B(x_0, t_-) B ˉ ( x 0 , t − ) に入っていることを確かめてください。
解答 (i) J ( x , y ) = ( 2 x 2 y 1 − 1 ) J(x, y) = \begin{pmatrix} 2x & 2y \\ 1 & -1 \end{pmatrix} J ( x , y ) = ( 2 x 1 2 y − 1 ) です。( 1 , 0 ) (1,0) ( 1 , 0 ) では J = ( 2 0 1 − 1 ) J = \begin{pmatrix} 2 & 0 \\ 1 & -1\end{pmatrix} J = ( 2 1 0 − 1 ) 、F = ( 0 , 1 ) F = (0, 1) F = ( 0 , 1 ) で、J d = − F J d = -F J d = − F を解くと d = ( 0 , 1 ) d = (0, 1) d = ( 0 , 1 ) 、次の点は ( 1 , 1 ) (1, 1) ( 1 , 1 ) です。そこでは J = ( 2 2 1 − 1 ) J = \begin{pmatrix} 2 & 2 \\ 1 & -1\end{pmatrix} J = ( 2 1 2 − 1 ) 、F = ( 1 , 0 ) F = (1, 0) F = ( 1 , 0 ) で、d = ( − 1 / 4 , − 1 / 4 ) d = (-1/4, -1/4) d = ( − 1/4 , − 1/4 ) 、次の点は ( 3 / 4 , 3 / 4 ) (3/4, 3/4) ( 3/4 , 3/4 ) です。
(ii) J ( p ) − J ( p ′ ) = ( 2 ( p 1 − p 1 ′ ) 2 ( p 2 − p 2 ′ ) 0 0 ) J(p) - J(p') = \begin{pmatrix} 2(p_1 - p'_1) & 2(p_2 - p'_2) \\ 0 & 0\end{pmatrix} J ( p ) − J ( p ′ ) = ( 2 ( p 1 − p 1 ′ ) 0 2 ( p 2 − p 2 ′ ) 0 ) の作用素ノルムは第 1 行のベクトルのノルム 2 ∥ p − p ′ ∥ 2\|p - p'\| 2∥ p − p ′ ∥ なので、γ = 2 \gamma = 2 γ = 2 (どこでも成り立つので R R R は任意)です。J 0 = ( 3 / 2 3 / 2 1 − 1 ) J_0 = \begin{pmatrix} 3/2 & 3/2 \\ 1 & -1\end{pmatrix} J 0 = ( 3/2 1 3/2 − 1 ) 、det J 0 = − 3 \det J_0 = -3 det J 0 = − 3 で、J 0 − 1 = ( 1 / 3 1 / 2 1 / 3 − 1 / 2 ) J_0^{-1} = \begin{pmatrix} 1/3 & 1/2 \\ 1/3 & -1/2\end{pmatrix} J 0 − 1 = ( 1/3 1/3 1/2 − 1/2 ) です。この行列の 2 つの列は直交し、ノルムは 2 / 3 \sqrt2/3 2 /3 と 2 / 2 \sqrt2/2 2 /2 なので、∥ J 0 − 1 v ∥ 2 = v 1 2 ⋅ 2 / 9 + v 2 2 ⋅ 1 / 2 \|J_0^{-1}v\|^2 = v_1^2 \cdot 2/9 + v_2^2 \cdot 1/2 ∥ J 0 − 1 v ∥ 2 = v 1 2 ⋅ 2/9 + v 2 2 ⋅ 1/2 から β = ∥ J 0 − 1 ∥ = 2 / 2 \beta = \|J_0^{-1}\| = \sqrt2/2 β = ∥ J 0 − 1 ∥ = 2 /2 です。F ( x 0 ) = ( 1 / 8 , 0 ) F(x_0) = (1/8, 0) F ( x 0 ) = ( 1/8 , 0 ) 、J 0 − 1 F ( x 0 ) = ( 1 / 24 , 1 / 24 ) J_0^{-1}F(x_0) = (1/24, 1/24) J 0 − 1 F ( x 0 ) = ( 1/24 , 1/24 ) なので η = 2 / 24 \eta = \sqrt2/24 η = 2 /24 です。L = β γ = 2 L = \beta\gamma = \sqrt2 L = β γ = 2 、h = L η = 1 / 12 < 1 / 2 h = L\eta = 1/12 \lt 1/2 h = L η = 1/12 < 1/2 で、
t − = 1 − 5 / 6 2 ≈ 0.06161 t_- = \frac{1 - \sqrt{5/6}}{\sqrt2} \approx 0.06161 t − = 2 1 − 5/6 ≈ 0.06161 です。零点との距離は 2 ( 3 / 4 − 1 / 2 ) ≈ 0.06066 < t − \sqrt2\,(3/4 - 1/\sqrt2) \approx 0.06066 \lt t_- 2 ( 3/4 − 1/ 2 ) ≈ 0.06066 < t − で、たしかに球に入っています。さらに 3 から、半径 1 / L ≈ 0.707 1/L \approx 0.707 1/ L ≈ 0.707 の開球には零点がこれ以外にありません(もう 1 つの零点 ( − 1 / 2 , − 1 / 2 ) (-1/\sqrt2, -1/\sqrt2) ( − 1/ 2 , − 1/ 2 ) との距離は約 2.06 2.06 2.06 です)。
山本哲朗『数値解析入門』(増訂版)サイエンス社、2003 — 非線形方程式の章(二分法・Newton 法・Kantorovich の定理)。
杉原正顯・室田一雄『数値計算法の数理』岩波書店、1994 — 非線形方程式の章。
K. E. Atkinson, An Introduction to Numerical Analysis , 2nd ed., Wiley, 1989 — Chapter 2(求根法、割線法の次数)。
J. M. Ortega, W. C. Rheinboldt, Iterative Solution of Nonlinear Equations in Several Variables , Academic Press, 1970(SIAM Classics in Applied Mathematics 30, 2000)— 収束の次数と Newton–Kantorovich の定理。
P. Deuflhard, Newton Methods for Nonlinear Problems: Affine Invariance and Adaptive Algorithms , Springer, 2004 — Chapter 2。
R. P. Brent, Algorithms for Minimization without Derivatives , Prentice-Hall, 1973 — Chapter 4(Brent 法)。
執筆: Claude Code(Opus 5.5)、Claude Code(Sonnet 5.5)・監査 2 回(Claude Code)・最終 2026-10-12