• bzoj3667: Rabin-Miller算法


    Description

     

    Input

    第一行:CAS,代表数据组数(不大于350),以下CAS行,每行一个数字,保证在64位长整形范围内,并且没有负数。你需要对于每个数字:第一,检验是否是质数,是质数就输出Prime
    第二,如果不是质数,输出它最大的质因子是哪个。

    Output

    第一行CAS(CAS<=350,代表测试数据的组数)
    以下CAS行:每行一个数字,保证是在64位长整形范围内的正数。
    对于每组测试数据:输出Prime,代表它是质数,或者输出它最大的质因子,代表它是和数

    Sample Input

    6
    2
    13
    134
    8897
    1234567654321
    1000000000000

    Sample Output

    Prime
    Prime
    67
    41
    4649
    5

    HINT

    数据范围:

    保证cas<=350,保证所有数字均在64位长整形范围内。

    Source

    题解:

    这是miller rabin和pollard rho的裸题,具体见:http://blog.csdn.net/thy_asdf/article/details/51347390

    code:

     1 #include<cstdio>
     2 #include<iostream>
     3 #include<cmath>
     4 #include<cstring>
     5 #include<algorithm>
     6 #define lowbit(x) ((x)&(-(x)))
     7 using namespace std;
     8 typedef long long int64;
     9 char ch;
    10 bool ok;
    11 void read(int &x){
    12     ok=0;
    13     for (ch=getchar();!isdigit(ch);ch=getchar()) if (ch=='-') ok=1;
    14     for (x=0;isdigit(ch);x=x*10+ch-'0',ch=getchar());
    15     if (ok) x=-x;
    16 }
    17 void read(int64 &x){
    18     ok=0;
    19     for (ch=getchar();!isdigit(ch);ch=getchar()) if (ch=='-') ok=1;
    20     for (x=0;isdigit(ch);x=x*10+ch-'0',ch=getchar());
    21     if (ok) x=-x;
    22 }
    23 int cases;
    24 int64 n;
    25 const int prime[10]={2,3,5,7,11,13,17,19,23,29};
    26 int64 mul(int64 a,int64 b,int64 mod){
    27     int64 d=((long double)a/mod*b+1E-8);
    28     int64 res=a*b-d*mod;
    29     if (res<0) res+=mod;
    30     return res;
    31     /*int64 t;
    32     for (t=0;b;b>>=1,a<<=1,a%=mod) if (b&1) t=(t+a)%mod;
    33     return t;*/
    34 }
    35 int64 ksm(int64 a,int64 b,int64 mod){
    36     int64 t;
    37     for (t=1;b;b>>=1,a=mul(a,a,mod)) if (b&1) t=mul(t,a,mod);
    38     return t;    
    39 }
    40 bool miller_rabin(int64 n){
    41     for (int i=0;i<10;i++) if (n==prime[i]) return true;
    42     if (!(n&1)) return false;
    43     int64 bit=lowbit(n-1),s=0,d;
    44     while (bit!=1) s++,bit>>=1;
    45     d=(n-1)>>s;
    46     for (int i=0;i<10;i++){
    47         int64 x=ksm(prime[i],d,n);
    48         for (int j=1;j<=s;j++){
    49             int64 xx=mul(x,x,n);
    50             if (xx==1&&x!=n-1&&x!=1) return false;
    51             x=xx;
    52         }
    53         if (x!=1) return false;
    54     }
    55     return true;
    56 }
    57 int cnt;
    58 int64 list[50];
    59 int64 random(int64 lim){return ((1LL*rand()<<31)+rand())%lim;}
    60 int64 f(int64 x,int64 mod,int64 c){return (mul(x,x,mod)+c+mod)%mod;}
    61 int64 pollard_rho(int64 n,int64 c){
    62     int64 x,y,d=1; x=random(n),y=f(x,n,c);
    63     while (d==1){
    64         d=__gcd(abs(x-y),n);
    65         x=f(x,n,c),y=f(f(y,n,c),n,c);
    66     }
    67     return d;
    68 }
    69 void work(int64 n){
    70     if (miller_rabin(n)){list[++cnt]=n;return;}
    71     int64 d=pollard_rho(n,random(n-1));
    72     while (d==n||d==1) d=pollard_rho(n,random(n));
    73     work(d),work(n/d);
    74 }
    75 void decompose(int64 n){
    76     cnt=0,work(n),sort(list+1,list+cnt+1);
    77     printf("%lld
    ",list[cnt]);
    78 }
    79 int main(){
    80     srand(19990617);
    81     for (read(cases);cases;cases--){
    82         read(n);
    83         if (miller_rabin(n)) puts("Prime");
    84         else decompose(n);
    85     }
    86     return 0;
    87 }
  • 相关阅读:
    17.10.13
    17.10.12
    17.10.11
    17.10.10
    17.10.05
    17.10.04
    17.10.03
    17.10.02
    17.10.01
    17.9.29
  • 原文地址:https://www.cnblogs.com/chenyushuo/p/5472155.html
Copyright © 2020-2023  润新知