素数筛法

1. 怎样判断一个数是素数

最直接的判断某个正整数是否是素数的方法是考察这个数是否有小于它本身且大于1的约数，实际上只需考察它是否有大于1且不大于它的算数平方根的约数即可。

素数还有一个特点，就是除2以外，所有的素数都是奇数。

于是有如下C语言代码来判断一个正整数是否为素数：

int is_prime(int n) {   // n > 0. if n is a prime number, returns 1; else returns 0
int i, lim;

if (n == 1)
return 0;
if (n == 2)
return 1;

lim = (int)sqrt(n) + 1;
for (i = 3; i < lim; i += 2) { // O(sqrt(n))
if (n % i == 0)
return 0;
}
return 1;
}

不考虑sqrt(n)的计算复杂度的情形下， 该算法的时间复杂度为O（N^(1/2)）。

若利用该方法求1,2, ... , N中所有的素数，则时间复杂度为O（N^(3/2)）。

2. 素数筛法： Sieve of Eratosthenes

对该筛法，Wiki上有非常详细的描述（网址： http://en.wikipedia.org/wiki/Sieve_of_Eratosthenes）。
算法的基本思想是一个素数的除它本身外的所有倍数均为合数。下图（来自wikipedia）非常直观的描述了该素数筛法的基本原理：



下面是找出0, 1, 2, ..., N - 1 中所有素数的C语言代码：

int is_prime
;

void find_primes() {
int i, j;
is_prime[0] = is_prime[1] = 0;
for (i = 2; i < N; i++) {
is_prime[i] = 1;
}
for (i = 2; i < N; i++) {
if (is_prime[i]) {
for (j = i * i; j < N; j += i) {
is_prime[j] = 0;
}
}
}
}

第 5~8行初始化is_prime， 第9~15行执行筛法。 i的每一次迭代中， 若is_prime[i]为1， 则i为素数（i若不是素数，i一定是比它小的某个素数的倍数，

is_prime[i]就会被标记为 0， 矛盾）， 然后将i的所有倍数均标记为合数, 即将对应的is_prime元素设置为0。

值得注意的是，第11行中j从i*i开始迭代, 原因是比i*i小的i的倍数已被比i小的素数标记过了。 即便做了上述优化，仍然存在重复标记的情况， 如18， 分别被2和3各标记一次。

上述算法的时间复杂度为O(N*loglogN)，空间复杂度为O（N）。

算法时间复杂度计算：



在通常的应用场景下，可以认为算法是线性的。 因为对于函数y = lnlnx（这里使用lnlnx只是为了方便地说明问题）, 随着x的增加， y增长的十分缓慢， 当x=5*10^8时，y还不超过3。

下面是从google上截图得到的y=lnlnx的函数图像：



3. 线性时间复杂度素数筛法

算法的基本思想是任何一个除1之外的数都可以表示成若干个素数的乘积（算数基本定理）， 一个合数必能够被比它小的素数整除。

算法的C语言代码：

int is_prime
;
int primes
;
int nps;

void find_primes()  {
int i, j;

is_prime[0] = is_prime[1] = 0;
for (i = 2; i < N; i++) {
is_prime[i] = 1;
}

for (nps = 0, i = 2; i < N; i++) {
if (is_prime[i])
primes[nps++] = i;
for(j = 0; j < nps && i * primes[j] < N; j++) {
is_prime[i * primes[j]] = 0;
if(i % primes[j] == 0)
break;
}
}
}

上述伪代码能在线性时间内找出1, 2, ..., N-1中的所有的素数, 即算法的时间复杂度为O（N）。 算法的空间复杂度为O(N)。

对伪代码作简要的解释：

首先说明一下全局变量所代表的含义， is_prime数组标记某个数是否为素数, 若是素数, 则对应的数组元素为1, 否则为0; primes数组记录已找到素数;

nps 表示找到的素数的数目。

第18~19行代码可能比较难以理解， 其实这行代码是根据算法的基本思想对算法所做的优化。 一个合数必能够被比它小的素数整除。 当i被某个素数整除时,

i的倍数必然能被该素数整除, 所以i的倍数必然能被整除i的最小素数所标记。 第18~19行就避免了合数被重复标记的情况, 使得每一个合数只被能整除它

的最小素数标记一次。

下面以求1,2,...,19中所有素数为例给出打印程序执行中数据的代码和编译执行后的输出结果。

【C语言代码】

#include <stdio.h>

#define N  20

int is_prime
;
int primes
;
int nps;

void print(int s[], int n) {
int i;
for (i = 0; i < n; i++) {
printf("%2d ", s[i]);
}
printf("\n");
}

void find_primes()  {
int i, j;

is_prime[0] = is_prime[1] = 0;
for (i = 2; i < N; i++) {
is_prime[i] = 1;
}

for (nps = 0, i = 2; i < N; i++) {
if (is_prime[i])
primes[nps++] = i;
printf("i = %2d : primes :", i);
print(primes, nps);
for(j = 0; j < nps && i * primes[j] < N; j++) {
is_prime[i * primes[j]] = 0;
printf("\tp = %2d : ", primes[j]);
print(is_prime, N);
if(i % primes[j] == 0)
break;
}
}
}

int main() {
find_primes();
return 0;
}

【程序输出】

i =  2 : primes : 2
p =  2 :  0  0  1  1  0  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1
i =  3 : primes : 2  3
p =  2 :  0  0  1  1  0  1  0  1  1  1  1  1  1  1  1  1  1  1  1  1
p =  3 :  0  0  1  1  0  1  0  1  1  0  1  1  1  1  1  1  1  1  1  1
i =  4 : primes : 2  3
p =  2 :  0  0  1  1  0  1  0  1  0  0  1  1  1  1  1  1  1  1  1  1
i =  5 : primes : 2  3  5
p =  2 :  0  0  1  1  0  1  0  1  0  0  0  1  1  1  1  1  1  1  1  1
p =  3 :  0  0  1  1  0  1  0  1  0  0  0  1  1  1  1  0  1  1  1  1
i =  6 : primes : 2  3  5
p =  2 :  0  0  1  1  0  1  0  1  0  0  0  1  0  1  1  0  1  1  1  1
i =  7 : primes : 2  3  5  7
p =  2 :  0  0  1  1  0  1  0  1  0  0  0  1  0  1  0  0  1  1  1  1
i =  8 : primes : 2  3  5  7
p =  2 :  0  0  1  1  0  1  0  1  0  0  0  1  0  1  0  0  0  1  1  1
i =  9 : primes : 2  3  5  7
p =  2 :  0  0  1  1  0  1  0  1  0  0  0  1  0  1  0  0  0  1  0  1
i = 10 : primes : 2  3  5  7
i = 11 : primes : 2  3  5  7 11
i = 12 : primes : 2  3  5  7 11
i = 13 : primes : 2  3  5  7 11 13
i = 14 : primes : 2  3  5  7 11 13
i = 15 : primes : 2  3  5  7 11 13
i = 16 : primes : 2  3  5  7 11 13
i = 17 : primes : 2  3  5  7 11 13 17
i = 18 : primes : 2  3  5  7 11 13 17
i = 19 : primes : 2  3  5  7 11 13 17 19