0%

质数筛

质数,也叫素数。作为自然界的神奇规律,在数学、物理、计算机等学科有着重要的地位。本篇文章将讨论计算机筛素数的三种经典算法。

朴素法

首先我们介绍朴素法。根据质数的定义,质数是指在大于 $1$ 的自然数中,除了1和它本身以外不再有其他因数的自然数。因此我们只需要一个一个判断,如果确实没有其他的因数,就是质数了。这就是判断质数的代码,十分暴力。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
#include <stdio.h>
bool isprime(int n) {
if(n <= 1) return false;
for(int i=2; i<=n; ++i) {
if(n % i == 0) return false;
}
return true;
}
int main() {
for(int i=1; i<=100; ++i) {
if(isprime(i)) printf("%d ", i);
}
return 0;
}

在这段代码中,有一个地方可以优化,那就是 $i$ 不一定非要从 $2$ 循环到 $n$ ,只需要循环到根号 $n$ 就行了。因为 $n$ 的因子除了根号 $n$ ,必然有一个大于根号 $n$ ,一个小于根号 $n$ ,如果没有发现小于根号 $n$ 的那个,也就必然不会再发现大于根号 $n$ 的那个了。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
#include <stdio.h>
bool isprime(int n) {
if(n <= 1) return false;
for(int i=2; i*i<=n; ++i) {
if(n % i == 0) return false;
}
return true;
}
int main() {
for(int i=1; i<=100; ++i) {
if(isprime(i)) printf("%d ", i);
}
return 0;
}

埃氏筛

接下来我们介绍埃氏筛法。那么我们既然已经有了朴素法,为什么还要学习埃氏筛法呢?那就是因为朴素法慢,求 $n$ 以内质数需要 $O(n\sqrt{n})$ 的时间复杂度。接下来我们学习一个时间复杂度为 $O(n\log{\log{n}})$ 的算法,埃氏筛算法。考虑合数的特征,一定存在两个除了1和它本身之外的因子,我们只需要从前往后把所有质数的倍数筛掉就行了。例如刚开始我们知道 $2$ 是质数,那么我们就把 $4, 6, 8, 10, 12 …$ 筛掉,然后我们知道 $3$ 是质数,就把 $6, 9, 12, …$ 筛掉,这样留下来没被筛掉的数就是质数了。这就是埃氏筛的代码,可以证明,它的时间复杂度是$O(n\log{\log{n}})$ ,当n很大时比朴素法有更大的优势。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
#include <stdio.h>
#define MAXN (1000005)
int prime[MAXN];
bool isprime[MAXN];
int sieve(int n) {
int p=0;
for(int i=0; i<=n; ++i) isprime[i]=true;
isprime[0]=isprime[1]=false;
for(int i=2; i<=n; ++i) {
if(isprime[i]) {
prime[p++]=i;
for(int j=2; j*i<=n; ++j) isprime[j*i]=false;
}
}
return p;
}
int main() {
int n = sieve(100);
for(int i=0; i<n; ++i) {
printf("%d ", prime[i]);
}
}

欧拉筛

最后,我们介绍欧拉筛,既然已经有这么优秀的埃氏筛,为什么还要学习欧拉筛呢?不妨观察,刚才筛 $2$ 的倍数和 $3$ 的倍数时,它们筛掉了很多共同的数,如 $6, 12, …$ 而这是没有必要的,每个数其实筛一次就够了。我们先写出来被筛掉的合数,观察一下,看看能不能发现什么规律。我们发现,从右向左看,每列的被乘数都是从 $2$ 开始的质数;从上至下看,每行的乘数都会逐渐增大,直到被乘数能整除乘数为止。所以算法的思路就很明确了,只需要用一个循环变量 $i$ ,如果 $i$ 没有被筛掉就把 $i$ 加入质数表,并筛掉质数表中现有的质数与 $i$ 的乘积,当被乘数整除乘数时, $i$ 自增 $1$ 。这就是欧拉筛的代码,虽然较为复杂,但算法效率高,复杂度是线性的。

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
#include <stdio.h>
#define MAXN (1000005)
int prime[MAXN];
int isprime[MAXN];
int sieve(int n) {
int p=0;
for(int i=0; i<=n; ++i) isprime[i]=true;
isprime[0]=isprime[1]=false;
for(int i=2; i<=n; ++i) {
if(isprime[i]) {
prime[p++] = i;
}
for(int j=0; j<p; ++j) {
if(prime[j]*i > n) break;
isprime[prime[j]*i] = false;
if(i % prime[j] == 0) break;
}
}
return p;
}
int main() {
int n = sieve(100);
for(int i=0; i<n; ++i) {
printf("%d ", prime[i]);
}
}