在Haskell的Prime Sieve

Question

我对Haskell很新,我只想找到前200万个素数的总和.我正试图用筛子生成素数(我想是Eratosthenes的筛子？),但它真的很慢,我不知道为什么.这是我的代码.

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

提前致谢.

Answer 1

它非常慢,因为算法是一个试验分区,不会停留在平方根.

如果仔细观察算法的作用,可以看出,对于每个素数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,因为复合材料通常很便宜)⁽¹⁾并且列表操作少得多.此外,通过跳过偶数(除了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

⁽¹⁾残酷的估计是

2000000
   ? ?k ? 4/3*?2*10^9
 k = 1

评估为两位有效数字.由于大多数数字都是具有较小素数因子的复合数 - 一半数字是偶数且只占一个分数 - 这大大高估了所需的分割数.

通过仅考虑质数,可以获得所需分区数的下限:

   ? ?p ? 2/3*N^1.5/log N
 p < N
p prime

其中,N = 2000000约为1.3*10 ⁸.这是正确的数量级,但低估了一个非平凡的因素(对于增长缓慢减少到1 N,从不超过2 N > 10).

除了素数之外,素数的平方和两个接近素数的乘积也要求试验除法达到(几乎)?k,因此如果有足够的数量,则对整体工作有显着贡献.

然而,处理半数所需的分割数量受到常数倍的限制

N^1.5/(log N)^2

因此,对于非常大的N,相对于处理素数的成本而言,它变得可以忽略不计.但是在试验部门完全可行的范围内,它们仍然有很大贡献.

Answer 2

这是Eratosthenes的一个筛子:

P = { 3,5,...}\⋃ {{ p ²,p ² + 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 

与奥古斯都的答案中的试验分区代码相比,它在产生200k质数时快1.9倍,在400k时快2.1倍,经验时间复杂度为O(n^1.12..1.15)vs O(n^1.4),在所述范围内.它产生1百万个素数的速度快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 == 95,或更不用5或7.

这是该代码中唯一最严重的问题,尽管在每篇人都提到的那篇文章的开头它都被认为是不相关的.通过推迟过滤器的创建来修复它可以实现显着的加速.