POJ_2429_GCD & LCM Inverse

PollardRho算法的应用.

输入a,b,处理为b/=a.

先使用Pollard算法求出b的所有的质因子,求出的信息中包含每个质因子在b中的个数.

设结果为x,y,

那么b=(x1^p1)*(x2^p2)*(x3^p3)*...,化简为b=X1*X2*X3*...*Xn,其中Xi=xi^pi.(因为x和y是互质的,所以x和y会独占Xi,这是化简的理由.)

通过打表可以推出n最大为15.

那么,那么位运算一下即可,(O(1<<n))

for(i=0;i<=(1<<n);i++)

{

　　LL tmp=1;

　　for(int j=0;j<n;j++)

　　{if(i&(1<<j)) tmp*=Xj;}

　　//下面是判断处理...求出最有结果即可...

}

最后x*=a,y*=a.

附上代码:

#include <map>
#include <iostream>
#include <fstream>
#include <algorithm>
#include <cmath>
#include <cstring>
#include <climits>
#include <cstdio>
#include <ctime>




using namespace std;

typedef long long  LL;

LL factor[1000];
LL cnt;


inline LL gcd(LL a,LL b)
{
    LL tmp;
    while(b)
    {
        tmp=a%b;
        a=b;
        b=tmp;
    }
    return a;
}

//求a*b%n,(a*b可能大于LLONG_MAX)
inline LL mulMod(LL a,LL b,LL n)
{
    LL res=0,e=a%n;
    while(b)
    {
        if(b&1)
        {
            res+=e;
            if(res>n)
            res-=n;
        }
        b>>=1;
        e<<=1;
        if(e>n)
        e-=n;
    }
    return  res;
}

//幂取模
LL expMod(LL a,LL b,LL n)
{
    LL res=1,e=a%n;
    while(b)
    {
        if(b&1)
        {
            res=mulMod(res,e,n);
        }
        b>>=1;
        e=mulMod(e,e,n);
    }

    return res;
}


bool millerRabin(LL n,LL times)
{
    if(n<2)
    return false;
    if(n==2)
    return true;
    if(!(n&1))
    return false;

    long long l=0,a,x,y;
    long long m=n-1;
    while(!(m&1))
    l++,m>>=1;
    srand(time(NULL));
    while(times--)
    {
        a=rand()%(n-1)+1;
        x=expMod(a,m,n);
        if(x==1)
        continue;
        for(int i=0;i<l;i++)
        {
            y=mulMod(x,x,n);
            if(y==1 && 1!=x && (n-1)!=x)
            return false;
            x=y;
        }
        if(y!=1)
        return false;
    }
    return true;
}

//质因子分解
LL pollardRho(long long n,long long c)
{
    LL x,y,d,i=1,k=2;
    srand(time(NULL));

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

//寻找所有的质因子
void findFactor(LL n,LL k)
{
    if(n==1)
    return;

    if(millerRabin(n,10))
    {
        factor[cnt++]=n;
        return;
    }
    LL p=n;
    while(p>=n)
    {
        p=pollardRho(p,k--);
    }
    findFactor(p,k);
    findFactor(n/p,k);
}

int main()
{

    long long a,b;
    long long nMin=LLONG_MAX;
    LL x,y;
    map<LL,LL> hash;
    while(scanf("%I64d %I64d",&a,&b)!=EOF)
    {
        if(a==b)
        {
            printf("%I64d %I64d\n",a,b);
            continue;
        }
        hash.clear();
        b/=a;
        factor[0]=1;
        cnt=0;
        findFactor(b,200);
        if(0==cnt)
        {
            printf("%I64d %I64d\n",a,a*b);
            continue;
        }
        for(int i=0;i<cnt;i++)
        {
            if(hash.find(factor[i])==hash.end())
            {
                hash[factor[i]]=factor[i];
            }
            else
            {
                hash[factor[i]]*=factor[i];
            }
        }
        cnt=0;
        for(map<LL,LL>::iterator itr=hash.begin();itr!=hash.end();itr++)
        {
            factor[cnt++]=itr->second;
        }
        LL tmp;
        nMin=LLONG_MAX;
        for(int i=0;i<(1<<cnt);i++)
        {
            tmp=1;
            for(int j=0;j<cnt;j++)
            {
                if(i&(1<<j))
                tmp*=factor[j];
            }
            //printf("tmp %I64d\n",tmp);
            if(nMin>tmp+b/tmp)
            {
                nMin=tmp+b/tmp;
                x=min(tmp,b/tmp);
                y=max(tmp,b/tmp);
            }
        }
        printf("%I64d %I64d\n",x*a,y*a);





    }


    return 0;
}