SICP Reading #13 [s1.3.4] 手続きを返す手続き, 抽象化

年末年始は、ありえないほどの寝正月っぷりだったので本日のSICP勉強会でいろいろ切り替えていく!

1章最後の問題を解いていく。

Newton法,不動点探索 な手法を使って、解を求めたりsqrt手続きを生成したり1章の締めくくりにふさわしい感じ。
最後らへんは、抽象化しすぎな感じもあったけど:-)

問題1.40

Newotn法を使って、

を満たすxを求める問題

Newton法を実装するために、不動点を求める手続き(fixed-point)や微分(近似:deriv)も全部自前で実装(すばらしい)

;; 不動点演算
(define (fixed-point f first-gueses)
  (define tolerance 0.000001)
  (define (close-enough? a b)
    (> tolerance (abs (- b a))))
  (define (try gueses)
    (let ((next (f gueses)))
      (if (close-enough? gueses next)
          gueses
          (try next))))
  (try first-gueses))

;;微分近似
(define (deriv g)
  (define dx 0.00001)
  (lambda (x)
    (/ (- (g (+ x dx)) (g x)) dx)))

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

;; Newton法の実装
(define (newton-transform g)
  (lambda (x)
    (- x (/ (g x) ((deriv g) x)))))

(define (newton-method g guess)
    (fixed-point (newton-transform g) guess))

(define (sqrt x)
  (newton-method (lambda (y) (- (square y) x))
                 1.0))

;; x^3+ax^2+bx+c
(define (cubic a b c)
  (lambda (x) (+ (cube x) (* (square x) a) (* x b) c)))

;; 零点を求める手続き
(define (cubic-zero a b c)
  (newton-method (cubic a b c)
                 1.0))

こんな感じで、全部実装。

で、テスト

;; x^3+23x^2-12x+823の零点を求める
(cubic-zero 23 -12 823)
;gosh> -24.819511855264825

;;実際に試してみる
((cubic 23 -12 823) (cubic-zero 23 -12 823))
;gosh> -2.368763716731337e-5

うーむ、ちゃんとできてる。凄い。

問題1.41

手続きを2回繰り返す手続きを作る問題。

そんなの速攻で実装

(define (double f)
  (lambda (x) (f (f x))))

(define (inc x) (+ x 1))

手続き引数にとり、その手続を2回繰り返す手続きを返す手続き(meta 手続きですね)

で、本当の問題は

(((double (double double)) inc) 5)

がどんな値を返すか?というもの
これ、最初はサクッとdouble3回使ってるからで5+8=13でしょ。
とか思ったりしたけど、実際に実行してみると違った。で、考えてみたらそのとおりなことに気づいた。

結果

(((double (double double)) inc) 5)
;gosh> 21

;;(double inc) -> ~を二回作用 (+= 2)
;;((double double) inc) -> 「~を2回作用」を2回作用 -> ~を4回作用 -> (+= 4)
;;((double (double double)) inc) -> 「~を4回作用」を4回作用 -> (+= 16)
;; 以下、2^n で指数的に増加

“「~を2回作用」を2回作用”
って、対角化ですよね、これ。ややこしい。

頭の中でちゃんと評価してなかったから最初は間違えたけど、面白い問題。

問題1.42

手続きを合成する手続きを実装する問題。
これって、手続きを引数にとって手続きを返す手続き(ややこしい)の典型的な例だ。

(define (square x) (* x x))
(define (inc x) (+ x 1))

;;合成手続きを返す手続き
(define (compose f g)
        (lambda (x) (f (g x))))

((compose square inc) 6)
;gosh> 49

問題1.43

先に実装したcomposeを用いて、手続きをn回作用させる手続きを返す手続きrepeatedを実装する問題。

(add-load-path ".")
(load "q-1-42.scm")

(define (repeated f count)
  (let iter ((p f) (count count))
    (if (< count 2) p
        (iter (compose f p) (- count 1)))))

で、square手続きを実際にrepeateしてみる。

;;(5^2)^2 = 5^4
((repeated square 2) 5)
;gosh> 625

これって、2^n回作用させる2^repeated手続きとかだったら、(double double)な繰り返しでいけるね。
repeated手続きに条件分岐でいれるのもあり(fast-expt的な)。

問題1.44

(n重)平滑化を実装する問題。平滑化って、画像処理で最近習ったのでちょっと得意げな感じががしたり。
ノイズを軽減したり、全体をぼかす処理に使われるそれ。

x-dx,x,x+dx の3点の平均をとって平滑化された値を返す手続きを実装。

(add-load-path ".")
(load "q-1-43.scm")

