【数论】Miller_Rabin 素性测试

· · 个人记录

本篇博客将对素性测试的方法进行讲解。

Part 1. 小素数判定

没什么好说的，根号复杂度筛质因子，可以跑单个 n\le 10^{12} 级别，跑多个用线性筛即可。

代码：

bool isprime(int x){
    if(x<2)return 0;
    for(int i=2;i*i<=x;i++)if(x%i==0)return 0;
    return 1;
}

Part 2. 大素数素性测试

素性测试是指在不对给定数进行质因子分解的情况下判断这个数是否为素数，通常包括确定性测试和概率性测试。

因为确定性测试十分复杂且比概率性测试慢很多，这里只讲解概率性测试。

Part 3. Fermat 素性测试

这种测试是使用费马小定理进行素数判定的方法。

设当前测试的数为 n，基本思想是不断地选取在 [2, n-1] 中的基 a，并检验是否每次都有 n 与 a 互质且 a^{n-1} \equiv 1 \pmod n。

证明显然，单次测试复杂度为 O(\log_2 n)，设测试 k 次，总时间复杂度即为 O(k\log_2 n)，对于单个 n\le 10^{18} 的数大约可以测试 1.5\times10^4 次。

然而这种方式会有较大的错误率，有一些合数对于一些基 a 也满足 a^{n-1}\equiv 1\pmod{p}，比如 n=341,a=2。更为不幸的是，有一些合数无论如何改变基都能通过 Fermat 素性测试，这类合数被称为 Carmichael 数，具体详见 OI wiki Link。

Part 4. Miller_Rabin 素性测试

这种测试是基于 Fermat 素性测试上优化后的素性测试。

先来证明一个定理：

证明：移项得 x^2-1\equiv 0 \pmod p，即 (x+1)(x-1)\equiv 0 \pmod p。

因为 p 为奇素数，所以有 x+1\equiv 0 \pmod p 或 x-1\equiv 0 \pmod p。

解得 x\equiv 1 \pmod p 或 x\equiv p-1 \pmod p。

得证。

具体步骤如下：

  1. 选取一个素数 p，判断 n 是否等于 p 或者为 p 的倍数并做出判断。

  2. 用费马小定理判断其是否极有可能是一个素数。

  3. 利用二次探测定理判断其是否是一个素数。具体做法如下：

将 n-1 拆分成 u\times 2^t 的形式，则费马小定理中的 p^{n-1}\equiv 1\pmod n 可化为 p^{u\times 2^t}\equiv 1\pmod n。

然后每次将左侧开方，第 i 次将左侧化为 p^{u\times2^{t-i}}，判断其是否等于 1 或 n-1，如果不满足则可直接判断 n 不是奇素数。

一个小细节：

如果在某个时刻 p^{u\times 2^s}\equiv n-1\pmod n，则可以直接判定 n 通过本轮测试，因为我们无法再使用二次探测定理继续测试。

注意：二次探测定理只能判奇素数，因此测试前就要特判素数 2。

核心代码如下：

int pp[20]={2,3,5,7,11,13,17,19,131,1145141};
int y=5;
ll qpow(ll a,ll b,ll md){
    a%=md;
    ll ans=1;
    while(b){
        if(b&1)ans=ans*a%md;
        a=a*a%md;
        b>>=1;
    }
    return ans;
}
struct Miller_Rabin{
    inline bool check(ll x,ll p){
        if(x%p==0||qpow(p,x-1,x)!=1)return 0;
        ll k=x-1;
        while(k%2==0){
            ll t=qpow(p,k>>=1,x);
            if(t!=1&&t!=x-1)return 0;
            if(t==x-1)return 1;
        }
        return 1;
    }
    inline bool isprime(ll x,ll y){
        if(x<2)return 0;
        if(x==2)return 1;
        for(ll i=0;i<y;i++){
            if(x==pp[i])return 1;
            if(check(x,pp[i])==0)return 0;
        }
        return 1;
    }
}mr;

模版题：P5285 2p部分。