問題3.52-3.61

ストリームの続きです。

3.52

accumを少し変更して、sumの値変化時に印字するようにします。

(define sum 0)
(define (accum x)
  (set! sum (+ x sum)) (display sum) (display "  accum") (newline) sum)

> (define seq (stream-map accum (stream-enumerate-interval 1 20)))
1 -- accum
> (define y (stream-filter even? seq))
3 -- accum
6 -- accum
> (define z (stream-filter (lambda (x) (= (remainder x 5) 0)) seq))
10 -- accum
> (stream-ref y 7)
15 -- accum
21 -- accum
28 -- accum
36 -- accum
45 -- accum
55 -- accum
66 -- accum
78 -- accum
91 -- accum
105 -- accum
120 -- accum
136 -- accum
136
> (display-stream z)
10
15
45
55
105
120
153 -- accum
171 -- accum
190 -- accum
190
210 -- accum
210
done

streamは遅延評価なので、(stream-ref y 7)ではsumは1+...+20=20*21/2=210ではなく、even?であるような7番目(0から数えて)の要素136が見つかった時点で値が返ります。
(display-stream z)はseqの5の剰余が0のものを全て印字します。
accumの印字を見れば分かる通り、sumの値の変更はseqの要素数と同じ20回しかされていません。これはdelayのメモ化がされていることを意味します。

次にメモ化のされていない場合を見てみます。昨日のエントリで作ったdelayとforceを使います。

> (define seq (stream-map accum (stream-enumerate-interval 1 20)))
1 -- accum
> (define y (stream-filter even? seq))
3 -- accum
6 -- accum
> (define z (stream-filter (lambda (x) (= (remainder x 5) 0)) seq))
8 -- accum
11 -- accum
15 -- accum
> (stream-ref y 7)
19 -- accum
24 -- accum
30 -- accum
37 -- accum
...()
145 -- accum
162 -- accum
162
> (display-stream z)
15
167 -- accum
173 -- accum
180 -- accum
180
188 -- accum
...()
305 -- accum
305
323 -- accum
342 -- accum
362 -- accum
done

メモ化がされていないので必要以上にsumの値を変更しています。

最後にdelayを特殊形式として定義せず、メモ化もしない場合を考えます。

> (define seq (stream-map accum (stream-enumerate-interval 1 20)))
1 -- accum
3 -- accum
6 -- accum
...()
190 -- accum
210 -- accum
> (define y (stream-filter even? seq))
> (define z (stream-filter (lambda (x) (= (remainder x 5) 0)) seq))
> (stream-ref y 7)
136
> (display-stream z)
10
15
45
55
105
120
190
210
done
> (define seq (stream-map accum (stream-enumerate-interval 1 20)))
211 -- accum
213 -- accum
216 -- accum
...()
400 -- accum
420 -- accum

ということで、seqに先に値を計算してから入れるようになり、遅延評価ではなくなります。

3.53

たぶん2のべき乗の並び(1, 2, 4, 8, ...)になる。

