前置知识:Millar Robin

之前已经学习了,试除法,线性筛等素数算法。

今天学习一个能在$O(n^{\frac{1}{4}})$时间内分解一个数的算法:Pollard-Rho。

首先介绍一下什么叫生日悖论。

有k个人,他们的生日均不相同的概率是多少?

按照理论上的计算的话,为$\prod\limits_{i=1}^k \frac{365-i+1}{365}$

生日有相同的情况就是,$1-\prod\limits_{i=1}^k \frac{365-i+1}{365}$

经过计算,当k=23时,有相同的概率就已经大于$\text{50%}$。

当$k=57$时,已经到达了$\text{99%}$。

但是在实际情况看,一个班内相同生日的同学往往少之又少。

由此就产生了生日悖论。

从随机角度来讲,生日悖论作为一个例子,很好的解释了只要我随机的够多,我就能骗分

对于一个数$x$,$x\le10^{16}$,要求出它的因子。

暴力枚举,枚中一个的几率是$\frac{1}{10^{16}}$

试除法,枚中一个的几率是$\frac{1}{10^{8}}$,还是不能接受。

不妨从枚举因子,变成$\gcd(b,N)>1$,这样几率又上升了,因为一个数,有可能含有$N$的因子,但不是$N$的因子。

选出这样一组$x_1,x_2,x_3,x_4,\cdots x_m$,若有$\gcd(|x_i-x_j|,N)>1$,则$\gcd(|x_i-x_j|,N)$是$N$的一个因子。

有早期论文中证明了大概要选出$N^{\frac{1}{4}}$大小的组即可选出一个因子,但是如何储存?这又是个问题,否则复杂度又会回弹。

于是有了一个随机序列$f(x)=x^2+c(\text{mod}~N)$,$c$是一个随机数。

但是这个序列,会有一个循环,像$\rho$一样。

因而有两种方式优化这个随机序列。

  1. Floyd判环,不常用,常数大,就是一个跳得慢,一个跳得快。
  2. 不判环,用小倍增来避免在环上绕太久。

主要讲解第二种方式,

由于$\gcd(a,b)=\gcd(b,a~\text{mod}~b)$,因此可以考虑把一组序列乘起来再进行一次$\gcd$

记录一段$[2^{k-1},2^k]$区间的贡献,与之前的区间贡献乘起来。

经过实验,当$2^k$达到128时,最优。

#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cmath>
#include<ctime>
#include<algorithm>
#define gc getchar()
#define ll long long
#define ull unsigned long long
using namespace std;
const int N=1e5+5,M=8555;
ll a[20]={2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 24251, 61};
ll inf=1ll<<62;
template<class o>
inline void qr(o &x)
{
    x=0;char c=gc;int f=1;
    while(c<'0'||c>'9'){if(c=='-')f=-1;c=gc;}
    while(c>='0'&&c<='9'){x=x*10+(c^48);c=gc;}
    x*=f;
}
template<class o>
void qw(o x)
{
    if(x<0)putchar('-'),x=-x;
    if(x/10)qw(x/10);
    putchar(x%10+48);
}
inline ull mul(ll a,ll b,ll p)
{
    ll x=(long double)a/p*b;
    ll c=1ull*a*b-1ull*x*p;
    return (c%p+p)%p;
}
inline ll fpow(ll a,ll b,ll p){ll ans=1;for(;b;b>>=1,a=mul(a,a,p))if(b&1)ans=mul(ans,a,p);return ans;}
ll gcd(ll a,ll b){return b?gcd(b,a%b):a;}
void exgcd(ll a,ll b,ll &x,ll &y){if(!b){x=1;y=0;return ;}exgcd(b,a%b,y,x);y-=(a/b)*x;}
inline ll f(ll x,ll c,ll p){return (mul(x,x,p)+c)%p;}
ll pr(ll n)
{
    ll s=0,t=0,c=1ll*rand()%(n-1)+1,val=1;
    for(int i=1;;i<<=1,s=t)
    {
        for(int j=1;j<=i;j++)
        {
            t=f(t,c,n);
            val=mul(val,abs(t-s),n);
            if(j%127==0)
            {
                ll d=gcd(val,n);
                if(d>1)return d;
            }
        }
        ll d=gcd(val,n);
        if(d>1)return d;
    }
}
bool mr(ll x,ll n)
{
    ll d=n-1;int cnt=0;
    while(d%2==0)d/=2,++cnt;
    ll ret=fpow(x,d,n);
    if(ret==1||ret==n-1)return 1;
    for(int i=1;i<=cnt;i++)
    {
        ret=mul(ret,ret,n);
        if(ret==n-1)return 1;
    }
    return 0;
}
bool is_prime(ll n)
{
    if(n<2||n==46856248255981ll)return 0;
    for(int i=0;i<12;i++)
        if(a[i]==n)return 1;
    return mr(2,n)&&mr(61,n);
}
ll mx;
void solve(ll n)
{
    if(n<=mx||n<2)return ;
    if(is_prime(n)){mx=mx<n?n:mx;return ;}
    ll p=n;
    while(p>=n)
        p=pr(n);
    while(n%p==0)n/=p;
    solve(n),solve(p);
}
int main()
{
    int T;qr(T);
    while(T--)    
    {
        srand(time(NULL));
        mx=0;ll x;qr(x);
        solve(x);
        if(mx==x)puts("Prime");
        else qw(mx),puts("");
    }
    return 0;
}