Miller Rabin & Pollard Rho
Miller_Rabin
用于检测大数素性($ \sqrt{n} \ge 1e8 $).
对于素数 $ P $ ,有费马小定理:
- 对于任意 $ a \in \lbrack 1,P) , a^{P-1} \equiv 1 : (mod : P) $
枚举 $ \lbrack 1,P ) $ 中的任意 $ a $ ,若均满足其逆定理 $ a^{P-1} \equiv 1 : ( mod : p ) $ ,则 $ P $ 大概率 为素数。
上述即费马素性检验。
存在卡迈克尔数/费马伪素数也满足上式 (如 $ 561=3 * 11 * 17 $) ,故需要优化费马素性检验。
二次探测定理:
- 若 $ P $ 为奇素数,则 $ x^2 \equiv 1 : (mod : P) $ 的解为 $ x \equiv \pm 1 : (mod : P) $ .
设当前选取底数为 $ a $ , $ u \cdot 2^t=n-1 : $ ($ u $ 为奇数),对于 $ x=u \cdot 2^y : (y \in \lbrack 0,t \rbrack) $ ,满足 $ a^{2x} \equiv 1 : (mod : n)$ 且 $ a^x \not\equiv \pm1 : (mod : n) $ ,则 $ n $ 为合数,否则大概率为质数。
bool check(ll a,ll n){
ll u=n-1,t=0;
while(!(u&1))u>>=1,t++;
ll x=qpow(a,u,n);
if(x==1)return false;
for(ll i=1,tmp;i<=t;i++,x=tmp){
tmp=qmul(x,x,n);//原题需要转__int128类型
if(x!=1&&x!=n-1&&tmp==1)
return true;
}
return x!=1;
}
bool Miller_Rabin(ll n){
if(n==2)return true;
if(n<2||!(n&1))return false;
for(ll i=1;i<=Test;i++){
if(check(rd()%(n-1)+1,n))//mt19937
return false;
}
return true;
}
Miller-Rabin的错误率是 $ 4^{-Test} $ ,为提高程序效率取了 \(10\) 。(参考 这里)
Pollard Rho
用 $ O(n^{\frac{1}{4}}) $ 的时间分解大数。
考虑随机函数 $ f(x)=(x^2+c) : mod : n $ ,初始值为 \(x_0\) ,令 $ x_1= f(x_0)$ ,然后再将 \(x_1\) 代入函数,数列轨迹形似 $ \rho $ 状,故称Pollard Rho.
对于 {\(x\)} ,该数列基本随机,取相邻两项取差作gcd.
基于Floyd算法的优化
数列有环存在,可以用一种Floyd的判圈方法:
设 $ t=r=0 ,2v_t=v_r$ ,则在第 \(i\) 时刻,\(x_t\) 对应 $ x_i $ , \(x_r\) 对应 \(x_{2i}\) ,两者相遇时结束。
基于Floyd判圈方法的Pollard Rho算法:
ll f(ll x,ll c,ll n){
return (qmul(x,x,n)+c)%n;
}
ll Pollard_Rho(ll n){
ll c=rd(n-1)+1;
ll t=f(0,c,n),r=f(t,c,n);
while(t!=r){
ll d=gcd(abs(t-r),n);
if(d>1)return d;
t=f(t,c,N),r=f(f(r,c,n),c,n);
}
return N;//多测保证正确性
}
路径倍增优化与乘法累积
上述优化瓶颈为 \(gcd\) 的 \(log\) 复杂度,由于当 \((a,b)>1\) 时,一定有 \((ac \: mod \: b ,b)>1\),可以通过乘法累积优化复杂度。
考虑倍增思想,在路径中取一段 \(\lbrack l=2^{k-1},r=2^k \rbrack\) ,以所有 \(|x_i-x_l| \: (l < i \le r)\) 为样本,以避免在环上过久或gcd的复杂度过大。
需要规定一个样本数限制相乘个数,这里取 \(127=2^7-1\)(据说经实验得出迭代7次作为上界为较优方案,具体不详。)
ll f(ll x,ll c,ll n){
return (qmul(x,x,n)+c)%n;
}
ll Pollard_Rho(ll n){
if(!(n&1))return 2;
ll x=rd()%(n-1)+1,y=x,c=rd()%(n-1)+1,ret;
for(ll i=1;;i<<=1,y=x,ret=1){
for(ll j=1;j<=i;j++){
x=f(x,c,n);
ret=qmul(ret,abs(x-y),n);
if(!(j&127)){
ll d=__gcd(ret,n);
if(d>1)return d;
}
}
ll d=__gcd(ret,n);
if(d>1)return d;
}
}
粗略大概就是这么多,放个模板代码。
#include<bits/stdc++.h>
#define ll long long
#define Test 10
using namespace std;
mt19937 rd(233);
ll read(){
ll x=0;
char ch=getchar();
while(!isdigit(ch))ch=getchar();
while(isdigit(ch))x=x*10+(ch^48),ch=getchar();
return x;
}
ll qmul(ll k,ll b,ll p){
return (__int128)k*b%p;
}
ll qpow(ll k,ll b,ll p){
ll val=1;k%=p;
while(b){
if(b&1)val=qmul(val,k,p);
k=qmul(k,k,p),b>>=1;
}
return val;
}
ll T,n,ans;
bool check(ll a,ll n){
ll u=n-1,t=0;
while(!(u&1))u>>=1,t++;
ll x=qpow(a,u,n);
if(x==1)return false;
for(ll i=1,tmp;i<=t;i++,x=tmp){
tmp=qmul(x,x,n);
if(x!=1&&x!=n-1&&tmp==1)
return true;
}
return x!=1;
}
bool Miller_Rabin(ll n){
if(n==2)return true;
if(n<2||!(n&1))return false;
for(ll i=1;i<=Test;i++){
if(check(rd()%(n-1)+1,n))
return false;
}
return true;
}
ll f(ll x,ll c,ll n){
return (qmul(x,x,n)+c)%n;
}
ll Pollard_Rho(ll n){
if(!(n&1))return 2;
ll x=rd()%(n-1)+1,y=x,c=rd()%(n-1)+1,ret;
for(ll i=1;;i<<=1,y=x,ret=1){
for(ll j=1;j<=i;j++){
x=f(x,c,n);
ret=qmul(ret,abs(x-y),n);
if(!(j&127)){
ll d=__gcd(ret,n);
if(d>1)return d;
}
}
ll d=__gcd(ret,n);
if(d>1)return d;
}
}
void find(ll n){
if(n<=ans)return;
if(Miller_Rabin(n)){
ans=max(ans,n);
return;
}
ll p=n;
while(p>=n)p=Pollard_Rho(n);
while(n%p==0)n/=p;
find(p),find(n);
}
int main(){
T=read();
while(T--){
n=read();
ans=1,find(n);
if(ans==n)printf("Prime\n");
else printf("%lld\n",ans);
}
return 0;
}