制作一个简单的筛子很容易：
for (int i=2; i<=N; i++){
if (sieve[i]==0){
cout << i << " is prime" << endl;
for (int j = i; j<=N; j+=i){
sieve[j]=1;
}
}
cout << i << " has " << sieve[i] << " distinct prime factors\n";
}
但是当 N 非常大并且我不能在内存中保存那种数组时呢？我查找了分段筛方法，它们似乎涉及在 sqrt(N) 之前找到素数，但我不明白它是如何工作的。如果 N 非常大（比如 10^18）怎么办？
分段筛的基本思想是选择小于n的平方根的筛分素数，选择一个相当大但仍适合内存的片段大小，然后依次筛选每个片段，从最小的片段开始。在第一段，计算段内每个筛素的最小倍数，然后以正常方式将筛素的倍数标记为合数；当所有的筛选素数都用完后，该段中剩余的未标记数是素数。然后，对于下一个片段，对于每个筛分素数，您已经知道当前片段中的第一个倍数（它是结束前片段中那个素数的筛分的倍数），所以您在每个筛分素数上进行筛分，依此类推直到你完成。
n的大小无关紧要，只是较大的n比较小的n需要更长的筛分时间；重要的是段的大小，它应该尽可能大（例如，机器上主内存缓存的大小）。
您可以在此处查看分段筛的简单实现。请注意，分段筛将比另一个答案中提到的 O'Neill 的优先队列筛快得多；如果你有兴趣，这里有一个实现。
编辑：我写这个是为了不同的目的，但我会在这里展示它，因为它可能有用：
尽管埃拉托色尼筛法非常快，但它需要 O(n) 空间。通过在连续的段中执行筛选，这可以减少到筛选素数的 O(sqrt(n)) 加上位数组的 O(1)。在第一段，计算段内每个筛素的最小倍数，然后以正常方式将筛素的倍数标记为合数；当所有的筛选素数都用完后，该段中剩余的未标记数是素数。然后，对于下一段，每个筛分素数的最小倍数是在前一段中结束筛分的倍数，因此筛分一直持续到完成。
考虑在 20 段中从 100 到 200 筛分的示例。五个筛分素数是 3、5、7、11 和 13。在从 100 到 120 的第一段中，位数组有十个槽，槽 0 对应于 101 ，slot k对应100+2k+1，slot 9对应119。segment中3的最小倍数是105，对应slot 2；槽 2+3=5 和 5+3=8 也是 3 的倍数。槽 2 处 5 的最小倍数是 105，槽 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 是一个链表，包含在相应筛选素数的当前段中的最小倍数的筛位数组中的offest。在每个段之后，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。除了上面显示的链接之外，如果您对使用素数进行编程感兴趣，我会在我的博客上谦虚地推荐这篇文章。
有一个基于优先级队列的 Sieve 版本，可以根据您的要求产生尽可能多的素数，而不是所有素数都达到上限。它在经典论文“The Genuine Sieve of Eratosthenes”中进行了讨论，谷歌搜索“eratosthenes 优先队列的筛子”在各种编程语言中出现了相当多的实现。
只是我们正在用我们拥有的筛子进行分段。基本思想是假设我们必须找出 85 到 100 之间的素数。我们必须应用传统的筛子，但采用如下所述的方式：
所以我们取第一个素数 2 ，将起始数字除以 2(85/2) 并取整到较小的数字我们得到 p=42，现在再次乘以 2 我们得到 p=84，从这里开始添加 2直到最后一个数字。所以我们所做的是我们已经删除了范围内的所有因素 2(86,88,90,92,94,96,98,100)。
我们取下一个素数 3，将起始数除以 3(85/3)，然后四舍五入到较小的数，得到 p=28，现在再次乘以 3，得到 p=84，从这里开始加 3 直到最后一个数字。所以我们所做的是我们已经删除了范围内的所有因素 3(87,90,93,96,99)。
取下一个素数=5 以此类推......继续做上述步骤。你可以得到素数（2,3,5,7, ...）通过使用传统的筛子达到 sqrt(n)。然后将其用于分段筛子。
如果有人想看 C++ 实现，这里是我的：
void sito_delta( int delta, std::vector<int> &res)
{
std::unique_ptr<int[]> results(new int[delta+1]);
for(int i = 0; i <= delta; ++i)
results[i] = 1;
int pierw = sqrt(delta);
for (int j = 2; j <= pierw; ++j)
{
if(results[j])
{
for (int k = 2*j; k <= delta; k+=j)
{
results[k]=0;
}
}
}
for (int m = 2; m <= delta; ++m)
if (results[m])
{
res.push_back(m);
std::cout<<","<<m;
}
};
void sito_segment(int n,std::vector<int> &fiPri)
{
int delta = sqrt(n);
if (delta>10)
{
sito_segment(delta,fiPri);
// COmpute using fiPri as primes
// n=n,prime = fiPri;
std::vector<int> prime=fiPri;
int offset = delta;
int low = offset;
int high = offset * 2;
while (low < n)
{
if (high >=n ) high = n;
int mark[offset+1];
for (int s=0;s<=offset;++s)
mark[s]=1;
for(int j = 0; j< prime.size(); ++j)
{
int lowMinimum = (low/prime[j]) * prime[j];
if(lowMinimum < low)
lowMinimum += prime[j];
for(int k = lowMinimum; k<=high;k+=prime[j])
mark[k-low]=0;
}
for(int i = low; i <= high; i++)
if(mark[i-low])
{
fiPri.push_back(i);
std::cout<<","<<i;
}
low=low+offset;
high=high+offset;
}
}
else
{
std::vector<int> prime;
sito_delta(delta, prime);
//
fiPri = prime;
//
int offset = delta;
int low = offset;
int high = offset * 2;
// Process segments one by one
while (low < n)
{
if (high >= n) high = n;
int mark[offset+1];
for (int s = 0; s <= offset; ++s)
mark[s] = 1;
for (int j = 0; j < prime.size(); ++j)
{
// find the minimum number in [low..high] that is
// multiple of prime[i] (divisible by prime[j])
int lowMinimum = (low/prime[j]) * prime[j];
if(lowMinimum < low)
lowMinimum += prime[j];
//Mark multiples of prime[i] in [low..high]
for (int k = lowMinimum; k <= high; k+=prime[j])
mark[k-low] = 0;
}
for (int i = low; i <= high; i++)
if(mark[i-low])
{
fiPri.push_back(i);
std::cout<<","<<i;
}
low = low + offset;
high = high + offset;
}
}
};
int main()
{
std::vector<int> fiPri;
sito_segment(1013,fiPri);
}
基于 Swapnil Kumar 的回答，我用 C 语言做了以下算法。它是用 mingw32-make.exe 构建的。
#include<math.h>
#include<stdio.h>
#include<stdlib.h>
int main()
{
const int MAX_PRIME_NUMBERS = 5000000;//The number of prime numbers we are looking for
long long *prime_numbers = malloc(sizeof(long long) * MAX_PRIME_NUMBERS);
prime_numbers[0] = 2;
prime_numbers[1] = 3;
prime_numbers[2] = 5;
prime_numbers[3] = 7;
prime_numbers[4] = 11;
prime_numbers[5] = 13;
prime_numbers[6] = 17;
prime_numbers[7] = 19;
prime_numbers[8] = 23;
prime_numbers[9] = 29;
const int BUFFER_POSSIBLE_PRIMES = 29 * 29;//Because the greatest prime number we have is 29 in the 10th position so I started with a block of 841 numbers
int qt_calculated_primes = 10;//10 because we initialized the array with the ten first primes
int possible_primes[BUFFER_POSSIBLE_PRIMES];//Will store the booleans to check valid primes
long long iteration = 0;//Used as multiplier to the range of the buffer possible_primes
int i;//Simple counter for loops
while(qt_calculated_primes < MAX_PRIME_NUMBERS)
{
for (i = 0; i < BUFFER_POSSIBLE_PRIMES; i++)
possible_primes[i] = 1;//set the number as prime
int biggest_possible_prime = sqrt((iteration + 1) * BUFFER_POSSIBLE_PRIMES);
int k = 0;
long long prime = prime_numbers[k];//First prime to be used in the check
while (prime <= biggest_possible_prime)//We don't need to check primes bigger than the square root
{
for (i = 0; i < BUFFER_POSSIBLE_PRIMES; i++)
if ((iteration * BUFFER_POSSIBLE_PRIMES + i) % prime == 0)
possible_primes[i] = 0;
if (++k == qt_calculated_primes)
break;
prime = prime_numbers[k];
}
for (i = 0; i < BUFFER_POSSIBLE_PRIMES; i++)
if (possible_primes[i])
{
if ((qt_calculated_primes < MAX_PRIME_NUMBERS) && ((iteration * BUFFER_POSSIBLE_PRIMES + i) != 1))
{
prime_numbers[qt_calculated_primes] = iteration * BUFFER_POSSIBLE_PRIMES + i;
printf("%d\n", prime_numbers[qt_calculated_primes]);
qt_calculated_primes++;
} else if (!(qt_calculated_primes < MAX_PRIME_NUMBERS))
break;
}
iteration++;
}
return 0;
}
它设置要找到的最大素数，然后用已知的素数（如 2、3、5...29）初始化一个数组。所以我们创建一个缓冲区来存储可能的素数段，这个缓冲区不能大于最大初始素数的幂，在这种情况下是 29。
我确信可以进行大量优化来提高性能，例如并行化段分析过程和跳过 2、3 和 5 的倍数的数字，但它可以作为低内存消耗的示例。
一个数是素数，如果没有一个较小的素数能整除它。由于我们按顺序遍历素数，我们已经将所有可以被至少一个素数整除的数字标记为可整除。因此，如果我们到达一个单元格并且它没有被标记，那么它就不能被任何更小的素数整除，因此它必须是素数。
请记住以下几点：-
// Generating all prime number up to R
// creating an array of size (R-L-1) set all elements to be true: prime && false: composite
#include<bits/stdc++.h>
using namespace std;
#define MAX 100001
vector<int>* sieve(){
bool isPrime[MAX];
for(int i=0;i<MAX;i++){
isPrime[i]=true;
}
for(int i=2;i*i<MAX;i++){
if(isPrime[i]){
for(int j=i*i;j<MAX;j+=i){
isPrime[j]=false;
}
}
}
vector<int>* primes = new vector<int>();
primes->push_back(2);
for(int i=3;i<MAX;i+=2){
if(isPrime[i]){
primes->push_back(i);
}
}
return primes;
}
void printPrimes(long long l, long long r, vector<int>*&primes){
bool isprimes[r-l+1];
for(int i=0;i<=r-l;i++){
isprimes[i]=true;
}
for(int i=0;primes->at(i)*(long long)primes->at(i)<=r;i++){
int currPrimes=primes->at(i);
//just smaller or equal value to l
long long base =(l/(currPrimes))*(currPrimes);
if(base<l){
base=base+currPrimes;
}
//mark all multiplies within L to R as false
for(long long j=base;j<=r;j+=currPrimes){
isprimes[j-l]=false;
}
//there may be a case where base is itself a prime number
if(base==currPrimes){
isprimes[base-l]= true;
}
}
for(int i=0;i<=r-l;i++){
if(isprimes[i]==true){
cout<<i+l<<endl;
}
}
}
int main(){
vector<int>* primes=sieve();
int t;
cin>>t;
while(t--){
long long l,r;
cin>>l>>r;
printPrimes(l,r,primes);
}
return 0;
}