题目地址:http://poj.org/problem?id=2429
解析:题目的关键就是把lcm/gcd的大数分解成质因子之积,用到了pollard—rho寻找因子和Miller_rabin素数检测。(可以用容斥原理求出所有的情况,然后选出来a+b最小的情况。)
#include <iostream>
#include <cmath>
#include <string>
#include <cstring>
#include <cstdlib>
#include <ctime>
#include <algorithm>
#include <cstdio>
using namespace std;
typedef long long ll;
#define MAX(a,b) a>b?a:b
#define MIN(a,b) a>b?b:a
#define INF ((ll)1<<59)
#define N 501
#define C 201 //pollard_rho里面加的常数C
#define Times 11
ll mini;
int ct;//记录因子的个数
ll r[N]; //存储N的质因子
ll gcd(ll a,ll b){return b?gcd(b,a%b):a;} //欧几里得
ll random(ll n) //生成随机数
{
return (ll)((double)rand()/RAND_MAX*n+0.5);
}
ll mul(ll a,ll b,ll MOD) //a*b%MOD
{
ll ans=0;
while(b>0)
{
if(b&1)
{
ans=(ans+a)%MOD;
b--;
}
b=b>>1;
a=(a<<1)%MOD;
}
return ans;
}
ll pow_mod_1(ll x,ll y,ll MOD) //快速幂mod(非递归实现)
{
ll ans=1;
x%=MOD;
while(y>0) //x,y的值可能会很大,因此选择这个可能会较好。
{
if(y&1)
{
ans=mul(ans,x,MOD);
y--;
}
y=y>>1;
x=mul(x,x,MOD);
}
return ans;
}
bool Witness(ll a,ll n) //Miller_rabin,素数测试
{
ll m=n-1;
int j=0;
while(!(m&1))
{
j++;
m=m>>1;
}
ll x=pow_mod_1(a,m,n);
if(x==1||x==n-1) return false; //
while(j--)
{
x=x*x%n;
if(x==n-1) return false; //若x==n-1,则x*x也为n-1
}
return true;
}
bool Miller_rabin(ll n)
{
if(n<2) return false;
if(n==2) return true;
if(!(n&1)) return false;
for(int i=1;i<=Times;i++) //测试了Times次,错误概率降到很低,几乎为0
{
ll a=random(n-2)+1;
if(Witness(a,n)) return false;
}
return true;
}
ll pollard_rho(ll n,int c) //寻找N的因子。
{
ll x,y,d,i=1,k=2;
x=random(n-1)+1;
y=x;
while(1)
{
i++;
x=(mul(x,x,n)+c)%n;
d=gcd(y-x,n);
if(d>1&&d<n) return d;
if(y==x) return n;
if(i==k)
{
y=x;
k=k<<1;
}
}
}
void find(ll n,int k) //找到因子,为素数则存储。
{
if(n==1) return ;
if(Miller_rabin(n))
{
r[ct++] =n;
if(n<mini) mini=n;
return ;
}
ll p=n;
while(p>=n) p=pollard_rho(p,k--); //p为n的因子
find(p,k);
find(n/p,k);
}
int temp[N][2];
int q; //记录含有不同质因子的个数,为容斥原理做铺垫。
void settle() //整理得到的质因子到temp中,temp[N][0]记录的是质因子,temp[N][1]记录的是该质因子的个数
{
sort(r,r+ct);
q=0;
temp[q][0]=r[0];
temp[q][1]=1;
for(ll i=1;i<ct;i++)
{
if(r[i]==r[i-1]) temp[q][1]++;
else
{
q++;
temp[q][0]=r[i];
temp[q][1]=1;
}
}
q++; //符合习惯
}
ll Pow(ll x,int y) //x^y
{
if(y==1) return x;
ll ans=Pow(x,y/2);
return y%2?(ans*ans)*x:(ans*ans);
}
void swap(ll &a,ll &b) //交换
{
ll t=a;
a=b;
b=t;
}
void resolve(ll lcm_gcd,ll gcdd)
{
ll mini=INF;
ll ed=1<<q;
ll temp_p;
ll ans_a,ans_b;
ll mid_l;
ll mid_r;
for(ll i=1;i<ed-1;i++) //容斥,见有人是dfs实现的,好像复杂度差不多。
{
temp_p=i;
mid_l=1;
int k=0;
while(temp_p)
{
if(temp_p&1) mid_l*=Pow(temp[k][0],temp[k][1]);
k++;
temp_p=temp_p>>1;
}
ll mid_r=lcm_gcd/mid_l+mid_l; //(a+b)两者之和
if(mid_r<mini)
{
mini=mid_r;
ans_a=mid_l;
ans_b=mid_r-mid_l;
}
}
if(ans_b<ans_a) swap(ans_b,ans_a);//调整下。
cout<<ans_a*gcdd<<" "<<ans_b*gcdd<<endl;
}
int main()
{
int i,j;
ll lcmm,gcdd;
ll lcm_gcd;
while(cin>>gcdd>>lcmm)
{
if(gcdd==lcmm)
{
cout<<gcdd<<" "<<lcmm<<endl;
continue;
}
ct=0;
lcm_gcd=lcmm/gcdd;
find(lcm_gcd,C);
settle();
resolve(lcm_gcd,gcdd);
}
return 0;
}