SICP四日目～ 問題 1.7
昨日書いたブログが意外と閲覧されてしまってて、最終的に魔法陣問題解いてないので非常に申し訳ない気持ちです…。煽りタイトルよくないですね。
順列の実装方法はこのあたりからカンニングできるので、この紹介を以て許してほしいと思います。
問題1.7
問題文をよく読みます。「実際の計算機では, 算術演算は殆んどの場合, 限られた精度で実行される」とあるので計算誤差の問題のようです。
「非常に大きい数の場合, どのようにテストが失敗するか」
以下の平方根計算プログラムを動かしてみます。
(define (sqrtrec guess x) (if (good-enough? guess x) guess (sqrtrec (improve guess x) x))) (define (good-enough? guess x) ;中身をのぞいてみます (print guess) (< (abs (- (square guess) x)) 0.001) ) (define (improve guess x) (average guess (/ x guess)) ) (define (average x y) (/ (+ x y) 2)) (define (sqrt x) (sqrtrec 1.0 x) ) ;まず、十分に大きい数から (sqrt 50000000000000000000000)
これを実行すると、
1.0 2.5e22 1.25e22 6.25e21 3.125e21 1.5625e21 7.8125e20 3.90625e20 （中略） 2.2690447253825894e11 2.236307608276733e11 2.2360679903385666e11 2.2360679774997894e11 2.2360679774997894e11 2.2360679774997894e11 2.2360679774997894e11 2.2360679774997894e11 2.2360679774997894e11 （以下無限ループ）
となります。
なんとなく数値を見るに桁あふれっぽいですよね。
improve関数の中の問題かな？ あたりをつけて、実際にREPLでループしている数値をimproveと同じ形式で実行してみる。
gosh> (average 2.2360679774997894e11 (/ 50000000000000000000000.0 2.2360679774997894e11） 2.2360679774997894e11
おお、数値が変わっていませんな。では、averageを分解してみます。ついでにwindows標準の電卓（関数電卓モード）の結果と比べてみます。
gosh> (/ 50000000000000000000000 2.2360679774997894e11) 2.2360679774997897e11 #関数電卓だと 2.2360679774997899928183473374626e11 → 桁あふれてる gosh> (+ 2.2360679774997894e11 2.2360679774997897e11) 4.472135954999579e11 #関数電卓だと 4.4721359549995791e11→桁あふれてる gosh> (/ 4.472135954999579e11 2) 2.2360679774997894e11 #最初のguessと同じ数字になってしまった
というわけで、十分に大きい値でsqrtを求めた場合に正しい答えが得られない原因は桁あふれのようです。
「非常に小さい数の場合, どのようにテストが失敗するか」
平方根プログラムを以下のように修正します。
(define (sqrtrec guess x) (if (good-enough? guess x) guess (sqrtrec (improve guess x) x))) (define (good-enough? guess x) ;中身をのぞいてみます (print guess) (< (abs (- (square guess) x)) 0.001) ) (define (improve guess x) (average guess (/ x guess)) ) (define (average x y) (/ (+ x y) 2)) (define (sqrt x) (sqrtrec 1.0 x) ) ;面倒なので確認も同時にしてしまう (square (sqrt 0.00000001))
実行結果は
gosh> 1.0 0.500000000000005 0.2500000000000125 0.12500000000002626 0.06250000000005312 0.03125000000010656 #↓これは計算した平方根を二乗したもの。答えの確認。 9.7656250000666e-4
これは確認するまでもなく、good-enough? 関数で使用している許容値の問題だとわかります。許容値は平方根を求めたい値よりも小さい値である必要があります。
「good-enough?を実装するもう一つの戦略」
good-enough?関数を次のように修正します。
(define (good-enough? guess x) (print guess) (if (< (abs (- (improve guess x) guess)) 0.0001) #t (< (abs (- (square guess) x)) 0.001) ))
で、大きい数から実行してみます。
gosh>(sqrt 50000000000000000000000.0) gosh> 1.0 2.5e22 1.25e22 6.25e21 3.125e21 1.5625e21 7.8125e20 3.90625e20 1.953125e20 9.765625e19 4.8828125e19 2.44140625e19 1.220703125e19 6.103515625000002e18 3.051757812500005e18 1.5258789062500106e18 7.629394531250217e17 (中略） 2.2360679774997894e11 (ここで演算終了）
おおー、上手くいきました。
小さい方も実行してみます。
gosh>(sqrt 0.00000001) gosh> sqrt gosh> 1.0 0.500000005 0.25000001249999987 0.12500002624999892 0.06250005312499106 0.03125010656242753 0.03125010656242753
許容値に問題があるので、今回の修正の影響はありません。
