Posted at

SICP読書録 #7 1.2.4 指数計算 - 1.2.5 最大公約数

More than 1 year has passed since last update.

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

練習問題等はGaucheで実行。

前回はこちら


1.2.4 指数計算

指数計算を行う手続きexptを、以下のように定義する。


1_2_4_0.scm

(define (expt b n)

(if (= n 0)
1
(* b (expt b (- n 1)))))

この手続きの、計算時間の増加オーダーは$\Theta(n)$。空間も$\Theta(n)$。

では、反復プロセス版を定義してみよう。


1_2_4_1.scm

(define (iter-expt b n)

(define (iter n product)
(if (= n 0)
product
(iter (- n 1) (* product b))))
(iter n 1))

これだと、計算時間の増加オーダーは$\Theta(n)$。空間は$\Theta(1)$となる。

さて、ここで $b^n=(b^{n/2})^2$ を利用して、高速化を図ろう。

$b$を必ず1ずつ減らしていくなんてまどろっこしいだろ?


1_2_4_2.scm

(define (fast-expt b n)

(define (even? n)
(= (remainder n 2) 0))
(cond ((= n 0) 1)
((even? n) (square (fast-expt b (/ n 2))))
(else (* b (fast-expt b (- n 1))))))

この手続きの計算時間の増加オーダーは $\Theta(\log{n})$ となる。$\Theta(n)$よりずっと効率が良い。

空間については、これは再帰プロセスなので $\Theta(n)$ となる。


練習問題 1.16

fast-expt手続きをさらに反復プロセス化せよ。


ex1_16.scm

(define (fast-iter-expt b n)

(define (even? n)
(= (remainder n 2) 0))
(define (iter n b a)
(cond ((= n 0) a)
((even? n) (iter (/ n 2) (square b) a))
(else (iter (- n 1) b (* a b)))))
(iter n b 1))


練習問題 1.17

世界から掛け算が失われたので、自らの手で作り出さなければならなくなった。

せめて、計算時間のオーダーが $\Theta(\log{n})$ ぐらいのものを作ってみせなければ示しがつかない。

幸い、足し算 +、整数を倍にするdouble、整数を2で割るhalveは残されている。


ex1_17.scm

(define (fast-* a b) ;本記事内での命名規則を統一するためこの名前に

(define (even? n)
(= (remainder n 2) 0))
(define (double a) (+ a a))
(define (halve a) (/ a 2))
(cond ((= b 0) 0)
((even? b) (double (fast-* a (halve b))))
(else (+ a (fast-* a (- b 1))))))


練習問題 1.18

fast-*が定義されたことで、無事世界に掛け算が戻ってきたが、その計算空間の増加オーダーは$\Theta(n)$。これじゃ大きな数には使えない。

反復プロセス版を定義しよう。先ほどの、指数を求める手続きと同じ考え方でいける。


ex1_18.scm

(define (fast-iter-* a b)

(define (even? n)
(= (remainder n 2) 0))
(define (double a) (+ a a))
(define (halve a) (/ a 2))
(define (iter a b product)
(cond ((= b 0) product)
((even? b) (iter (double a) (halve b) product))
(else (iter a (- b 1) (+ product a)))))
(iter a b 0))


練習問題 1.19

$a,b$を以下のように変換することを考える。

\begin{align}

a &\leftarrow bq+aq+ap \\
b &\leftarrow bp+aq
\end{align}

この変換を$T_{pq}$と名付けよう。

$T_{pq}$を2回適用したらどうなる?

\begin{align}

a &\leftarrow (bp+aq)q+(bq+aq+ap)q+(bq+aq+ap)p\\
b &\leftarrow (bp+aq)p+(bq+aq+ap)q
\end{align}

この式を変換すると、次のような形になった。

\begin{align}

a &\leftarrow b(2pq+q^2)+a(p^2+q^2)+a(p^2+q^2)\\
b &\leftarrow b(p^2+q^2) + a(2pq+q^2)
\end{align}