(define (smooth f)
  (define dx 0.000001)
  (lambda (x) (/ (+ (f (- x dx)) (f x) (f (+ x dx))) 3)))

で、n重平滑化関数は先に実装したrepeatedを用いれば

;; return n-fold smoothed function (arg f)
(define (n-fold-smoothed f n)
  ((repeated smooth n) f))

と簡単に実装。first-classな手続きすばらしい。

問題1.45

不動点探索,平均緩和法を用いて、n乗根を求める手続きを実装する問題。
全て自前で(簡単に)実装できてる点がすごい。

(add-load-path ".")
(load "q-1-43.scm")
(load "q-1-40.scm")

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

(define (average-damp f)
  (lambda (x) (average x (f x))))

(define (equation x n) ;;(expt y n)
  (lambda (y) (/ x ((repeated (lambda (a) (* y a)) n) 1))))

で、平方根,立方根はこんなふうに表現できる。

(define (sqrt x)
  (fixed-point (average-damp (equation x 1)) 1.0))

(define (cube-root x)
  (fixed-point (average-damp (equation x 2)) 1.0))

ここで問題なのは、n乗根を求めるとき、何回平均緩和法を適用させれば、解に収束するか?という問題。
ここは直感的&実験的にlog2(n)回平均緩和法を繰り返せば良いことが判明。

僕のバージョンのgaucheでは、デフォルトのlogは自然底なやつしか求められないっぽいので、自前で簡易なlog2を実装してしのぐ(最新版いれろという話ですが、エラーが。。)

;; return integer
(define (log2 x)
  (let iter ((x x) (e -1) (n 1))
    (if (> n x) e
        (iter x (+ e 1) (* n 2)))))

(define (n-fold-root x n)
  (fixed-point ((repeated average-damp (log2 n))
                (equation x (- n 1))) 1.0))

ちゃんと動くか、試してみる

(n-fold-root (expt 5 4) 4)
;gosh> 5.000000000004688
(n-fold-root (expt 2 16) 16)
;gosh> 5.000000000004688

うん、ちゃんと動く。すばらしい。

問題1.46

1章最後の問題。

本章に述べた数値計算法のいくつかは、反復改良法(iterative improvement)という非常に一般的な計算戦略の特殊化である。

とのこと(SICP)。ほとんどの手続きが反復的に値を改良してきてた。

ということで、さらに抽象化されたiterative-improvementを実装して、1.1.7節で実装したsqrt手続きと、1.3.3節のfixed-pointを書き直せとのこと。

(add-load-path ".")
(load "q-1-45.scm")

(define (iterative-improve close-enough? improve)
  (lambda (guess)
    (let recur ((guess guess))
      (if (close-enough? guess) guess
          (recur (improve guess))))))

(define (sqrt x)
  (define tolerance 0.000001)
  (define (square x) (* x x))
  (define (close-enough? guess)
    (> tolerance (abs (- x (square guess)))))
  (define (improve guess)
    (average guess (/ x guess)))
  ((iterative-improve close-enough? improve) x))

(exact->inexact (sqrt 256))
;gosh> 16.00000000000039

(define (fixed-point f first-gueses)
  (define tolerance 0.000001)
  (define (close-enough? guess)
    (> tolerance (abs (- guess (f guess)))))
  ((iterative-improve close-enough? f) first-gueses))

(define (sqrt n)
  (fixed-point ((repeated average-damp 1) (equation n 1))
               1.0))

(sqrt 16)
;gosh> 4.000000636692939

すばらしい!
だけど、ここまで抽象化されると。。逆になにしてるか分からかったり:-)

SC 合格した!

遅ればせながら、あけましておめでとうございます。

正月はまわりにひかれるほどの寝正月っぷりでした。今年は大学生活で最も重要な年になる気がする。。

タイトルの通り

