Haskellは分数がお得意(eの計算)


2012年 11月 22日

『関数プログラミング』3.5 有理数
Haskells.jpgのサムネール画像
foldr について、格好良い使い方を考える前に、ちょっとだけ前に戻ろう。

まず、1 / 2 + 1 / 3 という分数計算をやってみよう。

Prelude Data.Ratio> 1 / 2 + 1 / 3
0.8333333333333333

浮動小数点数になる。実際に型を確認すると

Prelude Data.Ratio> :t 1 / 2 + 1 / 3
1 / 2 + 1 / 3 :: Fractional a => a

ところで、3.5 例:有理数 には、有理数 Ratioanl の定義などが載っている。

しかし、最新のGHCiなら、何も自分で用意しなくたって、インポートするだけで使える。
ネット上には、自分で有理数を用意する方法がいっぱい転がっているが、清く正しく、GHCiに用意されているのを使ってみよう。
より横着な方法を考えるのが、プログラマのあるべき姿だ。

ということで、用意されているものを利用できるようにしよう。

Prelude> import Data.Ratio
Prelude Data.Ratio>

import の代わりに、 :m Data.Ratio でもよい。

すると、分数同士の演算が分数のまま可能になる。

さて、 / の代わりに % を使って同じことをしてみよう。

Prelude Data.Ratio> 1 % 2 + 1 % 3
5 % 6
Prelude Data.Ratio> :t 1 % 2 + 1 % 3
1 % 2 + 1 % 3 :: Integral a => Ratio a
Prelude Data.Ratio>

% で、割り算結果が分数のままになる。

ということで、分数のまま、ちょっと級数展開の計算をしてみよう。

簡単なところで、

*Main Data.Ratio> exp 1
2.718281828459045

を、級数展開の和を実際に計算して求めることにしよう。

e = Σ1/n!  (n=0..) で求まることになっているので、
最初のn項までの和を求めてみよう。
この級数は、結構収束が良いはずだ。

ずは、階乗の計算だが、今回は foldr を使って作ってみた。
foldr を使えば、再帰の必要も無いのだ。

fact :: Integer -> Integer
fact n = foldr (*) 1 [1..n]

*Main Data.Ratio> fact 10
3628800
*Main Data.Ratio> fact 20
2432902008176640000

となって、動いているようだ。
後は、だらだらと級数展開の通りに書いてみた。
右のmap で、階乗を求め、次のmap で逆数にし、foldr で全部加えた。

import Data.Ratio

expseq    :: Integer -> Rational
expseq  n =  foldr (+) 0 (map (1 %) (map (fact) [0..n]))

fact :: Integer -> Integer
fact n = foldr (*) 1 [1..n]

さて、実行してみよう。

*Main Data.Ratio> expseq 3
8 % 3
*Main Data.Ratio> expseq 4
65 % 24

分数だとよく分からないが、差をとると、1つの項が出るはずだ。

*Main Data.Ratio> (expseq 4) - (expseq 3)
1 % 24
*Main Data.Ratio> (expseq 5) - (expseq 4)
1 % 120

これで、動作確認もできた訳だが、俗人は浮動小数点が好きらしいので、直してみよう。

*Main Data.Ratio> fromRational (expseq 4)
2.7083333333333335
*Main Data.Ratio> fromRational (expseq 5)
2.716666666666667
*Main Data.Ratio> fromRational (expseq 6)
2.7180555555555554
*Main Data.Ratio> fromRational (expseq 10)
2.7182818011463845
*Main Data.Ratio> fromRational (expseq 15)
2.7182818284589945
*Main Data.Ratio> fromRational (expseq 20)
2.718281828459045
*Main Data.Ratio>

ということで、e が計算できた。
実際の計算では、精度のことも考えて、どこで打ち切るかを考えねばならぬだろうが、省略する。

有理数として計算をしているので、foldl, foldr のどちらを使っても結果に誤差は出ない。
まだ、foldrの持ち味がほとんど出ていない。
ちゃんと foldr が生きてくる例を考えてみよう。