ここで、以下のようにすると。

\begin{align}

p' &\leftarrow p^2+q^2\\
q' &\leftarrow 2pq+q^2
\end{align}

すてきなことに、こうなる。

\begin{align}

a &\leftarrow bq'+aq'+ap' \\
b &\leftarrow bp'+aq'
\end{align}

つまり、$T_{p'q'}$は、$T_{pq}$を2回適用したのと同じになるのだ。

ところで、$T_{pq}$にて、$(p=0,q=1)$とすると、以下のようになる。

\begin{align}

a &\leftarrow a+b \\
b &\leftarrow a
\end{align}

これは・・・どこかで見たような。そうだ。フィボナッチ数列を生成するための変換だ。

さっきの$T_{p'q'}$にこれを適用すれば、1ステップで2回のフィボナッチ変換を行える。計算時間の増加オーダーが$\Theta(\log{n})$となるような手続きを書けるぞ。


ex1_19.scm

(define (fast-fib n)

(define (even? n)
(= (remainder n 2) 0))
(define (iter a b p q n)
(cond ((= n 0) b)
((even? n)
(iter a
b
(+ (* p p) (* q q))
(+ (* 2 p q) (* q q))
(/ n 2)))
(else (iter (+ (* b q) (* a q) (* a p))
(+ (* b p) (* a q))
p
q
(- n 1)))))
(iter 1 0 0 1 n))

以下のように、普通に実装したバージョンも用意して。


ex1_19_org.scm

(define (fib n)

(define (fib-iter a b count)
(if (= count 0)
b
(fib-iter (+ a b) a (- count 1))))
(fib-iter 1 0 n))

性能比較しつつ、結果が一致するか確認してみる。


ex1_19_result.scm

gosh> (= (time (fib 100000)) (time (fast-fib 100000))

)
;(time (fib 100000))
; real 0.605
; user 1.440
; sys 0.070
;(time (fast-fib 100000))
; real 0.025
; user 0.020
; sys 0.000
#t


1.2.5 最大公約数

ユークリッドの互除法で、$a$と$b$の最大公約数(GCD)を求める。

$a$を$b$で割った余りが$r$とすると、$a$と$b$のGCDと$b$と$r$のGCDは等しい。

具体例をみたほうがわかりやすい。

\begin{align}

GCD(206,40) &= GCD(40, 6)\\
&= GCD(6, 4) \\
&= GCD(4, 2) \\
&= GCD(2, 0) \\
&= 2
\end{align}

実装してみる。


1_2_5.scm

(define (gcd a b)

(if (= b 0)
a
(gcd b (remainder a b))))

このアルゴリズムの計算時間(ステップ数)の増加オーダーは$\Theta(\log{n})$となる。


練習問題 1.20

(gcd 206 40)を正規順序評価に基づいて展開すると。remainder手続きは何回呼ばれるだろうか?

ifは特殊形式なので、まず条件句を評価し、その結果に応じてthen句とelse句のいずれかが評価されるという流れ。


ex1_20_1.scm

> (gcd 206 40)

> (if (= 40 0)
206
(gcd 40 (remainder 206 40)))

> (gcd 40 (remainder 206 40))
> (if (= (remainder 206 40) 0)
40
(gcd (remainder 206 40) (remainder 40 (remainder 206 40))))
; (remainder 206 40) -> 6 : 1回
> (if (= 6 0)
40
(gcd (remainder 206 40) (remainder 40 (remainder 206 40))))

> (gcd (remainder 206 40) (remainder 40 (remainder 206 40)))
> (if (= (remainder 40 (remainder 206 40)) 0)
(remainder 206 40)
(gcd (remainder 40 (remainder 206 40)) (remainder (remainder 206 40) (remainder 40 (remainder 206 40))))))
; (remainder 206 40) -> 6 : 1回
; (remainder 40 6) -> 4 : 1回
> (if (= 4 0)
(remainder 206 40)
(gcd (remainder 40 (remainder 206 40)) (remainder (remainder 206 40) (remainder 40 (remainder 206 40))))))

