Miller-Rabin和Pollard-Rho

\(Miller-Rabin\)算法是一种高效的质数判断方法，虽然是一种不确定的质数判断方法，但在选择多种底数的情况下，正确率是可以接受的。
费马素性检验：
有费马小定理：若\(p\)为质数，则对任意\(1\leq x<p\)，有\(x^{p-1}\equiv 1\pmod p\)。
考虑其逆否命题：若存在\(1\leq x<p\)，有\(x^{p-1}\not\equiv 1\pmod p\)。
于是可以在\([1,p-1]\)中随机若干数，若存在\(x^{p-1}\not\equiv 1\pmod p\)，则\(p\)一定不为质数，否则\(p\)有较大概率为质数。
不幸的是，仍然有合数\(p\)，同样满足对任意\(1\leq x<p\)，有\(x^{p-1}\equiv 1\pmod p\)。
考虑优化费马素性检验。
二次探测定理：若\(p\)为奇素数，则\(x^2\equiv 1\pmod p\)的解为\(x\equiv \pm 1\pmod p\)。
于是考虑如下优化算法的正确性：
设选取的底数为\(x\)，初始化指数\(y=p-1\)。
判断\(x^y\mod p\)是否为\(1\)，不为\(1\)则返回\(p\)为合数。
令\(y=\frac{y}{2}\)，令\(v=x^y\ mod\ p\)，若\(v=-1\)返回质数，若\(v\neq 1\)返回合数，否则重复进行。
若\(y\)至奇数仍未错误，返回质数。
在程序竞赛范畴内，若要求在\(2^{32}\)内判素数，选取底数\(2,7,61\)即可，若要求在\(2^{64}\)内判素数，选取底数\(2,325,9375,28178,450775,9780504,1795265022\)即可。

\(Pollard-Rho\)算法能在期望\(O(n^{\frac{1}{4}})\)内计算合数\(n\)的某个非平凡因子。
先介绍生日驳论：一个班上有\(k\)个学生，问存在两学生生日相同概率。
\(P(k)=1-\prod_{i=0}^{k-1}\frac{365-i}{365}\),若要求\(P(k)>=\frac{1}{2}\)，解得\(k>=23\)。
因为\(gcd(k,n)\mid n\)，所以若\(gcd(k,n)>1\)，则找到\(n\)的一个非\(1\)因子，不妨设法找到这样的\(k\)。
不妨选取\(x_1,x_2,…,x_n\)，若存在\(gcd(|x_i-x_j|,n)>1\)，则找到\(n\)的一个因子。
考虑随机生成一个伪随机序列，相邻两项做差求\(gcd\)。
定义\(f(x)=(x^2+c)\mod n\)，\(x_i=f(x_{i-1})\)。
但这样生成的序列，会存在环，若经过环还未找到答案，则会进入死循环。所以需要一个方法判环，及时停止循环。
不妨让\(x_1\)以每次一步走环，\(x_2\)以每次两步走环，每次判断\(gcd(|x_1-x_2|,n)\)是否大于\(1\)。若出现\(x_1=x_2\)，说明\(x_2\)已经对\(x_1\)实现套圈，直接停止循环，并改换参数即可。
操作中大量调用\(gcd\)函数，复杂度较高，考虑优化。
\(gcd(kx\ mod\ n,n)=gcd(kx,n)>gcd(x,n)\)，所以可以考虑用多次\(|x_i-x_j|\)做积，统一检测。
但一次检测的正确性又会降低，考虑倍增，每次检测\(2^k\)步的乘积，即设初始为\(s\)，检测\(val=\prod_{i=1}^{2^k}|x_{s+i}-x_s|\)。
为防止\(k\)过大时复杂度过高，对于每一次倍增块内，每隔\(127\)步，进行一次检测。
总而言之，\(Pollard_Rho\)的目的只是找到\(n\)的任意一个因子，由于目标范围的宽容，可以考虑随机化的方法。
我们先通过生日驳论得知通过随机组合达到目标要远比随机单体达到目标概率大，于是用简单函数生成随机数列，并逐项组合。同时为了防止进入死循环，巧妙利用追击问题去判环。再运用数论的基础性质，去统一检测，减少\(gcd\)的调用。最后运用倍增的方式，并设定检测阈值，去平衡随机过程的正确性和复杂度。
一道模板题
\(P4718\)
https://www.luogu.com.cn/problem/P4718

#include <bits/stdc++.h>
#define int __int128
#define LL long long

using namespace std;

struct Miller_Rabin{
	int qpow(int a,int b,int p){
		int res=1;
		while(b){
			if(b&1) res=res*a%p;
			a=a*a%p;
			b>>=1;
		}
		return res;
	}
	bool check(int x,int y,int p){
		if(x>=y) return 1;
		if(qpow(x,y,p)!=1) return 0;
		while(y%2==0){
			y/=2;
			int v=qpow(x,y,p);
			if(v==p-1) return 1;
			if(v!=1) return 0;
		}
		return 1;
	}
	bool is_prime(int x){
		if(x==2) return 1;
		int p[12]={2,3,5,7,11,13,17,19,23,29,31,37};
		for(int i=0; i<12; i++) if(!check(p[i],x-1,x)) return 0;
		return 1;
	}
}mr;

struct Pollard_Rho{
	int f(int x,int c,int p){
		return (x*x+c)%p;
	}
	int gcd(int a,int b){
		return b?gcd(b,a%b):a;
	}
	int find(int x){
		int c=rand()%(x-1)+1;
		int s=0,goal=1;
		while(1){
			int t=s,v=1;
			for(int step=1; step<=goal; step++){
				t=f(t,c,x);
				v=v*(t>s?t-s:s-t)%x;
				if(!v) return x;
				if(step%127==0){
					int d=gcd(v,x);
					if(d>1) return d;
				}
			}
			int d=gcd(v,x);
			if(d>1) return d;
			s=t,goal*=2;
		}
	}
}pr;

void solve(int x,LL &mx){
	if(x<=mx||x<2) return;
	if(mr.is_prime(x)){
		if(x>mx) mx=x;
		return;
	}
	int t=0;
	while(1){
		t=pr.find(x);
		if(t!=x) break;
	}
	while(x%t==0) x/=t;
	solve(x,mx),solve(t,mx);
}

signed main(){
	LL T;
	cin >> T;
	while(T--){
		LL n;
		cin >> n;
		LL mx=0;
		solve(n,mx);
		if(n==mx) cout << "Prime" << endl;
		else cout << mx << endl;
	}
	return 0;
}