求质数2 筛法 C实现

该问题出自《C语言名题精选百则技巧篇》

不用除法，用筛子把i的所有倍数去掉，求出2~N之间的所有质数。这个方法叫做筛法（Sieve Method）

举例：求2~40之间的质数，2是质数，所以2的倍数一定被删除，所以在一开始可以考虑就不把2的倍数考虑在内，只考虑奇数，3,5,7,9,11,...

下一个数是自然数3，把3的倍数都删除，下一个是5， 把5的倍数都删除，下一个是7，把7的倍数都删除，接下来是11，。。。。

（1）用一个数组x[0],x[1],...来储存3,5,7,9,...这些奇数，x[i] = 2*i+3, i为x[i]的下标，2*i+3为x[i]的值,用一个数组sieve[]来记录被删除的值的下标，例如sieve[i] = 1,代表下标为i的值被删除，即2*i+3这个值被删除。

（2）在3,5,7,9,...这些奇数中，从x[0]开始筛选，先筛选能被3整除的数，标记下标，设置所有能被3整除的的sieve[i] = 1,然后筛选能被5整除的数，标记下标sieve[i] = 1;筛选能被7整除的数，标记下标sieve[i] = 1,9所在位置的sieve已经被标记sieve[i] = 1，因此不再寻找能被9整除的数，...一直遍历到x[n],最后得出的结果是2~2*n+3之间的质数。

（3）我们要寻找能被3,5,7,...这些数整除的那些数的位置，即下标的值i，并标记它，x这个数组全部是奇数，因此，不存在3,5,7，这些数字的偶数倍，因此我们要寻找的是3,5,7,...这些数字的奇数倍的位置。x[i] = 2*i+3,x[i]的奇数倍为(2*n+1)(2*i+3),n=1,2,3,...

(2*n+1)(2*i+3) = 2n(2i+3)+(2i+3) = 2[n(2i+3)+i]+3 ,这些值的下标为n(2i+3)+i,这个式子可以理解为，以i为起点，间隔为2i+3位置的那些数。

i       0  1  2  3   4   5     6    7    8    9    10    11   12  13   14  15  16  17  18  19  20  21  22  23  24  25  26  27  28  29  30

x[i]   3  5  7  9  11  13  15  17  19  21   23    25   27  29   31  33   35  37  39  41  43  45  47  49  51  53  55  57  59  61  63

i=0时，x[i] = 3.  i为0+3*1，0+3*2， 0+3*3， ...  ，  0+3*n位置的数都是3的倍数

i=1时，x[1] = 5,i为1+5*1 ，1+5*2,  1+5*3，  ...  ，  1+5*n位置的数都是5的倍数

.......

因此我们可以以下标i+(2i+3)为初始值，2i+3步进的所有值的标记sieve[i+n(2i+3)] =1,知道所有x都被遍历完，再遍历下标标记，如果被标记为sieve[i] = 1的值，就是被删除的值，而sieve[i] = 0的那个下标i对应的x[i]值就是保留的值.

（4）求x[0]~x[MAXSIZE]（2*MAXSIZE+3）之间的所有质数，数组sieve记录被删除的数的，sieve[i] = 1，代表x[i]被删除，sieve[i] = 0,代表x[i]没被删除，x[i]是质数。在开始，sieve数组所有的元素都被初始化为0（#define KEPT 0）,从x[0]开始遍历x数组，一直到x[MAXSIZE](for(i=0;i<MAXSIZE;i++))，先判断x[i]是否已经被删除，如果没被删除（sieve[i] == KEPT），将x[i+k(2*i+3)]全部删除（设置sieve[i+k(2i+3)] = DELETED(#define DELETED 1 )）,

for(i=0;i<=MAXSIZE;i++)
         if(sieve[i] == KEPT){
             prime = i+i+3;
             count ++;
             for(k=prime+i;k<=MAXSIZE;k+=prime)
                 sieve[k] = DELETED;
         }

(5)x所有的数值已经遍历完以后，遍历sieve数组，sieve[i] == KEPT的数值为质数，这个质数是2i+3

for(i=0,k=2;i<=MAXSIZE;i++)
        if(sieve[i] == KEPT){
            if(k>10){
                printf("\n");
                k=1;         
            }
            printf("%6d",2*i+3);
            k++;
        }

所有的实现代码如下：

#include <stdio.h>  
  
#define MAXSIZE      200  
#define DELETED      1  
#define KEPT         0  
  
int main(int argc,char *argv[])  
{  
    int sieve[MAXSIZE+1];  
    int count =1;    
    int i,k,prime;  
    
    printf("\nEratosthenes Sieve Method");
    printf("\n=========================");
    printf("\n\nPrime Numbers between 2 and %d\n",MAXSIZE*2+3);
    
    for(i=0;i<=MAXSIZE;i++)
        sieve[i] = KEPT;
    
    for(i=0;i<=MAXSIZE;i++)
         if(sieve[i] == KEPT){
             prime = i+i+3;
             count ++;
             for(k=prime+i;k<=MAXSIZE;k+=prime)
                 sieve[k] = DELETED;
         } 
      
    printf("\n%6d",2);
    for(i=0,k=2;i<=MAXSIZE;i++)
        if(sieve[i] == KEPT){
            if(k>10){
                printf("\n");
                k=1;         
            }
            printf("%6d",2*i+3);
            k++;
        }
    printf("\n\nThere are %d primes in total.",count);
        
    getchar();  
    return 0;  
  }    

运行结果：