Pollard_rho

#include <bits/stdc++.h>
#define LL long long 
using namespace std;

  LL lis[23333],liscnt,T,n;

  LL rnd(LL mo){
      return(((rand()*rand()+rand())%mo+mo)%mo);
  }
  
  LL gcd(LL x,LL y){
      if (x%y==0) return(y);else return(gcd(y,x%y));
  }

  inline LL mult(LL a, LL b, LL p){
    LL tmp=(a*b-(LL)((long double)a/p*b+1e-8)*p);
    while (tmp<0) tmp+=p;
    while (tmp>p) tmp-=p;
    return(tmp);
  }

  LL qpow(LL bas,LL powe,LL mo){
      LL ret=1;
      for (;powe;bas=mult(bas,bas,mo)){
        if (powe&1) ret=mult(ret,bas,mo);
      powe>>=1;    
    }
    return(ret);
  }

  int check(LL p,LL n){
      if (n==2) return(1);
      if (n>2&&n%2==0) return(0);
      LL bas=n-1,tim=0;while (bas%2==0) bas/=2,tim++;
      LL y=qpow(p,bas,n),x=qpow(p,bas,n);
      for (int i=1;i<=tim;i++){
        x=mult(x,x,n);
        if (x==1&&y!=1&&y!=n-1) return(0);
        y=x;
    }
    return(x==1);
  }

  int isprime(LL n){
      for (int i=0;i<30;i++)
        if (!check(rnd(n-1)+1,n)) return(0);
      return(1);
  }

  LL pollard_rho(LL n){
      LL c=rnd(n-1)+1;
      LL stp=2,lim=2,x=rnd(n),y=x,d=1;
      while (d==1){
        x=(mult(x,x,n)+c)%n;
      d=gcd(abs(y-x),n);
      if (stp==lim) y=x,lim*=2;
      stp++;    
    }
    return(d);
  }

  void solve(LL n){
      if (isprime(n)){
        lis[++liscnt]=n;
      return;    
    }
    
    LL t=pollard_rho(n);
    while (t==n)
      t=pollard_rho(n); 
    solve(t);solve(n/t);
  }

  int main(){      
      scanf("%lld",&T);
      while (T--){
        scanf("%lld",&n);
        liscnt=0;
      solve(n);
      sort(lis+1,lis+liscnt+1);    
      if (liscnt==1) printf("Prime\n");
        else printf("%lld\n",lis[liscnt]);
    }
  }
 

2018.5.4修改了mult函数。数据范围过大时导致精度误差。

posted @ 2017-11-06 11:35  z1j1n1  阅读(210)  评论(0编辑  收藏  举报