说到筛法,一般我们见到的模式是使用素数去筛除所有该素数的倍数。因为每个合数显然很有可能是多个素数的倍数,所有会被筛很多次,筛的次数自然就显得冗余。
第一次看到线性筛法的时候真的感觉非常不可思议,因为它是严格O(n)的,并且具备非常优雅的性质,以至于可以求素数表以外的很多积性函数的筛表,而且,还可以做因式分解的表,而它们因此就被全部被优化到了O(n)。
线性筛法
下面给出一个C的代码。
1 | typedef long long typec; |
2 | /****************************************************************************** |
3 | * Variables using below |
4 | ****************************************************************************/ |
5 | const int N = 1000000; |
6 | int prime[N + 1]; |
7 | int fac_table[N + 1][2]; |
8 | int phi_table[N + 1]; |
9 | int facque[100]; |
10 | |
11 | |
12 | /****************************************************************************** |
13 | * Prime Table |
14 | * O(n) |
15 | * Need: int* prime |
16 | ****************************************************************************/ |
17 | int prime_table() |
18 | { |
19 | memset(prime, 0, sizeof(prime)); |
20 | for(int i = 2; i <= N; i++) { |
21 | if(!prime[i]) prime[++prime[0]] = i; |
22 | for(int j = 1; j <= prime[0] && prime[j] <= N / i; j++) { |
23 | prime[prime[j] * i] = 1; |
24 | if(i % prime[j] == 0) break; |
25 | } |
26 | } |
27 | return prime[0]; |
28 | } |
这个筛法的动机是什么呢?
前面有提到过朴素筛法的冗余主要就在每个合数被筛到多次,如果把次数优化到1,那自然效率就得到了提高(如果确实要筛到所有的合数,那么这种优化就达到了上限)。而线性的筛法确实做到了这点!问题的关键在于,合数在什么情况下被筛到?
合数仅被其最小的素因子筛到
这个线性筛法就是做到了以上的这点。下面给出代码的分析
- line 19。初始化,prime[0]是当前的素数表中的素数个数
- line 20。逐个判断范围内的数i
- line 21。如果i没有被筛到过,那么它当然是素数,更新加入筛表
- line 22。在这个循环内筛合数。注意for的第二参数,这些是为了防止溢出的,同时也是避免筛到范围外的合数
- line 23。把prime[j]*i所对应的合数标记上。后面会说明为什么prime[j]总是最小因子
- line 24。如果prime[j]能整除i,那么退出当前i的筛合数过程。因为prime[j]是逐个递增的,可以预见到prime[j]在之前的过程中不会大于i的最小素因子。这导致prime[j]是prime[j]*i的最小素因子。因为每个合数的这种拆分是唯一的,对于一个prime[j],仅能找到一个i(并且一定能找到)满足要求。所以,每个合数被筛到有且仅有的一次
注意这个筛法的核心:合数仅被其最小的素因子筛到。这导致它可以用来构造因式分解的筛表!!!
素因子分解
依然给出代码
1 | /****************************************************************************** |
2 | * factor Table |
3 | * O(n) |
4 | * Note: You can get a prime table from this function too |
5 | * Need: int* prime, facque |
6 | ****************************************************************************/ |
7 | void get_factors() |
8 | { |
9 | memset(fac_table, 0, sizeof(fac_table)); |
10 | memset(prime, 0, sizeof(prime)); |
11 | prime[0] = 0; |
12 | fac_table[0][0] = 0, fac_table[0][1] = 1; |
13 | for(int i = 2; i <= N; i++) { |
14 | if(!prime[i]) { |
15 | prime[++prime[0]] = i; |
16 | fac_table[i][0] = i; |
17 | fac_table[i][1] = 1; |
18 | } |
19 | for(int j = 1; j <= prime[0] && prime[j] <= N / i; j++) { |
20 | prime[i * prime[j]] = 1; |
21 | if(fac_table[i][0] != 1) { |
22 | fac_table[i * prime[j]][1] = i; |
23 | } else { |
24 | fac_table[i * prime[j]][1] = fac_table[i][1]; |
25 | } |
26 | if(i % prime[j]) { |
27 | fac_table[i * prime[j]][0] = prime[j]; |
28 | } else { |
29 | fac_table[i * prime[j]][0] = 1; |
30 | break; |
31 | } |
32 | } |
33 | } |
34 | } |
35 | /****************************************************************************** |
36 | * return factors of given n, using fac_table above, O(1) |
37 | * facque[0] is the number of factors, facque[1~facque[0]] are factors |
38 | * Note!!!:every factor is prime. |
39 | * eg. 2 is one factor of 12, factor would return 2 instead of 4 |
40 | * Need: int* facque |
41 | ******************************************************************************/ |
42 | void getfac(int n) |
43 | { |
44 | facque[0] = 0; |
45 | do { |
46 | if(fac_table[n][0] != 1) { |
47 | facque[++facque[0]] = fac_table[n][0]; |
48 | } |
49 | n = fac_table[n][1]; |
50 | } while(n != 1); |
51 | } |
- 这段代码的核心在于:对于给定的合数
(x是n最小的素因子),则n的其它素因子可以有y来提供,所以问题被归约为求y的素因子。所以,每个数n,仅仅需要维护自己的最小素因子和另外一个因子就可以了,所以对空间的需求是2n。而这个素因子和另一个因子也就是之前筛法的print[j]和i。n=x⋅y - 代码还给出了一个优化:因为对于某个素因子,如果存在重复的话,会使得链过长。比如2
,会有30个2挂链。所以我在上面的代码里面给出了优化,把链缩短了(每个因子在链中只出现一次)。当然,有些因子会用1替换(仅有一个),这个是需要额外特判的。30 - 取因子采用了另一个函数void getfac(int n)。结果在facque中,函数上方有注释说明细节。
欧拉函数的线性筛法
1 | /****************************************************************************** |
2 | * Phi Table |
3 | * O(n) |
4 | * Need: int* phi_table, prime |
5 | * Note: You can get a prime table from this function, too. |
6 | ****************************************************************************/ |
7 | void get_phis() |
8 | { |
9 | memset(phi_table, 0, sizeof(phi_table)); |
10 | prime[0] = 0; |
11 | for(int i = 2; i <= N; i++) { |
12 | if(!phi_table[i]) { |
13 | prime[++prime[0]] = i; |
14 | phi_table[i] = i - 1; |
15 | } |
16 | for(int j = 1; j <= prime[0] && prime[j] <= N / i; j++) { |
17 | if(i % prime[j]) { |
18 | phi_table[i * prime[j]] = phi_table[i] * (prime[j] - 1); |
19 | } else { |
20 | phi_table[i * prime[j]] = phi_table[i] * prime[j]; |
21 | break; |
22 | } |
23 | } |
24 | } |
25 | } |
——— 欧拉函数
在数论中,对正整数n,欧拉函数是小于或等于n的正整数中与n互质的数的数目。此函数以其首名研究者欧拉命名,它又称为φ函数、欧拉商数等。
它是积性的(multiplicative),满足以下的性质。
ϕ(1)=1 ϕ(p)=p−1 ϕ(pk)=pk−1(p−1) 对于m,n,如果gcd(m,n)=1,我们有ϕ(mn)=ϕ(m)⋅ϕ(n)
φ(n)可以通过n的因式分解完成求解。而素数的线性筛法确实非常优雅,它能够递推的维护φ。正如我们在因式分解中所看到的特性,φ同样可以在因式分解的链上传递,通过一个最小的因子,把问题归约给到更小的规模。而这次,我们不需要显式的产生因式分解了。
Read full article from 素数线性筛法 - 子曰博文
No comments:
Post a Comment