「OI」Miller Rabin 素性测试

Miller Rabin 素性测试 与 Pollard’s Rho 因数分解以它们的随机性和高效性著称,甚至 Linux 下的 factor 都使用了它们。


Miller Rabin 素性测试

简介

用于判断一个数是不是素数。算法的渐进时间复杂度是 $\Omicron(k\log^{3}n)$,用 FFT 或者其他科技可以优化到 $\Omicron(k\log^{2}n\log\log n\log\log\log n)$,详见Miller–Rabin primality test – Wikipedia

注:没有通过 Miller Rabin 素数测试的一定是合数,而通过的不一定是质数。

理论基础

费马小定理

定理

若 $p\in\mathbb{P},p\nmid a$,则有$a^{p-1}\equiv 1\pmod p$。

证明

我们取一个不为 $p$ 倍数的数 $a$,以保证 $\gcd(a,p)=1$。

构造一个序列:$A={1,2,3,\dots,n-1}$,这个序列有着这样一个性质:

$$\prod_{i=1}^{n}A_i\equiv\prod_{i=1}^{n}(A_i\times a)\pmod p$$

又因为每一个 $A_i\times a\pmod p$ 都是独一无二的,且$A_i\times a\pmod p < p$,所以每一个$A_i\times a$都对应了一个$A_i$。

设 $f=(p-1)!$,则 $f\equiv a\times A_1\times a\times A_2\times a \times A_3\dots\times A_{p-1}\pmod p$。

$$a^{p-1}\times f\equiv f\pmod p$$

$$a^{p-1}\equiv 1\pmod p$$

证毕。

推论

若 $ p $ 为质数,那么 $\forall a\in\mathbb{N_{+}},a\in[2,p-1]$,有 $a^{p-1}\equiv 1\pmod p$。

不足之处

当 $p$ 为合数时,依然有可能满足 $\forall a\in\mathbb{N_{+}},a\in[2,p-1],a^{p}\equiv a\pmod p$,例如 $561=3\times 11\times 17$ 和 $41041=7\times 11\times 13\times 14$ 以及 $1105$ 等,更多请参考卡迈克尔数 OEIS:A002997/卡迈克尔数_百度百科

你可以通过下面的程序生成 $10^9$ 以内的卡迈克尔数。

//file:fermat.cpp
#include<bits/stdc++.h>
using namespace std;
#define reg register
typedef long long ll;
typedef unsigned long long ull;
inline int read(void){
    reg bool f=false;
    reg char ch=getchar();
    reg int res=0;
    while(ch<'0'||'9'<ch)f|=(ch=='-'),ch=getchar();
    while('0'<=ch&&ch<='9')res=10*res+ch-'0',ch=getchar();
    return f?-res:res;
}

#define mul(a,b,mod) ( (a) * (b) % (mod) )

inline ll pow(reg ll x,reg ll exp,reg ll mod){
    reg ll res=1;
    while(exp){
        if(exp&1)
            res=mul(res,x,mod);
        x=mul(x,x,mod);
        exp>>=1;
    }
    return res;
}

inline bool Miller_Rabin(reg ll n){
    if(n<3)
        return n==2; //二次探测定理仅仅使用于奇素数
    const int TEST_TIMES=6; //测试次数,尽量不要小于 8
    reg ll u=n-1,t=0;
    while((u&1)==0)
        ++t,u>>=1;
    for(reg int i=1;i<=TEST_TIMES;++i){
        reg ll x=rand()%(n-2)+2,v=pow(x,u,n);
        if(v==1||v==n-1)
            continue;
        reg int j;
        for(j=0;j<t;++j){
            v=mul(v,v,n);
            if(v==n-1)
                break;
        }
        if(j>=t)
            return false;
    }
    return true;
}

inline bool Fermat(reg ll n){
    for(reg ll i=2;i<=n;++i)
        if(__gcd(i,n)==1&&pow(i,n-1,n)!=1)
            return false;
    return true;
}

int main(void){
    ll n=1e9;
    for(reg ll i=2;i<=n;++i)
        if(!Miller_Rabin(i)&&Fermat(i))
            printf("%lld\n",i);
    return 0;
}

二次探测定理

定理

如果 $p$ 是奇素数,则 $x^{2}\equiv 1\pmod p$ 的解为 $x=1$ 或者 $x=p-1\pmod p$。

证明

$$x^{2}\equiv 1\pmod p\Leftrightarrow p\mid(x^{2}-1)=(x+1)(x-1)$$

所以有

$$p\mid(x+1)$$ 或 $$p\mid(x-1)$$

继而

$$x=p-1$$ 或 $$x=p+1$$

所以 $x^{2}\equiv 1\pmod p$ 的解为 $x=1$ 或者 $x=p-1 \pmod p$。

得证。

推论

如果 $p$ 是奇素数,那么将 $p-1$ 分解为 $p-1=u\times 2^t$,对于一个 $a\in[2,p-1]$,考虑序列 $(a^u,a^{2u},\cdots,a^{u\times 2^t})$,在模 $p$ 意义下必定形如 $(p-1,p-1,1,\cdots,1,1)$(前面可以没有 $p-1$)。

理论实现

根据卡迈克尔数的性质,可知其一定不是 $p^e$。

不妨将费马小定理和二次探测定理结合起来使用:

将 $n-1$ 分解为 $n-1=u\times 2^{t}$,不断地对 $u$ 进行平方操作,若发现非平凡平方根时即可判断出其不是素数。

代码实现

inline bool Miller_Rabin(reg ll n){
    if(n<3)
        return n==2; //二次探测定理仅仅使用于奇素数
    const int TEST_TIMES=6; //测试次数,尽量不要小于 8
    reg ll u=n-1,t=0;
    while((u&1)==0)
        ++t,u>>=1;
    for(reg int i=1;i<=TEST_TIMES;++i){
        reg ll x=rand()%(n-2)+2,v=pow(x,u,n);
        if(v==1||v==n-1)
            continue;
        reg int j;
        for(j=0;j<t;++j){
            v=mul(v,v,n);
            if(v==n-1)
                break;
        }
        if(j>=t)
            return false;
    }
    return true;
}