poj 2429 GCD & LCM Inverse

Given two positive integers a and b, we can easily calculate the greatest common divisor (GCD) and the least common multiple (LCM) of a and b. But what about the inverse? That is: given GCD and LCM, finding a and b.

Input

The input contains multiple test cases, each of which contains two positive integers, the GCD and the LCM. You can assume that these two numbers are both less than 2^63.

Output

For each test case, output a and b in ascending order. If there are multiple solutions, output the pair with smallest a + b.

Sample Input

3 60

Sample Output

12 15

a / gcd * b = lcm
a / gcd * b / gcd = lcm / gcd
a / gcd 和 b / gcd 互质。
pollard_rho因数分解 lcm / gcd，来求答案，对于大数判断素数采用rabin_miller强伪素数测试，小数查表。
两个东西都不会，看到别人写的。
代码：

#include <iostream>
#include <cstdio>
#include <vector>
#include <map>
#include <cstdlib>
using namespace std;
#define PMAX 10000
typedef long long LL;
vector<int> primes;
vector<bool> is_prime;
inline LL mod_mult(LL a,LL b,LL m) {///乘法取模
    LL d = 0;
    a %= m;
    while(b) {
        if(b & 1) {
            d += a;
            if(d > m)d -= m;
        }
        a <<= 1;
        if(a > m)a -= m;
        b >>= 1;
    }
    return d;
}
LL mod_pow(LL a, LL b, LL m) {///幂取模
    LL d = 1;
    a %= m;
    while(b) {
        if(b & 1)d = mod_mult(d,a,m);
        a = mod_mult(a,a,m);
        b >>= 1;
    }
    return d;
}
LL gcd(LL a, LL b) {///求最小公倍数
    return b ? gcd(b,a % b) : a;
}
bool rabin_miller(LL n, LL times)///Rabin-Miller强伪素数测试
{
    if (n < 2) return false;
    if (n == 2) return true;
    if (!(n & 1)) return false;

    LL q = n - 1;
    int k = 0;
    while (q % 2 == 0) {
        k++;
        q >>= 1;
    }
    /// n - 1 = 2^k * q (q是奇素数)
    /// n是素数的话，一定满足下面条件
    /// (i) a^q ≡ 1 (mod n)
    /// (ii) a^q, a^2q,..., a^(k-1)q 中的某一个对n求模为-1
    ///
    /// 所以、当不满足(i)(ii)中的任何一个的时候，就有3/4的概率是合成数
    for (int i = 0; i < times; ++i) {
        LL a = rand() % (n - 1) + 1; /// 从1,..,n-1随机挑一个数
        LL x = mod_pow(a, q, n);
        /// 检查条件(i)
        if (x == 1) continue;
        /// 检查条件(ii)
        bool found = false;
        for (int j = 0; j < k; j++) {
            if (x == n - 1) {
                found = true;
                break;
            }
            x = mod_mult(x, x, n);
        }
        if (found) continue;
        return false;
    }
    return true;
}
LL pollard(LL n,int c)///Pollard 因数分解算法 c自己定
{
    LL x = rand() % (n - 1) + 1;
    LL y = x;
    LL d = 1;
    while (d == 1) {
        x = mod_mult(x, x, n) + c;
        y = mod_mult(y, y, n) + c;
        y = mod_mult(y, y, n) + c;
        d = gcd((x - y >= 0 ? x - y : y - x), n);
    }
    if (d == n) return pollard(n, c + 1);
    return d;
}
void init_primes() {///筛选素数并生成小范围素数表
    is_prime = vector<bool>(PMAX + 1, true);
    is_prime[0] = is_prime[1] = false;
    for (int i = 2; i <= PMAX; ++i) {
        if (is_prime[i]) {
            primes.push_back(i);
            for (int j = i * i; j <= PMAX; j += i) {
                is_prime[j] = false;
            }
        }
    }
}
/// 判断是否是素数，优先查表，如果n很大用Rabin-Miller强伪素数测试
bool isPrime(LL n) {
    if (n <= PMAX) return is_prime[n];
    else return rabin_miller(n,1);
}

/// 分解成素因子，小数用素数表，大数用Pollard 因数分解算法
void factorize(LL n, map<LL, int>& factors) {
    if (isPrime(n)) {
        factors[n]++;
    }
    else {
        for (int i = 0; i < primes.size(); ++i) {
            int p = primes[i];
            while (n % p == 0) {
                factors[p]++;
                n /= p;
            }
        }
        if (n != 1) {
            if (isPrime(n)) {
                factors[n]++;
            }
            else {
                LL d = pollard(n, 1);
                factorize(d, factors);
                factorize(n / d, factors);
            }
        }
    }
}

pair<LL, LL> solve(LL a, LL b) {
    LL c = b / a;///  lcm / gcd
    map<LL, int> factors;
    factorize(c, factors);

    vector<LL> mult_factors;/// 每个质因子的n次方，比如2^2和5^1
    for (map<LL, int>::iterator it = factors.begin();it != factors.end();it ++) {
        LL mul = 1;
        while (it->second) {
            mul *= it->first;
            it->second--;
        }
        mult_factors.push_back(mul);
    }
    LL best_sum = 1e18, best_x = 1, best_y = c;
    for (int mask = 0; mask < (1 << mult_factors.size()); ++mask) {/// 这是一个取数的过程，一共 2^size 种情况
        LL x = 1;
        for (int i = 0;i < mult_factors.size();i ++) {
            if (mask & (1 << i)) x *= mult_factors[i];
        }
        LL y = c / x;
        if (x < y && x + y < best_sum) {
            best_sum = x + y;
            best_x = x;
            best_y = y;
        }
    }
    return make_pair(best_x * a, best_y * a);
}

int main() {
    init_primes();
    LL a, b;
    while(~scanf("%lld%lld",&a,&b)) {
        pair<LL, LL> ans = solve(a, b);
        printf("%lld %lld
",ans.first,ans.second);
    }
    return 0;
}