「SICP (問題1.*)」の編集履歴(バックアップ)一覧に戻る

SICP (問題1.*) - (2007/01/11 (木) 07:54:56) のソース

「第1章 手続きによる抽象の構築」で提示されている問題を解いています。

* 目次
#contents(fromhere)

* 問題
** 1.1
省略。

** 1.2
> (/ (+ 5 4 (- 2 (- 3 (+ 6 (/ 4/5))))) (* 3 (- 6 2) (- 2 7)))

** 1.3
少し悩んだ。最大値を返す(max a b c)と中間値を返す(mid a b c)を作れば良いのかと思って悩んだけれど、
> (define (func a b c)
>  (define (square x) (* x x))
>  (if (> a b)
>      (+ (square a) (square (if (> b c) b c)))
>      (+ (square b) (square (if (> a c) a c)))))
で良いみたい。

** 1.4
この表現は驚いた。ifが返すのは値だけでなく、+や-の演算子でも良いのか。これまでに無い考え方だったので、慣れるには時間がかかりそう。

** 1.5
手こずった。先に解答を書いておくと、
「正規順序ではpはpのまま展開しておき、いざ評価する段階では先に(= x 0)が評価され、pには触れないですむので0が返る。しかし作用的順序ではpを置換し始めてしまうため無限ループに陥る」
とのこと。(define (p) (p))という定義の意味が、けっきょく評価ができないものと理解することができなかったため、実際に実行するまで分からなかったことが理由。

** 1.6
問題1.5の補足で書いてあったことが1.6の答えだった。
> 特殊形式ifの評価規則は、解釈系が正規順序と作用的順序のどちらを
> 使うかに無関係とする:述語式を最初に評価し、その結果が帰結式と
> 代替式のいずれを評価するかを決める。
の中にあるとおり、「述語式の評価結果によって、帰結式と代替式どちらかの式は評価されない」という点が合成関数new-ifとの違いになっている。
new-ifでは帰結式・代替式(に相当する式)の双方を置換し切ってしまおうとするため、無限ループに陥ってしまう。
(ごめんなさい、何故そうなるかが厳密にはまだ理解できてないです)

** 1.7
> (define (good-enough2? guess x)
>  (< (abs (- guess (improve guess x))) 0.0001))

** 1.8
> (define (cbrt x)
>   (define (cbrt-iter guess)
>     (if (good-enough? guess)
>         guess
>         (cbrt-iter (improve guess))))
>   (define (improve guess)
>     (define (average a b c) (/ (+ a b c) 3))
>     (average guess guess (/ x (* guess guess))))
>   (define (good-enough? a)
>     (< (abs (- (* a a a) x)) 0.0001))
>   (cbrt-iter 1.0))

** 1.9
前者が再帰的プロセスで、後者が反復的プロセス。

** 1.10
+ (f n) → 2n
+ (g n) → 2^n
+ (h n) → ...
自分の数学力の無さを痛感した問題。
-(h 1)=2
-(h 2)=2^2
-(h 3)=2^2^2
-(h 4)=2^2^2^2
というルールは分かったのだけれど、これをどうやって数式に書いたら良いんだろう・・・

追記:数列で表せばいいのか。
-F_0 = 0, F_1 = 2
-F_n = (F_(n-1))^2

** 1.11
> (define (f_rec n)
>   (cond ((< n 3) n)
>         (else (+ (f_rec (- n 1))
>                  (* 2 (f_rec (- n 2)))
>                  (* 3 (f_rec (- n 3)))))))
> 
> (define (f_rep n)
>   (define (f_itr a b c rep)
>     (if (= rep 0)
>         a
>         (f_itr (+ a (* 2 b) (* 3 c)) a b
>                (- rep 1))))
>   (if (< n 3) n
>       (f_itr 2 1 0 (- n 2))))

反復式プロセスの方はイテレータとn<3以下の処理を別に分けてみたけれど、もっといい方法があるんだろうか。
(あ、f_itrの中で分岐してもいいわけか。でも少し汚いソースになってしまうかな。このままでいいや)

** 1.12
> (define (pascal row col)
>   (cond ((or (< row 1) (> col row) (< col 1)) -1)
>         ((= col 1) 1)
>         ((= col row) 1)
>         (else (+ (pascal (- row 1) (- col 1))
>                          (pascal (- row 1) col)))))

