Posted at

SICP読書録 #13 1.3.3 汎用手法としての手続き - 練習問題1.36まで

More than 1 year has passed since last update.

hiroshi-manabe様の日本語訳を使わせてもらいます。

練習問題等はGaucheで実行。

前回はこちら


1.3.3 汎用手法としての手続き


区間二分法によって方程式の根を求める

関数の零点、つまり $f(x) = 0$ となる $x$ 値を見つける手続きを書こう。

$f(a) < 0 < f(b)$ となる 範囲 $a, b$ を与え、繰り返しによってそれを狭めていくのだ。


_.scm

(define (search f neg-point pos-point)

(let ((midpoint (average neg-point pos-point)))
(if (close-enough? neg-point pos-point)
midpoint
(let ((test-value (f midpoint)))
(cond ((positive? test-value)
(search f neg-point midpoint))
((negative? test-value)
(search f midpoint pos-point))
(else midpoint))))))

(define (average x y) (/ (+ x y) 2))
(define (close-enough? x y) (< (abs (- x y)) 0.001))


与えた $a, b$ が有効かを検証した上で検索を行う手続きも用意しよう。


_.scm

(define (half-interval-method f a b)

(let ((a-value (f a))
(b-value (f b)))
(cond ((and (negative? a-value) (positive? b-value))
(search f a b))
((and (negative? b-value) (positive? a-value))
(search f b a))
(else
(error "Values are not of opposite sign" a b)))))

以下、実行結果。


_.scm

gosh> (half-interval-method sin 2.0 4.0)

3.14111328125
gosh> (half-interval-method (lambda (x) (- (* x x x) (* 2 x) 3))
1.0 2.0)
1.89306640625


関数の不動点を求める

不動点とは $f(x) = x$ となる $x$ のこと。

これを求めるには、$f(f(f(x)))$ のように、$f$ を収束するまで繰り返し適用していく。


_.scm

(define tolerance 0.00001)

(define (fixed-point f first-guess)
(define (close-enough? v1 v2)
(< (abs (- v1 v2))
tolerance))
(define (try guess)
(let ((next (f guess)))
(if (close-enough? guess next)
next
(try next))))
(try first-guess))

以下、実行結果。


_.scm

gosh> (fixed-point cos 1.0)

0.7390822985224023
gosh> (fixed-point (lambda (y) (+ (sin y) (cos y)))
1.0)
1.2587315962971173

この方法、うまく収束しない場合もあるのだ。。。

例えば、この手続きを使って平方根を求めよう。

$y^2 = x$ となる $y$ を求めればよい。つまり $f(y) = x / y$ の不動点を求める。


_.scm

(define (sqrt x)

(fixed-point (lambda (y) (/ x y))
1.0))

これを実行すると、無限ループに陥る。

$y_3 = x / y_2 = x / (x / y_1) = y_1$ となり、$y_1, y_2$ の間を行ったり来たりするのだ。

これを防ぐために、次の推定値を現在の推定値との平均としてみよう。

これを、平均緩和法という。


_.scm

(define (sqrt x)

(fixed-point (lambda (y) (average y (/ x y)))
1.0))

これならうまくいく。


_.scm

gosh> (sqrt 25)

5.0
gosh> (sqrt 2)
1.4142135623746899


練習問題 1.35

黄金比 $\varphi^2=\varphi+1$ を求めよう。

$f(x) = 1 + 1 / x$ の不動点を求めればいい。


_.scm

(define (golden-ratio)

(fixed-point (lambda (x) (+ 1 (/ 1 x)))
1.0))

実行結果。


_.scm

gosh> (golden-ratio)

1.6180327868852458
gosh> (expt (golden-ratio) 2)
2.6180300994356354


練習問題 1.36

推定値がどのように更新されていくかが表示されるよう、 fixed-point 手続きを修正する。


_.scm

(define tolerance 0.00001)

(define (fixed-point f first-guess)
(define (close-enough? v1 v2)
(< (abs (- v1 v2))
tolerance))
(define (try guess)
(format #t "guess: ~D\n" guess)
(let ((next (f guess)))
(if (close-enough? guess next)
next
(try next))))
(try first-guess))

$x^x = 1000$ の解を求めよう。$f(x)=log(1000)/log(x)$ の不動点を求めればいい。

ちなみにここで、初期値を1にするのはまずい。$log(1)=0$ なのでゼロでの割り算が行われてしまう。


_.scm

gosh> (fixed-point (lambda (x) (/ (log 1000) (log x))) 2.0)

guess: 2.0
guess: 9.965784284662087
guess: 3.004472209841214
guess: 6.279195757507157
guess: 3.759850702401539
guess: 5.215843784925895
guess: 4.182207192401397
guess: 4.8277650983445906
guess: 4.387593384662677
guess: 4.671250085763899
guess: 4.481403616895052
guess: 4.6053657460929
guess: 4.5230849678718865
guess: 4.577114682047341
guess: 4.541382480151454
guess: 4.564903245230833
guess: 4.549372679303342
guess: 4.559606491913287
guess: 4.552853875788271
guess: 4.557305529748263
guess: 4.554369064436181
guess: 4.556305311532999
guess: 4.555028263573554
guess: 4.555870396702851
guess: 4.555315001192079
guess: 4.5556812635433275
guess: 4.555439715736846
guess: 4.555599009998291
guess: 4.555493957531389
guess: 4.555563237292884
guess: 4.555517548417651
guess: 4.555547679306398
guess: 4.555527808516254
guess: 4.555540912917957
4.555532270803653

続いて、平均緩和法を使って実行してみる。


_.scm

gosh> (fixed-point (lambda (x) (average x (/ (log 1000) (log x)))) 2.0)

guess: 2.0
guess: 5.9828921423310435
guess: 4.922168721308343
guess: 4.628224318195455
guess: 4.568346513136242
guess: 4.5577305909237005
guess: 4.555909809045131
guess: 4.555599411610624
guess: 4.5555465521473675
4.555537551999825

なぜだか速くなった。

それぞれ検算してみよう。


_.scm

gosh> (expt 4.555532270803653 4.555532270803653)

999.9913579312362
gosh>(expt 4.555537551999825 4.555537551999825)
1000.0046472054871