タグ : SICP

SICP Reading #19 [2.2.3] 公認インターフェースとしての並び

春休みはSICPガッツリ進めるぜ!と意気込んでいたけど、微妙に時間を割けてない感じで不甲斐ない。。
とりあえず溜まってる更新をサラッと済ませる。

全てはリスト処理

まず例として以下のコード実装を

木を引数にとり、奇数である葉の二乗を返す
(define (sum-odd-squares tree)
  (cond ((null? tree) 0)
        ((not (pair? tree))
         (if (odd? tree) (squares tree) 0))
        (else (+ (sum-odd-squares (car tree))
                 (sum-odd-squares (cdr tree))))))
引数n以下の偶数のFibonacci数のリストを返す
(define (even-fibs n)
  (defien (next k)
    (if (> k n)
        '()
        (let ((f (fib k)))
          (if (even? f)
              (cons f (next (+ k 1)))
              (next (+ k 1)))))))

2つの手続きの実装を見てみると、当然全く違う処理の流れに見える。
(どちらも再帰だけど、木の処理は二重再帰だしね)

だけど、どちらも抽象的に記述すると似通っていて、かつどちらも

  • enumerate (数え上げ)
  • accumulate (リスト演算)
  • map (リスト変換)
  • filter (リストのフィルタ)

という高階な処理で記述できるよね、という感じ。

  • 木の葉を数え上げる[enumerate]
  • 奇数でフィルタ[filter]
  • 選ばれた要素を二乗[map]
  • 要素を0初期値で加算[accumulate]

及び

  • 整数を0からnまで数え上げる[enumerate]
  • それぞれの要素でFibonacci数を計算[map]
  • 偶数でフィルタ[filter]
  • 要素を空リスト初期値でcons[accumulate]

よって

処理なんて全部(?)こんな感じのリスト処理に落とせるじゃん!明白じゃん!
“リスト”というデータ構造の処理に落とし込むことで、プログラミングの各ステージが明瞭だし、独立性が高くなるね(っていうことかな?)。

問題2.33

map, append, lengthの実装穴埋め問題

(define (accumulate op initial sequence)
  (if (null? sequence)
      initial
      (op (car sequence)
          (accumulate op initial (cdr sequence)))))

