2005/06/07 数値解析学と演習 福井 先生 # 先週は二分法、今週は、ニュートン法 (やったひとは解ると思うが..)二分法は、一度の反復で、1 ビットだけ精度が上る # 単精度 ( 24 bit ) ならば、20 回程度で収束する 前回の演習で利用する f(x) ( = \sin{x} - x ) は、それほど時間がからない => 20 回やろうが、10 回やろうが、あまり気にならない。 もし、一回に 30 分かかる計算だったら.. 10 回と 20 回ではえらく違う 遲いが、まあ、安定性があるので、価値がある # 安定性と、速度のトレードオフは数値計算の世界では良くあること == [ニュートン方法] # f(x) が微分可能だと仮定して.. (x_i,f(x_i)) の値が分れば、そこでの微分係数 f'(x_i) を計算し、 点(x_i,f(x_i)) を通り、傾きが f'(x_i) の直線 l_i を考え l_i と x 軸との交点を x_{i+1} とする # 図形的な意味を考えれば、新しい x_{i+1} が、元の x_i より、 # f(x)=0 # の根に近付いていることが解る ## ただし、この性質は何時も成立するわけがなく、例えば、 ## x_i が変曲点の近く ( f'(x_i_) = 0 ) だと、x_{i+1} が無限大になってしまって、上手く行かない !! # 結局.. x_{n+1} = x_n - \frac{f(x_n)}{f'(x_n)} という漸化式で、反復を行う。 # この公式と、テーラー展開の関係 ( 実は、多くの数値解析の公式は、テーラー展開と密接な関係がある ) # テーラ展開 f(x_k+h) = f(k_k)+\frac{h}{1!}f'(x_k)+\frac{h^2}{2!}f''(x_k) + .. # テーラ展開そのものは無限に続くことになる # ニュートン法は、2 次以降の項目を無視する ## 3 次以降の項目を無視すれば、二次のニュートン法を作ることができる ### 手間の割にはメリットが少いので、普通はあまり使われていない ## シャープの電卓に入っていたが、あまり役にたたなかった # ニュートン法の為に、2 次以降の項目を無視した式が以下の通り f(x_k+h) '=, f(k_k)+\frac{h}{1!}f'(x_k) これに対して f(x_k+h) = 0 として、h を解けば、 x_{n+1} = x_n - \frac{f(x_n)}{f'(x_n)} すなわち、ニュートン法の公式ができる。 == [ニュートン法の適用例] # 適当(適切:数学の世界では、これが困る..)なスタート点を取れば.. 先週の問題 f(x) = \sin{x} - x の場合 x_{n+1} = x_n - \frac{\sin{x}}{\cos{x}} = x_n - \tan{x} で、反復すれば良い。 # [注意] この式では、\tan{x} があるので x が 0 の付近で変なことが起きる # ニュートン法が非常に上手く行く例 \sqrt{a} の計算 [\sqrt{a} をニュートン法で解く] 関数 f(x)=x^2-a の零点は f(x)=0 より、 x = \pm \sqrt{a} となる。 これにニュートン法を適用すれば、 x_{n+1} = x_n - \frac{f(x_n)}{f'(x_n)} = x_n - \frac{x_n^2-a}{2x_n} とすれば、高速に収束する。 == [初期値の問題] 根が複数あるばあい、ニュートン法では、初期値をかえると、根が変る # 二分法とことなり、何処になるかは予想ができない。 求めたい根を出すには、初期値を選ぶ必要があるのだが大変.. # \sin{x}-a みたいなものは、根が無限にあるので、大変 ## まあ、\sin{x} は周期関数なので、一つわかればたもわかるが.. 初期値の問題は、非常に重要な問題 # 物理学的には、通常欲しい根の範囲が予めわかっていることが多いので、その場合は、あまり困らない。 == [収束速度] 反復回数を x 軸、得られる桁数を y 軸として対数グラフにすれば、反復法の違いにより、グラフに違いがしょうじる 二分法の場合は、直線 => 一次式になるので、一次収束 反復回数 n に対して n bit しか精度がでない ニュートン法の場合は、 => ニ次式になるので、ニ次収束 n 回数目に 1 ビットの精度がでれば n + 1 回数目に 2 ビットの精度 n + 2 回数目に 4 ビットの精度 n + 3 回数目に 8 ビットの精度 # 倍々ゲーム !! # 三次より高次の収束もあるが、あまり意味がない.. 24 bit ( 単精度 ) から、52 bit ( 倍精度 ) にするには.. ? 二分法の場合 52 - 24 = 28 回の追加の反復が必要 ニュートン法の場合 1 回で、24 * 2 = 48 2 回で、48 * 2 = 96 2 回の追加の反復で十分 # ニュートン法では、初期値が 1 bit の精度があれば、その後は高速に収束する # \sqrt{a} の場合は、a に対して、適当な精度の、初期値 table を保持し、その結果を利用している ## table の大きさをふやせば、高速になるが、メモりを必要とするようになる、 ## 速度と記憶容量のトレードオフ == [ニュートン法の収束速度] # 実際に、ニュートン法の収束速度が二次収束であることを示す.. ## テキストの誤差の定義は、ちょっと違うので、訂正する 現在の値 x_k と、真の値 x_T の間の誤差 e_k を考えて、次のように x_k を 表現する x_k = e_k + x_T # e_k は、現在の値の x_k と真の値 x_T との誤差 すると、ニュートン法の式 x_{n+1} = x_n - \frac{f(x_n)}{f'(x_n)} に、この式を代入すれば、 e_{n+1} = e_n - \frac{f(x_n)}{f'(x_n)} となる。ここで、テーラ展開 f(x_k+h) = f(x_k)+\frac{h}{1!}f'(x_k)+\frac{h^2}{2!}f''(x_k) + .. を使って、 f(x_k) = f(x_T+e_k) = f(x_T)+\frac{h}{1!}f'(x_T)+\frac{h^2}{2!}f''(x_T) + .. f'(x_k) = f'(x_T+e_k) = f'(x_T)+\frac{h}{1!}f''(x_T)+\frac{h^2}{2!}f'''(x_T) + .. ここで、f(x_T)=0 に注意しながら、さっきの式に、これを代入すると.. e_{k+1}=e_k - \frac{e_k f' + e_k^2 f'' + .. }{f' + e_k f'' + .. } 分数部分の分母分子を、f' で割れば.. =e_k - \frac{e_k + e_k^2 \frac{f''}{f'} + .. }{1+e_k \frac{f''}{f'} + e_k^2 \frac{f'''}{f'} + .. } # Text p.231 は誤っている.. # 誤 \frac{f'''}{6f'}^3e_k # 正 \frac{f'''}{6f'}e_k^3 分母は、1+\frac{f''}{f'}e_k+.. の形なので、分母分子に 1-\frac{f''}{f'}e_k をかけると、分母は、ほぼ 1 になり、分子は、 = e_k - (e_k + e_k^2 \frac{f''}{f'} + .. )(1-e_k\frac{f''}{f'}) = \frac{f''}{2f'}e_k^2 + ... となるので、最上位の次数が、e_k^2 となり、二次収束であることがわかる # 重根の場合は、一次収束になってしまう !! # f'(x) が 0 に近いので、なかなか根に近付かない ( 図形的な性質から ) # # 式に代入しても、同様な結論が得られる ## 一般に、n 重根の場合は、更に収束が遅くなる !! == [二次のニュートン法] テーラ展開の三次以降を無視して、二次のニュートン公式を作ることができる # h の二次の式なので、公式は、二つ出てくる ( 二次方程式の根になる ) 式が複雑になって、色々面倒 そのわりには、それほど精度の旨味がない 100bit 程度だと、二次でも三次でも似た様なもの.. # 1000000 bit 位の精度がほしければ、意味があるが.. ## 福井先生御自身も、平方根の問題を二次のニュートン法を解いてみたが、あまり、芳しい結果がでなかった.. # 二次収束の証明が、前の本と今の本で違う。 # 収束の条件の話は、Text にあるが、ここでは説明しない。 ( あまり意味がない ) == 本日の演習 課題番号 06/07 課題 関数 f(x)=x-\cos{x} に対して f(x)=0 の根をニュートン法で求める。 ただし、初期値 x_0 は、 x_0=5 とする。 [ポイント] 誤差の振舞いが、以前の二分法とどのように違うをかを視る。 # 一旦、精度が 1 桁でる ( 誤差に、0 が現れる ) と、 # 後は、( 誤差の 0 の ) 桁数が倍々になることが解る。 ==