2010-12-28

フィボナッチ数列は、自然界に現れる単純にして基本的な数列である。たとえば、ヒマワリの種の並びやミツバチの家系図はフィボナッチ数列を成す。

再帰

フィボナッチ数列の漸化式をそのまま Haskell で実装すると以下のようになる。

fib :: Int -> Integer
fib 0 = 0
fib 1 = 1
fib n = fib (n - 1) + fib (n - 2)

fib 0 = 0 の代わりに fib 2 = 1 とする流儀もある。フィボナッチ数列は、もともとウサギのつがいの話だから後者の方が正統だが、前者の方が数学的には美しい。

漸化式は、すなわちトップダウン的な表現であり、理解しやすい。しかし、この実装は同じ値を何度も計算するので、とても遅い。

反復

単純な再帰を反復(ループ)に変形してみよう。それには、関数が末尾再帰の形になればよい。末尾再帰へ変換する定石は、引数に蓄積変数を追加することである。今回は、蓄積変数が2つ必要になる。

fib :: Int -> Integer
fib n = iter 1 1 0
  where
    iter m a1 a2
      | n == m    = a1
      | otherwise = iter (m + 1) (a1 + a2) a1

この実装は速い。今回は m を増加するカウンターとして使っているが、減少するカウンターを使っても実装できる。興味があれば考えてほしい。

反復はボトムアップな解決方法である。

汎関数

漸化式を素直に実装した版を高速にするには、メモ化が使える。メモ化のために、まず不動点コンビネータとフィボナッチ数列の汎関数を用意しよう。

fix :: ((Int -> Integer) -> Int -> Integer) -> Int -> Integer
fix f x = f (fix f) x

fib :: (Int -> Integer) -> Int -> Integer
fib _ 0 = 0
fib _ 1 = 1
fib f n = f (n-1) + f (n-2)

fib は fib という名前では再帰していないことに注意。以下のように動く。

fix fib 10
→ 55

もちろん、このままだと遅い。

メモ化(1)

不動点コンビネータにメモ化の機能を入れる。ここでは、分かりやすさのために禁断の unsafePerformIO 用いて副作用を発生させる。

import Data.IORef
import qualified Data.Map as M
import System.IO.Unsafe

mfix :: ((Int -> Integer) -> Int -> Integer) -> Int -> Integer
mfix g y = fix g y
  where
    memo = unsafePerformIO (newIORef M.empty)
    fix f x = unsafePerformIO $ do 
        m <- readIORef memo
        case M.lookup x m of
            Just v -> return v
            Nothing -> do
                let v = f (fix f) x
                modifyIORef memo (M.insert x v)
                return v

fib :: (Int -> Integer) -> Int -> Integer
fib _ 0 = 0
fib _ 1 = 1
fib f n = f (n-1) + f (n-2)

mfix を使うと、メモ化が効いて高速になる。

mfix fib 10
→ 55

メモ化(2)

HackageDB には、メモ化のパッケージとして data-memocombinators が登録されている。これを使えば、フィボナッチ数列を以下のように実装できる。

import qualified Data.MemoCombinators as Memo

fib :: Int -> Integer
fib = Memo.integral fib'
  where
    fib' 0 = 0
    fib' 1 = 1
    fib' n = fib (n - 1) + fib (n - 2)

メモ化(3)

フィボナッチ数列に対しては、遅延評価を用いたメモ化の方法が知られている。この実装は、芸術的である。

fib :: Int -> Integer
fib n = fibs !! n
  where
    fibs = 0 : 1 : zipWith (+) fibs (tail fibs)

分割統治

フィボナッチ数列は、分割統治することで、対数時間で計算できることが知られている。以下のようなコードになる。

fib :: Int -> Integer
fib = iter 1 0 0 1
  where
    iter a b p q count
      | count == 0 = b
      | even count = iter a b (p*p+q*q) (2*p*q+q*q) (count `div` 2)
      | otherwise  = iter (b*q+a*q+a*p) (b*p+a*q) p q (count - 1)

詳細を知りたければ、SICP の問題1.19を参照のこと。

最速

フィボナッチ数列は、黄金比と密接な関係があり、黄金比からの単純な計算で求められる。

phy :: Double
phy = (1 + sqrt 5) / 2

phy' :: Double
phy' = (1 - sqrt 5) / 2

fib :: Int -> Double
fib n = (phy ^ n - phy' ^ n) / sqrt 5

この式はダニエル・ベルヌーイが1728年に発見した。一般的には、母関数を使うと導出できる。

この実装では、答えが Double になってしまい、欲しい関数とはいえない。幸運にも、n 番目のフィボナッチ数列が黄金比の n 乗を sprt 5 で割った数に近いという性質を利用して、Integer を求める式が知られている。

phy :: Double
phy = (1 + sqrt 5) / 2

fib :: Int -> Integer
fib n = truncate $ (phy ^ n) / sqrt 5 + 1 / 2

あどけない話

Internet technologies

フィボナッチ数列のまとめ

再帰

反復

汎関数

メモ化(1)

メモ化(2)

メモ化(3)

分割統治

最速

あわせて読みたい