POJ 2429 GCD & LCM Inverse(Miller_Rabin大素数测试+Pollar_rho素因子分解)

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

题意：

给定两个数的gcd和lcm求这两个数，并要求两个数是满足条件的且之和最小的答案

分析：

$\because lcm(a,b) = \frac{a\times b}{gcd(a,b)}$

$\therefore lcm \times gcd = a\times b$

$\therefore \frac{lcm}{gcd} = \frac{a}{gcd}\times \frac{b}{gcd}$

$\therefore \frac{a}{gcd}和 \frac{b}{gcd}$是互质的

我们令$n = \frac{lcm}{gcd},k = \frac{a}{gcd},\frac{n}{k}=\frac{b}{gcd}$

所以我们只需要对n进行素因子分解，然后组成两个互质的因子相乘得到n即可，因为要求和最小所以要尽量接近$\sqrt{n}$

比如n分解成了$n = p_1^{a_1}p_2^{a_2}p_3^{a_3}\cdots p_r^{a_r}$

这样我们只需要组合素质因子的幂这样组成的一个因子和剩下的一个因子必然互质，然后再保证接近$\sqrt{n}$就行。

又因为n比较大，所以只能采用大素数测试即Miller_Rabin素性测试和Pollar_rho素因子分解

注意不知道为什么直接使用rand函数就行了，使用了srand函数就RE了

code：

#include<stdio.h>
#include<string.h>
#include<stdlib.h>
#include<time.h>
#include<iostream>
#include<cmath>
#include<algorithm>
using namespace std;

const int S=20;

long long mult_mod(long long a,long long b,long long c)
{
    a%=c;
    b%=c;
    long long ret=0;
    while(b)
    {
        if(b&1){ret+=a;ret%=c;}
        a<<=1;
        if(a>=c)a%=c;
        b>>=1;
    }
    return ret;
}

long long pow_mod(long long x,long long n,long long mod)
{
    if(n==1)return x%mod;
    x%=mod;
    long long tmp=x;
    long long ret=1;
    while(n)
    {
        if(n&1) ret=mult_mod(ret,tmp,mod);
        tmp=mult_mod(tmp,tmp,mod);
        n>>=1;
    }
    return ret;
}

bool check(long long a,long long n,long long x,long long t)
{
    long long ret=pow_mod(a,x,n);
    long long last=ret;
    for(int i=1;i<=t;i++)
    {
        ret=mult_mod(ret,ret,n);
        if(ret==1&&last!=1&&last!=n-1) return true;
        last=ret;
    }
    if(ret!=1) return true;
    return false;
}

bool Miller_Rabin(long long n)
{
    if(n<2)return false;
    if(n==2)return true;
    if((n&1)==0) return false;
    long long x=n-1;
    long long t=0;
    while((x&1)==0){x>>=1;t++;}
    for(int i=0;i<S;i++)
    {
        long long a=rand()%(n-1)+1;
        if(check(a,n,x,t))
            return false;
    }
    return true;
}

long long factor[100];
int tol;

long long gcd(long long a,long long b)
{
    if(a==0)return 1;
    if(a<0) return gcd(-a,b);
    while(b)
    {
        long long t=a%b;
        a=b;
        b=t;
    }
    return a;
}

long long Pollard_rho(long long x,long long c)
{
    long long i=1,k=2;
    long long x0=rand()%x;
    long long y=x0;
    while(1)
    {
        i++;
        x0=(mult_mod(x0,x0,x)+c)%x;
        long long d=gcd(y-x0,x);
        if(d!=1&&d!=x) return d;
        if(y==x0) return x;
        if(i==k){y=x0;k+=k;}
    }
}

void findfac(long long n)
{
    if(Miller_Rabin(n))
    {
        factor[tol++]=n;
        return;
    }
    long long p=n;
    while(p>=n)p=Pollard_rho(p,rand()%(n-1)+1);
    findfac(p);
    findfac(n/p);
}
//以上全模板
long long r[100];
int num;
long long k;

void dfs(long long now,int x,long long n)
{
    if(now>sqrt(n)) return;//尽可能接近根号n
    k=max(k,now);
    for(int i=x;i<=num;i++) dfs(now*r[i],i+1,n);
}

int main()
{
    long long gcd,lcm,n;
    while(scanf("%lld%lld",&gcd,&lcm)!=EOF)
    {
        if(gcd==lcm)
        {
            printf("%lld %lld\n",gcd,lcm);
            continue;
        }
        tol=0;
        n=lcm/gcd;
        findfac(n); sort(factor,factor+tol);
        num=0; for(int i=0;i<=50;i++) r[i]=1;
        r[num]=factor[0];
        for(int i=1;i<tol;i++)
        {
            if(factor[i]==factor[i-1]) r[num]=r[num]*factor[i];//素因子的幂的形式
            else
            {
                num++;
                r[num]=factor[i];
            }
        }
        k=1; dfs(1,0,n);
        printf("%lld %lld\n",gcd*k,gcd*(n/k));
    }
    return 0;
}