Pollard Rho 质因数分解

暴力算法

这是一个暴力递归分解质因数的伪代码

void rp(int n) {if (n 是 质数) n 是 rp 的一个质因子, return;枚举找到 n 的一个非 1 因子 p;rp(p), rp (n / p);
}

其中判断质数可以用 Miller Rabin。该暴力复杂度的瓶颈是找因子，如果给定的 nnn 是一个大质数如 109+710^9+7109+7，那么枚举的复杂度就一飞冲天了。

生日悖论

假设 n=109+7n=10^9+7n=109+7，只有一个非 111 因子。

如果一个班上有 kkk 个人，那么有两个生日相同的人的概率是多少？

答案为 P(k)=1−∏i=1k365−i+1365P(k)= 1 - \prod_{i=1}^k \frac{365-i+1}{365}P(k)=1−∏i=1k365365−i+1。当 k=48k=48k=48 时，P(k)≈0.960598P(k) \approx 0.960598P(k)≈0.960598。

那么我们可以用生日悖论来进行找因子

我们在 [1,n][1,n][1,n] 中随便选一个数，那么他是因子的概率为 1n\frac{1}{n}n1。

我们在 [1,n][1,n][1,n] 中随便选两个数 aaa 和 bbb，那么 ∣a−b∣|a-b|∣a−b∣ 是因子的概率为 2n\frac{2}{n}n2。这相当于 a=b+na = b + na=b+n 和 a=b−na = b - na=b−n 两个情况。

我们在 [1,n][1,n][1,n] 中随便选 kkk 个数 a1∼aka_1 \sim a_ka1∼ak，对于任意两个数 aia_iai 和 aja_jaj，如果有 gcd⁡(∣ai−aj∣,n)≠1\gcd(|a_i-a_j|,n) \not =1gcd(∣ai−aj∣,n)=1，意味着我们找到了 nnn 的一个因子，它的概率为 2×k(k−1)2n\frac{2 \times \frac{k(k-1)}{2}}{n}n2×2k(k−1)。当 k=nk = \sqrt{n}k=n 时，概率已经接近百分之百了。但如果求这 kkk 个数的两两差，复杂度又回到了 O(n)O(n)O(n)，显然不可行。

伪随机函数

fi=(fi−12+c)modnf_i = (f_{i-1^2+c}) \bmod n fi=(fi−12+c)modn

其中 ccc 可以是 [1,n][1,n][1,n] 内的任意一个数，我们任取一个数作为 f1f_1f1，然后剩下的数列根据这个函数随机生成。可以发现在modp\bmod pmodp 下，一旦有两个数相同，那么这个数列会陷入循环节，陷入循环节前的长度期望是 n\sqrt nn 的。

那么这个函数的优点是什么呢？如果 gcd⁡(fi−fj,n)≠1\gcd(f_i-f_j,n) \not=1gcd(fi−fj,n)=1，那么有 fi+1−fj+1=(fi−fj)(fi+fj)f_{i+1}-f_{j+1}=(f_i-f_j)(f_i+f_j)fi+1−fj+1=(fi−fj)(fi+fj)，gcd⁡(fi+1,fj+1,p)≠1\gcd(f_{i+1},f_{j+1},p) \not = 1gcd(fi+1,fj+1,p)=1。

算法实现

int f(int a, int n, int c)
{return (1ll * a * a + c) % n;
}
void factor(int n) {if (miller_rabin(n)) {z[++cnt] = n; return ;}int g = n;while (g == n) g = pollard_rho(n); //若 g == n，则这次 pollard_rho 的 rp 不好，没有找到因子factor(g), factorn(n / g);
}

现在我们要实现 Pollard_Rho 的主体了。Pollard_Rho 的思想就是在随机数列上选一些数与 nnn 求最大公约数，如果遇到环则跳出。现在有一个效率低下的双指针法。

int pollard_rho(int n) {int c = rand() % (n - 1) + 1;int v1 = c;int v2 = f(v1,n,c);while (v1 != v2) {int g = gcd(n, abs(v1 - v2));if (g > 1) return g;v1 = f(v1,n,c); //v1 跳一次v2 = f(v2,n,c); v2 = f(v2,n,c); //v2 跳两次}return n;
}

我们知道求 gcd⁡\gcdgcd 的复杂度是 log⁡\loglog 的，求太多次会使得代码时间复杂度低下。我们不如用倍增的思想，让 v2v2v2 在一个范围中跳上好几次，把这几次的结果一起求一下 gcd⁡\gcdgcd。这样看上去比较玄学，但是对复杂度却如虎添翼。

int pollard_rho(int n) {int c = rand() % (n - 1) + 1;int v1 = 0;for (int s = 1, t = 2; s <<= 1, t <<= 1) {int v2 = v1, mul = 1;for (int i = s, step = 1; i < t; i++, step++) {v2 = f(v2, n, c);mul = 1ll * mul * abs(v1 - v2) % n;if (step == 128) {step = 0;int v = gcd(mul, n);if (v > 1) return v;}}int v = gcd(mul, n);if (v > 1) return v;v1 = v2;}
}

现在时间复杂度得到了大大提升，当然 Pollard_rho 还有许多更优的实现方式，甚至有许多比 Pollard_rho 快了不知道几十倍的质因数分解算法，但在算法竞赛中此法足矣。