> (gcd (remainder 40 (remainder 206 40)) (remainder (remainder 206 40) (remainder 40 (remainder 206 40))))
> (if (= (remainder (remainder 206 40) (remainder 40 (remainder 206 40))) 0)
(remainder 40 (remainder 206 40))
(gcd (remainder (remainder 206 40) (remainder 40 (remainder 206 40))) (remainder (remainder 40 (remainder 206 40)) (remainder (remainder 206 40) (remainder 40 (remainder 206 40))))))
; (remainder 206 40) -> 6 : 2回
; (remainder 40 6) -> 4 : 1回
; (remainder 6 4) -> 2 : 1回
> (if (= 2 0)
(remainder 40 (remainder 206 40))
(gcd (remainder (remainder 206 40) (remainder 40 (remainder 206 40))) (remainder (remainder 40 (remainder 206 40)) (remainder (remainder 206 40) (remainder 40 (remainder 206 40))))))

> (gcd (remainder (remainder 206 40) (remainder 40 (remainder 206 40))) (remainder (remainder 40 (remainder 206 40)) (remainder (remainder 206 40) (remainder 40 (remainder 206 40)))))
> (if (= (remainder (remainder 40 (remainder 206 40)) (remainder (remainder 206 40) (remainder 40 (remainder 206 40)))) 0)
(remainder (remainder 206 40) (remainder 40 (remainder 206 40)))
(gcd (remainder (remainder 40 (remainder 206 40)) (remainder (remainder 206 40) (remainder 40 (remainder 206 40)))) (remainder (remainder (remainder 206 40) (remainder 40 (remainder 206 40))) (remainder (remainder 40 (remainder 206 40)) (remainder (remainder 206 40) (remainder 40 (remainder 206 40)))))))
; (remainder 206 40) -> 6 : 3回
; (remainder 40 6) -> 4 : 2回
; (remainder 6 4) -> 2 : 1回
; (remainder 4 2) -> 0 : 1回
> (if (= 0 0)
(remainder (remainder 206 40) (remainder 40 (remainder 206 40)))
(gcd (remainder (remainder 40 (remainder 206 40)) (remainder (remainder 206 40) (remainder 40 (remainder 206 40)))) (remainder (remainder (remainder 206 40) (remainder 40 (remainder 206 40))) (remainder (remainder 40 (remainder 206 40)) (remainder (remainder 206 40) (remainder 40 (remainder 206 40)))))))

> (remainder (remainder 206 40) (remainder 40 (remainder 206 40)))
; (remainder 206 40) -> 6 : 2回
; (remainder 40 6) -> 4 : 1回
; (remainder 6 4) -> 2 : 1回

> 2


計18回。

適用順序評価ならばどうか。


ex1_20.scm

> (gcd 206 40)

> (if (= 40 0)
206
(gcd 40 (remainder 206 40)))

> (gcd 40 (remainder 206 40))
; (remainder 206 40) -> 6 : 1回
> (gcd 40 6)
> (if (= 6 0)
40
(gcd 6 (remainder 40 6)))

> (gcd 6 (remainder 40 6))
; (remainder 40 6) -> 4 : 1回
> (gcd 6 4)
> (if (= 4 0)
6
(gcd 4 (remainder 6 4)))

> (gcd 4 (remainder 6 4))
; (remainder 6 4) -> 2 : 1回
> (gcd 4 2)
> (if (= 2 0)
4
(gcd 2 (remainder 4 2)))

> (gcd 2 (remainder 4 2))
; (remainder 4 2) -> 0 : 1回
> (gcd 2 0)
> (if (= 0 0)
2
(gcd 0 (remainder 2 0)))

> 2


計4回。