Pollard-Rho 算法

Miller-Rabin 素性检测

部分内容摘自 题解 P4718/论 Miller-Rabin 算法的确定性化 - It's LUNATIC time!)

根据费马小定理，若 \(p\) 为素数，那么对于 \(1 \leq a < p\)，都有 \(a^{p-1} \equiv 1 \pmod p\)。

因此，若存在 \(1 \leq a < p\) 使得 \(a^{p-1} \not\equiv 1 \pmod p\)，那么 \(p\) 一定不是素数。

于是有费马素性检测：在 \([1,p)\) 中随机选取几个 \(a\)，然后使用快速幂来检测是否有 \(a^{p-1} \equiv 1 \pmod p\)。如果对于每个选取的 \(a\)，\(p\) 都通过了检测，那么费马素性检测认为 \(p\) 是素数。

但是有的合数 \(p\) 满足对于 \(1 \leq a < p\)，都有 \(a^{p-1} \equiv 1 \pmod p\)。这些数被称为卡迈克尔数。因此，我们需要检测 \(p\) 的更多性质。

二次探测定理：若 \(p\) 为奇素数，则 \(x^2 \equiv 1 \pmod p\) 的解为 \(x \equiv \pm 1 \pmod p\)。

若 \(p\) 为奇数，显然 \(p-1\) 是偶数。首先令 \(d = p-1\)，然后判断 \(a^d \bmod p\) 是否为 \(1\)，若不是则 \(p\) 为合数。接下来，若 \(d\) 为偶数则令 \(d\) 除以 \(2\)，然后重复上述判断，直到 \(d\) 为奇数或 \(a^d \not\equiv 1 \pmod p\)。

如果最后一步中 \(a^d \not \equiv \pm 1 \pmod p\)，则 \(p\) 一定为合数，否则 Miller-Rabin 素性检测认为它是素数。

事实上只需要选择前 \(12\) 个素数为底数就可以适用于 \(2^{78}\) 内的素性检测。

namespace miller{
	int pr[20]={2,3,5,7,11,13,17,19,23,29,31,37};

	inline bool isprime(ll x){
		if(x<=40){
			for(int i=0;i<12;++i) if(x==pr[i]) return true;
			return false;
		}
		
		ll d=x-1,k=0;
		while(!(d&1)) d>>=1,++k;
		
		for(int i=0;i<12;++i){
			ll a=pr[i],u=qpow(a,d,x),v;
			for(int j=0;j<k;++j){
				v=u,u=u*u%x;
				if(u==1&&v!=x-1&&v!=1) return false;
			}
			if(u!=1) return false;
		}
		return true;
	}
}

Pollard-Rho 算法

Pollard-Rho 算法用于在 \(O(n^{1/4})\) 的复杂度计算合数 \(n\) 的某个非平凡因子。

考虑一个随机函数 \(f(x) = (x^2+c) \bmod n\)，那么取 \(f\) 的很多项相邻作差求绝对值，很大概率能够取到一个数和 \(n\) 的 GCD 不是 \(1\)。

考虑连边 \(x \to f(x),f(x) \to f(f(x))\)，那么它最终必然成环，此时就没必要继续走下去了。我们考虑 Brent 判圈法，初始取 \(s = t =0\)。第 \(i\) 轮令 \(t\) 向后走 \(2^i\) 步，在此过程中 \(s\) 保持不动。如果过程中出现了 \(s = t\)，那么就成环了，并且 \(t\) 已经完整地遍历了整个 \(\rho\) 形，此时应当直接返回。否则令 \(s \leftarrow t\)，然后进入第 \(i+1\) 轮。

但是，在这个过程中，频繁地调用 \(\gcd(|f(x)-x|,n)\) 是不优的。我们只能够退而求其次。具体来说：

• 在每一轮，维护变量 \(val\)，初始为 \(1\)；
• \(t\) 每往后走一步，令 \(val\) 变为 \(val \times |t-s| \bmod n\)。
• \(t\) 每走 \(127\) 步，或者到达该轮末尾，调用 \(\gcd(val,n)\)。若该 \(\gcd > 1\) 则返回。

这是流传较广的一种实现。这样看上去会损失一些效率（例如：找到了一个非平凡因子后很快又转回了 \(t = s\) 处，这样会直接使得 \(val\) 变为 \(0\)，从而返回 \(n\)；找到了两个相乘为 \(n\) 的非平凡因子也会返回 \(0\)），但是根据经验，这样的整体运行速度是最快的。

# include <bits/stdc++.h>
// # define int long long
const int N=100010,INF=0x3f3f3f3f;

typedef __int128 ll;

std::mt19937 rng(114);

inline ll qpow(ll d,ll p,ll mod){
	ll ans=1;
	while(p){
		if(p&1) ans=ans*d%mod;
		d=d*d%mod,p>>=1;
	}
	return ans;
}

namespace miller{
	int pr[20]={2,3,5,7,11,13,17,19,23,29,31,37};

	inline bool isprime(ll x){
		if(x<=40){
			for(int i=0;i<12;++i) if(x==pr[i]) return true;
			return false;
		}
		
		ll d=x-1,k=0;
		while(!(d&1)) d>>=1,++k;
		
		for(int i=0;i<12;++i){
			ll u=qpow(pr[i],d,x),v;
			for(int j=0;j<k;++j){
				v=u,u=u*u%x;
				if(u==1&&v!=x-1&&v!=1) return false;
			}
			if(u!=1) return false;
		}
		return true;
	}
}

ll res;

inline ll read(void){
	ll res,f=1;
	char c;
	while((c=getchar())<'0'||c>'9')
		if(c=='-') f=-1;
	res=c-48;
	while((c=getchar())>='0'&&c<='9')
		res=res*10+c-48;
	return res*f;
}

inline ll f(ll x,ll c,ll n){
	return (x*x+c)%n;
}
inline ll myabs(ll x){
	return (x>0)?x:-x;
}
inline ll gcd(ll a,ll b){
	return (!b)?a:gcd(b,a%b);
}

inline ll pollard(ll x){ // 有的时候会有 s = t，此时成环了，不需要继续找下去
	ll s=0,t=0,val=1,c=rng()%(x-1)+1,d;
	for(int lim=1;;lim<<=1,s=t,val=1){
		for(int i=1;i<=lim;++i){
			t=f(t,c,x),val=myabs(s-t)*val%x;
			if(i%127==0){ 
				d=gcd(val,x);
				if(d>1) return d;
			}
		}
		d=gcd(val,x);
		if(d>1) return d;
	}
	return 0;
}
inline void print(ll x){
	if(x<0) putchar('-'),x=-x;
	if(x>9) print(x/10);
	putchar(x%10+'0');
	return;
}

inline void solve(ll x){
	if(x==1) return;
	if(miller::isprime(x)) return res=std::max(res,x),void();
	ll cur=x;
	while(cur>=x) cur=pollard(x);
	while(x%cur==0) x/=cur;
	solve(cur),solve(x);
	return;
}


signed main(void){
	int T=read();
	while(T--){
		ll n=read();
		res=0,solve(n);
		if(res==n) puts("Prime");
		else print(res),puts("");
	}
	return 0;
}