JZTXT

Pollard-Rho算法

发布时间 2023-11-02 00:14:54作者: carp_oier

prelogue

怎么感觉我这个人和随机化关系这么好。

鲤鱼我是从这篇博客中进行学习的。

Pollard-Rho 算法

Pollard-Rho 算法是一种求非 1 非自身的因子的高效算法。

main body

我们求素数平常是用的复杂度为 \(O(sqrt(n))\) 的试除法，如果 \(n\) 这个数字很大，我们这个算法就很低效，需要一些卡常技巧来帮助我们卡过去这个题目。

但是想从根本来解决问题，我们就得要换一种更加高效的算法 Pollard-Rho 算法。

这个算法是一种带有随机化思想的方法。有点类似于猴子排序（就是那个很 naive 很 trivial 的算法，你能说他效率低嘛，那只能说明你 rp 不好）

我们如果单纯的用这个随机数来进行匹配，效率极其低下，甚至可以说是 \(O(+\infty)\) 的。那我们就要想着怎么对它尽心优化。

在这之前，我们要先引入一个知识点：生日悖论

生日悖论

最后反正就是一顿操作猛如虎，然后就给硬生生将一个 \(O(+\infty)\) 的做法砍成了 \(O(n ^ {\frac{1}{4}})\)

代码如下：

code time

#include <bits/stdc++.h>
using namespace std;
#define ll long long
#define rl register ll
#define fom(i, a) for(rl i=a; i; -- i)
#define foa(i, a, b) for(rl i=a; i < b; ++ i)
#define fos(i, a, b) for(rl i=a; i <= b; ++ i)
#define fop(i, a, b) for(rl i=a; i >= b; -- i)

namespace IO
{
    int pp1=-1; const int pp2=(1<<21)-1; char buf[1<<21],*p1=buf,*p2=buf,buffer[1<<21];
    inline void flush() {fwrite(buffer,1,pp1+1,stdout); pp1=-1;}
    inline void pc(const char ch) {if(pp1==pp2) flush(); buffer[++pp1]=ch;}
    inline char gc(){ return p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++;}
    template <class T> inline void read(T &res){char ch=gc();bool f=0;res=0; for(;!isdigit(ch);ch=gc()) f|=ch=='-'; for(;isdigit(ch);ch=gc()) res=(res<<1)+(res<<3)+(ch^48); res=f?~res+1:res;}
    template <class T> inline void ww(T x) { if(!x) pc('0'); else { static int stk[21]; int top = 0; if(x<0) pc('-'),x=~x+1; while(x) stk[top++]=x%10, x/=10; while(top--) pc(stk[top]^48);}}
}

#define ws IO::pc(' ')
#define wl IO::pc('\n')
#define ww(x) IO::ww(x)
#define read(x) IO::read(x)
#define flush() IO::flush()

ll max_factor, n;

// mt19937 rand(time(0));

inline ll gcd(ll a, ll b) { return b ? gcd(b, a % b) : a; }

inline ll qmi(ll a, ll k, ll p)
{
    ll res = 1;

    while(k)
    {
        if(k & 1) res = ((__int128)res * a) % p;
        a = ((__int128)a * a) % p, k >>= 1;
    }

    return res;
}

inline bool Miller_Rabin(ll p)
{
    if(p < 2) return false;
    if(p == 2 || p == 3) return true;
    ll d = p - 1, r = 0;

    while(!(d & 1)) ++ r, d >>= 1;
    
    foa(k, 0, 10) 
    {
        ll a = rand() % (p - 2) + 2;
        ll x = qmi(a, d, p);
        if(x == 1 || x == p - 1) continue;
        foa(i, 0, r - 1)
        {
            x = (__int128)x * x % p;
            if(x == p - 1) break;
        }
        if(x != p - 1) return false;
    }
    return true;
}

inline ll Pollard_Rho(ll n)
{
    ll s = 0, t = 0;
    ll c = (ll)rand() % (n - 1) + 1;
    ll step = 0, goal = 1, val = 1;

    for(goal = 1;; goal *= 2, s = t, val = 1)
    {
        for(step = 1; step <= goal; ++ step)
        {
            t = ((__int128)t * t + c) % n;
            val = ((__int128)val * abs(t - s)) % n;
            if((step % 127) == 0) { ll d = gcd(val, n); if(d > 1) return d; }
        }
        ll d = gcd(val, n);
        if(d > 1) return d;
    }
}

inline void fac(ll x)
{
    if(x <= max_factor || x < 2) return ;
    if(Miller_Rabin(x)) { max_factor = max(max_factor, x); return ; }

    ll p = x;
    while(p >= x) p = Pollard_Rho(x);
    while((x % p) == 0) x /= p;
    fac(x), fac(p); 
}

inline void solve()
{
    srand(time(0));
    max_factor = 0; read(n);
    fac(n);
    if(max_factor == n) flush(), puts("Prime");
    else ww(max_factor), wl;
}

int main()
{
    // freopen("1.in", "r", stdin);
    ll T; read(T); while(T --) solve(); 
    flush(); return 0;
}