アドベントカレンダー(自称)vol.16
皆さん、こんにちは。yumyum116 です。SWE転職を目指す者です。
アドベントカレンダーに便乗して、1日1記事投稿に挑戦するとともに、執筆活動の習慣化にも挑戦してみます。
本記事は、頭休めになるような、それでいて、学びにもなるtipsを取り扱います。
今回のテーマは ニュートン法 です。
コーヒー片手にお付き合いいただけますと幸いです。
1. ニュートン法とは
関数 f(x) が与えられた時、その導関数 f'(x) を用いて、f(x) = 0 となる解 x を反復計算によって求める方法を、ニュートン法 といいます。
ある適当な値 x_0 から計算をはじめ、次の漸化式を繰り返し計算することによって、真の解 α の近似値を求めます。
x_{n + 1} = x_n - \frac{f(x_n)}{f'(x_n)}
ここで、$n = 0, 1, 2, ... $ とする。
この反復計算を実行し、決められた条件が満たされた時に計算を停止し、その時の値 x_n を真の解 α の近似解とする。
具体例として、\sqrt{5} の値を求めてみましょう。
\sqrt{5} の値を求める
\sqrt{5} は、
の解の一つであるので、関数 f(x) を
と定義し、方程式 f(x) = 0 の解をニュートン法で求める。
初期値を
とする。f(x) の微分は
であるので、第1ステップ後の値は、(1)(2)(3)から
x_1 = x_0 - \frac{f(x_0)}{f'(x_0)} \\
= 5 - \frac{5^2 - 5}{2 * 5} \\
= 3 ...(4)
となる。第2ステップの解は、(1)(3)(4)より、
x_2 = x_1 - \frac{f(x_1)}{f'(x_1)}
= 3 - \frac{3^2 - 5}{2 * 3}
= 2.333333333 ...(5)
となる。第3ステップの値は、(1)(3)(5)より、
x_3 = x_2 - \frac{f(x_2)}{f'(x_2)}
= 2.333333333 - \frac{(2.333333333)^2 - 5}{2 * 2.333333333}
= 2.238095238095238 ...(6)
となる。第4ステップの値は、(1)(3)(6)より、
x_4 = x_3 - \frac{f(x_3)}{f'(x_3)}
= 2.238095238095238 - \frac{(2.238095238095238)^2 - 5}{2 * 2.238095238095238}
= 2.2360688956433634 ...(7)
となる。第5ステップの値は、(1)(3)(7)より、
x_5 = x_4 - \frac{f(x_4)}{f'(x_4)}
= 2.2360688956433634 - \frac{(2.2360688956433634)^2 - 5}{2 * 2.2360688956433634}
= 2.236067977499978 ...(8)
となる。
真の値が
\sqrt{5} = 2.23606797749979
であるので、第5ステップの値は、真の値と13桁ほど一致する。
ここで、反復計算を停止する条件を、各ステップ間の変化量の絶対値が 10^{-5} 以下となった時とすると、
|x_5 - x_4|
= |2.236067977499978 - 2.2360688956433634|
= 9.18 * 10^{-7}
であるから、x_5 と x_4 の差の絶対値が 10^{-5} 以下になる。よって、計算を第5ステップで停止とする。
2. ニュートン法を実装してみる
前章で紹介した漸化式
x_{n + 1} = x_n - \frac{f(x_n)}{f'(x_n)}
を変形します。ある数値 k の平方根を求める場合、f(x) = x^2 - k となるので
x_{n + 1} = x_n - \frac{x_n^2 - k}{2*x_k}
= \frac{1}{2} * (x_n + \frac{k}{x_n})
となります。この式をベースに、ニュートン法を実装しましょう。
newton.py
x = 5.0
while True:
x2 = x - (x * x - 2) / (x * 2)
if abs(x2 - x) < 0.0001:
break
x = x2
print(x)
比較的シンプルな実装に抑えることができました。
記事内に誤謬等ございましたら、修正いたします。その際は、ご連絡いただけますと幸いです。
それでは、また。
Discussion