手続きに(行,列)を渡すと、対応する値を返してくれる。余裕がすこしあったので、エラー処理も入れてみた。

** 1.13
省略。

** 1.14
省略。

** 1.15
省略。

** 1.16
ヒントにある、「状態変数aを用意し、状態の移り変わりで積ab^nを一定となるよう定義する」の意味がよく分からず解答をネットで見ることに。
「偶数:b^n=(b^2)^(n/2)」、「奇数:b^n=b*b^(n-1)」と分解していく中で、そのときにaを交えた式ab^nがそのままの状態となるよう、式を変形していけばいいみたい。

-偶数 : {a,b,n} → {a, (b^2), (n/2)}
-奇数 : {a,b,n} → {ab, b, (n-1)}

というわけで、

> (define (fast-expt-rep b n)
>   (define (fast-itr b n a)
>     (define (square x) (* x x))
>     (cond ((= n 0) a)
>           ((even? n) (fast-itr (square b) (/ n 2) a))
>           (else (fast-itr b (- n 1) (* a b)))))
>   (fast-itr b n 1))

** 1.17
-偶数 : xy → (x/2)*2y
-奇数 : xy → y+(x-1)y

と変換できるので、それを利用したのが以下のとおり。

> (define (double x) (* x 2))
> (define (halve x) (/ x 2))
> (define (even? n) (= (remainder n 2) 0))
> 
> (define (multiply-rec x y)
>   (cond ((= x 1) y)
>         ((even? x) (multiply-rec (halve x) (double y)))
>         (else (+ y (multiply-rec (- x 1) y)))))
> (define (f a b) (multiply-rec a b))

** 1.18
少し手こずった。理由は、xyが状態の移り変わりでも一定となる値nを用意したときに、xyn=Cとなるような計算式を一生懸命に考えていたため。問題1.16の実装方針に影響されすぎた。

-偶数 : xy+n → (x/2)*2y + n
-奇数 : xy+n → (x-1)y + (n+y)

とすればよく、以下のプロセスを書いた。

> (define (multiply-rep x y)
>   (define (itr x y n)
>     (cond ((= x 0) n)
>           ((even? x) (itr (halve x) (double y) n))
>           (else (itr (- x 1) y (+ n y)))))
>   (itr x y 0))
> (define (g x y) (multiply-rep x y))

** 1.19
a ← bq+aq+ap, b ← bp+aqなので、それぞれのa,bにもう一度代入してあげる。

 a' ← (bp+aq)q+(bq+aq+ap)q+(bq+aq+ap)p
     = bpq+aq^2+bq^2+aq^2+apq+bpq+apq+ap^2
     = ap^2+2aq^2+2aqp+bq^2+2bpq
     = a(p^2+2q^2+2pq) + b(q^2+2pq)
 b' ← (bp+aq)p+(bq+aq+ap)q
     = bp^2+apq+bq^2+aq^2+apq
     = a(q^2+2pq) + b(p^2+q^2)

ここでp'=(p^2+q^2), q'=(q^2+2pq)と置き換えれば、ふたたびa' ← bq'+aq'+ap', b' ← bp'+aq'と表現できる。
(最初、この「置き換える」という視点が欠けていて先に進めなかった)

これにより、Fibonacci数を対数的ステップ数で求めるプロセスは次のようになる。
> (define (fib n)
>   (fib-itr 1 0 0 1 n))
> 
> (define (fib-itr a b p q count)
>   (define (even? x) (= (remainder x 2) 0))
>   (define (square x) (* x x))
>   (cond ((= count 0) b)
>         ((even? count)
>          (fib-itr a b
>                   (+ (square p) (square q))
>                   (+ (square q) (* 2 p q))
>                   (/ count 2)))
>         (else (fib-itr (+ (* b q) (* a q) (* a p))
>                        (+ (* b p) (* a q))
>                        p q (- count 1)))))

** 1.20
省略。
勉強にはDrSchemeを使用しているのだけれど、どうしたらトレースが取得・表示出来るようになるんだろう・・・。

追記:[言語]-[ティーチパックの追加...]でcalltrace.ssを追加すれば良いみたい。

** 1.21
199,1999,19999のうち、199と1999が素数で19999は素数ではない(7で割れる)。

