梅森素數 判定總結 - Lucas-Lehmer算法 & Miller-rabin算法

梅森素數

定義：

• if m是一個正整數 and 2^m-1是一個素數 then m是素數
• if m是一個正整數 and m是一個素數 then M(m)=2^m-1被稱為第m個梅森數
• if p是一個素數 and M(p)是一個素數 then M(p)被稱為梅森素數

Lucas-Lehmer判定法：判定一個梅森數是否是梅森素數

設p是素數，第p個梅森數為M(p)為2^p-1,r1 = 4,對於k >= 2

r(k) = r(k+1)^2-2(modM(p)), 0 <= r(k) <= M(p)

可以得到r(k)序列,則有M(p)是素數，當且僅當r(p-1) = 0(mod M(p))

推論:設p是素數，M(p)為第p個梅森數，則算法復雜度為O(n^3)

梅森素數 - nefu 120

思路：R.1 = 4;R.k = (R.k-1 ^ 2 - 2) % Mp;

如果R.p-1 == 0，則是梅森素數，否則不是。
特殊判斷:p == 2，即Mp = 3是梅森素數。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cstdlib>

using namespace std;
typedef long long ll;

ll multi(ll a, ll b, ll m)
{
    ll ret = 0;
    while(b>0)
    {
        if(b&1)
        {
            ret = (ret+a)%m;
        }
        b >>= 1;
        a = (a<<1)%m;
    }
    return ret;
}
int main()
{
    ll sum = 1, data[66], tmp;
    int n, p;
    data[1] = 4;
    cin >> n;
    while(n--)
    {
        sum = 1;
        cin >> p;
        sum <<= p;
        sum -= 1;
        for(int i = 2; i <= p-1; i++)
        {
            tmp = multi(data[i-1],data[i-1],sum);
            data[i] = (tmp-2)%sum;
        }
        if(p == 2)
            cout << "yes" << endl;
        else
        {
            if(data[p-1] == 0)
                cout << "yes" <<endl;
            else
                cout << "no" << endl;
        }
    }
    return 0;
}

模板：

long long multi(long long a, long long b, long long m){//實現a * b % m的操作,用2 * 3 = 6模擬一下就懂了
    long long ans = 0;
    while(b > 0){
        if(b & 1)  ans = (ans+a) % m;
        b >>= 1;
        a = (a<<1) % m;
    }
    return ans;
}
//判斷是否是梅森素數
bool is_msPrime(int p){
    long long r[70];
    long long m = 1;
    m <<= p;  m -=1;//求出Mp;
    r[1] = 4LL;
    if(p == 2)  return true;
    for(int i = 2; i <= p-1; ++i)
        r[i] = (multi(r[i-1],r[i-1],m)-2) % m;
    if(r[p-1] == 0)  return true;
    return false;
}

Miller-rabin 素數測試：直接判斷M(p)是不是素數

理論知識：

費馬小定理： 對於素數p和任意整數a，有ap ≡ a(mod p)（同余）。反過來，滿足ap ≡ a(mod p)，p也幾乎一定是素數。

偽素數： 如果n是一個正整數，如果存在和n互素的正整數a滿足 an-1 ≡ 1(mod n)，我們說n是基於a的偽素數。如果一個數是偽素數，那么它幾乎肯定是素數。

Miller-Rabin測試： 不斷選取不超過n-1的基b(s次)，計算是否每次都有bn-1 ≡ 1(mod n)，若每次都成立則n是素數，否則為合數。

還有一個定理，能提高Miller測試的效率：

二次探測定理： 如果p是奇素數，則 x2 ≡ 1(mod p)的解為 x = 1 || x = p - 1(mod p);

兩個高效求解a*b%m a^b%m的方法

// a * b % n
//例如: b = 1011101那么a * b mod n = (a * 1000000 mod n + a * 10000 mod n + a * 1000 mod n + a * 100 mod n + a * 1 mod n) mod n 

ll mod_mul(ll a, ll b, ll n) {
    ll res = 0;
    while(b) {
        if(b&1)    res = (res + a) % n;
        a = (a + a) % n;
        b >>= 1;
    }
    return res;
}

//a^b % n
//同理
ll mod_exp(ll a, ll b, ll n) {
    ll res = 1;
    while(b) {
        if(b&1)    res = mod_mul(res, a, n);
        a = mod_mul(a, a, n);
        b >>= 1;
    }
    return res;
}

代碼如下：

bool miller_rabin(ll n) {
    if(n == 2 || n == 3 || n == 5 || n == 7 || n == 11)    return true;
    if(n == 1 || !(n%2) || !(n%3) || !(n%5) || !(n%7) || !(n%11))    return false;

    ll x, pre, u;
    int i, j, k = 0;
    u = n - 1;    //要求x^u % n

    while(!(u&1)) {    //如果u為偶數則u右移，用k記錄移位數
        k++; u >>= 1;
    }

    srand((ll)time(0));
    for(i = 0; i < S; ++i) {    //進行S次測試
        x = rand()%(n-2) + 2;    //在[2, n)中取隨機數
        if((x%n) == 0)    continue;

        x = mod_exp(x, u, n);    //先計算(x^u) % n，
        pre = x;
        for(j = 0; j < k; ++j) {    //把移位減掉的量補上，並在這地方加上二次探測
            x = mod_mul(x, x, n);
            if(x == 1 && pre != 1 && pre != n-1)    return false;    //二次探測定理，這里如果x = 1則pre 必須等於 1，或則 n-1否則可以判斷不是素數
            pre = x;
        }
        if(x != 1)    return false;    //費馬小定理
    }
    return true;
}