(define (map p sequence)
  (accumulate (lambda (x y) (cons (p x) y)) '() sequence))

(map (lambda (x) (* x x)) '(1 2 3 4 5 6 7 8 9 10))
;gosh> (1 4 9 16 25 36 49 64 81 100)

(define (append seq1 seq2)
  (accumulate cons seq2 seq1))

(append '(1 2 3) '(4 5))
;gosh> (1 2 3 4 5)

(define (length sequence)
  (accumulate (lambda (x y) (+ y 1)) 0 sequence))

(length '(1 2 3 4 4 890))
;gosh> 6

問題2.34

hornerの方法を使って、xの多項式の演算をaccumulateで実装する問題

(add-load-path ".")
(load "q2-33")

(define (horner-eval x coefficient-sequence)
  (accumulate (lambda (this-coeff higher-terms)
                (+ this-coeff (* higher-terms x)))
              0
              coefficient-sequence))

(horner-eval 2 '(1 3 0 5 0 1))
;gosh> 79
;; 実質こういうこと
;; (let1 x 2 (+
;;            1
;;            (* 3 x)
;;            (* 5 (expt x 3))
;;            (expt x 5)))

問題2.35

2.2.2節のcount-leavesをaccumulateで。
提示された条件では、mapを使ってたけど最初に思い浮かんだのは以下な実装

(add-load-path ".")
(load "q2-33")

;;2.2.2の実装は
(define (count-leaves x)
  (cond ((null? x) 0)
        ((not (pair? x)) 1)
        (else (+ (count-leaves (car x))
                 (count-leaves (cdr x))))))

;;accumulate版. map使わない方が自然に思いついた
(define (count-leaves t)
  (accumulate (lambda (x y)
                (+ (if (pair? x) (count-leaves x) 1)
                   y))
              0
              t))

(count-leaves '(1 2 (3 4) (5 6) (7 8 (9 1))))
;gosh> 10

一応、指定通りmapを使ってみた版も

;;map使う場合(提示された条件)
(define (count-leaves t)
  (accumulate (lambda (x y) (+ x y))
              0
              (map (lambda (x) (if (pair? x)
                                   (count-leaves x)
                                   1)) t)))

accumulateに渡すのはリスト!!
ってのが自然だと思うかから、mapの中で再帰して最終的に各部分木の葉の数のリストに変換して一気に処理するほうが見通しがいいかもしれないと思った。

問題2.36

accumulateのmulti-list対応版

(add-load-path ".")
(load "q2-33.scm")

(define (accumulate-n op init seqs)
  (if (null? (car seqs))
      '()
      (cons (accumulate op init (map
                                 (lambda (x) (car x))
                                 seqs))
            (accumulate-n op init (map (lambda (x) (cdr x))
                                       seqs)))))

テストしてみないと不安..

(add-load-path ".")
(load "q2-37.scm")

(use gauche.test)
(test-start "accumulate-n")

(define s '((1 2 3) (4 5 6) (7 8 9) (10 11 12)))

(print "s = " s)
(test* "(accumulate-n + 0 s)"
       '(22 26 30)
       (accumulate-n + 0 s))
(test-end)
SHINYA% gosh q2-36.test.scm
Testing accumulate-n ...
s = ((1 2 3) (4 5 6) (7 8 9) (10 11 12))
test (accumulate-n + 0 s), expects (22 26 30) ==> ok
passed.

(テスト少な!!)まぁOK?

問題2.37

ベクトル演算。 ベクトルなんて所詮リストのリストだよね!
小さい手続きから実装していくと、意外に簡単に実装できた。 C++のオブジェクトで実装して結構苦戦した記憶が。。(まぁそっちでは逆行列とかも導出してたわけだけど)

(add-load-path ".")
(load "q2-36")

(define map (with-module gauche map))
(define v1 '(1 2 3))
(define v2 '(4 5 6))
(define m  '((1 2 3) (4 5 6) (7 8 9)))

(define (dot-product v w)
  (accumulate + 0 (map * v w)))

(dot-product v1 v2)
;gosh> 32

(define (matrix-*-vector m v)
  (map (lambda (x) (dot-product x v)) m))

(matrix-*-vector m v1)
;gosh> (14 32 50)

(define (transpose mat)
  (accumulate-n cons '() mat))

m
;gosh> ((1 2 3) (4 5 6) (7 8 9))
(transpose m)
;gosh> ((1 4 7) (2 5 8) (3 6 9))

(define (matrix-*-matrix m n)
  (let ((cols (transpose n)))
    (map (lambda (x)
           (map (lambda (y) (dot-product x y))
                cols)) m)))

(define m1 '((1 2) (3 4)))

(matrix-*-matrix m1 m1)
;((7 10) (15 22))

(matrix-*-matrix '((1 2 3)
                   (4 5 6)
                   (7 8 9))
                 '((1 0 0)
                   (0 1 0)
                   (0 0 1)))
;基本行列Eとの掛け算. E*A = A*E = A
;gosh> ((1 2 3) (4 5 6) (7 8 9))

問題2.34

accumulateは要素同士の演算を右から左に行なう。逆に演算したい場合は
再帰で評価を保留しておくべし。

(define (fold-right op initial sequence)
  (if (null? sequence)
      initial
      (op (car sequence)
          (fold-right op initial (cdr sequence)))))

(define (fold-left op initial sequence)
  (define (iter result rest)
    (if (null? rest)
        result
        (iter (op result (car rest))
              (cdr rest))))
  (iter initial sequence))

この二つのaccumulate実装で、結果が同じとなるような演算opの定義についての問題。
「op: 交換則がなりたつ演算の場合(加算,乗算) fold-right/leftの結果は同値」かな。

;;成り立つ
(fold-right + 0 (list 1 2 3 4 5 6))
;gosh> 21
(fold-left  + 0 (list 1 2 3 4 5 6))
;gosh> 21
(fold-right * 1 (list 1 2 3 4 5 6))
;gosh> 720
(fold-left  * 1 (list 1 2 3 4 5 6))
;gosh> 720

;;成り立たない
(fold-right / 1 (list 1 2 3))
;(/ 1 (/ 2 (/ 3 1))) -> 3/2
(fold-left  / 1 (list 1 2 3))
;(/ (/ (/ 1 1) 2) 3) -> 1/6
(fold-right list '() (list 1 2 3))
;gosh> (1 (2 (3 ())))
(fold-left list '() (list 1 2 3))
;gosh> (((() 1) 2) 3)

問題2.39

実装の穴埋め問題

(add-load-path ".")
(load "q2-38.scm")
;; with fold-right
(define (reverse-fr sequence)
  (fold-right (lambda (x y)
                (append y (list x))) '() sequence))
;(append (append (append '() '(3)) '(2)) '(1))

;; with fold-left
(define (reverse-fl sequence)
  (fold-left (lambda (x y)
               (cons y x)) '() sequence))
;(cons 3 (cons 2 (cons 1 '())))

appendはリストの末尾までたどってポインタを入れ替えるからO(n), consはセルで包むだけなのでO(1) (多分)。
かつfold-rightは再帰ので(null? sequence)という最終条件が来るまで式が展開されるので、再帰が深いとスタックをヒープにコピーする作業でガツンと遅くなる場合がある(らしい)。

(use srfi-1)
(time (begin (reverse-fr (iota 1200 0)) #t))
;gosh> ;(time (begin (reverse-fr (iota 1200 0)) #t))
; real   0.271
; user   0.270
; sys    0.000
;#t

(time (begin (reverse-fl (iota 500000 0)) #t))
;gosh> ;(time (begin (reverse-fl (iota 500000 0)) #t))
; real   0.236
; user   0.240
; sys    0.010
;#t

;;スケールするのは反復的プロセス!!

(time (begin (fold-right * 1 (iota 10000 0)) #t))
(time (begin (fold-left  * 1 (iota 10000 0)) #t))

;append -> O(n)
;cons   -> O(1)

まぁ、appendの遅延なんて取るに足らないと思うけど。

上記の件はid:yamanetoshiさんのブログでshiroさんがコメントしてくれてたり。

goshを-fcollect-starsオプションで起動するとスタックオーバーフローの回数やら時間やらが見れるらしい。使わせてもらいますとも。

SICP Reading #18 [2.2.2] 階層構造

前節で扱った非常にシンプルなデータ構造リストを組み合わせて、木などの階層構造も表現できるよね。という感じらしい。

問題2.24

(list 1 (list 2 (list 3 4)))

を評価した結果のデータ構造を図で書け、という問題。書いたけど見せれる物ではないのでここではスルー。

問題2.25

提示されたリストから、”7″をとりだす操作をcar と cdrで書け!という問題。

(define l1 '(1 3 (5 7) 9))
(car (cdr (car (cdr (cdr l1)))))
;gosh> 7
(car (cdaddr l1))
;gosh> 7

(define l2 '((7)))
(car (car l2))
;gosh> 7
(caar l2)
;gosh> 7

(define l3 '(1 (2 (3 (4 (5 (6 7)))))))
(car (cdr (car (cdr (car (cdr (car (cdr (car (cdr (car (cdr l3))))))))))))
;gosh> 7

問題2.26

リスト演算を評価した結果を予想する問題。
とりあえず結果から。

(define x (list 1 2 3))
(define y (list 4 5 6))

(append x y)
;gosh> (1 2 3 4 5 6)

(cons x y)
;gosh> ((1 2 3) 4 5 6)

(list x y)
;gosh> ((1 2 3) (4 5 6))

リスト は お尻がnilな対の連鎖。
appendは第一引数のリストの末尾のnilを台に引数に変更する命令。
なので第一引数はリストじゃなくてはならないけど、第二引数はリストじゃなくても良い。
でも、その場合は結果はリストじゃなくなる(末尾がnilでない)。

(append x 3)
;gosh> (1 2 3 . 3)
(list? (append x 3))
;gosh> #f

consした結果が((1 2 3) 4 5 6)ってのは一瞬アレ?って思ったけど、cdrした結果が(4 5 6)となることを考えれば当然。

listはx, yのリストを作ってx, yをそれぞれ置き換えればok.

問題

listの入れ子にも対応したreverse手続きdeep-reverseを定義する問題。
reverse対象のlistの要素がlistなら、再帰的にreverseすれば良いので、いかな定義に

(define (deep-reverse lst)
  (let iter ((lst lst) (rev '()))
    (if (null? lst) rev
        (iter (cdr lst) (cons
                         (if (pair? (car lst))
                             (deep-reverse (car lst))
                             (car lst))
                         rev)))))

テスト

僕のとは気合い違うテスト。
実はid:yamanetoshiさんが書いてたり :-)

(add-load-path ".")
(load "q2-27")

(use gauche.test)
(test-start "deep-reverse")
(test-section "deep-reverse")
(test* "'() to '()"
      '()
      (deep-reverse '()))

(test* "'(1) to '(1)"
      '(1)
      (deep-reverse '(1)))

(test* "'(1 2) to '(2 1)"
      '(2 1)
      (deep-reverse '(1 2)))

(test* "'((1 2) 3 4) to '(4 3 (2 1))"
      '(4 3 (2 1))
      (deep-reverse '((1 2) 3 4)))

(test* "'((1 2) (3 4)) to '((4 3) (2 1))"
      '((4 3) (2 1))
      (deep-reverse '((1 2) (3 4))))

(test* "'(1 2 3 (4 5) 6 (7 8)) to '((8 7) 6 (5 4) 3 2 1)"
      '((8 7) 6 (5 4) 3 2 1)
      (deep-reverse '(1 2 3 (4 5) 6 (7 8))))

(test* "'(1 (2 3 (4 5)) 6 (7 (8))) to '(((8) 7) 6 ((5 4) 3 2) 1)"
      '(((8) 7) 6 ((5 4) 3 2) 1)
      (deep-reverse '(1 (2 3 (4 5)) 6 (7 (8)))))
(test-end)
SHINYA% gosh q2-27.test.scm
Testing deep-reverse ...
<deep-reverse>-----------------------------------------------------------------
test '() to '(), expects () ==> ok
test '(1) to '(1), expects (1) ==> ok
test '(1 2) to '(2 1), expects (2 1) ==> ok
test '((1 2) 3 4) to '(4 3 (2 1)), expects (4 3 (2 1)) ==> ok
test '((1 2) (3 4)) to '((4 3) (2 1)), expects ((4 3) (2 1)) ==> ok
test '(1 2 3 (4 5) 6 (7 8)) to '((8 7) 6 (5 4) 3 2 1), expects ((8 7) 6 (5 4) 3 2 1) ==> ok
test '(1 (2 3 (4 5)) 6 (7 (8))) to '(((8) 7) 6 ((5 4) 3 2) 1), expects (((8) 7) 6 ((5 4) 3 2) 1) ==> ok
passed.

問題2.28

木構造(入れ子リスト)のデータを受け取り、木のすべての葉を左から右に順に並べたリストを返す手続きを実装する問題
要は、tree->list

(define x (list (list 1 2) (list 3 4)))

(define (fringer tree)
  (cond ((null? tree) '())
        ((pair? tree)
         (append (fringer (car tree))
                 (fringer (cdr tree))))
        (else (list tree))))

なんで”fringer”って名前なのか微妙に謎

以下な感じで動いてる。

(fringer x)
;gosh> (1 2 3 4)

(fringer (list x x))
;gosh> (1 2 3 4 1 2 3 4)
(fringer (cons x x))
;gosh> (1 2 3 4 1 2 3 4)
(fringer (append x x))
;gosh> (1 2 3 4 1 2 3 4)

問題2.29

木構造を用いて2進モービルを表現し、色々な演算を定義していく。
最初2進モービルってなんぞ?って感じだったけど、まぁ天秤みたいなモノのことらしい。
天秤に天秤をぶらさげることもできるので、まぁcompositeパターン。

定義は以下。

;;2進モービル(天秤)
(define (make-mobile left right)
  (list left right))

(define (make-branch length structure)
  (list length structure))

選択子(インタフェース)の定義

consではなくlistなので、2番目の要素が取りたい場合はcadr(cdrするとlistが帰る)。

(define (left-branch mobile)
  (car mobile))

(define (right-branch mobile)
  (cadr mobile))

(define (branch-length branch)
  (car branch))

(define (branch-structure branch)
  (cadr branch))

モービルの全重量を返す手続き

ブランチの先がモービルの場合と、重りの場合で処理を分けてモービルなら再帰すれば良い.

(define (mobile-weight mobile)
  (+ (branch-weight (left-branch mobile))
     (branch-weight (right-branch mobile))))

(define (branch-weight branch)
  (if (pair? (branch-structure branch))
      (mobile-weight (branch-structure branch))
      (branch-structure branch)))

ブランチが釣り合ってるか判定する手続き

rootのブランチだけでなく、全ての部分ブランチが釣り合ってる場合にのみ#tを返すという制約が。

(define (branched? mobile)
  (if (not (pair? mobile)) #t
      (let ((lb (left-branch mobile))
            (rb (right-branch mobile)))
        (if (= (* (branch-weight lb) (branch-length lb))
               (*   (branch-weight rb) (branch-length rb)))
            (and (branched? (branch-structure lb))
                 (branched? (branch-structure rb)))
            #f))))

この実装だと、部分木の重さを調べるmobile-weightが余分に呼び出されてしまうけど、memoizeすればいいじゃん。ってことでスルー:-)

ここでテスト

(add-load-path ".")
(load "q2-29.scm")

(use gauche.test)
(test-start "mobile test")

(define b1 (make-branch 4 8))
(define b2 (make-branch 8 4))
(define m1 (make-mobile b1 b2))

(test* "mobile test 1"
       #t
       (branched? m1))

(define b3 (make-branch 3 6))
(define b4 (make-branch 9 2))
(define m2 (make-mobile b3 b4))

(test* "mobile test 2"
       #t
       (branched? m2))

(define b5 (make-branch 2 m1))
(define b6 (make-branch 3 m2))
(define m3 (make-mobile b5 b6))

(test* "mobile test 3"
       #t
       (branched? m3))

(test-end)
SHINYA% gosh q2-29.test.scm
Testing mobile test ...
test mobile test 1, expects #t ==> ok
test mobile test 2, expects #t ==> ok
test mobile test 3, expects #t ==> ok
passed.

あまり気合いの入ってないテストだけど、とりあえずOKぽい!

mobileのlistをconsにした場合

mobileのデータ構造をlistからconsの入れ子に変更した場合どうなる?という問題。
データ構造と処理を、インタフェースで区切っていたおかげでインタフェースの変更のみでOK

;;list を cons にした場合
(define (make-mobile left right)
  (cons left right))

(define (make-branch length structure)
  (cons length structure))

;;インターフェースのみを変更 cadr -> cdr
(define (right-branch mobile)
  (cdr mobile))

(define (branch-structure branch)
  (cdr branch))

2進の場合はconsが最適。 n進(そんなのあるのか?)だとlistな必要が。

問題2.30

問題2.21で実装したsquare-listの木構造版を実装する問題。

(define (square x) (* x x))

(define (square-tree tree)
  (if (null? tree) '()
      (cons (if (pair? (car tree))
                (square-tree (car tree))
                (square (car tree)))
            (square-tree (cdr tree)))))

mapを使うと以下な感じで実装できる。

(define (map proc items)
  (if (null? items) '()
      (cons (proc (car items))
            (map proc (cdr items)))))

(define (mapped-square-tree tree)
  (map (lambda (t) (if (pair? t)
                       (mapped-square-tree t)
                       (square t))) tree))

テスト

一応書いたのでテスト

(add-load-path ".")
(load "q2-30")

(use gauche.test)
(test-start "square-tree")
;;(test-section "make-rat")
(test* "square-tree"
       '(1 (4 (9 16) 25) (36 49))
       (square-tree
        (list 1
              (list 2 (list 3 4) 5)
              (list 6 7))))

(test* "square-tree (use map)"
       '(1 (4 (9 16) 25) (36 49))
       (mapped-square-tree
        (list 1
              (list 2 (list 3 4) 5)
              (list 6 7))))

(test-end)
SHINYA% gosh q2-30.test.scm
Testing square-tree ...
test square-tree, expects (1 (4 (9 16) 25) (36 49)) ==> ok
test square-tree (use map), expects (1 (4 (9 16) 25) (36 49)) ==> ok
passed.

問題2.31

2.30の実装を抽象化し、mapのtree版を実装する問題

(add-load-path ".")
(load "q2-30")

(define (tree-map proc tree)
  (if (null? tree) '()
      (cons (if (pair? (car tree))
                (tree-map proc (car tree))
                (proc (car tree)))
            (tree-map proc (cdr tree)))))

tree-mapを用いると、square-treeが

(define (square-tree proc tree) (treemap square tree))

と定義できる。これは嬉しい。

問題2.32

最後に若干質の違う問題。
冪集合の問題。
Xの冪集合は、Xの先頭要素xを除いた集合の冪集合とその各要素にxを加えた集合の和集合。
日本語でも再帰的な定義。以下のコードは上の定義をそのままコードに落とした感じ。

(define (subsets s)
  (if (null? s)
      (list '())
      (let ((rest (subsets (cdr s))))
        (append rest (map (lambda (x)
                            (cons (car s) x)) rest)))))

ちなみに、gaucheでsubsetsを実行すると

(subsets '(1 2 3))
;gosh> (() #0=(3) #1=(2) #2=(2 . #0#) (1) (1 . #0#) (1 . #1#) (1 . #2#))

となっているけど、これはどうやら同値なオブジェクトを共有しているということらしい。

(define s (subsets '(1 2 3)))
(eq? (list-ref s 3) (cdr (list-ref s 7)))
;gosh> #t

手続きeq?は
“ふたつとも同じ型で、ドット対かベクター、または、文字列でメモリの同じ場所にあるとき”
なので、メモリ節約のソレかな。実装はハッシュとか?

SICP Reading #16 [2.1.4] 区間算術演算

大学の講義がようやく一段落。この春休みでSICPとGCC読みを大幅に進めていきたいな。
(実はエンベデッドシステムスペシャリスト試験も気になってたりするけど、、まだ未定! )

SICPに関してはおよそ2週間ぐらい離れてたりした。しかも、その間C/C++でがりがりプログラミングをしてたので”lambdaの精神”を若干わすれてたり:-)

区間算術演算

誤差を考慮した計算の問題。
例として電気的な抵抗R1とR2の並列にした等価抵抗(高校[もしや中学?]の時にやった)を求める課題。

問題2.7

区間算術演算の構成子のインターフェースを定義

;cdr
(define (upper-bound i)
  (max (car i) (cdr i)))

;car
(define (lower-bound i)
  (min (car i) (cdr i)))

構成子intervalはただのdotted pairなので それぞれcar,cdrでOK。
でも、本の例だとmake-interval は引数の順番でそのままcons している。
なのでupperで上限,lowerで下限が取れるようにmin,max (ほんとはmake-intervalで順序づけたら良いと思うけど)

問題2.8

intervalの減算。きちんと上限,下限となるように構成子を組み合わせること。

(add-load-path ".")
(load "q2-7.scm")

(define (sub-interval x y)
  (make-interval (- (upper-bound x) (lower-bound y))
                 (- (lower-bound x) (upper-bound y))))

問題2.9

区間の幅とそれぞれの演算結果との関係をとく。

区間の幅(誤差の幅)は、以下な定義。

(add-load-path ".")
(load "q2-8.scm")

(define (interval-width x)
  (- (upper-bound x) (lower-bound x)))

加算(減算)の場合

(define r1 (make-interval-with-tolerance 6.8 0.1))
(define r2 (make-interval-with-tolerance 4.7 0.05))

(+ (interval-width r1) (interval-width r2))
;gosh> 1.830000000000001
(interval-width (add-interval r1 r2))
;gosh> 1.83
(interval-width (sub-interval r1 r2))
;gosh> 1.830000000000001

加減の場合は、結果の幅は引数の幅の和(関数)となる。
これは簡単。

乗算(除算)

(define s1 (make-interval 3 5))
(define s2 (make-interval -4 -2))
(define s3 (make-interval -1 1))

(interval-width s1)
;gosh> 2
(interval-width s2)
;gosh> 2
(interval-width s3)
;gosh> 2

(interval-width (mul-interval s1 s2))
;14
(interval-width (mul-interval s2 s3))
;8
(interval-width (mul-interval s1 s3))
;10

乗除の場合は、引数の幅がおなじであっても結果の幅が異なる場合がある。
乗除の場合幅よりも実際の値が関数となる、という感じで。

問題2.10

零除算時のエラー処理を追加せよという問題。
もともとの除算手続きdiv-intervalは

(div-interval (make-interval 2 3)
              (make-interval 0 2))
;gosh> (1.0 . +inf.0)

となる。最もgaucheの場合はerrorは投げられずinf.0 となるだけ。
(これはSchemeの仕様?)

gaucheでerror手続きはメッセージを出力してstack traceの表示。(例外を投げるのはraise か? いつか試す)

(define (div-interval x y)
  (let ((uy (upper-bound y))
        (ly (lower-bound y)))
    (if (= 0 (* uy ly))
        (error "divided by 0 - " y)
        (mul-interval x (make-interval
                         (/ 1 (upper-bound y))
                         (/ 1 (lower-bound y)))))))

(div-interval (make-interval 2 3)
              (make-interval 0 2))

;gosh> *** ERROR: divided by 0 -  (0 . 2)
;Stack Trace:
;_______________________________________

問題2.11

これは意外と考え込んだ。というのも、上限と加減を持つ区間2つの演算で、9つという場合の数がピンと来なかったから。

(ちょっと時間がかかったけど)”端点の符号を調べる”というBenの発言で、1つの区間は3つの場合に分けられることに気付いた。

  1. 正区間:上限,加減がともに正
  2. 零区間:区間に零を含む
  3. 負区間:上限,加減がともに負

ここまで気付けば、Benの発言の意味が全て分かる:-)

プログラムでは、条件分岐で区間判定の条件式がならぶのがいやだったので関数でまとめた。

(add-load-path ".")
(load "q2-8.scm")

;;正区間
(define (positive-interval? x)
  (> (lower-bound x) 0))

;;負区間
(define (negative-interval? x)
  (< (upper-bound x) 0))

;;零区間(0を含む)
(define (0-included-interval? x)
  (not (> (* (lower-bound x) (upper-bound x)) 0)))

で、実装。

;;posi,nega,eleseでもちろんOK
;;でも、分かりすさのためposi,nega,0-includedで場合分け
(define (mul-interval x y)
  (let ((lx (lower-bound x))
        (ux (upper-bound x))
        (ly (lower-bound y))
        (uy (upper-bound y)))
    (cond ((positive-interval? x)
           (cond ((positive-interval? y)
                  (make-interval (* lx ly) (* ux uy)))
                 ((negative-interval? y)
                  (make-interval (* ux ly) (* lx uy)))
                 ((0-included-interval? y)
                  (make-interval (* ux ly) (* ux uy)))))
          ((negative-interval? x)
           (cond ((positive-interval? y)
                  (make-interval (* lx uy) (* ux ly)))
                 ((negative-interval? y)
                  (make-interval (* ux uy) (* lx ly)))
                 ((0-included-interval? y)
                  (make-interval (* lx uy) (* lx ly)))))
          ((0-included-interval? x)
           (cond ((positive-interval? y)
                  (make-interval (* lx uy) (* ux uy)))
                 ((negative-interval? y)
                  (make-interval (* ux ly) (* lx ly)))
                 ((0-included-interval? y)
                  (make-interval (min (* ux ly) (* lx uy))
                                 (max (* ux uy) (* lx ly)))))))))

長い。。書いてる間頭の中で区間がイメージできてないとかなりやっかいな問題
(図を書けば余裕だけど)

テスト

テストケースを生成するのが面倒くさかったので乱数列で生成して計算することに。
そのさい3つの区間が均等に生成されるようにちょっと工夫。

(add-load-path ".")
(load "q2-11.scm")
(use math.mt-random)

(use gauche.test)

(define mt (make <mersenne-twister> :seed (sys-time)))
(define (_mul-interval x y)
  (let ((p1 (* (lower-bound x) (lower-bound y)))
        (p2 (* (lower-bound x) (upper-bound y)))
        (p3 (* (upper-bound x) (lower-bound y)))
        (p4 (* (upper-bound x) (upper-bound y))))
    (make-interval (min p1 p2 p3 p4)
                   (max p1 p2 p3 p4))))

(define (mt-random-num mt range sign)
  (* (expt -1 sign)
     (mt-random-integer mt range)))

(define (make-random-interval mt range sign-pattern)
  (let1 pattern (remainder sign-pattern 3)
        (cond ((= pattern 0) ;;正区間ペア
               (make-interval
                (mt-random-num mt range 0)
                (mt-random-num mt range 0)))
              ((= pattern 1) ;;負-生区間
               (make-interval
                (mt-random-num mt range 1)
                (mt-random-num mt range 0)))
              ((= pattern 2) ;;負区間ペア
               (make-interval
                (mt-random-num mt range 1)
                (mt-random-num mt range 1))))))

(test-start "mul-interavl")
(let iter ((count 20)
           (interval1 (make-random-interval mt 10 20))
           (interval2 (make-random-interval mt 10 20)))
  (cond ((> count 0)
         (test* #`",interval1 * ,interval2"
                (_mul-interval interval1 interval2)
                (mul-interval interval1 interval2))
         (iter (- count 1)
               (make-random-interval mt 10 (- count 1))
               (make-random-interval mt 10 (- count 1))))
        (else (test-end))))

で、実行

SHINYA% gosh q2-11.test.scm                                                                                 [~/ln/sicp/q/s2]Testing mul-interavl ...
test (-1 . -5) * (-8 . -1), expects (1 . 40) ==> oktest (-6 . 7) * (-3 . 4), expects (-24 . 28) ==> ok
test (3 . 4) * (2 . 7), expects (6 . 28) ==> oktest (-8 . -7) * (-4 . -4), expects (28 . 32) ==> ok
test (-5 . 4) * (-8 . 5), expects (-32 . 40) ==> ok
test (3 . 6) * (9 . 3), expects (9 . 54) ==> ok
test (-6 . -5) * (0 . -4), expects (0 . 24) ==> ok
test (-2 . 7) * (0 . 5), expects (-10 . 35) ==> ok
test (9 . 5) * (5 . 0), expects (0 . 45) ==> ok
test (-3 . -8) * (0 . -6), expects (0 . 48) ==> ok
test (-9 . 7) * (-5 . 0), expects (-35 . 45) ==> ok
test (6 . 6) * (9 . 8), expects (48 . 54) ==> ok
test (-7 . -4) * (-3 . 0), expects (0 . 21) ==> ok
test (-7 . 0) * (-9 . 1), expects (-7 . 63) ==> ok
test (6 . 6) * (1 . 8), expects (6 . 48) ==> ok
test (-4 . -2) * (-4 . -8), expects (8 . 32) ==> ok
test (-4 . 8) * (-1 . 9), expects (-36 . 72) ==> ok
test (9 . 7) * (0 . 1), expects (0 . 9) ==> oktest (-6 . -7) * (-4 . -4), expects (24 . 28) ==> ok
test (-7 . 8) * (-9 . 7), expects (-72 . 63) ==> ok
passed.

テストケースを生成して走らせるのは面白いし大量の結果がみれて楽。
素晴らしい。

問題2.12

下限と上限ではなく、中央値と誤差率で区間を元に区間を生成する手続きの作成。
(じつは、問題の前に即実装してたり)

(add-load-path ".")
(load "q2-8.scm")

;;q2-7で定義したのが答え
(define (make-interval-with-tolerance val tol)
  (make-interval (* val (+ 1 tol)) (* val (- 1 tol))))

;;一応、定義しなおす
(define (make-center-percent center percent)
  (make-interval
   (* center (- 1 percent))
   (* center (+ 1 percent))))

(define (center i)
  (/ (+ (lower-bound i) (upper-bound i)) 2))

(define (width i)
  (/ (- (upper-bound i) (lower-bound i)) 2))

(define (percent i)
  (/ (width i) (abs (center i))))

問題2.13

パーセント相対許容誤差が小さい場合の、mul-intervalの近似式を実装せよという問題。

誤差が小さいのなら、中心値どうしを直接掛けて、結果をそのまま中心値として区間を生成すればよいと判断。
で、誤差は和を取ってみた(ちょっと安易?)。

(add-load-path ".")
(load "q2-12.scm")

;;mul-intervalの近似式(正区間のみと仮定)
(define (rough-mul-interval x y)
  (make-center-percent (* (center x) (center y))
                       (+ (percent x) (percent y))))

テスト

(add-load-path ".")
(load "q2-13.scm")

(define s1 (make-center-percent 7 0.1))
(define s2 (make-center-percent 5 0.1))

(mul-interval s1 s2)
;gosh> (28.349999999999998 . 42.35000000000001)
(rough-mul-interval s1 s2)
;gosh> (27.999999999999996 . 42.00000000000001)

う〜む。一応近似できてると思う。。次!

問題2.14

この問題の重要なところ、等価な計算でも誤差の幅に変化があるという問題

並列抵抗の式

par1と

par2は等価な式というのは自明だけど、プログラムで計算すると結果が異なる。

(add-load-path ".")
(load "q2-12.scm")

(define (par1 r1 r2)
  (div-interval (mul-interval r1 r2)
                (add-interval r1 r2)))

(define (par2 r1 r2)
  (let1 one (make-interval 1 1)
        (div-interval one
                      (add-interval
                       (div-interval one r1)
                       (div-interval one r2)))))
(add-load-path ".")
(load "q2-14.scm")

(define r1 (make-center-percent 6.8 0.1))
(define r2 (make-center-percent 4.7 0.05))

(define p1 (par1 r1 r2))
p1 ;gosh> (2.201031010873943 . 3.4873689182805863)
(define p2 (par2 r1 r2))
p2 ;gosh> (2.581558809636278 . 2.97332259363673)

(center p1)
;gosh> 2.844199964577265
(center p2)
;gosh> 2.777440701636504
(percent p1)
;gosh> 0.22613352145193358
(percent p2)
;gosh> 0.0705260392723452

par1よりもpar2のほうが誤差が少ないことがわかる(精度が3倍ほど良い?)。par2の方が演算は多いけど。

考察

これは、幅が中央値に比べて小さいパーセントの区間の場合加減では誤差がへり、乗除では誤差が増えるから(と実験して確認)。

(define s1 (make-center-percent 2.3 0.1))
(define s2 (make-center-percent 1.8 0.1))

(percent (add-interval s1 s2))
;gosh> 0.09999999999999999
(percent (sub-interval s1 s2))
;gosh> 0.8200000000000006
(percent (mul-interval s1 s2))
;gosh> 0.19801980198019803
(percent (div-interval s1 s2))
;gosh> 0.19801980198019808

だけど、除算で区間の逆元を取る場合は、誤差に変化はでない。

(define one (make-interval 1 1))

;;(div-interval one s1)
;;(div-interval one s2)
;;の計算は誤差が増えない

(percent (div-interval one s1))
;gosh> 0.09999999999999998
(percent (div-interval one s2))
;gosh> 0.10000000000000006

問題2.15

par1よりpar2が優れてる理由は、par2では誤差が増える乗算及び誤差を持ったintervalどうしの除算を行なっているけど、par2では誤差の減る加算及び誤差の変わらない逆元の演算のみで構成されているところ。

問題2.16

等価な式でも異なる答えが出る理由は、ある一つの式は別の演算構成による等価な式に変換可能で、また演算によって誤差がことなるため。
この欠点を補うには、式自体を誤差のすくない演算で構成された式に最適化すれば良い(かも)。
先の例で思いつくのは、誤差を持ったinterval同士の乗除をなるべく行なわず、なるべく加減,逆元を求める演算に置き換えるといったところかな。

SICP Reading #15 [2.1.3] データとは何か

データとは何か、ということでconsを実装してみるやらという問題。
言語の最も基本的なデータ構造を自前で実装できるのかよ!とか、一瞬思ったけどなんてことはなく

問題2.4

(define (cons x y)
  (lambda (m) (m x y)))

(define (car z)
  (z (lambda (p q) p)))

こんな感じ。
つまり、レキシカルクロージャの中に変数を閉じ込めておという手法。
で、carの定義も素晴らしい。p,qを引数にとってp(つまりcar)を返す手続きをconsセルに適用。
lambdaさえあればなんでもできんな。

で、cdrな実装は下のようにcdr部を返す手続きを。

(define (cdr z)
  (z (lambda (p q) q)))

問題2.5

こんどは、consセルに保存する値が非負の正数に限定する場合、a,bのペアを

で表現してみろという問題。

(define (cons x y)
  (* (expt 2 x) (expt 3 y)))

(define (car x)
  (let iter ((e 0) (t x))
    (if (= 0 (remainder t 2))
        (iter (+ e 1) (/ t 2))
        e)))

(define (cdr x)
  (let iter ((e 0) (t x))
    (if (= 0 (remainder t 3))
        (iter (+ e 1) (/ t 3))
        e)))

これって、2,3の乗数の積じゃなくても、任意の素数の乗数の積でも定義できるね。

(cons 0 4)
;gosh> 81
(car (cons 0 4))
;gosh> 0
(cdr (cons 0 4))
;gosh> 4

(cons 3 40)
;gosh> 97261323672455430408
(car (cons 3 40))
;gosh> 3
(cdr (cons 3 40))
;gosh> 40

問題2.6

Church数に関数問題。

最初にChurch数での0と、1を足す演算の定義が

(define zero (lambda (f) (lambda (x) x)))

(define (add-1 n)
  (lambda (f) (lambda (x) (f ((n f) x)))))

で与えられている。

oneとtwoを直接定義せよという問題。

とりあえず、(add-1 zero)を置き換えていくことに。

(define one (add-1 zero))

(define one
  (lambda (f)
    (lambda (x)
      (f (((lambda (f)
             (lambda (x)
               x)) f) x)))))

う〜む、複雑すぎる。。とかちょっと悩んでた。
そもそもChurch数の定義が全く分からなかったので、改めて0の定義を見直すことに。

(define zero (lambda (f) (lambda (x) x)))

((zero (lambda ())) 0)
;gosh> 0

fを引数にうけとって「xを受け取りxを返す手続き」を返す手続き。
ん?fの使って無いじゃんとか思ったり。
で、Church数について調べてみると、手続きfを受け取りn回xにfを適用する手続きで自然数を表しているとのこと。

なるほどね、と。タネが分かればadd-1手続きの書き下しも楽に。(簡約の形が明確に見えた)

(define one
  (lambda (f)
    (lambda (x)
      (f (((lambda (f)
             (lambda (x)
               x)) f) x)))))
;;は、
;;(((lambda (f)
;;  (lambda (x)
;;    x)) f) x)
;;=== x (恒等写像)
;;で簡約できるので
(define one (lambda (f) (lambda (x) (f x))))

なので、(形は見えてるけど)twoもadd-1の書き下しから簡約し直接定義いてみる。

(define two
  (lambda (f)
    (lambda (x)
      (f (((lambda (f)
             (lambda (x)
               (f (((lambda (f)
                      (lambda (x)
                        x)) f) x)))) f) x)))))

;;ここも、恒等写像((zero f) x)の部分を評価すれば
(define two
  (lambda (f)
    (lambda (x)
      (f (((lambda (f)
             (lambda (x)
               (f x))) f) x)))))

;;さらに評価
;;(((lambda (f)
;;           (lambda (x)
;;             (f x))) f) x)
;;=== (f x)

(define two
  (lambda (f)
    (lambda (x)
      (f (f x)))))

Church数の意味さえ分かってしまえば楽勝。
で、追加問題としてChurch数の加算手続きを定義せよという問題が。

m回作用とn回作用の関数をうけとり、m+n回作用させる関数を返せば良いので

(define (+ a b)
  (lambda (f)
    (lambda (x)
      ((a f) ((b f) x)))))

とサクッと実装。

実装が正しいか確かめるために、手続きとしてchurch-to-num手続きを実装(インクリメント手続き)
これを0に何回適用させるかで、そのままChurch数から自然数に変換できるはず。

(define (church-to-num c)
  ((c (lambda (x) ((with-module gauche +) x 1))) 0))

(define three (+ two one))

(church-to-num zero)
;gosh> 0
(church-to-num one)
;gosh> 1
(church-to-num two)
;gosh> 2
(church-to-num three)
;gosh> 3

それぞれが正しく変換されてることから、実装が正しいと分かりました。

SICP Reading #14 [2.1.1 ~ 2.1.2] データ抽象

講義の合間にちょくちょく進めてようやく1章を抜けました。
ここにきて、ようやくcons手続きが! なんか嬉しかったり。

データ抽象の入門ということで、consを使えばデータがペアで扱える。
で、consを組み合わせることでさらに複雑なデータを表現できるという感じで進んでいく。
重要なのはデータ表現は、データを使うプログラムから独立に定義すべきで、さらにインターフェースを提供して抽象データを実装するとのこと。

問題2.1

正負両方の引数から分数を作る手続きmake-ratの実装。
分母は正に正規化するというしばりつき。

分母,分子がそれぞれ正負の場合で4つの場合分けがどうのだとめんどくさいなーとか考えてたけど、
要は分母が負の場合だけ分母,分子に負をかければよろしいだけと気づく。

;; gcd is natural number
(define (gcd x y)
  (let iter ((x (abs x)) (y (abs y)))
    (if (= 0 (remainder x y)) y
        (iter y (remainder x y)))))

(define (make-rat n d)
  (let* ((d-sign (if (> d 0) 1 -1))
        (g (gcd n d)))
  (cons (* (/ n g) d-sign) (* (/ d g) d-sign))))

(define (numer x) (car x))
(define (denom x) (cdr x))

(define (print-rat x)
  (print (numer x) "/" (denom x)))

分母が正か負かで場合分けして演算を2つ記述する実装でもよかったけど、符号をそのまま演算に組み込むことで演算は1つで統一(結局,符号抜き出すときにif使ってるけど)

テスト

ちゃんとtestしてみる。

(add-load-path ".")
(load "q2-1.scm")

(use gauche.test)
(test-start "make-rat")
(test-section "make-rat")
(test* "24/3"
       (cons 8 1)
       (make-rat 24 3))

(test* "-5/2"
       (cons -5 2)
       (make-rat -5 2))

(test* "-90/-3"
       (cons 30 1)
       (make-rat -90 -3))

(test* "6/-8"
       (cons -3 4)
       (make-rat 6 -8))

(test-end)
SHINYA% gosh q2-1.test.scm
Testing make-rat ...
<make-rat>---------------------------------------------------------------------
test 24/3, expects (8 . 1) ==> ok
test -5/2, expects (-5 . 2) ==> ok
test -90/-3, expects (30 . 1) ==> ok
test 6/-8, expects (-3 . 4) ==> ok
passed.

test が通るのは気持ち良い!

問題2.2

線分をデータ抽象して、さらに中間点を求める手続きを実装する問題。

(define (make-point x y)
  (cons x y))

(define (x-point p)
  (car p))

(define (y-point p)
  (cdr p))

(define (print-point p)
  (print "\n(" (x-point p) ", " (y-point p) ")"))

(define (make-segment start end)
  (cons start end))

(define (start-segment segment)
  (car segment))

(define (end-segment segment)
  (cdr segment))

(define (average x y)
  (/ (+ x y) 2))

(define (midpoint-segment segment)
  (let ((sx (x-point (start-segment segment)))
        (sy (y-point (start-segment segment)))
        (ex (x-point (end-segment segment)))
        (ey (y-point (end-segment segment))))
    (make-point (average sx ex) (average sy ey))))

x,y座標を抜き出すインターフェースであるx-point,y-pointの実装は、carとcdrの再実装。
なら(define x-point car)で良いじゃんとか思ったけど、しっかりコメントに説明があった
(define x-point car)とすると、make-ratがconsを呼びだすのではなくmake-ratがconsということになるので、トレースやブレークなどのデバッグ作業に支障がでるとのこと。
なるほどねなー、という感じ。*本質的に同じようなことを最近研究室えも教えられたり。

テスト

テストはprint-pointで事足りる問題だけど、testを書くクセを付ける。

(add-load-path ".")
(load "q2-2.scm")

(use gauche.test)
(test-start "segment")
(test-section "make-segment")
(test* "segment (0 0) to (10 10)"
       (cons (cons 0 0) (cons 10 10))
       (make-segment (make-point 0 0)
                     (make-point 10 10)))

(test-section "midpoint-segment")
(test* "midpoint of segment (0 0) to (10 10)"
       (make-point 5 5)
       (midpoint-segment
        (make-segment (make-point 0 0)
                      (make-point 10 10))))
(test-end)
SHINYA% gosh q2-2.test.scm
Testing segment ...
<make-segment>-----------------------------------------------------------------
test segment (0 0) to (10 10), expects ((0 . 0) 10 . 10) ==> ok
<midpoint-segment>-------------------------------------------------------------
test midpoint of segment (0 0) to (10 10), expects (5 . 5) ==> ok
passed.

passed.

問題2.3

線分というデータ表現から、長方形のデータ表現とちょっとレベルアップした問題。(平面上の〜ということで、どっちも2次元ですが)
これ、ヒントに「問題2.2のmidpoint-segment手続きが使いたくなるだろう」みたいなことが書かれてあるけど、長方形の表現で中間点を使う表現はわかりにくいだろうということでボツ:-)

結局、[原点(左下座標)と幅,高さ], [2つの対角点]から長方形を表現する感じで2つ実装してみた。

(add-load-path ".")
(load "q2-2.scm")

;;基本となる長方形の表現 原点座標と幅/高さ
;;幅/高さは負でもOK(正規化する)
(define (make-rectangle origin width height)
  ;;引数を正規化して再帰
  (cond ((and (< width 0) (< height 0))
         (make-rectangle
          (make-point (+ (x-point origin) width)
                      (+ (y-point origin) height))
          (- width) (- height)))
        ((< width 0)
         (make-rectangle
          (make-point (+ (x-point origin) width)
                      (y-point origin))
          (- width) height))
        ((< height 0)
         (make-rectangle
          (make-point (x-point origin)
                      (+ (y-point origin) height))
          width (- height)))
        (else (cons origin (cons width height)))))

(define (width rectangle)
  (cadr rectangle)) ;; (cadr hoge) == (car (cdr hoge))

(define (height rectangle)
  (cddr rectangle)) ;; (cddr hoge) == (cdr (cdr hoge))

(define hoge (make-rectangle (make-point 0 0) 10 20))

(define (perimeter rectangle)
  (* (+ (width rectangle) (height rectangle)) 2))

(define (area rectangle)
  (* (width rectangle) (height rectangle)))

;;rectangle の別表現
;;2対角点での表現(組み合わせ、順番 don't care)
(define (make-rectangle-from-2-opposite-corner p1 p2)
  (make-rectangle p1
                  (- (x-point p2) (x-point p1))
                  (- (y-point p2) (y-point p1))))

一番シンプル(と思う)な、原点と幅,高さでの長方形を基本として実装。
ここで、幅と高さは負のあたいでも動くように実装。

(make-rectangle (make-point 0 0) 5 5)
;gosh> ((0 . 0) 5 . 5)
(make-rectangle (make-point 0 5) 5 -5)
;gosh> ((0 . 0) 5 . 5)
(make-rectangle (make-point 5 0) -5 5)
;gosh> ((0 . 0) 5 . 5)
(make-rectangle (make-point 5 5) -5 -5)
;gosh> ((0 . 0) 5 . 5)

こんな感じで、引数が負でも(つまり、原点が左下の点じゃなくても)幅,高さを正に正規化する感じで実装。
こうすることで、引数の順番を気にしないで2つの対角点から長方形が実装できるように(組み合わせも)。

テスト

testするぜー、超testするぜー。

(add-load-path ".")
(load "q2-3.scm")

(use gauche.test)
(test-start "make-rectangle")
(test-section "rectangle from origin,width,height")
(define square5x5 (make-rectangle
                   (make-point 0 0) 5 5))
(test* "perimeter of square5x5"
       20
       (perimeter square5x5))

(test* "area of square5x5"
       25
       (area square5x5))

(test-section "rectangle from 2 opposite corner")
(define rectangle3x8
  (make-rectangle-from-2-opposite-corner
   (make-point 0 0)
   (make-point 3 8)))

(test* "perimeter of rectangle3x8"
       22
       (perimeter rectangle3x8))

(test* "area of rectangle3x8"
       24
       (area rectangle3x8))
(test-end)
SHINYA% gosh q2-3.test.scm
<rectangle from origin,width,height>-------------------------------------------
test perimeter of square5x5, expects 20 ==> ok
test area of square5x5, expects 25 ==> ok
<rectangle from 2 opposite corner>---------------------------------------------
test perimeter of rectangle3x8, expects 22 ==> ok
test area of rectangle3x8, expects 24 ==> ok
passed.

よしよし。

今回(から)は

id:yamanetoshiさん曰く「テストを書くくせをつける。テストから先に書くのもあり(要求は問題で指定されてるからね)。」とのことなので、testちゃんとしてみた。
2.1章の残り、は問題が大量にあるっぽいので若干 test手抜きする可能性ありかも:-)

 
Better Tag Cloud