1. 程式人生 > >poj-2429 Miller_Rabin強偽素數測試 +pollard_rho質因數分解

poj-2429 Miller_Rabin強偽素數測試 +pollard_rho質因數分解

題目:

GCD & LCM Inverse
Time Limit: 2000MSMemory Limit: 65536K
Total Submissions: 17978Accepted: 3331

Description

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

題意:輸入兩個數的最大公約數和最小公倍數,求這兩個數。如果有多組,輸出和最小的一對。

思路:

設這兩個數分別為:x , y。

那麼  必然有 gcd(x,y)*k1 = x ; gcd(x,y)*k2 = y。

由 :x*y / gcd(x,y)= lcm(x,y)    變形得=>>     x*y = lcm*gcd;

則:(gcd*k1) * (gcd*k2) = lcm*gcd.

化簡得:k1*k2 = lcm/gcd

若我們找到一組k1,k2,滿足:k1*k2 = lcm/gcd 的關係,那麼k1*gcd就是可能的x,k2*gcd就是可能的y。

我們還應該注意到由於 k1 = x/gcd , k2 = y/gcd 。那麼一定有gcd(k1,k2)=1 。

所以我們可以將lcm/gcd進行質因數分解,之後將質因數分別再湊為兩堆,一堆相乘結果就是x ,令一堆就是y。特別要注意的是,由於gcd(k1,k2)=1 ,那麼為了保證這個性質,分解出的質因數若存在多個,要再乘回去,否則gcd(k1,k2)就會等於那個質因數了。

題設還有最後一個條件:輸出和最小的一對。由於x*y =gcd*lcm , 那麼求和最小的一對也就是求差最小的一對解。而要使,差最小,無非就是使分解出的兩個數都接近sqrt(lcm/gcd),所以我們將所有質因數的組合進行dfs,找到最優解即可。

程式碼:

#include<cstdio>  
#include<cstdlib>  
#include<cmath>  
#include<cstring>  
#include<iostream>  
#include<algorithm>  
#include<queue>  
#include<map>  
#include<stack> 
#include<set>

typedef long long ll;

using namespace std;

const int INF = 0x3f3f3f3f;
const int N = 5500;
const int Times = 10; 

ll ct,cnt;
ll fac[N];
ll f[N];
ll ans;

ll gcd(ll a , ll b){
	return b==0? a:gcd(b,a%b);
}

ll multi(ll a , ll b , ll m){
	ll ans = 0;
	a %= m;
	while(b){
		if(b&1){
			ans = (ans+a)%m;
			b--;
		}
		b >>= 1;
		a = (a+a)%m;
	}
	return ans;
}

ll quick_mod(ll a , ll b , ll m){
	ll ans = 1;
	a %= m;
	while(b){
		if(b&1){
			ans = multi(ans , a , m);
			b--;
		}
		b>>=1;
		a = multi(a,a,m);
	}
	return ans;
}

bool Miller_Rabin(ll n){
	if(n==2)	return true;
	if(n<2 || !(n&1)) return false;
	ll m = n-1;
	int k = 0;
	while((m&1)==0){
		k++;
		m>>=1;
	}
	for(int i = 0 ; i<Times ; i++){
		ll a = rand()%(n-1)+1;
		ll x = quick_mod(a,m,n);
		ll y = 0;
		for(int j = 0 ; j<k ; j++){
			y = multi(x,x,n);
			if(y==1 && x!=1 && x!=n-1)	return false;
			x = y;
		}
		if(y!=1)	return false;
	}
	return true;
}

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

void find(ll n , int c){
	if(n == 1)	return;
	if(Miller_Rabin(n)){
		fac[ct++] = n;
		return;
	}
	ll p = n;
	ll k = c;
	while(p >= n) p = pollard_rho(p , c--);
	find(p , k);
	find(n/p , k);
}

void findx(ll i , ll x , ll q){
	if(i==cnt)	return;
	if(x>ans && x<=q) ans = x;
	findx(i+1 , x , q);
	
	x *= f[i];
	if(x>ans && x<=q) ans = x;
	findx(i+1 , x , q);
}

int main(){
	ll a,b;
	while(scanf("%lld %lld",&a,&b)!=EOF){
		if(a == b){
			printf("%lld %lld\n",a,b);
			continue;
		}
		
		ct = 0;
		ll n = b/a;
		find(n,120);
		sort(fac,fac+ct);
		
		f[0] = fac[0]; cnt = 0;
		for(int i = 1 ; i<ct ; i++){
			if(fac[i] == fac[i-1]){
				f[cnt] = f[cnt] * fac[i];
			}
			else{
				cnt++;
				f[cnt] = fac[i];
			}
		}
		cnt++;
		
		ll goal = (ll)sqrt(n*1.0);
		ans = 1;
		findx(0,1,goal);
		printf("%lld %lld\n" , ans*a , n/ans*a);
	}
	return 0;
}