10

私は Haskell の初心者で、最初の 200 万個の素数の合計を見つけようとしています。ふるいを使って素数を生成しようとしていますが (エラトステネスのふるいだと思いますか?)、実際には非常に遅く、その理由がわかりません。これが私のコードです。

sieve (x:xs) = x:(sieve $ filter (\a -> a `mod` x /= 0) xs)
ans = sum $ takeWhile (<2000000) (sieve [2..])

前もって感謝します。

4

5 に答える 5

16

アルゴリズムが平方根にとどまらない試行分割であるため、非常に遅いです。

アルゴリズムの動作をよく見ると、素数ごとpに、それより小さい素約数を持たないその倍数が候補のリストから削除されていることがわかります (小さい素数約数を持つ倍数は以前に削除されました)。

したがって、各数値は、最小の素因数の倍数として削除されるか、素数の場合は残りの候補のリストの先頭に表示されるまで、すべての素数で除算されます。

合成数の場合、ほとんどの合成数には小さな素因数があり、最悪の場合、 の最小の素因数nは を超えないため、これは特に悪いことではありません√n。

しかし、素数はすべての小さい素数で除算されるため、k番目の素数が素数であることが判明するまで、それはすべてのk-1小さい素数で除算されています。mlimit 未満の素数がある場合n、それらすべての素数を見つけるために必要な作業は次のとおりです。

(1-1) + (2-1) + (3-1) + ... + (m-1) = m*(m-1)/2

部門。素数定理により、以下の素数の数nは漸近的ですn / log n(ここでlogは自然対数を表します)。複合体を除去する作業は、n * √n分割によって大まかに制限される可能性があるためn、素数に費やされる作業と比較して無視できるほど小さいわけではありません。

200 万までの素数の場合、ターナーふるいは約 10 10の分割を必要とします。さらに、多くのリスト セルを分解して再構築する必要があります。

平方根で止まる試行分割、

isPrime n = go 2
  where
    go d
      | d*d > n        = True
      | n `rem` d == 0 = False
      | otherwise      = go (d+1)

primes = filter isPrime [2 .. ]

必要な分割数は 1.9*10 9未満 (すべての isPrime nチェックが行われた場合の残忍な見積もり√n- 実際には、コンポジットは一般に安価であるため、179492732 しかかからない) (1)であり、リスト操作ははるかに少なくなります。さらに、この試行分割は、除数候補として偶数 ( を除く2) をスキップすることで簡単に改善できます。これにより、必要な分割数が半分になります。

エラトステネスのふるいは分割を必要とせず、O(n * log (log n))操作のみを使用します。これはかなり高速です。

primeSum.hs:

module Main (main) where

import System.Environment (getArgs)
import Math.NumberTheory.Primes

main :: IO ()
main = do
    args <- getArgs
    let lim = case args of
                (a:_) -> read a
                _     -> 1000000
    print . sum $ takeWhile (<= lim) primes

そして、1000万の制限でそれを実行します:

$ ghc -O2 primeSum && time ./primeSum 10000000
[1 of 1] Compiling Main             ( primeSum.hs, primeSum.o )
Linking primeSum ...
3203324994356

real    0m0.085s
user    0m0.084s
sys     0m0.000s

試行分割を 100 万まで実行します (型を に固定Int):

$ ghc -O2 tdprimeSum && time ./tdprimeSum 1000000
[1 of 1] Compiling Main             ( tdprimeSum.hs, tdprimeSum.o )
Linking tdprimeSum ...
37550402023

real    0m0.768s
user    0m0.765s
sys     0m0.002s

そして、ターナーふるいは 100000 までしかありません。

$ ghc -O2 tuprimeSum && time ./tuprimeSum 100000
[1 of 1] Compiling Main             ( tuprimeSum.hs, tuprimeSum.o )
Linking tuprimeSum ...
454396537

real    0m2.712s
user    0m2.703s
sys     0m0.005s

(1)残忍な見積もりは

2000000
   ∑ √k ≈ 4/3*√2*10^9
 k = 1

2 桁の有効数字に評価されます。ほとんどの数は小さな素因数を持つ合成であるため (半分の数は偶数であり、1 つの除算しか必要としません)、必要な除算の数を大幅に過大評価します。

必要な分割数の下限は、素数のみを考慮することによって得られます。

   ∑ √p ≈ 2/3*N^1.5/log N
 p < N
p prime

これにより、N = 2000000おおよそ 1.3*10 8が得られます。これは正しい大きさのオーダーですが、自明でない要因によって過小評価されています (成長する では 1 までゆっくりと減少しN、 では 2 を超えることはありませんN > 10)。

素数に加えて、素数の 2 乗と 2 つの近接した素数の積も (ほぼ) まで試行分割を行う必要があるため、√k十分な数がある場合は全体の作業に大きく貢献します。

ただし、半素数を処理するために必要な除算の数は、次の定数倍によって制限されます。

N^1.5/(log N)^2

そのため、非常に大きいN場合、素数を処理するコストに比べて無視できます。しかし、試行分割がまったく可能な範囲では、それらは依然として大きく貢献しています。

