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}\) 内的素性检测。

struct MillerRabin{
	int prime[20]={2,3,5,7,11,13,17,19,23,29,31,37};
	inline bool test(int v,ll p){
		ll d=p-1,cur=qpow(v,d,p);
		if(cur!=1) return false;
		while(!(d&1)){
			d>>=1,cur=qpow(v,d,p);
			if(cur==p-1) return true;
			else if(cur!=1) return false;
		}
		return true;
	}
	inline bool check(ll x){
		if(x<=40){
			for(int i=0;i<12;++i) if(x==prime[i]) return true;
			return false;
		}
		for(int i=0;i<12;++i) if(!test(prime[i],x)) return false;
		return true;
	}
}T;

Pollard-Rho 算法

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

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

我们倍增一下。每次取 \([2^k,2^{k+1})\) 这一段，设段头的 \(f(x)\) 值为 \(s\)，这一段里面其它的值为很多个 \(t\)，那么每次把 \(|s - t|\) 乘起来对 \(t\) 取模，如果乘了 \(127\) 次就取一下模，然后求一遍 GCD。

反正这是玄学，它能找就当它能找吧。

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

typedef __int128 ll;

std::mt19937 rng(time(0));

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;
}

struct MillerRabin{
	int prime[20]={2,3,5,7,11,13,17,19,23,29,31,37};
	inline bool test(int v,ll p){
		ll d=p-1,cur=qpow(v,d,p);
		if(cur!=1) return false;
		while(!(d&1)){
			d>>=1,cur=qpow(v,d,p);
			if(cur==p-1) return true;
			else if(cur!=1) return false;
		}
		return true;
	}
	inline bool check(ll x){
		if(x<=40){
			for(int i=0;i<12;++i) if(x==prime[i]) return true;
			return false;
		}
		for(int i=0;i<12;++i) if(!test(prime[i],x)) return false;
		return true;
	}
}T;
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 PollardRho(ll x){ // 有的时候会有 s = t，此时成环了，不需要继续找下去
	ll s=0,t=0,val=1,c=rng()%(x-1)+1,d;
	for(int ga=1;;ga<<=1,s=t,val=1){
		for(int cur=1;cur<=ga;++cur){
			t=f(t,c,x),val=myabs(s-t)*val%x;
			if(cur==127){ 
				d=gcd(val,x);
				if(d>1) return d;
			}
		}
		d=gcd(val,x);
		if(d>1) return d;
	}
	assert(false);
	return 114514;
}
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(T.check(x)){
		res=std::max(res,x);
		return;
	}
	ll cur=x;
	while(cur>=x) cur=PollardRho(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;
}