算法详解 - 素数

定义

素数，又被称为质数。我们一般定义质数为除了 \(1\) 和自身之外没有其他因数的数，特别的 $1 $ 不是质数。

判断

试除法

我们考虑如何判断质数，一种很显然的思路就是枚举这个数的因数。假设这个数是 \(x\)，则我们可以枚举 \([2,\,\sqrt{x}]\) 之间的整数，判断是否可以整除，如果没有找到一个可以整除的数，那么可以认为这是一个质数。时间复杂度是 \(O(\sqrt{x})\)​​。显然这并不是很快。

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

Fermat 素数测试

我们都知道，对于素数 \(p\)，根据费马小定理： \[ a^{p-1} \equiv 1 \pmod p \] 那么我们猜想是否这个逆定理也成立，也就是说，如果 \(a^{p - 1} \equiv 1 \pmod p\)，那么 \(p\) 是质数。很显然，你可以举出很多反例，大悲。

于是我们再次猜想多枚举几个 \(a\) 是否可以判断掉更多的伪质数？显然是可以的，但是有一些数即使你遍历完 \([2,\,p - 1]\) 的所有数，也无法判断出来。这些数被称为卡迈克尔数（Carmichael Number），也被称为强伪素数。举个例子，\(561\) 就是最小的卡迈克尔数。

很不幸的是，卡迈克尔数有无数个，但是很幸运的是 \(10^8\) 以内的卡迈克尔数只有 \(255\) 个，特判一下即可。显然复杂度更加劣了，正确性也堪忧。

bool fermat(int n) {
    if (n < 3) return n == 2;
    // test_time 为测试次数,建议设为不小于 8
    // 的整数以保证正确率,但也不宜过大,否则会影响效率
    for (int i = 1; i <= test_time; ++i) {
        int a = rand() % (n - 2) + 2;
        if (quickPow(a, n - 1, n) != 1) return 0;
    }
    return 1;
}

Miller–Rabin 素性测试

Miller-Rabin 是一种极其快速的新型素数判定方法，其本质是对 Fermat 素数测试的一些优化和新的判断方式。首先我们需要知道一些定理。

• 二次探测定理。我们令 \(p\) 为任意一个奇素数，则若 \(x^2 \equiv 1 \pmod p\) 的解为 \(x \equiv 1 \pmod p\) 或 \(x \equiv p - 1 \pmod p\)。

于是我们可以将二次探测定理和费马素数测试结合起来，得到一个较稳稳定的算法。

我们考虑将费马素数判断中的式子 \(a^{p - 1} \equiv 1 \pmod p\) 改写成 \(\left(a^{\frac{p - 1}{2}}\right)^2 \equiv 1 \pmod p\)。为什么我们可以这样改呢？很显然 ， \(p\) 是质数当且仅当 \(p\) 是奇数（\(2\) 除外），因此 \(p - 1\) 必然是偶数。于是我们可以用二次探测定理来将 \(\left(a^{\frac{p - 1}{2}}\right)^2 \equiv 1 \pmod p\) 化为 \(a^{\frac{p - 1}{2}} \equiv 1 \pmod p\) 或 \(a^{\frac{p - 1}{2}} \equiv p - 1 \pmod p\)，判断是否符合费马素数测试。如果符合，那么继续分解测试下去，，直到变成奇数。注意不要产生 \(a^{p - 1} \equiv 0 \pmod p\) （\(p\,|\,a^{p -1 }\)） 的情况。

然后就是底数的选择了。经过数千万 oier 的验证，在 int 范围内，使用 \(\{2,\,7,\,61\}\) 这三个底数是足够的。在 long long 的范围内，使用 \(\{2,\,325,\,9375,\,28178,\,450775,\,9780504,\,1795265022\}\) 这 \(7\) 个底数也是完全足够的。时间复杂度为 \(\mathcal{O}(\log_3n)\) 带一个 \(7\) 或 \(3\)​​ 倍的常数，非常优秀！

下面给出一种正确性较高的实现方式（源自 oi-wiki）

```cpp bool millerRabin(int n) { if (n < 3 || n % 2 == 0) return n == 2; int u = n - 1, t = 0; while (u % 2 == 0) u /= 2, ++t; // test_time 为测试次数，建议设为不小于 8 // 的整数以保证正确率，但也不宜过大，否则会影响效率 for (int i = 0; i < test_time; ++i) { int a = rand() % (n - 2) + 2, v = quickPow(a, u, n); if (v == 1) continue; int s; for (s = 0; s < t; ++s) { if (v == n - 1) brea