情報セキュリティスペシャリスト試験に合格しました。(すごい嬉しい

得点

得点

ちょ、午後2!!
まぁ、ギリギリとはいえ合格しました。

本試験は

2009春に午後2が3点差で落ちて、2回目の挑戦でした。
春の時点では、ネットワーク関連の知識がほぼ0で、サーバもいじったことがなく落ちて当然な感じが。
2回目の今回は、試験勉強は全っくしてない(2回目ということで、モチベーションが。。)。それでも受かったってのは、大学のサーバ班での活動とかがつながったんだろうなぁ、としんみり。

研究室も決まったし、今年もがんばるぞ!!

SICP Reading #12 [s1.3.3] 零点,不動点,無限連分数

この章では関数の零点・不動点を求める手法などを例に手続きを引数として受け取る手続きで実装いていった。
勉強会で解いた問題も含めて、問題の回答を。

問題1.35

x |-> 1 + 1/x の不動点となるxを求める問題。
ここで、不動点は黄金比 φ となる。

ちなみに、上の定義は黄金比の元の定義の変形

不動点を求める手続きは、f(x) = x を誤差が無視できるほど小さくなるまで繰り返していくという単純な実装で以下のように

;;search fixed-point
(define (fixed-point f first-gueses)
  (define (close-enough? a b)
    (> tolerance (abs (- b a))))
  (define (try gueses)
    (let ((next (f gueses)))
      (if (close-enough? gueses next)
          gueses
          (try next))))
  (try first-gueses))

(define tolerance 0.000001)

;;calculate golden-ratio
(fixed-point (lambda (x) (+ 1 (/ 1 x))) 1.0)
;gosh> 1.618034447821682

(define golden-ratio (/ (+ 1 (sqrt 5)) 2))
golden-ratio
;gosh> 1.618034447821682

きちんと黄金比が出力されてますね!

問題1.36

fixed-pointに出力をかまして、平均緩和法を用いた場合と用いない場合のx^x=1000の解を求める。
出力を噛ましたfixed-pointは

;;search fixed-point
(define (fixed-point f first-gueses)
  (define (close-enough? a b)
    (> tolerance (abs (- b a))))
  (define (try gueses)
    (let ((next (f gueses)))
      (newline)
      (display gueses)
      (if (close-enough? gueses next)
          gueses
          (try next))))
  (try first-gueses))

(define tolerance 0.000001)

な感じ。

平均緩和法を使わないで、式をそのまま使用した場合。

;;search fixed-point x^x = 1000
(fixed-point (lambda (x) (/ (log 1000) (log x))) 2)
;; gosh>
;; 2
;; 9.965784284662087
;; 3.004472209841214
;; 6.279195757507157
;; 3.759850702401539
;; 5.215843784925895
;; 4.182207192401397
;; 4.8277650983445906
;; 4.387593384662677
;; 4.671250085763899
;; 4.481403616895052
;; 4.6053657460929
;; 4.5230849678718865
;; 4.577114682047341
;; 4.541382480151454
;; 4.564903245230833
;; 4.549372679303342
;; 4.559606491913287
;; 4.552853875788271
;; 4.557305529748263
;; 4.554369064436181
;; 4.556305311532999
;; 4.555028263573554
;; 4.555870396702851
;; 4.555315001192079
;; 4.5556812635433275
;; 4.555439715736846
;; 4.555599009998291
;; 4.555493957531389
;; 4.555563237292884
;; 4.555517548417651
;; 4.555547679306398
;; 4.555527808516254
;; 4.555540912917957
;; 4.555532270803653
;; 4.555537970114198
;; 4.555534211524127
;; 4.555536690243655
;; 4.555535055574168
;; 4.55553613360814

;;"try" was called 40-times!

平均緩和を使った場合。

;;search fixed-point x^x = 1000, with average dumping
(fixed-point
 (lambda (x)
   (* (/ 1 2) (+ x (/ (log 1000)(log x)))))
             2)
;; gosh>
;; 2
;; 5.9828921423310435
;; 4.922168721308343
;; 4.628224318195455
;; 4.568346513136242
;; 4.5577305909237005
;; 4.555909809045131
;; 4.555599411610624
;; 4.5555465521473675
;; 4.555537551999825
;; 4.555536019631145

;;"try" was called 11-times!!!

(expt 4.555536019631145 4.555536019631145)
;gosh> 1000.000791229235

平均緩和法を使う、使わないで試行回数は40 : 11と大幅に差が出た。
どちらの式もx^x=1000を表す全く等価な式だけど、ほんの少し形を工夫する事で解の探索でとても差がでるなー。
(ここら辺は授業でやってたり)

問題1.37

無限連分数を一般化した手続きを作成し、黄金比を求める問題。

Recursive implementation

再帰的な実装。

(define (cont-frac ni di k)
  (let recur ((i 0))
    (if (> i k) 0
        (/ (ni i) (+ (di i) (recur (+ i 1)))))))

(define (calc-golden-ratio k)
  (/ 1 (cont-frac (lambda (i) 1.0)
                  (lambda (i) 1.0)
                  k)))


(define golden-ratio (/ (+ 1 (sqrt 5)) 2))

;;4桁の精度を求めるには k = 10
(calc-golden-ratio 10)
;gosh> 1.6179775280898876
golden-ratio
;gosh> 1.618033988749895

Iterative implementation

で、若干はまってしまったのが反復的な実装(頭が固まってたなー)。
自分の頭の中では連分数の浅い方(分母の方)から先に計算していくイメージでしか考えてなかったので、再帰的な実装しか思い浮かばなかった。
id:yamanetoshiさんは難なく反復的な実装を考えついたり。言われるまで気付かなかったりした。用は、深い方の分数から先に計算していく実装だと、途中でも値を計算できるということで。

とりあえず、実装を見れば分かる感じで、

(define (cont-frac ni di k)
  (let iter ((i (- k 1)) (rslt (/ (ni k) (di k))))
    (if (< i 0) rslt
        (iter (+ i 1) (/ (ni i) (+ (di i) rslt))))))

うーん、すばらしい。

問題1.38

連分数で自然底eの近似値を求める問題。
分母数列(di)を簡単にでっちあげて実装。

(define (di i)
  (let1 _i (modulo (+ i 2) 3)
        (if (= 0 _i) (* (/ (+ i 2) 3) 2)
            1)))

(map di '(0 1 2 3 4 5 6 7 8 9 10))
;gosh> (1 2 1 1 4 1 1 6 1 1 8)

(define (euler k)
  (cont-frac (lambda (i) 1.0)
             di
             k))

(euler 20)
;gosh> 0.7182818284590453

(use math.const)
(- (euler 10) (- e 2))
;gosh> 6.746947445179785e-9

1週間schemeから離れていただけでファイルの読み込みをちょっと忘れてたり(use math.const)。
とりあえず、eの近似が実装できた。

問題1.39

tan(x)の連分数展開な問題。
タンジェントは分子の第一項がx,第二項からx^2で分母が奇数列な無限連分数で実装できるとのこと。

以下が実装、分子は全てx^2を返して最後でxで割っているところがミソ

(define (tan-cf x k)
  (/ (cont-frac (lambda (i) (* x x))
                (lambda (i) (+ (* i 2) 1))
                k) x))

(use math.const)
(exact->inexact (tan-cf (/ pi 3) 20))
;gosh> 1.7320508075688767
(tan (/ pi 3))
;gosh> 1.7320508075688767

すばらしい!

SICP勉強会 #5

宜野湾Macで久しぶりにid:yamanetoshiさんとの勉強会。
Android/スマートフォン なお話とか聞けたりした。

で、勉強会の方ですが、、、やっぱり時間とって人と話しながら問題の解いたりしてると頭回転しますね!
勉強会までの間に進めておいて、勉強会では確認とか相談とかでガンガン進めて行く感じが良いと思う(今は自分が進めるのが遅かったり…)。

ようやっと1章が終わりそうな感じなので、年末も頑張っていきたいです!

回答の投入は後で。

SICP Reading #11 [1.3.2] lambda

研究室配属やら発表やらでSICPをⅠ週間ほど放置ぎみに。いかん、年末までにいろいろ進めるぞ!

とりあえずリハビリで復習から。

lambda

手続きを作る特殊形式lambda

lambdaで生成した手続きの中の変数は? -> 引数は手続きの内部での変数(局所変数)として、それ以外の変数(自由変数)は手続きが定義された環境の変数を参照。

Goのチュートリアルを英語で読み始めたり、徐々に英語と友達になれる気がするのでen.wikipediaの記事を読んでみたり。

Closure(computer science)

In computer science, a closure is a first-class function with free variables that are bound in the lexical environment. Such a function is said to be “closed over” its free variables

first-class functionというのは、関数を引数として/戻り値として使える. 変数に格納できるなど自由に扱えるといった意味合い(まじ?)。

Goではfunction literalとかいってたな。つまり自由変数を静的に解決するfirst-classなソレがclosureということ?

局所変数を使いたい

(let ((x 2) (y 10))
     (expt x y))

((lambda (x y)
         (expt x y)) 2 10)

として評価される(syntax-sugar)。

使い捨て関数とか,局所変数とかはlambda作用でなんとかしろという話ですか。

問題1.34

(define (square n) (* n n))
(define (f g)
  (g 2))

(f square)
;gosh> 4

(f (lambda (z) (* z (+ z 1))))
;gosh> 6

(f f)

を実行したとき、どうなるか?という問題

手続きfは、引数に手続きgを受け取りそれに(g 2)を実行するという手続き。
それに自分自身を与えると

(f f)
;gosh> *** ERROR: invalid application: (2 2)
;Stack Trace:

;;(f f) #=> ((lambda (g) (g 2)) (lambda (g) (g 2)))
;;#=> ((lambda (g 2)) 2) #=> (2 2)

こんな感じになる。

(f f)で不動点オペレーターを思い出してみたり。あれは本当に頭の体操。

 
Better Tag Cloud