Miller-Rabin素數測試算法

\(Miller-Rabin\)​素數測試

用途

判斷整數\(n\)是不是質數，在\(n\)較小的狀況下，可使用試除法，時間複雜度爲\(O(\sqrt n)\)。但當\(n\)的值較大的時候，樸素的試除法已經不能在規定時間內解決問題。此時，咱們能夠用\(Miller-Rabin\)素數測試算法，時間複雜度能夠下降至\(O(\log_2n)\)。

引理

費馬小定理

若\(a,p \in \mathbb{Z}\)，\(p\)爲質數，則
\[ a^{p-1} \equiv 1(mod\;p) \]
在此不給出證實。

二次探測定理

描述

若\(a,p \in \mathbb{Z}\)，\(a^{2} \equiv 1(mod\;p)\)，\(p\)爲質數，則\(a \equiv 1(mod\;p)\)或\(a \equiv p-1(mod\;p)\)。

證實

\[ \begin{aligned} &\because a^{2} \equiv 1(mod\;p)\\ &\therefore p \mid (a^{2}-1)\\ &\therefore p \mid (a+1)(a-1)\\ &\because p爲質數\\ &\therefore p \mid (a+1) 或(a-1)\\ &\therefore a+1 \equiv 0(mod\;p)或a-1 \equiv 0(mod\;p)\\ &\therefore a \equiv 1 (mod\;p)或a \equiv p-1 (mod\;p)\\ \end{aligned} \]

過程

根據費馬小定理，咱們能夠獲得一個真命題：若\(p\)爲質數，則\(a^{p-1} \equiv 1(mod\;p)\)。咱們考慮這一命題的逆命題：若\(a^{p-1} \equiv 1\)，則\(p\)爲質數。咱們會驚訝地發現，這一逆命題在大多數狀況下居然成立。也就是說，咱們獲得了一種有效地判斷質數的方法，即取一個底數\(a\)，判斷它與所需判斷的數\(p\)是否知足這一等式。儘管有時可能出錯，但這一算法的效率相比起樸素算法來講有了很大的提高。

接下來咱們要作的就是提升這一算法的正確性。首先想到的天然是取多個\(a\)值，在常見的題目中，取\([2,29]\)大概就能經過測試，固然也能夠隨機生成，注意\(a\)的值應該小於\(p\)。第二個優化是基於二次探測定理的。設\(p=2^nm+1\)，則可先算出\(a^m\)，而後再平方\(n\)次，求得\(a^{p-1}\)。在這一過程當中，若某次平方後所得的結果爲\(1\)但上次平方後的結果不等於\(p-1\)或\(1\)，就出現了矛盾，從而就不知足\(p\)爲質數這一前提。最後再次判斷是否知足等式便可。

注意乘法可能越界，應拆成相似快速冪的算法。

代碼

const int prime[10]={2,3,5,7,11,13,17,19,23,29};
long long multi(long long a,long long b,long long p)
{
    long long t=0;
    while(b)
    {
        if(b&1)
            t=(t+a)%p;
        a=(a<<1)%p;
        b>>=1;
    }
    return t;
}
long long power(long long a,long long b,long long p)
{
    long long t=1;
    while(b)
    {
        if(b&1)
            t=multi(t,a,p);
        a=multi(a,a,p);
        b>>=1;
    }
    return t;
}
bool Miller_Rabin(long long x)
{
    if(x==2)
        return true;
    if(!(x&1)||x<2)
        return false;
    long long t=x-1,exponent=0;
    while(!(t&1))
    {
        t>>=1;
        ++exponent;
    }
    for(int i=0;i<10&&prime[i]<x;++i)
    {
        long long m=power(prime[i],t,x);
        for(int j=0;j<exponent;++j)
        {
            long long n=multi(m,m,x);
            if(n==1&&m!=1&&m!=x-1)
                return false;
            m=n;
        }
        if(m!=1)
            return false;
    }
    return true;
}