用筛选法求质数

作者：cywater2000

日期：2005-4-5

出处：http://blog.csdn.net/cywater2000/

 

一．起因：

很早以前曾在CSDN某位网友的blog中看到过有关“筛选法求质数”的报道。由于当时没有具体研究，所以并没有看懂那段代码。觉得有些不可思议！

…

前几天正好有机会在某门课上学了一下相关知识，于是拿出来大家分享一下。

在讲之前，大家可以先看看下面这篇帖子：

http://www.math.utah.edu/classes/216/assignment-07.html

 

是不是觉得有些神奇？如果只是拿来用的话，也没多少好说的。但正所谓“知其然，知其所以然也”。我们来看看算法是怎么得来的。

 

二．算法描述

我们先设保存整数0~N的数组为sieve[N+1]，用“逐步求精”的方法来得出算法：

 

[定理1]若比素数P小的所有素数的倍数均已从sieve中删去，且P*P > N,则sieve = primes(2…N)。 其中：primes(2..N)表示2..N之间所有的素数。

 

因此算法的框架应该是：

//////////////////////////////////////////////////////////////////////////////////

P = 3; //(P=2不考虑的原因是因为偶数一开始就删除了)

while(P * P <= N)

{

       从sivev中删去P的倍数;

　  P =下一个素数；

}

//此时sieve中留下的就是素数了

//////////////////////////////////////////////////////////////////////////////////

 

那么问题1：如求P的倍数呢？

由于P是质数，所以P肯定是奇数，因此：奇数＋2＝奇数；

所以P的倍数可以表示为：

…P*(P-4),  P*(P-2),  P*P,  P*(P+2),  P*(P+4)…

而“比素数P小的所有素数的倍数均已从sieve中删去”，所以只考虑:

P*P,  P*(P+2),  P*(P+4)…即可：

//////////////////////////////////////////////////////////////////////////////////

//“从sivev中删去P的倍数;”

  t = P * P;

s = 2 * P;

while(t <= N)

{

  从sieve中删去t;

t = t +s;

}

//////////////////////////////////////////////////////////////////////////////////

 

问题2：如何得到下一个素数？

[定理2]　令P`是sieve中大于P的最小元素，则P`是下一素数。

 

//////////////////////////////////////////////////////////////////////////////////

//“P =下一个素数；”

 P = P+1;

while(P不在sieve中) P = P+1;

//////////////////////////////////////////////////////////////////////////////////

 

最后一个问题3：算法是否会终止呢？即在P和P*P之间是否一定能找到一个素数？

[定理3]　已知素数P，在P和P*P之间至少存在一个素数。

 

 

最后，根据上面的算法，我们用一个简单的C++程序来描述就是：

 

/************************************************************

 *   使用改进的The Sieve of Eratosthenes（爱拉托逊斯筛选法）求质数

 *  cywater2000

 *  2005-4-5

 /************************************************************/

#include <iostream.h>

 

#define N  100

 

void main()

{

 

       int sieve[N + 1]; //N以内的数

 

       //step 1:

       //初始化（sieve[i] = 0 表示不在筛中，即不是质数；1表示在筛中）

      for(int i = 2; i <= N; i++) //从2开始，因为0和1不是质数

        sieve[i] = 1;

 

       //step 2:

       //偶数(2的倍数)肯定不是质数，所以应该先筛除

       for(i = 2; i <= N / 2; i++)

        sieve[i * 2] = 0;

 

       int p = 2; //第一个质数是2

 

       //step 3:从sieve中删去P的倍数

       while(p * p <= N)

       {

              //选下一个p

              p = p + 1;

              while(sieve[p] == 0)

                     p++;

 

              int t = p * p;

              int s = 2 * p;

              while(t <= N)

              {

                     sieve[t] = 0;  //删除

                     t = t + s;

              }

       }

      

 

  //step 4:　输出结果

       for(i = 2; i <= N; i++)

       {    

              if(sieve[i] != 0)

              {

                     cout<<i<<", ";

              }

       }