BZOJ.3667.Rabin-Miller算法(MillerRabin PollardRho)

题目链接

Pollard_Rho:http://blog.csdn.net/thy_asdf/article/details/51347390

#include<cstdio>

#include<cctype>

#include<algorithm>

#define gc() getchar()

const int p[]={2,3,5,7,11,13,17,19};

typedef long long LL;

LL Ans;

inline LL read()

{

	LL now=0,f=1;register char c=gc();

	for(;!isdigit(c);c=gc()) if(c=='-') f=-1;

	for(;isdigit(c);now=now*10+c-'0',c=gc());

	return now*f;

}

inline LL Mult(LL a,LL b,LL p)//O(1)快速乘

{

	LL tmp=a*b-(LL)((long double)a/p*b+1e-8)*p;

	return tmp<0?tmp+p:tmp;

}

LL Fast_Pow(LL n,LL k,LL p)

{

	LL t=1;

	for(;k;k>>=1,n=n*n%p)

		if(k&1) t=t*n%p;

	return t;

}

bool Miller_Rabin(LL n)

{

	if(n==2) return 1;

	if(!(n&1)||n==1) return 0;

	for(int i=0;i<8;++i)

		if(n==p[i]) return 1;

		else if(!(n%p[i])) return 0;

	LL u=n-1,now,las; int t=0;

	while(!(u&1)) u>>=1,++t;

	for(int i=0;i<8;++i)

	{

		now=Fast_Pow(p[i],u,n);

		for(int j=1;j<=t;++j)

		{

			las=now, now=Mult(now,now,n);

			if(now==1&&las!=1&&las!=n-1) return 0;

		}

		if(now!=1) return 0;

	}

	return 1;

}

LL gcd(LL x,LL y)

{

	return y?gcd(y,x%y):x;

}

LL Rho(LL n,LL delta)

{//现要分解n，有两个随机数x,y，若p=gcd(x-y,n)!=1&&p!=n，那么p为n的一个约数...省略

	LL x=rand()%n,y=x,p=1; int k=2;//设定k为此次路径长

	for(int i=1;p==1;++i)

	{

		x=(Mult(x,x,n)+delta)%n;//随机函数f(x)=x*x+d

		p=gcd(std::abs(x-y),n);//多次生成随机数，直至找到p是n的一个因子

		if(i==k) y=x,k<<=1;//达到k次后把y赋值为x。路径每次倍长

	}

	return p;

}

void Find(LL n)

{

	if(n==1) return;

	if(Miller_Rabin(n)) {Ans=std::max(Ans,n);/*fac[++cnt]=n;*/ return;}

	LL t=n;

	while(t==n) t=Rho(n,rand()%(n-1)+1);

	//t=n说明这个随机函数会导致走到n的环上，再换一个重试即可

	Find(t), Find(n/t);

}

int main()

{

#ifndef ONLINE_JUDGE

	freopen("3667.in","r",stdin);

#endif

	int t=read();LL n;

	while(t--)

		n=read(),Ans=0,Find(n),Ans==n?puts("Prime"):printf("%lld\n",Ans);

	return 0;

}

巴特西

BZOJ.3667.Rabin-Miller算法(MillerRabin PollardRho)

最新文章

热门文章