Loading

POJ - 2429 GCD & LCM Inverse

GCD & LCM Inverse

大整数质因数分解

主要是一个 \(Pollard_Rho\) 板子想放上来

\(GCD\)\(LCM\) 质因数分解之后,对于每个质因数,有相应的两个幂 \(a^{p_1}\)\(a^{p_2}\),且 \(a^{p_1} <> a^{p_2}\),对于每个这样的质因数,就考虑如何分配到 \(a\)\(b\) 之中是最优的

这个推不出来,所以就直接二进制枚举,找到最优答案

!!!这题有个巨大坑,他输出两个数,要求第一个是比较小的,第二个是比较大的

#include <iostream>
#include <cstdio>
#include <ctime>
#include <cstdlib>
#include <map>
#include <vector>
using namespace std;
typedef long long ll;

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

inline ll ksc(ll x, ll y, ll mod)
{
    return (x * y - (ll) ( (long double) x / mod*y )*mod + mod ) % mod;     
}

// 快速幂
ll qpow(ll x, ll n, ll m)
{
    ll ans = 1;
    x %= m;
    while (n)
    {
        if (n & 1)
            ans = ksc(x, ans, m);
        x = ksc(x, x, m);
        n >>= 1;
    }
    return ans;
}

// 判断大质数
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;
    for (ll k = 0; k < 10; k++)
    {
        ll a = rand() % (p - 2) + 2;
        ll x = qpow(a, d, p);
        if (x == 1 || x == p - 1)
            continue;
        for (int i = 0; i < r - 1; i++)
        {
            x = ksc(x, x, p);
            if (x == p - 1)
                break;
        }
        if (x != p - 1)
            return 0;
    }
    return 1;
}

// 大质数质因数分解,分治
ll Pollard_Rho(ll x)
{
    ll s = 0, t = 0;
    ll c = (ll)rand() % (x - 1) + 1;
    int step = 0, goal = 1;
    ll val = 1;
    for (goal = 1;; goal <<= 1, s = t, val = 1)
    {
        for (step = 1; step <= goal; step++)
        {
            t = (ksc(t, t, x) + c) % x;
            val = ksc(val, abs(t- s), x);
            if ((step % 127) == 0)
            {
                ll d = gcd(val, x);
                if (d > 1)
                    return d;
            }
        }
        ll d = gcd(val, x);
        if (d > 1)
            return d;
    }
}

// 多次试探,结果储存在中
map<ll, int> mp;
void fac(ll x)
{
    if(x == 1) return; // 切记
    if (Miller_Rabin(x))
    {
        mp[x]++;
        return;
    }
    ll p = x;
    while (p >= x)
        p = Pollard_Rho(x);
    // while ((x % p) == 0)
    //     x /= p;
    fac(x / p);
    fac(p);
}

ll ax[110], px[110], py[110];
int main()
{
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);
    srand(time(0)); // 种子
    ll a, b;
    while(cin >> a >> b)
    {
        int tp = 0;
        mp.clear();
        fac(b);
        for(map<ll, int>::iterator it = mp.begin(); it != mp.end(); it++)
        {
            ax[tp] = it->first;
            px[tp] = 1;
            for(int i=0; i<it->second; i++)
                px[tp] *= ax[tp];
            tp++;
        }
        mp.clear();
        fac(a);
        int st = 0;
        for(int i=0; i<tp; i++) py[i] = 1;
        for(map<ll, int>::iterator it = mp.begin(); it != mp.end(); it++)
        {
            while(ax[st] != it->first) st++;
            for(int i=0; i<it->second; i++)
                py[st] *= ax[st];
        }
        ll ans = a + b, ans_a = a, ans_b = b;
        for(int i=(1<<tp)-1; i>=0; i--)
        {
            ll xa = 1, xb = 1;
            for(int j=0; j<tp; j++)
            {
                if(i >> j & 1)
                {
                    xa *= px[j];
                    xb *= py[j];
                }
                else
                {
                    xb *= px[j];
                    xa *= py[j];
                }
            }
            if(xa + xb < ans)
            {
                ans = xa + xb;
                ans_a = xa;
                ans_b = xb;
            }
        }
        if(ans_a > ans_b) swap(ans_a, ans_b);
        cout << ans_a << " " << ans_b << endl;
    }
    return 0;
}
posted @ 2022-10-11 14:54  dgsvygd  阅读(20)  评论(0编辑  收藏  举报