Pollard-Rho

一种复杂度大概在$O(n^{\frac 1 4} \log n)$的分解质因数方法。

Miller-Rabin

给定一个$10^{18}$范围的数，判断质数

由于费马小定理，我们知道当$a\bmod p \neq 0$且$p$为 质数 时，$a^{p-1} \equiv 1 \pmod p$，而如果$p$不为质数，这东西有一定几率成立。

所以我们有一个想法，随机一些$a$，对每个$a$用这个式子判断一下，如果有一个错误了它就直接是合数了。

但是有些合数更特殊，对于$1\leq a < p$都成立这个式子，比如$561$。

于是有了个二次探测定理。如果$x^2 \equiv 1 \pmod p$则有$x\equiv \pm 1 \pmod p$，其中$p$为质数。

证明很显然，$p | x^2 - 1 \to p | (x+1)(x-1)$，所以$p | (x+1)$或$p | (x-1)$。

然而仍然，对于合数这个东西有一定几率不成立。

这个算法的做法就是，检测$a^{p-1}$模$p$下的值，如果不是 1 或者 -1 那么必然不是合数。如果算出来是 1 并且$2 | \frac{p-1} 2$，那么递归下去算$\frac{p-1} 2$，如果算出来是$p-1$直接返回为质数。

这样做并不是一定正确的，所以我们得多拿几个$a$来测。

复杂度是$O(k\log^2 n)$的。

具体而言：

-$n \leq 4\times 10^9$时，$2,7,61$即可保证正确
-$n \leq 10^{18}$时，$2,3,5,7,11,13,17,19,23$可以保证正确
-$n \leq 10^{19}$时，$2,3,5,7,11,13,17,19,23,29,31,37$可以保证正确。

虽然多数时候这个复杂度不会成为瓶颈，但是它仍然可以被优化。可以预处理出$p - 1$的分解中 2 的个数$k$然后先算出$a^{\frac{p-1}{2^k}}$再一步步平方回去，这样复杂度就可以优化掉 ksm 变成$O(k\log n)$

Pollard-Rho

现在已经完全明白了如何判质数，那么怎么质因数分解呢？

首先用 Miller-Rabin 判一下它是不是质数，是就直接 return 了。

否则，它一定有一个$\leq \sqrt n$的因子$p$

于是有一个乱搞思路，随机选出$n^{\frac 1 4}$个数，两两做差再求 gcd ，期望情况下可以分解出来，但是这样并不优秀。

所以考虑设函数$f(i) = (f(i-1)^2 + a) \pmod n$。其中$a$是一个常数。我们发现，这样的话$f$必然会循环。它最后的形状很类似一个$\rho$，先经历尾巴再进入一个循环。（同时这也是为啥名字叫 Pollard-Rho）。

显然，即使是在$\bmod p$的意义下它仍然是循环的，我们可以认为这个 函数 的数字出现是随机的，所以$\bmod p$意义下的最终环的长度是$O(p^\frac 1 2) = O(n^ \frac 1 4)$的。

具体应该怎么做呢？ Pollard-Rho 有两种实现方法：

• Floyd 判环

  一种比较好写好理解的方法，但是貌似不如第二种优秀。

  做法是，每次让$s$走一步，$t$走两步，然后求一下$y-x$和$n$的$\gcd$。但是如果$x = y$了就说明进入了$\bmod n$的循环，就挂掉了，分解失败，需要再次随机继续来。

  如果我们可以在某一时刻发现$\gcd(y-x , n) \neq 1$这就意味着我们找到了一个$\bmod p$意义下的环，并且没有找到$\bmod n$意义下的环。成功进行了一次分解！

好现在可以总结一下，Pollard-Rho 的本质是在环上跑，目的是在进入$\bmod n$的环之前进入$\bmod p$的环。如果做到了，那就可以分解出最小质因数$p$的某个倍数（并且也是$n$的约数），否则（比如同时进入这两个环）就需要调参重来。

然后是一种更优秀的判环方法：

• Brent 判环

  每次$x$走一步，当它走了$2^k$时让$y = x$，然后$x$继续走。同时在$x$走了$128$的倍数步的时候，我们 check 一次。

  我也不知道为啥$128$这个数字很优秀，但是这样做确实很优秀。

#include "iostream"
#include "algorithm"
#include "cstring"
#include "cstdio"
#include "cmath"
#include "vector"
using namespace std;
typedef long long ll;

ll mul( ll a , ll b , ll md ) {
    return ( a * b - (ll)( (long double)a / md * b + 0.5 ) * md + md ) % md;
}
ll Pow( ll x , ll a , ll md ) {
    ll cur = x % md , ans = 1;
    while( a ) {
        if( a & 1 ) ans = mul( ans , cur , md );
        cur = mul( cur , cur , md ) , a >>= 1;
    }
    return ans;
}

const int ck[] = {2,3,5,7,11,13,17,19,23} , _l = 8;
bool miller( ll n ) {
    if( n == 1 ) return false;
    ll t = n - 1; int cn = 0;
    while( !( t & 1 ) ) t >>= 1 , ++ cn;
    for( int i = 0 ; i < _l ; ++ i ) {
        if( n == ck[i] ) return true;
        ll a = Pow( ck[i] , t , n ) , nex = a;
        for( int j = 1 ; j <= cn ; ++ j ) {
            nex = mul( a , a , n );
            if( nex == 1 && a != 1 && a != n - 1 ) return false;
            a = nex;
        }
        if( a != 1 ) return false;
    }
    return true;
}

inline ll f( ll x , ll c , ll md ) {
    return ( mul( x , x , md ) + c ) % md;
}

inline ll _rand(  ) {
    return (ll) rand() << 32 | rand();
}
inline ll _randw() {
    return (ll)rand() << 48 | (ll)rand() << 32 | rand() << 16 | rand();
}
inline ll _abs( ll x ) {
    return x > 0 ? x : -x;
}
inline ll gcd( ll a , ll b ) {
    return !b ? a : gcd( b , a % b );
}

inline ll pollard_rho( ll n ) {
    ll s = 0 , t = 0 , c = _rand() % ( n - 1 ) + 1 , val = 1;
    for( int cir = 1 ; ; cir <<= 1 , s = t , val = 1 ) {
        for( int i = 0 ; i < cir ; ++ i ) {
            t = f( t , c , n ) , val = mul( val , _abs( t - s ) , n );
            if( i % 127 == 0 ) {
                ll g = gcd( val , n );
                if( g != 1 ) return g;
            }
        }
        ll g = gcd( val , n );
        if( g != 1 ) return g;
    }
}

vector<ll> divs;
inline void analyze( ll n ) {
    if( n == 1 ) return;
    if( miller( n ) ) { divs.push_back( n ); return; }
    ll d = n;
    while( d == n ) d = pollard_rho( n );
    n /= d;
    analyze( n ) , analyze( d );
}

int main() {
    srand( time(0) );
    int T;cin >> T;
    ll n;
    while( T-- ) {
        scanf("%lld",&n);
        if( miller( n ) ) { puts("Prime"); continue; }
        divs.clear();
        analyze( n );
        printf("%lld\n",*max_element( divs.begin() , divs.end() ) );
    }
}