円周率の計算

30分プログラム、その573。円周率を計算してみる。
SICP曰く
\frac{1}{1 \cdot 3} + \frac{1}{5 \cdot 7} + \frac{1}{9 \cdot 11} + \cdot\cdot\cdot
の無限級数π/8が計算できるらしい。
無限ストリームを使って、だんだんと精度が上っていくπのストリームを作ってみた。

使い方

gosh> (stream-for-each print pi)
0.0
2.6666666666666665
2.895238095238095
2.976046176046176
3.017071817071817
3.0418396189294024
3.058402765927332
3.0702546177791836
3.079153394197426
3.0860798011238333
3.0916238066678385
3.096161526463641
3.099944032373807
3.1031453128860114
3.105889738271946
3.1082685666989462
3.110350273698686
3.1121872426998345
3.113820229023573
3.1152814162381857
3.1165965567938323

ソースコード

#! /opt/local/bin/gosh
;; -*- mode:scheme; coding:utf-8 -*-
;;
;; pi.scm -
;;
;; Copyright(C) 2009 by mzp
;; Author: MIZUNO Hiroki / mzpppp at gmail dot com
;; http://howdyworld.org
;;
;; Timestamp: 2009/05/01 21:26:42
;;
;; This program is free software; you can redistribute it and/or
;; modify it under MIT Lincence.
;;

(use util.stream)

(define pi/8
  (stream-unfoldn (lambda (seed)
		    (let* [(n   (car seed))
			   (sum (cdr seed))
			   (a   (+ 1 (* 4 n)))]
		      (values (cons (+ 1 n)
				    (+ sum (/ 1 (* a (+ a 2)))))
			      (list sum))))
		  '(0 . 0) 2))
(define pi
  (stream-map (cut *. 8 <>)
	      pi/8))