円周率を連分数展開で求める

仕事にはまったく役に立ちませんが、円周率を連分数展開で求めてみます。

円周率の連分数展開

https://www.timedia.co.jp/tech/wikipediahaskell-haskell13-let-ratsqrt2-n-foldr-x-y/
書くのが面倒なので他人任せ。

Scheme での実装例

参考文献　計算機プログラムの構造と解釈 問題1.37 1.38 1.39
https://sicp.iijlab.net/fulltext/x133.html

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

(define (cont-frac n d k)
  (define (iter i result)
    (if (= i 0)
        result
        (iter (- i 1) (/ (n i) (+ (d i) result)))))
  (iter k 0))

(define (pi k)
  (define (n i)
    (print i)
    (square (- (* 2.0 i) 1)))
  (define (d i) 6.0)
  (+ 3 (cont-frac n d k)))

cont-frac 手続きは分数になっている部分を求めますから、最初の 3 は計算結果に別途加える必要があります。

実行結果

gosh> (pi 10000)
3.1415926535895435