Pollard Rho 质因数分解

暴力算法

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

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

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

生日悖论

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

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

答案为 $P(k)=1-{\prod }_{i=1}^k\frac{365-i+1}{365}$。当 $k=48$ 时，$P(k)\approx 0.960598$。

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

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

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

我们在 $[1,n]$ 中随便选 $k$ 个数 $a_1\sim a_k$，对于任意两个数 $a_i$ 和 $a_j$，如果有 $\gcd (|a_i-a_j|,n)\ne 1$，意味着我们找到了 $n$ 的一个因子，它的概率为 $\frac{2\times \frac{k(k-1)}{2}}{n}$。当 $k=\sqrt{n}$ 时，概率已经接近百分之百了。但如果求这 $k$ 个数的两两差，复杂度又回到了 $O(n)$，显然不可行。

伪随机函数

$$f_i=(f_{i-1^2+c})\mathrm{mod}n$$

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

那么这个函数的优点是什么呢？如果 $\gcd (f_i-f_j,n)\ne 1$，那么有 $f_{i+1}-f_{j+1}=(f_i-f_j)(f_i+f_j)$，$\gcd (f_{i+1},f_{j+1},p)\ne 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 的思想就是在随机数列上选一些数与 $n$ 求最大公约数，如果遇到环则跳出。现在有一个效率低下的双指针法。

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$ 的复杂度是 $\log$ 的，求太多次会使得代码时间复杂度低下。我们不如用倍增的思想，让 $v2$ 在一个范围中跳上好几次，把这几次的结果一起求一下 $\gcd$。这样看上去比较玄学，但是对复杂度的提升却很大。

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 快了不知道几十倍的质因数分解算法，但在算法竞赛中此法足矣。