> (define s (cons-stream 1 (add-stream s s)))
> s
(1 . #<struct:promise>)
> (stream-ref s 1)
2
> (stream-ref s 2)
4
> (stream-ref s 3)
8
> (stream-ref s 10)
1024

なってました。

3.54

(define (mul-streams s1 s2) (stream-map * s1 s2))
(define factorials (cons-stream 1 
                                (mul-streams factorials
                                             (stream-cdr integers))))

> factorials
(1 . #<struct:promise>)
> (stream-ref factorials 1)
2
> (stream-ref factorials 2)
6
> (stream-ref factorials 3)
24
> (stream-ref factorials 10)
39916800

3.55

(define (partial-sums s)
  (define ps (cons-stream (stream-car s)
                          (add-stream ps (stream-cdr s))))
  ps)

> (define s (partial-sums integers))
> (stream-ref s 1)
3
> (stream-ref s 2)
6
> (stream-ref s 3)
10
> (stream-ref s 9)
55
> (define int (partial-sums ones))
> (stream-ref int 0)
1
> (stream-ref int 1)
2
> (stream-ref int 9)
10

3.56

確認用にstreamの先頭からn個の要素をリストにして返すdisp-stream-nを作る。
Sは以下のようになる。

(define (disp-stream-n s n)
  (if (= n 0)
      '()
      (cons (stream-car s)
            (disp-stream-n (stream-cdr s) (- n 1)))))

> (define S (cons-stream 1
                         (merge (scale-stream S 2)
                                (merge (scale-stream S 3)
                                       (scale-stream S 5)))))
> (disp-stream-n S 30)
(1 2 3 4 5 6 8 9 10 12 15 16 18 20 24 25 27 30 32 36 40 45 48 50 54 60 64 72 75 80)

3.57

メモ化ありの時。

> (time (stream-ref fibs 100))
cpu time: 0 real time: 1 gc time: 0
354224848179261915075

メモ化なしの時。

> (time (stream-ref fibs 20))
cpu time: 260 real time: 261 gc time: 0
6765
> (time (stream-ref fibs 21))
cpu time: 425 real time: 426 gc time: 0
10946
> (time (stream-ref fibs 22))
cpu time: 836 real time: 837 gc time: 152
17711
> (time (stream-ref fibs 23))
cpu time: 1108 real time: 1106 gc time: 0
28657
> (time (stream-ref fibs 24))
cpu time: 1916 real time: 1915 gc time: 112
46368

遅いです。fibs,add-streamsの定義は以下のようになってます。

(define fibs (cons-stream 0
                          (cons-stream 1
                                       (add-streams fibs
                                                    (stream-cdr fibs)))))
(define (add-streams s1 s2) (stream-map + s1 s2))

メモ化delayでは、一回forceされたものを覚えておくので1ステップでアクセスすることができます。なのでadd-streamsに渡されているストリームがfibs自身であることから、第n項のfibonacci数の計算に必要な加算の回数はn-1回となります。
一方、非メモ化delayではforceしたものを覚えておく事が出来ないので、第n項のfibonacci数の計算に必要な加算回数をa(n) (n≧2)としたとき、

a(0) = 0
a(1) = 0
a(n) = a(n-1) + a(n-2) + 1

となります。
これが指数関数的に増大することを確認したいのですが、上の漸化式をそのまま実装すると結局同じくらい加算回数が増えるので確認できません。そこでメモ化手続きa_nを使って確認します。

(define (make-table) (list '**table**))
(define (lookup key table)
  (let ((record (assoc key (cdr table))))
    (if record
        (cdr record)
        false)))
(define (insert! key value table)
  (let ((record (assoc key (cdr table))))
    (if record
        (set-cdr! record value)
        (set-cdr! table (cons (cons key value) (cdr table))))))
(define (assoc key records)
  (cond ((null? records) false)
        ((equal? key (caar records)) (car records))
        (else (assoc key (cdr records)))))
(define (memoize f)
  (let ((table (make-table)))
    (lambda (x)
      (let ((previously-computed-result (lookup x table)))
        (or previously-computed-result
            (let ((result (f x)))
              (insert! x result table)
              result))))))
(define a_n
  (memoize
   (let ((table (make-table)))
     (lambda (n)
       (cond ((= n 0) 0)
             ((= n 1) 0)
             (else (+ (a_n (- n 1)) (a_n (- n 2)) 1)))))))

> (a_n 10)
88
> (a_n 20)
10945
> (a_n 30)
1346268
> (a_n 50)
20365011073
> (a_n 100)
573147844013817084100

ということで、加算回数が爆発的に増えることが分かります。
ところで、a_nは階差をとるとfibonacci数列になるので、一般項も求められます。

 a_n=a_0+\sum_{i=0}^{n-1}{fib(i)}

 fib(n) = \frac{1}{\sqrt{5}} ((\frac{1+\sqrt{5}}{2})^n - (\frac{1-\sqrt{5}}{2})^n)

この式からも、a_nが爆発的に増えることが分かります。

3.58

(define (expand num den radix)
  (cons-stream (quotient (* num radix) den)
               (expand (remainder (* num radix) den) den radix)))

expandが返すストリームは、第一要素が小数点第一位、第二要素が小数点第二位、...となるnum/denのradix進数での値。ただし、0 < num/den < 1でないとexpandの挙動は保証されない感じ。
手続きがちょうど筆算の手順みたいになっていておもしろいです。

10進数

> (define s (expand 1 3 10))
> (disp-stream-n s 10)
(3 3 3 3 3 3 3 3 3 3)            => 0.33333333...
> (define s (expand 11 12 10))
> (disp-stream-n s 10)
(9 1 6 6 6 6 6 6 6 6)            => 0.91666666...
> (define s (expand 1 7 10))

> (disp-stream-n s 20)
(1 4 2 8 5 7 1 4 2 8 5 7 1 4 2 8 5 7 1 4)  => 0.142837142857......
> (define s (expand 3 8 10))
> (disp-stream-n s 20)
(3 7 5 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0)  => 0.375000000000......

2進数

> (define s (expand 5 10 2))
> (disp-stream-n s 10)
(1 0 0 0 0 0 0 0 0 0)            => 0.10000000... = 1/2
> (define s (expand 7 8 2))
> (disp-stream-n s 10)
(1 1 1 0 0 0 0 0 0 0)            => 0.11100000... = 1/2 + 1/4 + 1/8
> (define s (expand 9 16 2))
> (disp-stream-n s 10)
(1 0 0 1 0 0 0 0 0 0)            => 0.1001000000... = 1/2 + 1/16

3.59

a.
(define (integrate-series s)
  (define frac-integers (stream-map (lambda (x) (/ 1 x)) integers))
  (mul-streams s frac-integers))

> (disp-stream-n (integrate-series integers) 20)
(1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1)

integers便利。

b.
(define cosine-series
  (cons-stream 1 (stream-map (lambda (x) (- 0 x)) (integrate-series sine-series))))
(define sine-series
  (cons-stream 0 (integrate-series cosine-series)))

> (disp-stream-n cosine-series 10)
(1 0 -1/2 0 1/24 0 -1/720 0 1/40320 0)
> (disp-stream-n sine-series 10)
(0 1 0 -1/6 0 1/120 0 -1/5040 0 1/362880)

3.60

べき級数の乗算は以下の式で表されます。
(a_0+a_1x+...+a_nx^n+...)(b_0+b_1x+...+b_mx^m+...)
\hspace{30pt}=a_0b_0+(a_0b_1+a_1b_0)x+(a_0b_2+a_1b_1+a_2b_0)x^2+...
\hspace{140pt}+(a_0b_{n+m}+a_1b_{n+m-1}+...+a_{n+m}b_0)x^{n+m}+...
これをそのまま実装すればいいです。

(define (mul-series s1 s2)
  (cons-stream (* (stream-car s1) (stream-car s2))
               (add-streams (scale-stream (stream-cdr s2) (stream-car s1))
                            (mul-series (stream-cdr s1) s2))))

> (define s
    (add-streams (mul-series sine-series sine-series)
                 (mul-series cosine-series cosine-series)))
> (disp-stream-n s 100)
(1
 0
 0
 0
 ...()
 0
 0)

ちゃんと\sin^2(x) + \cos^2(x) = 1になってます。

3.61

問題文に書いてある通りに実装します。

(define (invert-unit-series s)
  (cons-stream 1
               (mul-series (stream-cdr s)
                           (stream-map (lambda (x) (- 0 x)) (invert-unit-series s)))))

> (define s (invert-unit-series ones))
> (disp-stream-n s 20)
(1 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0)
> (define s (invert-unit-series integers))
> (disp-stream-n s 20)
(1 -2 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0)

ということで、ちゃんとS*X=1を満たすようになっています。

残り文字数が足りないので今日はここまで。