理系学生日記

おまえはいつまで学生気分なのか

問題 3-65 (3.5.3 Exploiting the Stream Paradigm)

今度は,\ln 2 = \sum_{i=0}^{\infty}(-1)^i \frac 1n を使って,\ln 2 の真の値を求めるよー.といってもほとんど SICP にやりかたが書いてあるので,ぼくは Scheme を書くだけだ.

まずはこれが上記の総和をとってるやつの各項をストリームにしたもの.なるほどとても直感的だ.

(define (ln2-summands n)
  (cons-stream (/ 1.0 n)
               (stream-map - (ln2-summands (+ n 1)))))

次に部分和を取る.

(define ln2-stream
  (partial-sums (ln2-summands 1)))

今 ln2-summands の第 i 要素を l_i とすると,ストリーム ln2-stream の第 n 要素は \sum_{i=0}^{n-1}l_i になってるよー.ということは n\rightarrow \infty とすると \ln 2 の値が求められているはず…!!

dump してみる.

(dump-stream ln2-stream 10)
; 1.0, 0.5, 0.8333333333333333, 0.5833333333333333, 0.7833333333333332, 0.6166666666666666, 0.7595238095238095, 0.6345238095238095, 0.7456349206349207, 0.6456349206349207, done

このように,意外と収束がおそい.しかし世界はよくできていて,オイラーが,連続する二項の符号が逆になっている数列の部分和を早く計算できることを示しているらしくて,こんな風にしたら早い!!
S_{n+1}-\frac{\left(S_{n+1}-S_n\right)^2}{S_{n-1}-2S_n+S_{n+1}}
これは SICP 中で euler-transform って関数として定義されている.

(define (euler-transform s)
  (let ((s0 (stream-ref s 0))           ; Sn-1
        (s1 (stream-ref s 1))           ; Sn
        (s2 (stream-ref s 2)))          ; Sn+1
    (cons-stream (- s2 (/ (square (- s2 s1))
                          (+ s0 (* -2 s1) s2)))
                 (euler-transform (stream-cdr s)))))

dump してみるよー

(dump-stream (euler-transform ln2-stream) 10)
; 0.7, 0.6904761904761905, 0.6944444444444444, 0.6924242424242424, 0.6935897435897436, 0.6928571428571428, 0.6933473389355742, 0.6930033416875522, 0.6932539682539683, 0.6930657506744464, done

収束はやい!!

なんか次は tableau だか,そういうデータ構造が紹介されている.このへんよくわからん!!

(define (make-tableau transform s)
  (cons-stream s
               (make-tableau transform
                             (transform s))))
(define (accelerated-sequence transform s)
  (stream-map stream-car
              (make-tableau transform s)))

これつかったらもっと収束早い!!

(dump-stream (accelerated-sequence euler-transform ln2-stream) 10)
; 1.0, 0.7, 0.6932773109243697, 0.6931488693329254, 0.6931471960735491, 0.6931471806635636, 0.6931471805604039, 0.6931471805599445, 0.6931471805599427, 0.6931471805599454, done