筛法求素数(转载)

TAG 素数 数论

素数总是一个比较常涉及到的内容,掌握求素数的方法是一项基本功。

基本原则就是题目如果只需要判断少量数字是否为素数,直接枚举因子2 。。N^(0.5) ,看看能否整除N。

如果需要判断的次数较多,则先用下面介绍的办法预处理。

一般的线性筛法

首先先介绍一般的线性筛法求素数

[cpp] view plaincopy

  1. void make_prime() {
  2. memset(prime, 1, sizeof(prime));
  3. prime[0]=false;
  4. prime[1]=false;
  5. int N=31700;
  6. for (int i=2; i
  7. if (prime[i]) {
  8. primes[++cnt ]=i;
  9. for (int k=i*i; k
  10. prime[k]=false;
  11. }
  12. return;
  13. }

这种方法比较好理解,初始时,假设全部都是素数,当找到一个素数时,显然这个素数乘上另外一个数之后都是合数(注意上面的 ii , 比 i2 要快点 ),把这些合数都筛掉,即算法名字的由来。

但仔细分析能发现,这种方法会造成重复筛除合数,影响效率。比如10,在i=2的时候,k=215筛了一次;在i=5,k=56 的时候又筛了一次。所以,也就有了快速线性筛法。

快速线性筛法

快速线性筛法没有冗余,不会重复筛除一个数,所以"几乎"是线性的,虽然从代码上分析,时间复杂度并不是O(n)。先上代码

[cpp] view plaincopy

  1. #include
  2. using namespace std;
  3. const long N = 200000;
  4. long prime[N] = {0},num_prime = 0;
  5. int isNotPrime[N] = {1, 1};
  6. int main()
  7. {
  8. for(long i = 2 ; i < N ; i ++)
  9. {
  10. if(! isNotPrime[i])
  11. prime[num_prime ++]=i;
  12. //关键处1
  13. for(long j = 0 ; j < num_prime && i * prime[j] < N ; j ++)
  14. {
  15. isNotPrime[i * prime[j]] = 1;
  16. if( !(i % prime[j] ) ) //关键处2
  17. break;
  18. }
  19. }
  20. return 0;
  21. }

首先,先明确一个条件,任何合数都能表示成一系列素数的积。

不管 i 是否是素数,都会执行到"关键处1",

①如果 i 都是是素数的话,那简单,一个大的素数 i 乘以不大于 i 的素数,这样筛除的数跟之前的是不会重复的。筛出的数都是 N=p1*p2的形式, p1,p2之间不相等

②如果 i 是合数,此时 i 可以表示成递增素数相乘 i=p1p2...*pn, pi都是素数(2<=i<=n), pi<=pj ( i<=j )

p1是最小的系数。

根据"关键处2"的定义,当p1==prime[j] 的时候,筛除就终止了,也就是说,只能筛出不大于p1的质数*i。

我们可以直观地举个例子。i=235

此时能筛除 2i ,不能筛除 3i

如果能筛除3i 的话,当 i' 等于 i'=335 时,筛除2i' 就和前面重复了。

需要证明的东西:

  1. 一个数会不会被重复筛除。
  2. 合数肯定会被干掉。

根据上面红字的条件,现在分析一个数会不会被重复筛除。

设这个数为 x=p1p2...*pn, pi都是素数(1<=i<=n) , pi<=pj ( i<=j )

当 i = 2 时,就是上面①的情况,

当 i >2 时, 就是上面②的情况, 对于 i ,第一个能满足筛除 x 的数 y 必然为 y=p2*p3...*pn(p2可以与p1相等或不等),而且满足条件的 y 有且只有一个。所以不会重复删除。

证明合数肯定会被干掉? 用归纳法吧。

类比一个模型,比如说我们要找出 n 中2个不同的数的所有组合 { i , j } ,1<=i<=n, 1<=j<=n,

我们会这么写

for (i=1; i

for (j=i+1; j<=n; ++j)

{

/////

}

我们取 j=i+1 便能保证组合不会重复。快速筛法大概也是这个道理,不过这里比较难理解,没那么直观。

1楼提供的方法,我整理下

//偶数显然不行,所以先去掉偶数。可以看作上面第一种的优化吧。

//不过这种方法不太直观,不太好理解。

[cpp] view plaincopy

1. 我推荐这个算法! 易于理解。 只算奇数部分,时空效率都还不错! 
2. half=SIZE/2; 
3. int sn = (int) sqrt(SIZE); 
4. for (i = 0; i < half; i++) 
5. p[i] = true;// 初始化全部奇数为素数。p[0]对应3,即p[i]对应2*i+3 
6. for (i = 0; i < sn; i++) { 
7. if(p[i])//如果 i+i+3 是素数 
8. { 
9. for(k=i+i+3, j=k*i+k+i; j < half; j+=k) 
10. // 筛法起点是 p[i]所对应素数的平方 k^2 
11. // k^2在 p 中的位置是 k*i+k+i 
12. // 下标 i k*i+k+i 
13. //对应数值 k=i+i+3 k^2 
14. p[j]=false; 
15. } 
16. } 
17. //素数都存放在 p 数组中,p[i]=true代表 i+i+2 是素数。 
18. //举例,3是素数,按3*3,3*5,3*7...的次序筛选,因为只保存奇数,所以不用删3*4,3*6.... 
19.