2日目

1.2.6 素数

最小の除数を求める。この通りにやるより、除数のリストを作るほうが簡単かな。

divisers :: Int -> [Int]
divisers n = [ x | x <- [2 .. limit], n `mod` x == 0]
             where limit = floor (sqrt (fromIntegral n))

isPrime :: Int -> Bool
isPrime n = case (divisers n) of
                (x:_) -> False
                _     -> True

そうでもないか。isPrime では divisers の1つ目の要素だけが評価されるようにしたつもりだけど、どうだろう。
次は確率的手法だけど、予想通り乱数を使うとモナドが出てくるみたいで、うむむ。まずこれを読もう -> All About Monads

理解しないまま System.Random を説明する

モナドのすべてのI部を読んだ。うーん。とりあえずやってみよう。
Haskell で標準に用意されている乱数生成機構は System.Random。乱数を得るには生成機 RandomGen が必要で、グローバルな生成機は getStdGen :: IO StdGen で取得できる。StdGen は RandomGen クラスのインスタンス。それを使って例えばこうする。

import System.Random
randomIOInt :: IO Int
randomIOInt = do gen <- getStdGen
                 let (r, g) = random gen
                 setStdGen g
                 return r

で、randomIOInt >>= putStr.show とかやると乱数が表示される。関数 random は Random クラスに書いてあって、Int や Float などがインスタンスになってる。IO (Int) と指定しているので r が Int と判断してくれているらしい。賢いね。あと、setStdGen g という風にグローバルな生成機をセットしなおさないと状態が変わらなくて延々と同じ数字が出るみたい。で実は、

getStdRandom :: (StdGen -> (a, StdGen)) -> IO a

っていう関数を使うと明示的に StdGen を使わなくても同じことができる。そうすると

randomIOInt = getStdRandom random

だけで済む。もっと言うと Random クラスの関数 randomIO はこれそのものなので、自分で関数を書かなくても (randomIO :: IO Int) >>= putStr.show とやれば乱数を出力できる。
乱数の無限配列を使うといいとGoogle様のお告げにあったので、SICP の問題用に範囲つきの整数乱数列を作ると

randomRInts :: (Int, Int) -> IO [Int]
randomRInts range = do gen <- getStdGen
                       return (randomRs range gen)

何度やっても同じ数列がでるけどまあいいや。これを使って p29 の Fermat テストを書くと、

expmod :: Int -> Int -> Int -> Int
expmod base exp m
       | exp == 0 = 1
       | even exp = (expmod base (exp `div` 2) m)^2 `mod` m
       | otherwise = (expmod base (exp - 1) m) * base `mod` m

fermatTest :: Int -> Int -> Bool
fermatTest n test = expmod test n n == test

isPrimeFast :: Int -> Int -> IO Bool
isPrimeFast n times = do rands <- randomRInts (1, n - 1)
                         return (all (fermatTest n) (take times rands))

最初 fermatTest の

expmod test n n

test^n `mod` n

と書いたら、オーバーフローするらしくしばらくはまった。やってることは同じに見えるけどなぁ。なぞ。構造をみてみると、乱数生成は基本のアルゴリズムに内包できないのでこうしたんだけど、(1, n - 1) の範囲を isPrimeFast の方で指定したりして、処理が分散しているように思える。うーん。

問題1.25

あー expmod と base^exp `mod` m は全然違うね。一目で気付け。