Loading

Pollard-Rho算法 学习笔记

Miller Rabin 素性测试

假设我们要判定一个数 \(m\) 是否是素数。

Fermat素性测试

可以随便找到一个数 \(a\),然后判定 \(a^{m - 1} \bmod m\) 是否为 \(1\)

如果我们找到 \(k\)\(a\) 都满足 \(a^{m - 1} \equiv 1\),那么我们就把 \(m\) 判定为素数。

这样子的判定事实上是不对的,因为为有很多 (无限个) 合数满足对于任意 \(a\)\(a^{m - 1} \equiv 1 \pmod m\)

二次探测定理

对于任意奇素数 \(p\)\(x^2 \equiv 1 \pmod p\) 的解只有 \(1\)\(p - 1\)

证明:\(x^2 - 1 \equiv 0 \pmod p\)\((x + 1) (x - 1) \pmod p\),因此要么 \(p | x + 1\),要么 \(p | x - 1\),因此 \(x\)\(1\)\(p - 1\)


我们把二次探测定理和费马小定理配合起来使用。

\(n - 1\) 分解成 \(t \times 2^k\)。我们判一下 \(k = a^t \pmod m\)。如果不是 \(1\) 且不是 \(m - 1\),那么继续判断。

如果 \(m\) 是素数,那么平方之后为 \(1\) 且原先不是 \(1\) ,那么一定是由 \(m - 1\) 平方而来。所以我们只要平方,然后如果没有遇到 \(m - 1\) 就认为他不是素数。

一次成功判定的概率是 \(\frac{3}{4}\)

(代码见下面 Pollard-Rho 算法 的)

Pollard-Rho 算法

生日悖论

随机找 \(23\) 个人,他们中有人在同一天生日的概率是 \(50 \%\)

推广:随机写在 \([1, n]\) 中的数,它们中有两个数相同的期望次数是 \(\sqrt n\)

算法

(假设我们要找 \(n\) 的因子)

我们生成一个序列 \(F\)\(f_i = F(f_{i - 1})\)\(F(x)\) 是一个函数,这里 \(F(x) = ( x^2 + seed ) \bmod n\)\(seed\) 是随机选择的一个数。

\(F(x)\) 是几乎随机的,\(F(x) \bmod q\)\(F(x) \bmod n\) 也几乎随机。

假设 \(n\) 有一个因数 \(q\) ,那么 \(F(x) \bmod q\) 必然会重复,因此有循环节,且循环节长度期望为 \(\sqrt q\)

如果找到循环节,那么用循环节中重复的两个数 \(val1, val2\)\(val1 \equiv val2 \pmod q\),我们可以用 \(\gcd(n, |val1 - val2|)\) 来找到 \(q\)

如何找循环节?

可以设两个变量 \(a\)\(b\) (初值相等),每次让 \(a = F(a), b = F(F(b))\)。最后他们会相遇,这样子我们就找到了环。

找环的期望步数是 \(\sqrt q\) 的。\(q\) 可以取到 \(n\) 的最小因数,所以步数是期望 \(n^{\frac{1}{4}}\)

但是我们找环的时候每次都要做一遍 \(\gcd\),时间复杂度 \(\Theta(n^{\frac{1}{4}} \log n)\)

但是我们可以优化:分段处理。每 \(lim\) 个数处理一次,处理出 \(mul = \prod\limits_{i = 1}^{lim} p_i \bmod n\),然后判一下 \(\gcd(mul, n)\) 是否为 \(1\)

如果 \(\gcd\) 不为 \(1\) 就判一下找到环的位置。时间复杂度变成了 \(\Theta (n^{\frac{1}{4}})\)

代码:

#include<bits/stdc++.h>
#define L(i, j, k) for(int i = j, i##E = k; i <= i##E; i++)
#define R(i, j, k) for(int i = j, i##E = k; i >= i##E; i--)
#define ll long long
#define ull unsigned long long
#define db double
#define pii pair<int, int>
#define mkp make_pair
using namespace std;
mt19937_64 orz(233);
const int N = 1e5 + 7;
ll Sum(ll aa, ll bb, ll mod) {
	aa += bb;
	return aa >= mod ? aa - mod : aa;
}
ll Mul(ll a, ll b, ll mod) {
	ll sav = (a * b - (ll)((long db) a * b / mod + 0.1) * mod) % mod; 
	return sav < 0 ? sav + mod : sav;
} 
ll mypow(ll x, ll y, ll mod) {
	ll res = 1;
	for(; y; x = Mul(x, x, mod), y >>= 1) if(y & 1) res = Mul(res, x, mod);
	return res;
} 
const int testp[8] = {2, 3, 5, 7, 13, 19, 61, 233};
bool Miller_Rabin(ll x) {
	if(x < 3) return x == 2;
	if(x <= testp[7]) {
		L(i, 2, sqrt(x)) if(x % i == 0) return 0;
		return 1;
	}
	int k, j; ll a;
	for(k = 0, a = x - 1; ! (a & 1); k ++) a >>= 1;
	L(i, 0, 7) {
		ll v = mypow(testp[i], a, x);
		if(v == 1 || v == x - 1) continue;
		for(j = 0; j < k; j ++) {
			v = Mul(v, v, x);
			if(v == x - 1) break;
		}
		if(j == k) return 0;
	}
	return 1;
}
int total;
ll f[N], seed;
ll F(ll x, ll mod) {
	return Sum(Mul(x, x, mod), seed, mod);
}
ll Pollard_Rho(ll x) {
	total = 0, seed = orz() % (x - 1) + 1;
	ll val1 = rand() % (x - 1) + 1, val2 = val1, mul = 1, sav;
	while(1) {
		val1 = F(val1, x), val2 = F(F(val2, x), x);
		f[++total] = abs(val1 - val2);
		mul = Mul(mul, f[total], x);
		if(total == 127) {
			if(__gcd(mul, x) > 1) {
				L(i, 1, total) if((sav = __gcd(f[i], x)) > 1) return sav;
				return x; 
			}
			total = 0;
		}
	}
}
int tot;
ll ans[N];
void edsolve(ll x) {
	if(x < 2) return;
	if(Miller_Rabin(x)) return ans[++tot] = x, void();
	ll sav = Pollard_Rho(x);
	while(sav == x) sav = Pollard_Rho(x);
	edsolve(sav), edsolve(x / sav);
}
void solve(ll x) {
	tot = 0;
	L(i, 2, testp[7]) while(x % i == 0) ans[++tot] = i, x /= i;
	edsolve(x), sort(ans + 1, ans + tot + 1);
}
ll n;
void Main() {
	scanf("%lld", &n), tot = 0;
	if(Miller_Rabin(n)) puts("Prime");
	else solve(n), printf("%lld\n", ans[tot]);
}
int main() {
	int T; scanf("%d", &T);
	while(T--) Main();
	return 0;
}

参考资料:command_block Mr素性测试+Prho分解小记OI-wiki 素数

posted @ 2021-01-19 19:11  zhoukangyang  阅读(220)  评论(0编辑  收藏  举报