c - Eratosthenes 筛法 - 按位优化问题

Question

想必大家都遇到过埃拉托色尼筛按位运算的优化代码吧。我正在努力解决这个问题，但我对这个实现中的一个操作有疑问。这是来自 GeeksforGeeks 的代码:

bool ifnotPrime(int prime[], int x) { 
    // checking whether the value of element 
    // is set or not. Using prime[x/64], we find 
    // the slot in prime array. To find the bit 
    // number, we divide x by 2 and take its mod 
    // with 32. 
    return (prime[x / 64] & (1 << ((x >> 1) & 31))); 
} 

// Marks x composite in prime[] 
bool makeComposite(int prime[], int x) { 
    // Set a bit corresponding to given element. 
    // Using prime[x/64], we find the slot in prime  
    // array. To find the bit number, we divide x 
    // by 2 and take its mod with 32. 
    prime[x / 64] |= (1 << ((x >> 1) & 31)); 
} 

// Prints all prime numbers smaller than n. 
void bitWiseSieve(int n) { 
    // Assuming that n takes 32 bits, we reduce 
    // size to n/64 from n/2. 
    int prime[n / 64]; 

    // Initializing values to 0 . 
    memset(prime, 0, sizeof(prime)); 

    // 2 is the only even prime so we can ignore that 
    // loop starts from 3 as we have used in sieve of 
    // Eratosthenes . 
    for (int i = 3; i * i <= n; i += 2) { 

        // If i is prime, mark all its multiples as 
        // composite 
        if (!ifnotPrime(prime, i)) 
            for (int j = i * i, k = i << 1; j < n; j += k) 
                makeComposite(prime, j); 
    } 

    // writing 2 separately 
    printf("2 "); 

    // Printing other primes 
    for (int i = 3; i <= n; i += 2) 
        if (!ifnotPrime(prime, i)) 
            printf("%d ", i); 
} 

// Driver code 
int main() { 
    int n = 30; 
    bitWiseSieve(n); 
    return 0; 
} 

所以我的问题是:

1. (prime[x/64] & (1 << ((x >> 1) & 31)) 是什么意思？更具体地说 (1 << ((x >> 1) & 31));
2. 在prime[x/64]为什么我们除以 64而不是 32 ，当我们使用 32 位整数时；
3. 是int prime[n/64]如果 n < 64 正确？

Answer 1

代码中存在多个问题:

• makeComposite()应该有返回类型 void .
• (1 << ((x >> 1) & 31))有未定义的行为如果 x == 63因为1 << 31溢出类型 int 的范围.你应该使用 1U或者最好是 1UL以确保 32 位。
• 移动表格条目并屏蔽最后一位会更简单:return (prime[x / 64] >> ((x >> 1) & 31)) & 1;
• 而不是假定类型 int有 32 位，你应该使用 uint32_t作为位数组的类型。
• 数组int prime[n / 64];如果 n 太短不是 64 的倍数。使用 uint32_t prime[n / 64 + 1];反而。这对您的示例来说是一个问题 n = 30 , 所以创建的数组长度为 0 .
• ifnotPrime(n)仅返回奇数的有效结果。更改此函数并将其命名为 isOddPrime() 可能会更好且更具可读性.

关于您的问题:

what is the meaning of (prime[x/64] & (1 << ((x >> 1) & 31)) and more specifically (1 << ((x >> 1) & 31))?

x先除以 2(右移一位)，因为只有奇数在数组中有位，然后结果用 31 掩码以保留低 5 位作为字中的位号。对于任何无符号值 x , x & 31相当于x % 32 . x / 64是测试该位的字号。

如上所述，1是一个 int因此不应向左移动 31 个位置。使用 1UL确保该类型至少有 32 位，并且可以移动 31 个位置。

in prime[x/64] why do we divide by 64 and not with 32, when we work with 32-bit integer?

数组中的位对应奇数，所以一个 32 位的字包含 64 个数的素数信息:32 个已知为合数的偶数和 32 个奇数，如果该数为合成的。

Is int prime[n/64] correct if n < 64?

不是不是，如果n是不正确的不是 64 的倍数: 大小表达式应为 (n + 63) / 64 , 或更好 int prime[n/64 + 1]

---

这是修改后的版本，您可以在其中传递命令行参数:

#include <stdio.h>
#include <stdint.h>
#include <stdlib.h>
#include <stdbool.h>
#include <string.h>

bool isOddPrime(const uint32_t prime[], unsigned x) {
    // checking whether the value of element
    // is set or not. Using prime[x/64], we find
    // the slot in prime array. To find the bit
    // number, we divide x by 2 and take its mod
    // with 32.
    return 1 ^ ((prime[x / 64] >> ((x >> 1) & 31)) & 1);
}

// Marks x composite in prime[]
void makeComposite(uint32_t prime[], unsigned x) {
    // Set a bit corresponding to given element.
    // Using prime[x/64], we find the slot in prime
    // array. To find the bit number, we divide x
    // by 2 and take its mod with 32.
    prime[x / 64] |= (1UL << ((x >> 1) & 31));
}

// Prints all prime numbers smaller than n.
void bitWiseSieve(unsigned n) {
    // Assuming that n takes 32 bits, we reduce
    // size to n/64 from n/2.
    uint32_t prime[n / 64 + 1];

    // Initializing values to 0 .
    memset(prime, 0, sizeof(prime));

    // 2 is the only even prime so we can ignore that
    // loop starts from 3 as we have used in sieve of
    // Eratosthenes .
    for (unsigned i = 3; i * i <= n; i += 2) {
        // If i is prime, mark all its multiples as composite
        if (isOddPrime(prime, i)) {
            for (unsigned j = i * i, k = i << 1; j < n; j += k)
                makeComposite(prime, j);
        }
    }

    // writing 2 separately
    if (n >= 2)
        printf("2\n");

    // Printing other primes
    for (unsigned i = 3; i <= n; i += 2) {
        if (isOddPrime(prime, i))
            printf("%u\n", i);
    }
}

// Driver code
int main(int argc, char *argv[]) {
    unsigned n = argc > 1 ? strtol(argv[1], NULL, 0) : 1000;
    bitWiseSieve(n);
    return 0;
}