c - 使用位掩码使用 Sieve 方法列出素数

Question

我编写了以下代码，使用 Sieve 的方法列出了 20 亿以内的所有质数。我使用位掩码来标记目的。虽然我能够正确地得到素数，但每次都会丢失一些开头的素数。请帮我找出程序中的错误。

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <stdbool.h>

#define MAX 2000000000

char* listPrimes(){
int block = sqrt(MAX);
char* mark = calloc((MAX/8),sizeof(char));
int i = 2;
int j;
char mask[8];
for(j=0;j<8;j++)
    mask[j] = 0;
mask[7] = 1;
mask[6] |= mask[7] << 1;
mask[5] |= mask[7] << 2;
mask[4] |= mask[7] << 3;
mask[3] |= mask[7] << 4;
mask[2] |= mask[7] << 5;
mask[1] |= mask[7] << 6;
mask[0] |= mask[7] << 7;

for(j=0;j<8;j++)
    printf("%d ",mask[j]);
mark[0] |= mask[0];
mark[0] |= mask[1];

while (i < block){

        for (j = 2; i*j <= block; j++)
                mark[(i*j) / 8] |= mask[((i*j) % 8 )];
        i++;
    }
printf("\n");
printf("The block size is\t:\t%d\n",block);


j = 2;
while(j<=block){
    if((mark[j / 8] & mask[j]) == 0 ){
        for(i = 2;i <= MAX; i++){
            if((i%j) == 0){
                mark[i / 8] |= mask[(i % 8)];
            }
        }
    }
while((mark[++j / 8] & mask[j % 8]) != 0);
}


for(j=0;j<=MAX;j++)
        if((mark[j / 8] & mask[(j % 8)]) == 0)
            printf("%d\n", ((8*(j / 8)) + (j % 8)));

return mark;
}   

int main(int argc,char* argv[]){

listPrimes();

return 0;
}

Answer 1

作为ArunMK也就是说，在第二个 while 循环中，您将质数 j 本身标记为 j 的倍数。而作为 Lee Meador说，您需要对 mask 索引取 j 模 8，否则您将越界访问并调用未定义的行为。

进一步调用未定义行为的地方是

while((mark[++j / 8] & mask[j % 8]) != 0);

在没有干预序列点的情况下使用和修改j。您可以通过编写来避免这种情况

do {
    ++j;
}while((mark[j/8] & mask[j%8]) != 0);

或者，如果您坚持使用空主体的 while 循环

while(++j, (mark[j/8] & mask[j%8]) != 0);

您可以使用逗号运算符。

更多未定义的行为通过访问 mark[MAX/8] 未分配

for(i = 2;i <= MAX; i++){

和

for(j=0;j<=MAX;j++)

另外，如果 char 是有符号的并且是八位宽的，

mask[0] |= mask[7] << 7;

是实现定义的(并且可能引发实现​​定义的信号)因为

的结果

mask[0] | (mask[7] << 7)

(int 128)不能表示为 char。

但为什么要将每个数字除以所有不超过第二个 while 循环中边界平方根的素数？

    for(i = 2;i <= MAX; i++){
        if((i%j) == 0){

这使您的算法不是埃拉托色尼筛法，而是试验除法。

为什么不在那里也使用第一个 while 循环中的技术？ (然后，为什么要有两个循环？)

while (i <= block){
    if ((mark[i/8] & mask[i%8]) == 0) {
        for (j = 2; i*j < MAX; j++) {
            mark[(i*j) / 8] |= mask[((i*j) % 8 )];
        }
    }
    i++;
}

不会溢出(对于 MAX 的给定值，如果它可以表示为 int)，并更快地产生正确的输出数量级。