您的位置:首页 > 其它

Miller_Rabin和Pollard_rho素数模板——POJ-2429的题解

2018-03-14 16:38 483 查看
网络上有很多讲解,不再赘述。此处记载此题题解以供复习之用。#include<cstdio>
#include<iostream>
#include<cstdlib>
#include<cmath>
#include<set>
#include<algorithm>
using namespace std;
typedef long long ll;
ll temp,pr[101],f[101],t,ans;
bool book[2000001];
ll random(ll n)//随机数
{
return ((double)rand()/RAND_MAX*n+0.5);
}
ll q_mul(ll a,ll b,ll mod)//快速乘
{
ll ans=0;
while(b)
{
if(b&1)
ans=(ans+a)%mod;
b>>=1;
a=(a+a)%mod;
}
return ans;
}
ll q_pow(ll a,ll b,ll mod)//快速幂
{
ll res=1;
while(b)
{
if(b&1)
res=q_mul(res,a,mod);
b>>=1;
a=q_mul(a,a,mod);
}
return res;
}
bool go(ll a,ll n)//Miller-rabin的核心算法
{
ll r=n-1,j=0;
while(r&1==0)
{
r>>=1;
j++;
}
ll x=q_pow(a,r,n);
if(x==1||x==n-1)
return true;
while(j--)
{
x=q_mul(x,x,n);
if(x==n-1)
return true;
}
return false;
}
bool miller_rabin(ll n)
{
if(n==2)
return true;
if(n<2||n&1==0)
return false;
for(int i=1;i<=20;i++)
{
ll a=random(n-2)+1;
if(!go(a,n))
return false;
}
return true;
}
ll gcd(ll a,ll b)
{
return b==0?a:gcd(b,a%b);
}
ll pollard_rho(ll p,ll c)
{
ll i=1,k=2,x=random(p-1)+1,y,d;
y=x;
while(1)
{
i++;
x=(q_mul(x,x,p)+c)%p;
d=gcd(y-x,p);
if(d>1&&d<p)//找到
return d;
if(y==x)//循环了没找到
return p;
if(i==k)//优化功能
y=x,k<<=1;
}
}
void insert(ll n)//pr数组记录约数值,f数组记录个数
{
int i;
for(i=0;i<t;i++)
if(pr[i]==n)
break;
if(i<t)
f[i]++;
else
pr[t]=n,f[t++]=1;
}
void find(ll n,ll c)//将辨别素数和查找合数约数的两个算法放在一个函数里了
{
if(n==1)
return;
//这里大的数字用Miller-rabin查的,小的素数用筛子
if(n>=2000000)
{
if(miller_rabin(n))
{
insert(n);
return;
}
}
else
{
if(!book
)
{
insert(n);
return;
}
}

ll p=n;
while(p>=n)//直到随机到他的因数
p=pollard_rho(p,c--);
find(p,c);
find(n/p,c);
}
void search(ll i,ll x,ll q)//查找最接近sqrt的那个
{
if(i==t)return;
if(x>ans&&x<=q) ans=x;
search(i+1,x,q);
x*=f[i];
if(x>ans&&x<=q) ans=x;
search(i+1,x,q);
}
int main()
{
ll a,b,n,m;

book[0]=book[1]=true;//筛小素数
for(int i=2;i*i<=2000000;i++)
{
if(!book[i])
for(int j=2*i;j<=2000000;j+=i)
book[j]=true;
}

while(~scanf("%lld%lld",&n,&m))
{
t=0;
ll multiple=m/n;
find(multiple,21373425);//常数c是随意定的
temp=floor(sqrt(multiple)+0.5);
for(int i=0;i<t;i++)
{
ll k=1;
for(int j=1;j<=f[i];j++)
k*=pr[i];
f[i]=k;//把同样的质数乘成一个整体以保证因数互质,且保存
}
ans=1;
search(0,ans,temp);//查找最接近sqrt的那个数
b=n*multiple/ans,a=n*ans;
printf("%lld %lld\n",a,b);
}
return 0;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
相关文章推荐