于 2012-08-02T01:31:18.293 に答える
6

これがエラトステネスのふるいです：

P = { 3、5 、 ... }\ ⋃ {{ p2、p 2 + 2p、...} | p in P }

（2なし）。:)または「機能的」、つまりリストベースのHaskellでは、

primes = 2 : g (fix g)  where
   g xs = 3 : (gaps 5 $ unionAll [[p*p, p*p+2*p..] | p <- xs])

unionAll ((x:xs):t) = x : union xs (unionAll $ pairs t)  where
  pairs ((x:xs):ys:t) = (x : union xs ys) : pairs t 
fix g = xs where xs = g xs
union (x:xs) (y:ys) = case compare x y of LT -> x : union  xs (y:ys)
                                          EQ -> x : union  xs    ys 
                                          GT -> y : union (x:xs) ys
gaps k s@(x:xs) | k<x  = k:gaps (k+2) s    
                | True =   gaps (k+2) xs

augustssによる回答の試行割り法コードと比較すると、200k素数の生成で1.9倍、400kで2.1倍高速であり、経験的な時間計算量はO(n^1.12..1.15)vsO(n^1.4)です。1mlnの素数を生成する場合は2.6倍高速です。

ターナーふるいがとても遅い理由

なぜなら、それは倍数を開くからです-各プライムのストリームをフィルタリングするのが早すぎて、それらの数が多すぎてしまいます。入力にその正方形が表示されるまで、素数でフィルタリングする必要はありません。

ストリーム処理パラダイムの下で見sieve (x:xs) = x:sieve [y|y<-xs, rem y p/=0]られるように、それが機能しているときに、それ自体の背後にストリームトランスデューサーのパイプラインを作成していると見なすことができます。

[2..] ==> sieve --> 2
[3..] ==> nomult 2 ==> sieve --> 3
[4..] ==> nomult 2 ==> nomult 3 ==> sieve 
[5..] ==> nomult 2 ==> nomult 3 ==> sieve --> 5
[6..] ==> nomult 2 ==> nomult 3 ==> nomult 5 ==> sieve 
[7..] ==> nomult 2 ==> nomult 3 ==> nomult 5 ==> sieve --> 7
[8..] ==> nomult 2 ==> nomult 3 ==> nomult 5 ==> nomult 7 ==> sieve

ここでnomult p = filter (\y->rem y p/=0)。3^2 == 9ただし、8は、5または7は言うまでもなく、よりも小さいため、まだ3で除算可能かどうかをチェックする必要はありません。

これは、そのコードで最も深刻な問題の1つですが、誰もが言及している記事の冒頭で無関係として却下されています。フィルタの作成を延期して修正すると、劇的なスピードアップが実現します。

于 2012-08-02T16:05:30.630 に答える
3

あなたがしたことは、エラトステネスのふるいではありません。それは試行分割です (mod オペレーターに注意してください)。エラトステネスのふるいの私のバージョンは次のとおりです。

import Control.Monad (forM_, when)
import Control.Monad.ST
import Data.Array.ST
import Data.Array.Unboxed

sieve :: Int -> UArray Int Bool
sieve n = runSTUArray $ do
    let m = (n-1) `div` 2
        r = floor . sqrt $ fromIntegral n
    bits <- newArray (0, m-1) True
    forM_ [0 .. r `div` 2 - 1] $ \i -> do
        isPrime <- readArray bits i
        when isPrime $ do
            forM_ [2*i*i+6*i+3, 2*i*i+8*i+6 .. (m-1)] $ \j -> do
                writeArray bits j False
    return bits

primes :: Int -> [Int]
primes n = 2 : [2*i+3 | (i, True) <- assocs $ sieve n]

http://ideone.com/mu1RNで実行できます。

于 2012-08-01T23:57:19.607 に答える
2

個人的には、素数を生成するこの方法が好きです

primes :: [Integer]
primes = 2 : filter (isPrime primes) [3,5..]
  where isPrime (p:ps) n = p*p > n || n `rem` p /= 0 && isPrime ps n

また、ここで提案されている他のいくつかの方法に比べて非常に高速です。まだ試行分割ですが、素数でのみテストします。(ただし、このコードの終了証明は少しトリッキーです。)

于 2012-08-02T09:00:57.660 に答える
1

使用しているアルゴリズムはまったくふるいではありません。したがって、遅いという点では、試行分割を使用することを期待する必要があります。

素数は、平方根関数の頻度で大まかに発生しています...つまり、 1 とnの間に球場n/log(n)素数があります。したがって、最初の 200 万の素数では、最大で約 3200 万が必要になります。しかし、これらの素数が通過する 200 万要素のデータ構造を構築しています。これで、なぜこれが非常に遅かったのかがわかります。実際は O(n^2) です。O(n*(log n)*log(log n)) に削減できます

これは、さまざまな治療法に関するページで、それを少し減らす方法を説明しています. http://en.literateprograms.org/Sieve_of_Eratosthenes_(Haskell)

于 2012-08-02T00:21:48.450 に答える