POJ 2429 GCD & LCM Inverse（Miller-Rabbin素性测试，Pollard rho质因子分解）

x = lcm/gcd，假设答案为a,b，那么a*b = x且gcd(a,b) = 1，因为均值不等式所以当a越接近sqrt(x)，a+b越小。

x的范围是int64的，所以要用Pollard_rho算法去分解因子。因为a,b互质，所以我们把相同因子一起处理。

最多16个不同的因子：2,3,5,7,11,13,17,19,23,29,31,37,41,43,47, 乘积为 614889782588491410， 乘上下一个质数53会爆int64范围。

所以剩下暴力枚举一下就好。

#include<cstdio>
#include<iostream>
#include<string>
#include<cstring>
#include<queue>
#include<vector>
#include<stack>
#include<vector>
#include<map>
#include<set>
#include<algorithm>
#include<cmath>
#include<ctime>
//#include<bits/stdc++.h>
using namespace std;

typedef long long ll;

ll mulMod(ll a, ll b, ll m) // a %= m
{
    ll re = 0;
    while(b){
        if(b&1){
            re = (re+a); if(re >= m) re -= m;
        }
        a <<= 1;
        if(a >= m) a -= m;
        b >>= 1;
    }
    return re;
}

ll powMod(ll a,ll q,ll m)
{
    ll re = 1;
    while(q){
        if(q&1){
            re = mulMod(re,a,m);
        }
        a = mulMod(a,a,m);
        q >>= 1;
    }
    return re;
}


bool witness(ll a,ll d,int t,ll n)
{
    ll x = powMod(a,d,n), y;
    while(t--){
        y = mulMod(x,x,n);
        if(y == 1 && x != 1 && x != n-1) return true;
        x = y;
    }
    return x != 1;
}

bool RabinMiller(ll n,int k = 10)
{
    if(n == 2) return true;
    if(n < 2 || ((n&1)^1) ) return false;
    int t = 0;
    ll d = n-1;
    while((d&1)^1) { d>>=1; t++; }
    while(k--){
        if(witness(rand()%(n-1)+1,d,t,n)) return false;
    }
    return true;
}

ll gcd(ll a,ll b)
{
    while(b){
        ll t = a%b;
        a = b;
        b = t;
    }
    return a < 0? -a: a;
}

ll rho(ll n)
{
    ll x = rand()%(n-1) + 1;
    ll y = x, d;
    int i = 1, k = 2;
    ll c = 1+rand()%(n-1);
    while(true){
        i++;
        x = (mulMod(x,x,n) + c);
        if(x >= n) x -= n;
        d = gcd(y-x,n);
        if(d != 1 && d != n) return d;
        if(y == x) return n;
        if(i == k){
            y = x;
            k <<= 1;
        }
    }
}

ll Prms[150], stk[150];
int tot;

void pollard(ll n)
{
    tot = 0;
    int top = 0;
    stk[++top] = n;
    while(top){
        n = stk[top--];
        if(RabinMiller(n)){
            Prms[tot++] = n;
        }else{
            ll p = n;
            while(p >= n) p = rho(p);
            stk[++top] = p;
            stk[++top] = n/p;
        }
    }
}

ll mpow(ll a,int q)
{
    ll re = 1;
    while(q){
        if(q&1) re *= a;
        a *= a;
        q>>=1;
    }
    return re;
}

//#define LOCAL
int main()
{
#ifdef LOCAL
    freopen("in.txt","r",stdin);
#endif
    ll g,l;
    while(~scanf("%I64d%I64d",&g,&l)){
        if(l == g){
            printf("%I64d %I64d
",g,l);
            continue;
        }
        ll x = l/g;
        pollard(x);
        sort(Prms,Prms+tot);
        Prms[tot] = -1;
        int k = 0;
        for(int j = 0,i = 1, cnt = 1; i <= tot; i++) {
            if(Prms[i] == Prms[j]) cnt++;
            else {
                Prms[k++] = mpow(Prms[j],cnt);
                cnt = 1; j = i;
            }
        }
        ll sqr = floor(sqrt(x)), best = 1;
        for(int S = 1<<(tot = k); --S; ){
            ll cur = 1;
            for(int i = 0; i < tot; i++){
                if(S>>i&1) cur *= Prms[i];
                if(cur > sqr) break;
            }
            if(cur <= sqr && cur > best) best = cur;
        }
        printf("%I64d %I64d
",g*best,l/best);
    }

    return 0;
}