読者です 読者をやめる 読者になる 読者になる

Haskellでエラトステネスの篩

優先順位つきキューを用いたエラトステネスの篩の変形がある。
優先順位つきキューの削除と挿入のコストをO(log n)とすると、エラトステネス全体のオーダーはO(n log^2 n)である。

追記:以下は整数値オーバーフローのため値がおかしい。かつ、修正すると低速になる可能性が高い。

(奇数だけ走査するなどすればさらに速くなるかもしれない)

Eratosthenes' sieve in Haskell and Data.Map.Map — Gist

import qualified Data.Map

sieve :: (Data.Map.Map Int [Int], Int) -> [Int]
sieve (m,i) =
        let v = Data.Map.lookup i m
        in case v of
              Nothing -> i : sieve (Data.Map.insert (i*i) [i] m,i+1)
              Just l -> sieve (foldl f m l, i+1)
        where
              f mm p = Data.Map.insertWith' (++) (i+p) [p] m

primes = sieve (Data.Map.empty, 2)

main :: IO ()
main = do
        putStrLn $ show $ primes !! 9999999

これは以下のC++コードを、Data.Map.Mapを用いてHaskellに移植したものである。
(Priority Queueは標準では用意されてないようなので、Mapで代用した)

#include <cstdio>
#include <queue>
#include <algorithm>
using namespace std;

int main(int argc, char **argv) {
    priority_queue<pair<int,int> > pq;
    printf("%d\n", 2);
    pq.push(make_pair(-4, 2));
    for(int i = 2; i < 100; ) {
        int x = -pq.top().first;
        int q = pq.top().second;
        pq.pop();
        if(x-i==2) {
            int p = i+1;
            printf("%d\n", p);
            pq.push(make_pair(-p*p, p));
        }
        pq.push(make_pair(-x-q, q));
        i = x;
    }
    return 0;
}

なお、上記のHaskellコードを手元で実行したところ、204秒ほどの時間がかかった。参考程度に。