二分法と勾配法
近似解へ近づく反復と、計算式の形が誤差へ与える影響を調べます。
説明
数値計算では、方程式の解や最小値を反復によって近似します。停止条件と誤差の見積もりがなければ、表示された桁がどこまで信頼できるか分かりません。また、数学的に等しい式でも浮動小数点計算では丸めの影響が異なります。
x²−2=0の正の解を二分法で挟み、同じ解を最小点に持つ関数で勾配更新を試します。更新幅が大きすぎる場合を比較し、最後に差の計算を別の式へ変えて桁落ちを確かめます。
1解を挟む区間
f(x)=x²−2は連続です。1では負、2では正なので、その間に少なくとも一つの解があります。
import math比較用の平方根; f = lambda x: x*x-2解きたい方程式; lo,hi = 1.0,2.0初期区間; f(lo),f(hi),math.sqrt(2)端点の符号と参照値二分法には連続性と端点の符号の条件を使います。比較用に平方根関数の値も表示していますが、反復の区間更新には使いません。
2中点で半分に絞る
中点1.5ではfが正なので、負の端1と中点1.5の間へ絞ります。
mid = (lo+hi)/2区間の中点; f(mid), mid中点の関数値; hi = mid正側の端を置き換える; lo,hi,hi-lo区間と幅新区間は[1,1.5]、幅0.5です。解を含む条件を保ちながら幅を半分にしています。最終式の出力は区間と幅です。
3二分法の反復
区間幅が指定値以下になったら止めます。反復ごとの端点と中点を記録します。
def bisect_root(lo,hi,tolerance):区間と停止幅
trace = []反復の記録
while hi-lo > tolerance:区間幅で停止
mid = (lo+hi)/2中点
trace.append((lo,hi,mid,f(mid)))区間と関数値を保存
if f(mid) > 0: hi = mid正なら右端を更新
else: lo = midそれ以外は左端を更新
return (lo+hi)/2, (hi-lo)/2, trace近似値・誤差上限・記録この関数は正側の根を挟み、左が負・右が正の今回の区間を想定しています。戻り値の第2項は、最後の中点から真の解までの距離の上限です。
4反復数と誤差上限
初期区間[1,2]から幅0.001以下まで繰り返します。最後の区間を参照値と照合します。
root,bound,trace = bisect_root(1.0,2.0,0.001)反復を実行; trace[:3], len(trace), root, bound初期記録・回数・近似・上限, abs(root-math.sqrt(2)) <= bound参照値との誤差を照合10回で幅1/1024となり、中点の誤差上限は1/2048です。反復回数だけでなく、区間に解を含む条件が保たれていることによって誤差を保証しています。
5最小化する関数と勾配
g(x)=x⁴/4−x²は、正の範囲で√2に最小点を持ちます。勾配はx³−2xです。xから学習率×勾配を引いて更新します。
def descent(rate,steps):学習率と更新回数
x = 1.0初期値
history = [x]履歴
for _ in range(steps):指定回数だけ更新
x -= rate*(x**3-2*x)勾配と逆向きへ更新
history.append(x)新しい値を保存
return history履歴を返す勾配が負ならxを増やし、正なら減らします。これは方程式の符号で区間を保つ二分法とは別の手順です。更新幅を決める学習率によって挙動が変わります。
6小さい更新幅と大きい更新幅
学習率0.1と1.0で、それぞれ6回更新します。
descent(0.1,6)小さい学習率, descent(1.0,6)大きい学習率0.1では正の最小点へ近づきます。1.0では1→2→−2→2と振動し、目的の解へ収束しません。勾配更新を行ったという事実だけでは収束を保証できません。
7近い値の引き算
ε=10⁻¹⁶として√(1+ε)−1を計算します。数学的には小さな正の値ですが、浮動小数点の丸めが影響します。
epsilon = 1e-16小さい増分; direct = math.sqrt(1+epsilon)-1元の式; 1+epsilon, direct加算の丸めと差この環境では1+εが1.0に丸められ、差は0になります。近い値の差では有効な桁が失われやすく、この例では引き算より前の加算でも情報が消えています。
8等しい式へ変形する
分子を有理化すると、ε/(√(1+ε)+1)になります。小さいεを分子へ残して計算します。
stable = epsilon/(math.sqrt(1+epsilon)+1)差を避けた式; stable, epsilon/2, direct == stable一次近似と元の計算を比較約5×10⁻¹⁷が得られます。数学的に同じ式でも、有限桁の演算では誤差が異なります。近似法の誤差と浮動小数点演算の誤差を分けて考える必要があります。
まとめ
二分法は解を含む区間を保ち、区間幅から近似誤差を見積もります。
勾配更新の収束は関数・初期値・学習率に依存します。
等価な式でも丸め誤差は同じとは限らず、計算の順序や形を工夫します。
操作履歴と解説を振り返れます。
Python 3
import math; f = lambda x: x*x-2; lo,hi = 1.0,2.0; f(lo),f(hi),math.sqrt(2)(-1.0, 2.0, 1.4142135623730951)
mid = (lo+hi)/2; f(mid), mid; hi = mid; lo,hi,hi-lo(0.25, 1.5) (1.0, 1.5, 0.5)
def bisect_root(lo,hi,tolerance):
trace = []
while hi-lo > tolerance:
mid = (lo+hi)/2
trace.append((lo,hi,mid,f(mid)))
if f(mid) > 0: hi = mid
else: lo = mid
return (lo+hi)/2, (hi-lo)/2, traceroot,bound,trace = bisect_root(1.0,2.0,0.001); trace[:3], len(trace), root, bound, abs(root-math.sqrt(2)) <= bound([(1.0, 2.0, 1.5, 0.25), (1.0, 1.5, 1.25, -0.4375), (1.25, 1.5, 1.375, -0.109375)], 10, 1.41455078125, 0.00048828125, True)
def descent(rate,steps):
x = 1.0
history = [x]
for _ in range(steps):
x -= rate*(x**3-2*x)
history.append(x)
return historydescent(0.1,6), descent(1.0,6)([1.0, 1.1, 1.1869, 1.2570776452091001, 1.309844207706923, 1.3470841462563339, 1.3720542776288311], [1.0, 2.0, -2.0, 2.0, -2.0, 2.0, -2.0])
epsilon = 1e-16; direct = math.sqrt(1+epsilon)-1; 1+epsilon, direct(1.0, 0.0)
stable = epsilon/(math.sqrt(1+epsilon)+1); stable, epsilon/2, direct == stable(5e-17, 5e-17, False)
演習完了