Haskellで黄金比を求めよう


2012年 11月 24日

RenbunsuunoFushigi.jpgHaskellで、foldrについて考えていたら、なぜか BLUE BACKS の『連分数の不思議』が机の上にあった。
表紙(帯)には、黄金比を連分数で示していた。

黄金比Φ= 1+1/(1+1/(1+1/(1+1/(1+1/(1+…..

実際には、これをちゃんと分数の形に書くと、とても美しい連分数になる。

この連分数を見ていると、右の方から畳み込んでいることが分かる。
つまり、これって、foldr そのものではないだろうか。

ここで、連分数を、もうすこし一般的な、正則連分数について考えることにしよう。

正則連分数 = a0+1/(a1+1/(a2+1/(a3+1/(a4+1/(a5+…..

という連分数だ。分数のまま書くのは大変なので、[a0;a1,a2,a3,a4,…]という表記方法がある。

黄金比Φ の場合、 [1;1,1,1,…]  となる。

さて、[a0;a1,a2,a3,a4,…]をリストとして作り、これからfoldrを使って連分数を計算してみよう。

連分数の分数1回分の演算を★で表すと、

x ★ y ->  x + 1 / y

と書けるのだ。詳しくは、自分で確認しよう。

これを1つの補助関数として定義しても良いのだが、この程度は1つの関数の中に書き込んでしまいたい。
そんなことをするために、名前をつけなくても使える無名(ラムダ)関数というものがRealWorldHaskell.jpgある。
このあたりの知識は、オライリーの ”Real World Haskell” から得たものだ。
こちらの本も、おなじ訳者だ。

名前がない関数なので、名前を書かなくても良い場所に書かねばならない。
foldr の最初の引数は関数なので、ここに書いてしまえば良いはずだ。

foldr (\x y -> x+1/y)

さて、この演算を畳み込むときの最初の値 e は何が妥当であろうか。
0にしてしまうと、最初に無名関数のyのところに使われるのでまずい。
ということで、1にしてみよう。

foldr (\x y -> x+1/y) 1

あとは、1が延々と続くリストを書けば完成だ。

*Main Data.Ratio> foldr (\x y -> x+1/y) 1 [1,1..]
  C-c C-cInterrupted.
*Main Data.Ratio>

しかし、いつまで経っても帰ってこないので殺した。

[1,1…] という無限リストを与えたのがまずかったようだ。
リストの長さを短くして様子を見よう。

*Main Data.Ratio> foldr (\x y -> x+1/y) 1 [1,1]
1.5
*Main Data.Ratio> foldr (\x y -> x+1/y) 1 [1,1,1]
1.6666666666666665
*Main Data.Ratio>

連分数なのに、分数ではなく浮動小数点という誤差を伴う計算に成り下がってしまった。
有理数(Rational)で計算しなくてはダメだ。
どうすれば出来るだろうか?