** 1.22
まず、DrSchemeでは(runtime)という基本手続き(primitive)が存在しないので、代わりに(current-milliseconds)を使わないといけない。

> (define (runtime) (current-milliseconds)

そして奇整数を対象にtimed-prime-testを実行していく以下の手続きを作った。

> (define (search-for-primes from to)
>   (cond ((> from to) 0)
>         ((even? from) (search-for-primes (+ from 1) to))
>         (else ((timed-prime-test from)
>                (search-for-primes (+ from 2) to)))))

ところが、これだと最後の値を実行した後に「procedure application: expected procedure, given: #<void>; arguments were: 0」というエラーが出てしまう。
原因は2行目でfrom>toのときに0を返しているけれど、その返戻値を5行目で何も取り扱っていないから。頭の中で再帰を単なるgotoのようなものとして考えてしまっている自分に気がついた。
関数型言語なのに返戻値を無視するのはマズイ気がするのだけれど、どうしたらいいだろう・・・。
→結局、他のサイトの回答を参考に次のような手続きを書いた。

> (define (search-for-primes from end)
>   (define (itr n)
>     (timed-prime-test n)
>     (if (< n end) (itr (+ n 2))))
>   (if (even? from) (itr (+ from 1)) (itr from)))

でも、itrの中で「timed-prime-testを実行して、その次に終了条件を見て・・・」というように手続き的に書いてしまっているのがいやだなあ。ちなみに、DrSchemeでの実行速度が速すぎて、Θ(√n)の増加はよく分かりませんでした。

** 1.23
次の数を奇数にして返すnext手続きのみ作成。実行速度が速いので、比較してもあまり有為な結果が得られないし。

> (define (next n)
>   (if (even? n) (+ n 1) (+ n 2)))

** 1.24
これも1.23と同じ理由で、省略。

** 1.25
理屈は設問のAlyssaが言っていることで正しいので、悩んだ。
いろいろインターネットを漁り(こればっか)、計算の過程で、Alyssaの手続きでは値が非常に大きな値を持ってしまうことが分かった。それに対し、もともとの方法では各計算で剰余を求めることで値の増加を抑えつつ計算しているので、最終的な計算時間が短くてすむ。

「本当にそうなの?っていうか、p.28の計算って正しいの?」と思って同手続きをステップ実行してみたけれど、確かに計算は合っているし、値もmod未満の値に抑えられている。剰余の計算アルゴリズムとして覚えておいた方がいいみたい。

** 1.26
Louisの作った手続きは、square手続きでまとめていたexpmodを*によって展開した状態で記述している。これだと、squareの評価に入る前に1回だけ評価していたexpmodをどんなときも2回評価してしまうから。

当初「何でただ2乗するだけの手続きを別にsquareとして切り出していたんだろう?」と思っていたけれど、モジュール化のみならず評価の回数を少なくする効果があったんですね。
なるほど。

** 1.27
素数かどうかチェックするための手続きを以下のように作成。
> (define (is-prime? n)
>   (define (itr a)
>     (cond ((= a 0) #t)
>           (else (if (= (expmod a n n) a)
>                     (itr (- a 1)) #f))))
>   (if (itr (- n 1))
>       (display "I'ts prime.")
>       (display "I'ts not prime.")))
そして、これを元にCarmichael数(561,1105,1729,2465,2821,6601 -- p29の脚注から抽出)をチェック。

ふーん、確かに全部素数として判定されているなあ。

** 1.28
考えたけど、ギブアップ。解答から写経することにした。

> (define (mod a b) (remainder a b))
> 
> (define (expmod2 base exp m)
>   (cond ((= exp 0) 1)
>         ((even? exp)
>          (let ((tmp (mod (square (expmod2 base (/ exp 2) m)) m)))
>            (if (and (< 1 tmp) (< tmp m) (= (mod (square tmp) m) 1))
>                (display (list base exp tmp (mod (square tmp) m)))) tmp))
>         (else (mod (* base (expmod2 base (- exp 1) m)) m))))
> 
> (define (miller-rabin-test n)
>   (define (tt a i)
>     (cond ((= a n) i)
>           ((> (expmod2 a (- n 1) n) 1) (tt (+ a 1) (+ i 1)))
>           (else (tt (+ a 1) i))))
>   (tt 1 0))
> 
> (define (fermat-test2 n)
>   (define (try-it a)
>     (= (expmod2 a n n) a))
>     (try-it (+ 1 (random(- n 1)))))
> (define (fast-prime2? n times)
>   (cond ((= times 0) #t)
>         ((fermat-test2 n) (fast-prime2? n (- times 1)))
>         (else #f)))

??letって何?listは?
まだ一度も使ったことのない文法が出てきている。これらが実際にテキストで出てきたら、もう一度この問題を見直してみよう。

** 1.29
ここからEmacs+Gaucheを使用することに。でもデバッグはステップ実行のできるDrSchemeを使用してます。
> (define (simpson-integral f a b n)
>   (define h (/ (- b a) n))
>   (define (coef x)
>     (cond ((or (= x 0) (= x n)) 1)
>         ((even? x) 2)
>         (else 4)))
>   (define (itr k)
>     (cond ((< k 0) 0)
>         (else (+ (* (/ h 3)
>                     (* (coef k) (f (+ a (* k h)))))
>                  (itr (- k 1))))))
>   (itr n))
> 
> (define (f func a b n)
>   (if (even? n)
>       (simpson-integral func a b n)
>       #f))

この状態で、インターネットの解答例と比較。
しまった、ここではsum手続きを使わないといけないんだ。

> (define (simpson-integral-2 f a b n)
>   (define h (/ (- b a) n))
>   (define next inc)
>   (define (coef x)
>     (cond ((or (= x 0) (= x n)) 1)
>           ((even? x) 2)
>           (else 4)))
>   (define (term k) (* (coef k) (f (+ a (* k h)))))
>   (* (/ h 3) (sum term 0 next n)))

これでいいかな。
でも気になるのが、外にあるsum手続きの中から評価されている各種手続きで、なぜsimpson-integral-2内でのみ使われている値(a,h,nなど)が使用できているんだろう??
大事なことのような気がするけど、説明出来ない・・・。

** 1.30
> (define (sum term a next b)
>   (define (iter a result)
>     (if (> a b)
>         result
>         (iter (next a) (+ (term a) result))))
>   (iter a 0))

最初、a>bのとき(つまり、最後の計算の時)にresultではなく0を返してしまっていた。反復手続きの場合は最後に状態を返さないと、再帰的手続きとは違ってもう何も計算が行われずスタックトレースを遡るだけなので、これまでの計算がふいになってしまうことに注意しないと。

** 1.31
正月休みも終わって、昼は仕事をしているので進行ペースが落ちている。やばい?

> ; 指定した範囲の積を返す手続きproduct
> ;  1.反復的手続き版
> (define (product term a next b)
>   (define (iter a result)
>     (if (> a b)
>         result
>         (iter (next a) (* (term a) result))))
>   (iter a 1))
> ;  2.再帰的手続き版
> (define (product-rec term a next b)
>   (if (> a b)
>       1
>       (* (term a)
>          (product-rec term (next a) next b))))

で、これを使ったfactorial(n!のこと)は次のとおり。抽象化によって、前よりずっと簡単に手続きが定義できている。

> (define (factorial n)
>   (define (term x) x)
>   (product term 1 inc n))

で、次の問題で解答を見てしまった。なぜかというと、John Wallisが見つけたというπの近似式:

 π/4 = (2*4*4*6*6*8*8*...)/(3*3*5*5*7*7*...)

について、どうやったら(next a)に相当する手続きを書けるかが分からなかったため。
素直(というか手続き言語的に)に考えて、カウンタに内部状態を持たせて、1つおきに値を更新するような手続きを書こうとしていたけれど、ここまで学んだSchemeの知識では無理だった。

ギブアップして解答を見たら、式を(2*4)/(3*3), (4*6)/(5*5)というようにまとめて計算していた。

うーん、確かにそれで良いって言えば良いんだけど、これだと「第(2n-1)項までを計算せよ」という問題は解けないんじゃないの?
→ 調べたら、[[そもそもそういう公式>http://ja.wikipedia.org/wiki/円周率の歴史]]でした。なんだ、それなら問題ないか。

> (define (pi-wallis times)
>   (define (term n)
>     (/ (* (* 2.0 n) (* 2.0 (+ n 1.0)))
>        (square (+ (* 2.0 n) 1.0))))
>   (product term 1 inc times))
記事メニュー
目安箱バナー