埃拉托色尼钻头阵列筛

Question

我试图找到素数使用厄拉多塞筛位数组，但我使用的是无符号整数数组。我需要能够产生多达2，147，483，647个素数。我的代码工作正常，可以生成大约10，000，000个，但是当我增加数组的大小以容纳更大的数字时，它失败了。有人能指导我如何用c语言(不是c语言)使用位向量吗？谢谢

这是我的代码：

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

#define MAXBYTES 2000000
#define MAX 50000000
#define BITSIZE 32

void ClearBit(unsigned int [], unsigned int);
void SetBit(unsigned int [], unsigned int);
int  BitVal(unsigned int [], unsigned int);
void PrintBitStream(unsigned int [], unsigned long);
void PrintBitStreamData(unsigned int[], unsigned long);
int  Sieve(unsigned int[], unsigned int, unsigned int);

int main(int argc, char ** argv) {
    unsigned int maxsize = MAX;
    unsigned int i;
    //Set Bit Array
    unsigned int BitArray[MAXBYTES] = {0};
    SetBit(BitArray, 0);
    SetBit(BitArray, 1);
    i = 2;
    for (;i < maxsize;i++){

        if(Sieve(BitArray, i, maxsize)==0)
            break;
    }
    PrintBitStreamData(BitArray, maxsize-1);
    return EXIT_SUCCESS;
}

void PrintBitStreamData(unsigned int BitArray[], unsigned long maxsize) {
    unsigned int i;
    for (i = 0; i < maxsize; i++)
        if (!BitVal(BitArray, i))
            printf("%ld ", i);
    printf("\n");
}

void PrintBitStream(unsigned int BitArray[], unsigned long maxsize) {
    unsigned int i;
    for (i = 2; i < maxsize; i+=2)
        printf("%d", BitVal(BitArray, i));
    printf("\n");
}

void SetBit(unsigned int BitArray[], unsigned int pos) {
    BitArray[pos / BITSIZE] |= 1 << (pos % BITSIZE);
}

void ClearBit(unsigned int BitArray[], unsigned int pos) {
    BitArray[pos / BITSIZE] &= ~(1 << (pos % BITSIZE));
}

int BitVal(unsigned int BitArray[], unsigned int pos) {
    return ((BitArray[pos / BITSIZE] & (1 << (pos % BITSIZE))) != 0);
}

int Sieve(unsigned int BitArray[], unsigned int p, unsigned int maxsize) {
    unsigned int i;
    unsigned int j;
    j = 0;
    for (i = 2 * p; i < maxsize; i += p) {
        SetBit(BitArray, i);
        j++;
    }
    return j;
}

Answer 1

使用位访问整数的示例< br >注意< code>GetBit()和< code>SetBit()。< br >使用2的幂时，优化编译器会使< code>/和< code>%更快。

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

#define ubitsize (sizeof(unsigned)*8)

unsigned GetBit(const unsigned *List, unsigned long index) {
  return !!(List[index / ubitsize] & (1u << (index % ubitsize)));
  }

void SetBit(unsigned *List, unsigned long index) {
  List[index / ubitsize] |= (1u << (index % ubitsize));
  }

void Sieve_of_Eratosthenes_via_bit_array(unsigned long MaxCandidatePrime) {
  unsigned long uByteSize = MaxCandidatePrime/ubitsize + 1;
  unsigned *List = calloc(uByteSize, sizeof *List);
  if (List == 0) return;

  unsigned long PrimeCount = 0;
  unsigned long Index = 0;
  for (Index = 2; Index <= MaxCandidatePrime; Index++) {
    // If found a prime ...
    if (GetBit(List, Index) == 0) {
      PrimeCount++;
      // let's see the progress
      if ((PrimeCount % (256LU*1024)) == 0) printf("%lu\n", Index);

      // Mark subsequent multiples as not--a-prime
      unsigned long Index2 = Index*2;
      while (Index2 <= MaxCandidatePrime) {
        SetBit(List, Index2);
        Index2 += Index;
      }
    }
  }
  printf("X %lu\n", Index);
  free(List);
}

