为什么这个主要测试如此缓慢？

Question

此代码取自"Haskell之路逻辑,数学和编程"一书.它实现了eratosthenes算法的筛选并解决了Project Euler Problem 10.

sieve :: [Integer] -> [Integer]
sieve (0 : xs) = sieve xs
sieve (n : xs) = n : sieve (mark xs 1 n)
  where
    mark :: [Integer] -> Integer -> Integer -> [Integer]
    mark (y:ys) k m | k == m = 0 : (mark ys 1 m)
                    | otherwise = y : (mark ys (k+1) m)

primes :: [Integer]
primes = sieve [2..]

-- Project Euler #10
main = print $ sum $ takeWhile (< 2000000) primes

实际上它比天真的素数测试运行得更慢.有人可以解释这种行为吗？

我怀疑它与迭代标记函数中列表中的每个元素有关.

谢谢.

Answer 1

您正在使用此算法构建二次数的未评估的thunk.该算法如此严重依赖于懒惰,这也是它无法扩展的原因.

让我们来看看它是如何工作的,这有望使问题显而易见.简单地说,让我们说我们想要无限制print的元素primes,即我们想要一个接一个地评估列表中的每个单元格.primes定义为:

primes :: [Integer]
primes = sieve [2..]

因为2不为0,的第二个定义sieve适用,和2加入到素数的列表,并且该列表的其余部分是一个未计算的形实转换(I使用tail模式匹配,而不是n : xs在sieve对xs,因此tail是不实际被调用,并没有在下面的代码中添加任何开销; mark实际上是唯一的thunked函数):

primes = 2 : sieve (mark (tail [2..]) 1 2)

现在我们想要第二个primes元素.因此,我们遍历代码(为读者练习)并最终得到:

primes = 2 : 3 : sieve (mark (tail (mark (tail [2..]) 1 2)) 1 3)

同样的程序,我们想评估下一个素数......

primes = 2 : 3 : 5 : sieve (mark (tail (tail (mark (tail (mark (tail [2..]) 1 2)) 1 3))) 1 5)

这开始看起来像LISP,但我离题了...开始看到问题？对于primes列表中的每个元素,mark必须评估越来越多的应用程序堆栈.换句话说,对于列表中的每个元素,必须通过评估mark堆栈中的每个应用程序来检查该元素是否由任何前面的素数标记.因此,对于n~=2000000,Haskell运行时必须调用函数,导致调用堆栈的深度约为......我不知道,137900(let n = 2e6 in n / log n给出下限)？这样的事情.这可能是导致减速的原因; 也许vacuum可以告诉你更多(我现在没有配备Haskell和GUI的计算机).

Eratosthenes筛子在C语言中工作的原因是:

1. 您没有使用无限列表.
2. 由于(1),您可以在继续下一个之前标记整个数组n,从而导致根本没有调用堆栈开销.

Answer 2

不仅是让它变得非常慢的thunk,如果在有限位数组中用C实现,该算法也会非常慢.

sieve :: [Integer] -> [Integer]
sieve (0 : xs) = sieve xs
sieve (n : xs) = n : sieve (mark xs 1 n)
  where
    mark :: [Integer] -> Integer -> Integer -> [Integer]
    mark (y:ys) k m | k == m = 0 : (mark ys 1 m)
                    | otherwise = y : (mark ys (k+1) m)

对于每个素数p,此算法检查所有数字p+1是否为极限,是否为倍数p.它并没有像特纳的筛子那样划分,而是通过将计数器与素数进行比较.现在,比较两个数字比分割要快得多,但是为此付出的代价n是现在每个数字都要检查每个数字< n,而不仅仅是素数到达n最小的素数因子.

结果是该算法的复杂性为特纳筛的O(N ^ 2/log N)与O((N/log N)^ 2)(和O(N*log(log N))筛子Eratosthenes).

该 嵌套¹dflemstr提到的thunk堆叠加剧了问题²,但即便没有,算法也会比Turner更差.我同时感到震惊和着迷.

¹"嵌套"可能不正确.虽然每个markthunk只能通过它上面的那个来访问,但它们不会引用封闭thunk范围内的任何内容.

²虽然thunk的大小或深度没有任何二次方,但是thunk的表现相当不错.为了说明,我们假装mark是用反向参数顺序定义的.然后,当发现7是素数时,情况是

sieve (mark 5 2 (mark 3 1 (mark 2 1 [7 .. ])))
~> sieve (mark 5 2 (mark 3 1 (7 : mark 2 2 [8 .. ])))
~> sieve (mark 5 2 (7 : mark 3 2 (mark 2 2 [8 .. ])))
~> sieve (7 : mark 5 3 (mark 3 2 (mark 2 2 [8 .. ])))
~> 7 : sieve (mark 7 1 (mark 5 3 (mark 3 2 (mark 2 2 [8 .. ]))))

