1. 程式人生 > >POJ——2429(數論之大整數分解)

POJ——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;
}