poj 2429

令x=lcm/gcd
将x分解质因数->Pollard rho方法
枚举将x分解为两个互素的数相乘

代码:

#include<ctime>
#include<iostream>
#include<fstream>

using namespace std;
long long factor[1000],fac_top = -1;

//计算两个数的gcd
long long gcd(long long a,long long b)
{
	if(a==0)
		return b;
	long long c;
	while(b!=0)
	{
		c=b;
		b=a%b;
		a=c;
	}
	return a;
}
//ret = (a*b)%n (n<2^62)
long long muti_mod(long long a,long long b,long long n)
{
    long long  exp = a%n, res = 0;
    while(b)
	{
        if(b&1)
		{
            res += exp;
            if(res>n) res -= n;
        }
        exp <<= 1;
        if(exp>n) 
			exp -= n;
       
        b>>=1;
    }
    return res;
}
// ret = (a^b)%n
long long mod_exp(long long a,long long p,long long m)
{
    long long exp=a%m, res=1; //  
    while(p>1)
    {
        if(p&1)//
            res=muti_mod(res,exp,m);
        exp = muti_mod(exp,exp,m);
        p>>=1;
    }
    return muti_mod(res,exp,m);
}
//miller-rabin法测试素数, time 测试次数
bool miller_rabin(long long n, long long times)
{
    if(n==2)return 1;
    if(n<2||!(n&1))return 0;
   
    long long a, u=n-1, x, y;
    int t=0;
    while(u%2==0){
        t++;
        u/=2;
    }
    srand(time(0));
    for(int i=0;i<times;i++)
	{
        a = rand() % (n-1) + 1;
        x = mod_exp(a, u, n);
        for(int j=0;j<t;j++)
		{
            y = muti_mod(x, x, n);
            if ( y == 1 && x != 1 && x != n-1 )
                return false; //must not
            x = y;
        }
        if( y!=1) return false;
    }
    return true;
}

long long pollard_rho(long long n,int c)
{//找出一个因子
    long long x,y,d,i = 1,k = 2;
    srand(time(0));
    x = rand()%(n-1)+1;
    y = x;
    while(true) {
        i++;
         x = (muti_mod(x,x,n) + c) % n;
        d = gcd(y-x, n);
        if(1 < d && d < n)return d;
         if( y == x)       return n;
         if(i == k) {
            y = x;
             k <<= 1;
         }
    }
}
void findFactor(long long n,int k)
{//二分找出所有质因子,存入factor
    if(n==1)return;
    if(miller_rabin(n, 10))
	{
		factor[++fac_top] = n;
		return;
	}
	long long p = n;
	while(p >= n)
		p = pollard_rho(p,k--);//k值变化,防止死循环
	findFactor(p,k);
	findFactor(n/p,k);
}


long long a[1000];
long long m,minx,ans;


int cmp(const void *a,const void *b){
	return *(long long *)a-*(long long *)b;
}

void dfs(long long s,long long num,long long t){
	
	if(s==m+1)
	{
		if(minx==-1||(num+t/num<minx))
		{
			minx=num+t/num;
			ans=num;
		}
		return;
	}

	dfs(s+1,num*a[s],t);
	dfs(s+1,num,t);

}

void solve(){
	qsort(factor,fac_top+1,sizeof(long long),cmp);
	a[0]=factor[0];
	long long i;
	for(i=0;i<fac_top;i++)
		if(factor[i]==factor[i+1])
		{	
			a[m]*=factor[i+1];
		}
		else
		{
			m++;
			a[m]=factor[i+1];
		}
}

int main()
{
    long long s,t,n;
 //   ifstream cin("in.txt");

	while(cin>>s>>t)
	{
		n=t/s;
		if(s==t)
		{
			cout<<s<<' '<<t<<endl;
			continue;
		}
		minx=-1;
		m=0;
        fac_top = -1;
        findFactor(n,107);
		solve();
		dfs(0,1,n);
		if(ans>n/ans) ans=n/ans;
		cout<<s*ans<<' '<<s*(n/ans)<<endl;
    }
    return 0;
}

转载于:https://www.cnblogs.com/zhaozhe/archive/2011/04/12/2013979.html

  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值