2180198 Sample Output3 183 22 AuthorHIT Source2012 Multi-University Training Contest 5 Recommendzhuyuanchen520 | We have carefully selected several similar problems for you: 4348 4347 4346 4345 4343题解:Miller Rabbin + Pollard rho
Miller Rabbin 可以在O(s(logn)^3)的时间复杂度内判断一个数是否为素数,有2^(-s)的概率出错
Pollard rho 是基于Miller Rabbin的一种快速分解质因数的做法。
该算法的大致流程是先判断当前数是否为素数(用Miller Rabbin),如果是直接记录下该质数,直接返回。如果不是就试图找到一个因子(可以不是质因子),然后对于当前因子p,和n/p分别递归寻找质因子。
对于质因数的寻找,我们采用一种随机化的算法。我们假设要找到的质因子为p,先随机取一个x,然后用x构造y,使p=gcd(x-y,n)如果p不等于1那么就找到了一个质因子,如果找不到我们就不断的调整x,使x=x*x+c(c一般可以随机),直到x==y出现循环,则选取失败。重新选取x,重复上述过程。
#include<iostream>#include<cstdio>#include<cstring>#include<algorithm>#include<cmath>#define LL long long #define N 10000003using namespace std;LL n,mx,cnt,num[N],prime[N],c[N];LL mul(LL a,LL b,LL p){ LL ans=0; LL base=a%p; while (b) { if (b&1) ans=(ans+base)%p; b>>=1; base=(base+base)%p; } return ans;}LL quickpow(LL num,LL x,LL p){ LL ans=1; LL base=num%p; while (x) { if (x&1) ans=mul(ans,base,p); x>>=1; base=mul(base,base,p); } return ans;}bool miller_rabbin(LL n){ if (n==2) return true; if (!(n&1)) return false; LL t=0,a,x,y,u=n-1; while (!(u&1)) t++,u>>=1; for (int i=0;i<=10;i++) { a=rand()*rand()%(n-1)+1; x=quickpow(a,u,n); for (int j=0;j<t;j++) { y=mul(x,x,n); if (y==1&&x!=1&&x!=n-1) return false; x=y; } if (x!=1) return false; } return true;}LL gcd(LL x,LL y){ LL r; while (y) { r=x%y; x=y; y=r; } return x;}LL pollard_rho(LL n,LL c){ LL i=1,k=2; LL x=rand()%(n-1)+1,y=x,p=1; while (p==1) { i++; x=(mul(x,x,n)+c)%n; p=gcd((y-x+n)%n,n); if (y==x) return n; if (i==k) y=x,k<<=1; } return p;}void solve(LL n){ if (n==1) return; if (miller_rabbin(n)) { num[++cnt]=n; return; } LL p=n; while (p==n) p=pollard_rho(p,rand()%(n-1)+1); solve(p); solve(n/p);}int main(){ freopen("a.in","r",stdin); freopen("my.out","w",stdout); srand(2000001001); int T; scanf("%d",&T); while (T--) { scanf("%I64d",&n); //cout<<n<<endl; mx=0; cnt=0; solve(n); int k=0; sort(num+1,num+cnt+1); prime[++k]=num[1]; c[k]=1; for (int i=2;i<=cnt;i++) if (num[i]==num[i-1]) prime[k]*=num[i]; else prime[++k]=num[i]; printf("%d ",k); LL sum=0; for (int i=1;i<=k;i++) sum+=prime[i]; if (k==1) sum/=num[1]; printf("%I64d/n",sum); }}
新闻热点
疑难解答