我对 Haskell 很陌生，我只是想找到前 200 万个素数的总和。我正在尝试使用筛子生成素数（我认为是 Eratosthenes 的筛子？），但它真的很慢，我不知道为什么。这是我的代码。
sieve (x:xs) = x:(sieve $ filter (\a -> a `mod` x /= 0) xs)
ans = sum $ takeWhile (<2000000) (sieve [2..])
提前致谢。
问问题
10913 次
5 回答
16
它非常慢，因为该算法是一个不停留在平方根的试除法。
如果您仔细查看算法的作用，您会看到对于每个 prime p，其不具有较小素因数的倍数将从候选列表中删除（先前已删除具有较小素因数的倍数）。
因此，每个数字都被所有素数除，直到它作为其最小素数除数的倍数被删除，或者如果它是素数，它出现在剩余候选者列表的开头。
对于合数来说，这并不是特别糟糕，因为大多数合数都有小的素数除数，在最坏的情况下，最小的素数除数n不会超过√n.
但是素数被所有更小的素数整除，所以直到第 k个素数被发现是素数，它已经被所有k-1更小的素数整除。如果存在m低于极限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）作为候选除数，可以轻松改进此试验除法，这将所需除法的数量减半。
Eratosthenes 的筛子不需要任何除法，只使用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
评估为两位有效数字。由于大多数数字是具有小素因数的复合数 - 一半的数字是偶数并且只需要一个除法 - 这大大高估了所需的除法数。
通过仅考虑素数可以获得所需除法数的下限：
∑ √p ≈ 2/3*N^1.5/log N
p < N
p prime
其中， forN = 2000000大约给出 1.3*10 8。这是正确的数量级，但低估了一个重要的因素（增长缓慢下降到 1 N，并且永远不会大于 2 N > 10）。
除了素数之外，素数的平方和两个相近素数的乘积也需要试除法上升到（几乎）√k，因此如果有足够多的话，对整体工作的贡献很大。
然而，处理半素数所需的除法数受以下常数倍数的限制
N^1.5/(log N)^2
所以对于非常大N的情况，相对于处理素数的成本来说，它变得可以忽略不计。但在试验划分完全可行的范围内，它们仍然有很大的贡献。
于 2012-08-02T01:31:18.293 回答
6
这是埃拉托色尼的筛子：
P = { 3 , 5 , ...} \ ⋃ {{ p 2 , p 2 +2p , ...} | p中的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)vs O(n^1.4)，在所述范围内。生成 100 万个素数时速度快 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)。但是 8 还不需要检查是否可以被 3 整除，因为它小于3^2 == 9，更不用说被 5 或 7 整除了。
这是该代码中最严重的问题，尽管在每个人都提到的那篇文章的开头，它被认为是无关紧要的。通过推迟过滤器的创建来修复它可以实现显着的加速。
于 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 回答