Miller-Rabin素数测试算法

前提

\(Miller-Rabin\) 素数测试算法需要如下两个引理：

费马小定理

设 \(p\) 是素数，\(a\) 为整数，且 \((a,p)=1\)，则 \(a^{p-1}\equiv 1 \pmod p\)

Q:求证?

A:
考虑 \(1,2,3,\dots,(p-1)\) 这 \((p-1)\) 个数字，同时乘 \(a\) 的 \(a,2a,3a,\dots，(p-1)a\)。

\(∵ a\equiv b\pmod{p},(c,p)=1\)

\(∴ac\equiv bc\pmod{p}\)

\(∴1\times2\times\dots\times(p-1) \equiv 1\times2\times\dots\times(p-1)\pmod{p}\)

\(∴(p-1)! \equiv (p-1)!\times a^{p-1}\pmod{p}\)

\(∵((p-1),p)=1\)

\(∴a^(p-1)\equiv 1 \pmod{p}\)

二次探测定理

设 \(p\) 是一个素数，且 \(0<x<p\)，则方程 \(x^2\equiv 1 \pmod{p}\) 的解为 \(x=1/p-1\)

Q:求证。

A:
容易知道 \(x^2-1\equiv 0 \pmod{p}\)

\(∴(x+1)(x-1)\equiv 0 \pmod{p} \leftarrow因式分解\)

\(∴p|(x+1)(x-1)\)

\(∵p是素数\)

\(∴\begin{cases}x=1\\x=p-1&\end{cases}\)

算法介绍

首先，费马小定理只是判定 \(p\) 为质数的必要条件。即费马小定理不成立，\(n\) 一定是合数；费马小定理成立，\(n\) 也有可能不是质数。

下面来看 \(Miller-Rabin\) 算法的分析过程：

假设 \(n\) 为素数且 \(n>2\)，则 \(n-1\) 为偶数。令 \(n-1=2^q\cdot m\)

随机选取整数 \(a(0<a<n)\)。

由费马小定理可得 \((a^{2^q \cdot m}=a^{n-1})\equiv \pmod{n}\)

由二次探测定理可得 \(a^{2^{q-1}\cdot m}\equiv1\pmod{n}\) 或 \(a^{2^{q-1}\cdot m}\equiv n-1\pmod{n}\)。

若\(a^{2^{q-1}\cdot m}\equiv1\pmod{n}\) 成立，那么再由二次探测定理可得\(a^{2^{q-2}\cdot m}\equiv 1\pmod{n}\) 或 \(a^{2^{q-2}\cdot m}\equiv n-1\pmod{n}\ldots\)

反复使用二次探测定理拆解，直到拆解到 \(a^m \equiv 1 \pmod{n} 或 a^m \equiv n-1 \pmod{n}\)。

总结一下： \(a^m \equiv 1 \pmod{n}\) 或存在 \(0\le r < q\) 使得 \(a^{2^{r}\cdot m}\equiv n-1 \pmod{n}\)，则称 \(n\) 通过了 \(Miller-Rabin\) 素数测试算法。

但是还是存在一种伪素数，他们可以通过 \(Miller-Rabin\) 素数测试算法但不是素数。

但是经过证明\(Miller-Rabin\) 素数测试算法的错误率小于等于 \(\dfrac{1}{4}\)。若用不同的素数作为 \(a\) 测验 \(k\) 次，错误率可降至 \(4^{-k}\)。

这里还有一个小技巧。如果被测数小于\(4759123141\) 那么只需要测试 \(3\) 个底数 \(2,7,61\) 即可。当然，你测试的越多，范围也就越大。如果你用前 \(7\) 个素数 \((2,3,5,7,11,13,17)\) 测试，不超过 \(341550071728320\) 都是正确的。如果你选用 \(2,3,7,61,24251\) 作为底数，那么 \(10^{16}\)内唯一的伪素数是 \(46856248255981\)，是 \(4840261\) 的倍数。对于 \(OI\) 选手来说，已经非常够用了。

\(Miller-Rabin\) 素数测试算法代码\(\downarrow\)

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll prime[10] = {2,3,5,7,11,13,17,19,23,29};
ll Quick_Multiply(ll a, ll b, ll c)
{
    ll ans = 0, res = a;
    while(b)
    {
        if(b & 1)ans = (ans + res) % c;
        res = (res + res) % c;
        b >>= 1;
    }
    return ans;
}
ll Quick_Power(ll a, ll b, ll m)
{
    ll ans = 1;
    while(b)
    {
        if(b & 1)
        {
            ans = ans * a % m;
        }
        a = a * a % m;
        b >>= 1;
    }
    return ans;
}
bool Miller_Rabin(ll x)    //判断素数
{
    ll i, j, k;
    ll s = 0, t = x-1;
    if(x == 2) return true; //2是素数
    if(x < 2 || !(x & 1)) return false;  //如果x是偶数或者是0，1，那它不是素数
    while(!(t & 1)) //将x-1分解成(2^s)*t的样子
    {
        s++;
        t >>= 1;
    }
    for(i = 0; i < 10 && prime[i] < x; i++) //随便选一个素数进行测试
    {
        ll a = prime[i];
        ll b = Quick_Power(a, t, x);   //先算出a^t
        for(j = 1; j <= s; j++) //然后进行s次平方
        {
            k = Quick_Multiply(b,b,x);  //求b的平方
            if(k == 1 && b != 1 && b != x-1)
                return false;
            b = k;
        }
        if(b != 1) return false;    //用费马小定理判断
    }
    return true;    //如果进行多次测试都是对的，那么x就很有可能是素数
}

int main()
{
    ll x;
    scanf("%lld",&x);
    if(Miller_Rabin(x)) printf("Yes");
    else printf("No");
    return 0;
}

若文章有误，欢迎在下方评论区指出，谢谢！