分段筛法的基本思想是选择小于根号下n的筛素数,选择一个合理的、可以容纳于内存中的大分段大小,然后依次筛选每个分段,从最小的开始。在第一个分段中,计算位于该分段内的每个筛素数的最小倍数,然后按照通常的方式标记筛素数的倍数为合数;当使用完所有的筛素数时,分段中未标记的数字就是质数。然后对于下一个分段,对于每个筛素数,你已经知道当前分段中的第一个倍数(它是该素数在前一个分段中结束筛选的倍数),因此你对每个筛素数进行筛选,以此类推,直到完成。
n的大小并不重要,除了较大的n需要比较长的筛选时间以外;关键的大小是分段的大小,应该尽可能大(例如机器上的主存缓存大小)。
你可以在
这里看到分段筛法的简单实现。请注意,相对于 O'Neill 的优先队列筛法,分段筛法会快得多;如果你感兴趣,可以在
这里找到一个实现。
编辑:我写这个的初衷不同,但是我想在这里展示一下,因为它可能会有用:
虽然埃拉托斯特尼筛法非常快,但它需要 O(n)的空间。通过在连续的段中执行筛选,可以将其减少到 O(sqrt(n))的筛素数加上 O(1)的位数组。在第一个段中,计算每个筛素数在该段内的最小倍数,然后按照普通方式标记筛选的倍数为合数;当所有筛素数都被使用时,该段剩余未标记的数字即为质数。然后,对于下一个段,每个筛素数的最小倍数是结束上一段筛选的倍数,因此筛选将继续进行,直到完成。
考虑一个筛子的例子,从100到200分段为20个一组。五个筛选质数是3,5,7,11和13。在100到120的第一段中,位数组有十个插槽,插槽0对应于101,插槽k对应于100+2k+1,插槽9对应于119。该段中3的最小倍数是105,对应于插槽2;插槽2+3=5和5+3=8也是3的倍数。5的最小倍数是105,在插槽2上,插槽2+5=7也是5的倍数。7的最小倍数是105,在插槽2上,插槽2+7=9也是7的倍数。以此类推。
函数primesRange接收参数lo、hi和delta;lo和hi必须为偶数,且lo < hi,lo必须大于sqrt(hi)。段大小为两倍的delta。Ps是包含小于sqrt(hi)的筛素的链接列表,在其中2被删除,因为偶数被忽略。Qs是包含相应筛选质数在当前段中最小多个数的偏移量的链接列表。每个段结束后,lo会增加两倍的delta,因此对应于筛选位数组的索引i的数字为lo + 2i + 1。
function primesRange(lo, hi, delta)
function qInit(p)
return (-1/2 * (lo + p + 1)) % p
function qReset(p, q)
return (q - delta) % p
sieve := makeArray(0..delta-1)
ps := tail(primes(sqrt(hi)))
qs := map(qInit, ps)
while lo < hi
for i from 0 to delta-1
sieve[i] := True
for p,q in ps,qs
for i from q to delta step p
sieve[i] := False
qs := map(qReset, ps, qs)
for i,t from 0,lo+1 to delta-1,hi step 1,2
if sieve[i]
output t
lo := lo + 2 * delta
当调用primesRange(100, 200, 10)时,筛法素数ps为[3, 5, 7, 11, 13];qs最初为[2, 2, 2, 10, 8],对应于最小倍数105、105、105、121和117,并在第二段重置为[1, 2, 6, 0, 11],对应于最小倍数123、125、133、121和143。
您可以在
http://ideone.com/iHYr1f上看到此程序的运行情况。除了上面显示的链接之外,如果您对使用质数进行编程感兴趣,我谦虚地推荐一下我的博客中的这篇
文章。