void test(void) {
  Sieve_of_Eratosthenes_via_bit_array(200LU*1000*1000);
}

重写可以采用通常的建议，不保存偶数，将2视为特殊情况。这很有帮助，但我认为这是一个练习。我可以通过节省每30倍编码1个字节来节省大约4倍，因为每30个整数最多有8个素数。还有其他方案。

Answer 2

我肯定不会使用位数组，而是使用原生int (64位或32位，取决于体系结构)数组，并使用按位< code>|和< code >包装一个函数，将普通数字重新映射到正确的位置和位

还要考虑省略偶数，几乎没有一个是质数。因此，您可以将前128个数字存储在第一个64位数字中，然后将下一个128个数字存储在第二个数字中，等等。

这听起来有点复杂，但让它工作有点有趣！

欧拉项目似乎已经产生了一些非常好的解决方案。

好的一面是：要筛分，您不需要重新计算偶数转移，但您可以取消每隔三分之一位的筛分3，每隔5位的筛分5等。

如果您想快速的java解决方案作为详细参考，请进入聊天。

EDIT4：已更正工作代码，但速度较慢。备注：记住使用calloc！

#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include <limits.h>
#include <time.h>

typedef unsigned long long number;

number lookFor = 2147483648ULL;
number max = 2147483648ULL*10ULL; // hopefully more then every 10th uneven number is prime

unsigned long * isComposite;

number bitslong = CHAR_BIT*sizeof(long);

time_t rawtime;
struct tm * timeinfo;
char buffer[80];

// is not called during sieve, only once per sieving prime 
// and needed for reading if a number is prime
inline long getParts(number pos, number *slot, unsigned int *bit){
    *slot = pos / bitslong;
    *bit = (unsigned int)(pos % bitslong);
}

int isPrime(number n){
    if(n == 1){
        return 0;
    }

    if(n < 4){
        return 1;
    }

    if((n%2) == 0){
        return 0;
    }

    number slot=0;
    unsigned int bit=0;
    number pos = (number)(n-3)/2;
    getParts(pos, &slot, &bit);
    // printf("? n=%lld  pos = %lld slot = %lld bit = %lu ret %d \n", n, pos, slot, bit, !(isComposite[slot] & (1<<bit)));
    return !(isComposite[slot] & (1UL<<bit));
}

// start sieving at offset (internal position) offset with width step 
int doSieve(number offset, number step){
    rawtime = time(0);
    time (&rawtime);
    timeinfo = localtime (&rawtime);

    strftime(buffer, 80, "%Y-%m-%d %H:%I:%S", timeinfo);
    unsigned int bit=0;
    number slot=0;
    getParts(offset, &slot, &bit);
    printf("doSieve %s  %lld %lld  %lu \n", buffer, offset, step, isComposite[slot]);

    while(slot < max/bitslong){
        slot += (step + bit)/bitslong;
        bit = (step + bit) % bitslong;
        isComposite[slot] |= (1UL << bit);
    } 
    return 1;
}

int sieve(){
    number spot;
    spot=1;
    number pos;
    pos = 0;
    while(spot < 1 + sqrt((float)max)){
        spot+=2;
        if(! isPrime(spot)){
            pos++;
            continue;
        }
        doSieve(pos, spot);
        pos++;
    }
}

void main(int argc, char *argv[]){
    if(argc > 1){
        char *tp = malloc(sizeof(char*));
        max = strtol(argv[1], &tp, 10);
    }
    printf("max %lld , sq %ld, malloc: %lld\n", max, (long)(1 + sqrt((float)max)), 1+max/bitslong);
    isComposite = calloc((2+max/bitslong), sizeof(unsigned long)) ;
    if(! isComposite){
        printf("no mem\n");
        exit(5);
    }
    sieve();
    number i;
    number found = 0;
    for(i = 1; i<max && found < lookFor; i++){
        if(isPrime(i)){
            found++;
            // printf(" %30lld %30lld \n", found, i);
            if(found % 10000 == 0 ){
                printf("%30lld %30lld \n", found, i);
            }
        }
        /*
        if(i % 1000 == 17){
            printf("%5lld %5lld \n", i, found);
        }
        */
    }
}