快速素数测试(Miller-Robin)

typedef long long ll;
typedef __int128 lll;
const int S = 20;//随机算法判定次数,S越大,判错概率越小,一般5~20
inline ll fp(ll b, ll p, ll mod) {
    ll ans = 1;
    while (p) {
        if (p & 1)
            ans = (lll)ans * b % mod;
        p >>= 1;
        b = (lll)b * b % mod;
    }
    return ans;
}
bool check(ll a, ll n, ll x, ll t) {
    ll ret = fp(a, x, n);
    ll last = ret;
    for (int i = 1; i <= t; i++) {
        ret = (lll)ret * ret % n;
        if (ret == 1 && last != 1 && last != n - 1)
            return 1;//合数
        last = ret;
    }
    if (ret != 1)
        return 1;
    return 0;
}
bool Miller_Rabin(ll n) {
    if (n < 2) return 0;
    if (n == 2) return 1;
    if ((n & 1) == 0) return 0;//偶数
    ll x = n - 1;
    ll t = 0;
    while ((x & 1) == 0) {
        x >>= 1;
        t++;
    }
    for (int i = 0; i < S; i++) {
        ll a = rand() % (n - 1) + 1;
        if (check(a, n, x, t))
            return 0;//合数
    }
    return 1;
}
/*
只能测试3e9以内的数,如果要测试long long范围内的数,需要将快速幂中的乘法换成快速乘
本地测试估计:1e9内,时间约为O(100),1e18内,时间约为O(15000)
*/
inline ll fm(ll x,ll y,ll mod){
    ll ans=0;
    while(y>0){
        if(y&1)ans=(ans+x)%mod;
        y>>=1;
        x=(x+x)%mod;
    }
    return ans;
}
inline ll fp(ll b,ll p,ll mod){
    ll ans=1;
    while(p){
        if(p&1)ans=fm(ans,b,mod);
        p>>=1;
        b=b*b%mod;//在判断ll范围内的数时,要用下面这个代替
        //b=fm(b,b,mod);
    }
    return ans;
}
int strong_pseudo_primetest(ll n, int base) {
    ll n2 = n - 1, res;
    int s = 0;
    while (n2 % 2 == 0) n2 >>= 1, s++;
    res = fp(base, n2, n);
    if ((res == 1) || (res == n - 1)) return 1;
    s--;
    while (s >= 0) {
        res = fm(res, res, n);
        if (res == n - 1) return 1;
        s--;
    }
    return 0;// n is not a strong pseudo prime
}
int isprime(ll n) {
    static ll testNum[] = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37};
    static ll lim[] = {4, 0, 1373653ll, 25326001ll,25000000000ll, 2152302898747ll, 3474749660383ll, 341550071728321ll,0, 0, 0, 0};
    if (n < 2 || n == 3215031751ll) return 0;
    for (int i = 0; i < 12; ++i) {
        if (n < lim[i]) return 1;
        if (strong_pseudo_primetest(n, testNum[i]) == 0) return 0;
    }
    return 1;
}
posted @ 2020-12-02 18:43  UCPRER  阅读(293)  评论(0编辑  收藏  举报