算法分析:用筛选法求质数
作者: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<<", ";
}
}
}
三.改进
除了可以考虑将*2和/2改为<<1和>>1之外,还可以将第一步和第二步合在一起
//////////////////////////////////////////////////////////////////////////////////
//step 1&2:
for(int i = 3; i <= N; i++)
sieve[i] = i % 2;
sieve[2] = 1; //2可是质数!
//////////////////////////////////////////////////////////////////////////////////
当然还有其它地方可以改进的,这里就不讨论了。