下一个模式匹配sieve强制mark 7 1thunk,强制mark 5 3thunk,强制mark 3 2thunk,强制mark 2 2thunk,强制[8 .. ]thunk并用0替换头,并将尾巴包裹在mark 2 1thunk中.这会冒泡到sieve,它会丢弃0,然后强制下一堆thunk.

因此,对于从每一个号码p_k + 1到p_(k+1)(含),所述模式匹配中sieve的堆叠/链力k的形式的的thunk mark p r.这些中的每一个都(y:ys)从封闭的thunk([y ..]最里面的mark 2 r)中获取并将尾部包裹ys在一个新的thunk中,要么保持y不变,要么将其替换为0,从而在尾部构建一个新的堆栈/ thunk串.列表到达sieve.

对于每个找到的素数,在顶部sieve添加另一个mark p rthunk,所以最后,当找到大于2000000的第一个素数并且takeWhile (< 2000000)完成时,将有148933个级别的thunk.

这里堆叠的thunk不会影响复杂性,它只会影响常数因素.在我们正在处理的情况下,一个懒惰生成的无限不可变列表,没有太多可以做的事情来减少将控制从一个thunk转移到下一个thunk所花费的时间.如果我们处理一个有限的可变列表或不是懒惰生成的数组,就像在C或Java这样的语言中那样,那么让每个mark p完成它的完整工作会更好(这将是一个简单的for循环,开销更少在检查下一个数字之前,调用/转移控制的功能,所以永远不会有多个标记激活且控制传递较少.

Answer 3

好吧,你肯定是对的,它比天真的实现慢.我从维基百科那里拿了这个,然后用GHCI将它与你的代码进行比较:

-- from Wikipedia
sieveW [] = [] 
sieveW (x:xs) = x : sieveW remaining 
  where 
    remaining = [y | y <- xs, y `mod` x /= 0]

-- your code
sieve :: [Integer] -> [Integer]
sieve (0 : xs) = sieve xs
sieve (n : xs) = n : sieve (mark xs 1 n)
  where
    mark :: [Integer] -> Integer -> Integer -> [Integer]
    mark (y:ys) k m | k == m = 0 : (mark ys 1 m)
                    | otherwise = y : (mark ys (k+1) m)

跑步给

[1 of 1] Compiling Main             ( prime.hs, interpreted )
Ok, modules loaded: Main.
*Main> :set +s
*Main> sum $ take 2000 (sieveW [2..])
16274627
(1.54 secs, 351594604 bytes)
*Main> sum $ take 2000 (sieve [2..])
16274627
(12.33 secs, 2903337856 bytes)

为了尝试了解正在发生的事情以及mark代码的工作原理,我尝试手动扩展代码:

  sieve [2..]
= sieve 2 : [3..]
= 2 : sieve (mark [3..] 1 2)
= 2 : sieve (3 : (mark [4..] 2 2))
= 2 : 3 : sieve (mark (mark [4..] 2 2) 1 3)
= 2 : 3 : sieve (mark (0 : (mark [5..] 1 2)) 1 3)
= 2 : 3 : sieve (0 : (mark (mark [5..] 1 2) 1 3))
= 2 : 3 : sieve (mark (mark [5..] 1 2) 1 3)
= 2 : 3 : sieve (mark (5 : (mark [6..] 2 2)) 1 3)
= 2 : 3 : sieve (5 : mark (mark [6..] 2 2) 2 3)
= 2 : 3 : 5 : sieve (mark (mark (mark [6..] 2 2) 2 3) 1 5)
= 2 : 3 : 5 : sieve (mark (mark (0 : (mark [7..] 1 2)) 2 3) 1 5)
= 2 : 3 : 5 : sieve (mark (0 : (mark (mark [7..] 1 2) 3 3)) 1 5)
= 2 : 3 : 5 : sieve (0 : (mark (mark (mark [7..] 1 2) 3 3)) 2 5))
= 2 : 3 : 5 : sieve (mark (mark (mark [7..] 1 2) 3 3)) 2 5)
= 2 : 3 : 5 : sieve (mark (mark (7 : (mark [8..] 2 2)) 3 3)) 2 5)

我想我可能在那里结束时犯了一个小错误,因为7看起来它将变成0并被删除,但机制很清楚.此代码仅创建一组计数器,计数到每个素数,在正确的时刻发出下一个素数并将其传递到列表中.这相当于仅仅在初始实现中检查每个先前素数的除法,以及在thunk之间传递0或素数的额外开销.

这里可能还有一些我想念的细微之处.有埃拉托色尼的Haskell中筛的各种最佳化沿着一个非常详细的治疗在这里.