2013-team4/code/rabin-miller

从 Trac 迁移的文章

这是从旧校内 Wiki 迁移的文章，可能存在一些样式问题，您可以向 memset0 反馈。

原文章内容如下：

{{{
#include <ctime>
#include <cstdlib>
using namespace std;
typedef long long LL;

inline LL add(LL a, LL b, LL mod){
    LL res=a+b;
    res%=mod;
    return res;
}

//复杂度log2(min(a,b)),避免溢出,用于替换原来的mul
inline LL mul(LL a, LL b, LL mod){
    a=a%mod; b=b%mod;
    LL res=0; if(a<b)swap(a,b);
    for(;b;b>>=1,a=add(a,a,mod))
        if(b&1)res=add(res,a,mod);
    return res;
}

//普通mul
/*
inline LL mul(LL a, LL b, LL mod){
    return (long long)a*b%mod;
}
*/

LL modExp(LL a, LL n, LL mod){
    LL res=1;
    a=a%mod;
    for(;n;n>>=1,a=mul(a,a,mod))
        if(n&1)res=mul(res,a,mod);
    return res;
}

bool mr_test(LL a, LL n, LL d, int k){
    LL r=modExp(a,d,n);
    if(r==1)return true;
    while(k--){
        if(r==n-1)return true;
        r=mul(r,r,n);
    }
    return false;
}

bool easy_prime(LL n){
    switch(n){
        case 2: case 3: case 5: case 7: case 11: case 13: case 17: return true;
        default: return false;
    }
}

bool easy_comp(LL n){
    if(n%2==0||n%3==0||n%5==0||n%7==0||n%11==0||n%13==0||n%17==0)return true;
    return false;
}

//!!if n < 1,373,653, it is enough to test a = 2 and 3;         //verified
//!!if n < 9,080,191, it is enough to test a = 31 and 73;       //verified
//!!if n < 4,759,123,141, it is enough to test a = 2, 7, and 61;//mostly verified
//!!if n < 1,122,004,669,633, it is enough to test a = 2, 13, 23, and 1662803;
//!!if n < 2,152,302,898,747, it is enough to test a = 2, 3, 5, 7, and 11;
//!!if n < 3,474,749,660,383, it is enough to test a = 2, 3, 5, 7, 11, and 13;
//!!if n < 341,550,071,728,321, it is enough to test a = 2, 3, 5, 7, 11, 13, and 17.

const int prime_num=4, p[]={2,13,23,1662803};
bool miller_rabin_fixed(LL n){
    if(n==0||n==1)return false;
    if(easy_prime(n))return true;
    if(easy_comp(n))return false;
    int k=0; LL d=n-1;
    while(~d&1)d>>=1,k++;
    for(int i=0;i<prime_num;i++){
        if(p[i]==n)return true;
        if(!mr_test(p[i]%n,n,d,k))return false;
    }
    return true;
}

bool miller_rabin_random(LL n, int time){
    if(n==0||n==1)return false;
    if(easy_prime(n))return true;
    if(easy_comp(n))return false;
    int k=0; LL d=n-1;
    while(~d&1)d>>=1,k++;
    for(int i=0;i<time;i++){
        LL p=rand()%(n-1)+1;
        if(!mr_test(p%n,n,d,k))return false;
    }
    return true; //probability=1-(1/4)^time
}
}}}

#include <ctime>
#include <cstdlib>
using namespace std;
typedef long long LL;
inline LL add(LL a, LL b, LL mod){
    LL res=a+b;
    res%=mod;
    return res;
}
//复杂度log2(min(a,b)),避免溢出,用于替换原来的mul
inline LL mul(LL a, LL b, LL mod){
    a=a%mod; b=b%mod;
    LL res=0; if(a<b)swap(a,b);
    for(;b;b>>=1,a=add(a,a,mod))
        if(b&1)res=add(res,a,mod);
    return res;
}
//普通mul
/*
inline LL mul(LL a, LL b, LL mod){
    return (long long)a*b%mod;
}
*/
LL modExp(LL a, LL n, LL mod){
    LL res=1;
    a=a%mod;
    for(;n;n>>=1,a=mul(a,a,mod))
        if(n&1)res=mul(res,a,mod);
    return res;
}
bool mr_test(LL a, LL n, LL d, int k){
    LL r=modExp(a,d,n);
    if(r==1)return true;
    while(k--){
        if(r==n-1)return true;
        r=mul(r,r,n);
    }
    return false;
}
bool easy_prime(LL n){
    switch(n){
        case 2: case 3: case 5: case 7: case 11: case 13: case 17: return true;
        default: return false;
    }
}
bool easy_comp(LL n){
    if(n%2==0||n%3==0||n%5==0||n%7==0||n%11==0||n%13==0||n%17==0)return true;
    return false;
}
//!!if n < 1,373,653, it is enough to test a = 2 and 3;         //verified
//!!if n < 9,080,191, it is enough to test a = 31 and 73;       //verified
//!!if n < 4,759,123,141, it is enough to test a = 2, 7, and 61;//mostly verified
//!!if n < 1,122,004,669,633, it is enough to test a = 2, 13, 23, and 1662803;
//!!if n < 2,152,302,898,747, it is enough to test a = 2, 3, 5, 7, and 11;
//!!if n < 3,474,749,660,383, it is enough to test a = 2, 3, 5, 7, 11, and 13;
//!!if n < 341,550,071,728,321, it is enough to test a = 2, 3, 5, 7, 11, 13, and 17.
const int prime_num=4, p[]={2,13,23,1662803};
bool miller_rabin_fixed(LL n){
    if(n==0||n==1)return false;
    if(easy_prime(n))return true;
    if(easy_comp(n))return false;
    int k=0; LL d=n-1;
    while(~d&1)d>>=1,k++;
    for(int i=0;i<prime_num;i++){
        if(p[i]==n)return true;
        if(!mr_test(p[i]%n,n,d,k))return false;
    }
    return true;
}
bool miller_rabin_random(LL n, int time){
    if(n==0||n==1)return false;
    if(easy_prime(n))return true;
    if(easy_comp(n))return false;
    int k=0; LL d=n-1;
    while(~d&1)d>>=1,k++;
    for(int i=0;i<time;i++){
        LL p=rand()%(n-1)+1;
        if(!mr_test(p%n,n,d,k))return false;
    }
    return true; //probability=1-(1/4)^time
}