みずぴー日記

円周率を計算しよう

Haskell 30分プログラム

30分プログラム、その611。http://gauc.no-ip.org/awk-users-jp/blis.cgi/DoukakuAWK_111にインスパイアされたので、円周率を計算してみよう。
元の記事と同様にwikipedia:ガウス＝ルジャンドルのアルゴリズムを使うのだけれども、このアルゴリズムはパラメータが4つあって、それを各ステップで受け渡す必要がある。4つもあると、タプルで受け渡すのは厳しいので、パラメータ用のレコードを定義して、それを使うようにしてます。

使い方

*Main> take 10 piList
[2.9142137,3.1405795,3.1415927,3.1415927,3.1415927,3.1415927,3.1415927,3.1415927,3.1415927,3.1415927]

ソースコード

data Param = Param { a :: Float, b :: Float, t :: Float, p :: Float} deriving Show

toPI :: Param -> Float
initial :: Param
improve :: Param -> Param
piList :: [Float]

toPI param = (a param + b param) ** 2 / (4 * t param)
initial = Param {a = 1, b = sqrt 0.5, t = 0.25, p = 1}
improve param = Param {a = x,
                       b = y,
                       t = t param - (p param)*(a param - x)**2,
                       p = 2 * p param}
                where x = (a param + b param) / 2
                      y = sqrt $ a param * b param

piList = map toPI $ iterate improve initial

参考