有的时候，我们不想求前 $\pi(n)$ 个正素数，我们只想知道 $n$ 是不是素数，此时可以使用素性检验.

素性检验分为两种：

• 确定性检验.
• 概率性检验.

后者通常比前者快很多，但是有很小的概率错误地将合数识别为素数（反之则不会），因此，通过概率素性测试的数字被称为可能素数，直到它们的素数可以被确定性地证明.

而通过测试但实际上是合数的数字则被称为伪素数，有许多特定类型的伪素数，最常见的是 Fermat 伪素数，它们是满足 Fermat 小定理的合数.

这篇笔记将重点叙述一种确定性检验：试除法和两种概率性检验：Fermat 素性测试和 Miller-Rabin 素性测试.

试除法

这种方法是枚举从小到大的每个小于 $n$ 且不为 $1$ 的数，并观察它是否能整除 $n$：

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

优化

这种做法十分稳妥，但是我们容易注意到如下事实：若 $x$ 是 $n$ 的约数，则 $\dfrac nx$ 也是 $n$ 的约数.

这个事实告诉我们，我们可以只考察 $x$ 与 $\dfrac nx$ 的其中一项. 方便起见，我们考察其中较小的一项，由此，我们只需要考察 $\left[1, \sqrt n\right]$ 中的整数即可，而 $1$ 是任何数的约数，故不必考察 $1$.

另一种较为朴素的优化思想是，我们只考虑对奇数进行判断，因为除了 $2$ 以外的偶数一定有约数 $2$. 此外，$3$ 的倍数也一定有约数 $3$. 故综合以上叙述后，可考虑如下优化：

int isPrime(int n) {
    if (n < 3 || !(n % 2)) {
        return n == 2;
    } else if (!(n % 3)) {
        return n == 3;
    }
    for (int i = 2; (long long)i * i <= n; i ++) {
        if (n % i == 0) {
            return 0;
        }
    }
    return 1;
}

以上是关于试除法的叙述，在进行概率性检验的叙述前，需要先引入一条数论知识.

Fermat 小定理

该定理的内容如下： 设 $p$ 是素数，对于任意整数 $a$，若有 $p\nmid a$，则有 $a^{p-1}\equiv1\left(\mathrm{mod}\space p\right)$

另一种等价叙述是： 设 $p$ 是素数，对于任意整数 $a$，都有 $a^p\equiv a\left(\mathrm{mod}\space p\right)$

以下是直接证明：

$\mathbf{Proof.}$

设 $p$ 是素数且 $p\nmid a$，首先证明：对于 $i=1,2,\cdots,p-1$，余数 $ia\space\mathrm{mod}\space p$ 各不相同. 考虑反证法：如果存在 $1\leqslant i<j<p$ 使得

$$ia\space\mathrm{mod}\space p=ja\space\mathrm{mod}\space p $$

可得 $\left(i-j\right)a\equiv0\left(\mathrm{mod}\space p\right)$，而 $a$ 与 $i-j$ 显然都不整除 $p$，因此矛盾. 所以对于 $i=1,2,\cdots,p-1$，余数 $ia\space\mathrm{mod}\space p$ 各不相同，也即 $\{ia\space\mathrm{mod}\space p\}$ 是 $\{1,2,\cdots,p-1\}$ 的一个排列，因此

$$\prod_{i=1}^{p-1}i=\prod_{i=1}^{p-1}ia\space\mathrm{mod}\space p\equiv\prod_{i=1}^{p-1}ia=a^{p-1}\prod_{i=1}^{p-1}i\left(\mathrm{mod}\space p\right) $$

这说明

$$\left(a^{p-1}-1\right)\prod_{i=1}^{p-1}i\equiv0\left(\mathrm{mod}\space p\right) $$

也就是说，$\displaystyle p\mid\left(a^{p-1}-1\right)\prod_{i=1}^{p-1}i$，而 $i=1,2,\cdots,p-1$ 都不是 $p$ 的倍数，所以 $p\mid\left(a^{p-1}-1\right)$，即Fermat 小定理成立.

$$\mathfrak{Q.E.D}$$

Fermat 小定理的逆命题并不成立，即使对于所有与 $n$ 互素的 $a$，都有 $a^{n-1}\equiv1\left(\mathrm{mod}\space n\right)$，那么 $n$ 也未必是素数. 这也是 Fermat 素性测试中错误概率的来源.

Fermat 素性测试

Fermat 素性测试（Fermat primality test）是最简单的概率性素性检验.

我们可以根据 Fermat 小定理得出一种检验素数的思路，基本思想是不断地选取在 $\left[2, n-1\right]$ 中的基底 $a$，并检验是否每次都有 $a^{n-1}\equiv1\left(\mathrm{mod}\space n\right)$：

int Fermat(int n) {
    if (n < 3 || !(n % 2)) {
        return n == 2;
    } else if (!(n % 3)) {
        return n == 3;
    }
    for (int i = 1; i <= testTime; i ++) {
        int a = rand() % (n - 2) + 2;
        if (quickPow(a, n - 1, n) != 1) {
            return 0;
        }
    }
    return 1;
}

上述代码中，参数 testTime 是测试次数，一般建议设置为不小于 $8$ 的整数，在保证效率的同时也保证正确率（不为 $100\%$）.

如果 $a^{n-1}\equiv1\left(\mathrm{mod}\space n\right)$ 但 $n$ 不是素数，则称 $n$ 为以 $a$ 为底的 Fermat 伪素数.

在实践中我们观察到，如果 $a^{n-1}\equiv1\left(\mathrm{mod}\space n\right)$，那么 $n$ 通常是素数，但其实存在反例：对于 $n=341$ 且 $a=2$，虽然有 $2^{340}\equiv1\left(\mathrm{mod}\space 341\right)$，但是 $341=11\times31$ 是合数.

事实上，对于任何固定的基底 $a$，这样的反例都有无穷多个.

既然对于单个基底，Fermat 素性测试无法保证正确性，一个自然的想法就是多检查几组基底. 但是，即使检查了所有可能的与 $n$ 互素的基底 $a$，依然无法保证 $n$ 是素数. 在上述条件下无法保证是素数的数称为 Carmichael 数，它也有无穷多个.

这迫使我们寻找更为严格的素性测试.

Miller-Rabin 素性测试

Miller-Rabin 素性测试（Miller-Rabin primality test）是一种更好的素数判定方法.

在之前我们已经提到过，要确保是素数，需要使用确定性算法，而确定性算法通常要慢得多. 然而，实际上没有已知的实际上是合数的数字通过了 Miller-Rabin 测试等高级概率性测试，因此我们可以放心使用.

在不考虑乘法的复杂度时，对数 $n$ 进行 $k$ 轮测试的时间复杂度是 $O\left(k\log n\right)$.

Miller-Rabin 素性测试常用于对高精度数进行测试，此时时间复杂度是 $O\left(k\log ^3n\right)$，利用 FFT 等技术可以优化到 $O(k\log^2 n\log\log n\log\log\log n)$.

为了解决 Carmichael 数带来的挑战，Miller-Rabin 素性测试进一步考虑了素数的如下性质：

二次探测定理：

若 $p$ 是奇素数，则同余方程 $x^2\equiv1\left(\mathrm{mod}\space p\right)$ 的解只有 $x\equiv1\left(\mathrm{mod}\space p\right)$ 或者 $x\equiv p-1\left(\mathrm{mod}\space p\right)$.

$\mathbf{Proof.}$

容易验证， $p$ 为奇素数时， $x\equiv1\left(\mathrm{mod}\space p\right)$ 与 $x\equiv p-1\left(\mathrm{mod}\space p\right)$ 都可以使上式成立. 由 Lagrange 定理可知，这就是该方程的所有解.

$$\mathfrak{Q.E.D}$$

将 Fermat 小定理和二次探测定理结合起来使用，就能得到 Miller-Rabin 素性测试：

1. 将 $a^{n-1}\equiv1\left(\mathrm{mod}\space n\right)$ 中的指数 $n-1$ 分解为 $n-1=u\times2^t$.
2. 在每轮测试中对随机出来的 $a$ 先求出 $v=a^u\space\mathrm{mod}\space n$，之后对这个值执行最多 $t$ 次平方操作.
3. 在整个过程中，如果发现 $1$ 的非平凡平方根（即除了 $\pm1$ 之外的其他根），就可以判断该数不是素数.
4. 否则，再使用 Fermat 素性测试判断.

此处还有一些实现上的小细节：

• 对于一轮测试，如果某一时刻 $a^{u\times2^s}\equiv n-1\left(\mathrm{mod}\space n\right)$，则之后的平方操作全都会得到 $1$，则可以直接通过本轮测试.
• 如果找出了一个非平凡平方根 $a^{u\times2^s}\not\equiv n-1\left(\mathrm{mod}\space n\right)$，则之后的平方操作全都会得到 $1$. 可以选择直接返回 false（0），也可以放到 $t$ 次平方操作后再返回 false（0）.

由此一来，我们便得到了较正确的 Miller-Rabin 素性测试：

int millerRabin(int n) {
    if (n < 3 || !(n % 2)) {
        return n == 2;
    } else if (!(n % 3)) {
        return n == 3;
    }
    int u = n - 1, t = 0;
    while (u % 2 == 0) {
        t ++;
        u /= 2;
    }
    for (int i = 1; i <= testTime; i ++) {
        int a = rand() % (n - 3) + 2, v = quickPow(a, u, n)
        if (v == 1) {
            continue;
        }
        int s;
        for (s = 0; s < t; s ++) {
            if (v == n - 1) {
                break;
            }
            v = (long long)v * v % n;
        }
        if (s == t) {
            return 0;
        }
    }
    return 1;
}

以上代码中的参数 testTime 的要求与 Fermat 素性测试中的相同.

可以证明，奇合数 $n>9$ 通过随机选取的一个基底 $a$ 的 Miller-Rabin 素性测试的概率至多为四分之一. 因此，随机选取 $k$ 个基底后，仍将合数误判为素数的概率不超过 $\dfrac1{4^k}$.

在 OI 范围内，通常都是对 $\left[1,2^{64}\right)$ 内的数作素性检验：

• 对于 $\left[1,2^{32}\right)$ 内的数，选取 $\{2,7,61\}$ 三个数作为基底进行 Miller-Rabin 素性检验就可以确定素性.
• 对于 $\left[1,2^{64}\right)$ 内的数，选取 $\left\{2,325,9375,28178,450775,9780504,1795265022\right\}$ 七个数作为基底进行 Miller-Rabin 素性检验就可以确定素性.

注：如果要使用上面的数列中的数作为基底 $a$ 判断 $n$ 的素性：

• 所有的数都要取一遍.
• 需要把 $a$ 换成 $a\space\mathrm{mod}\space n$.
• 如果 $a\equiv0\left(\mathrm{mod}\space n\right)$ 或 $a\equiv\pm1\left(\mathrm{mod}\space n\right)$，则直接通过该轮测试.

以下是一个对于 $\left[1,2^{32}\right)$ 内的数 $100\%$ 正确的 Miller-Rabin 素性测试：

int millerRabin(int n) {
    if (n < 3 || !(n % 2)) {
        return n == 2;
    } else if (!(n % 3)) {
        return n == 3;
    }
    int u = n - 1, t = 0;
    while (u % 2 == 0) {
        t ++;
        u /= 2;
    }
    int tester[3] = {2, 7, 61};
    for (int i = 0; i < 3; i ++) {
        int a = tester[i] % n, v = quickPow(a, u, n)
        if (a == 0 || a == n - 1 || a == 1) {
            continue;
        }
        if (v == 1) {
            continue;
        }
        int s;
        for (s = 0; s < t; s ++) {
            if (v == n - 1) {
                break;
            }
            v = (long long)v * v % n;
        }
        if (s == t) {
            return 0;
        }
    }
    return 1;
}

这里是一条笔记，用于记录与素数相关的数论算法.

本题相当于求解前 $\pi(100)$ 个正素数，一个自然的想法是对每个小于等于 $n$ 的数进行一次质数检验，但这种暴力的做法显然不能达到最优复杂度.

下面介绍几种比暴力更优的素数筛法.

Eratosthenes 筛法

我们很容易发现一个事实：对于任意一个大于 $1$ 的正整数 $n$，它的 $x$ 倍是合数（$x>1$）.

利用这个事实，我们可以避免很多次不必要的检测.

如果我们从小到大考虑每个数，然后同时把当前这个数的所有（比自己大的）倍数记为合数，那么运行结束的时候没有被标记的数就是素数了.

并且，十分显然的是，要找到直到 $n$ 为止的所有素数，仅对不超过 $\sqrt n$ 的素数进行筛选就足够了.

其代码实现如下：

int isPrime[N];
void Eratosthenes(int n) {
    isPrime[0] = isPrime[1] = 0;
    for (int i = 2; i <= n; i ++) {
        isPrime[i] = 1;
    }
    for (int i = 2; (long long)i * i <= n; i ++) {
        if (isPrime[i]) {
            for (int j = i * i; j <= n; j ++) {
                isPrime[j] = 0;
            }
        }
    }
    return;
}

以上就是 Eratosthenes 筛法（埃拉托斯特尼筛法，简称埃氏筛法）的内容，其时间复杂度是 $O\left(n\log\log n\right)$.

以下是对于其时间复杂度的证明：

$\mathbf{Proof.}$

如果每一次对数组的操作花费 $1$ 个单位时间，则时间复杂度为：

$$O\left(\sum_{k=1}^{\pi\left(\lfloor\sqrt n\rfloor\right)}\dfrac n{p_k}\right)=O\left(n\sum_{k=1}^{\pi\left(\lfloor\sqrt n\rfloor\right)}\dfrac1{p_k}\right) $$

其中 $p_k$ 表示第 $k$ 小的素数，$\sum_{k=1}^{\pi\left(\lfloor\sqrt n\rfloor\right)}$ 表示第一层 for 循环，$\dfrac n{p_k}$ 表示第二层 for 循环的执行次数. 根据 Mertens 第二定理，存在常数 $B_1$，使得：

$$\sum_{k=1}^{\pi\left(\lfloor\sqrt n\rfloor\right)}\dfrac1{p_k}=\log\log\lfloor\sqrt n\rfloor+B_1+O\left(\dfrac1{\log\lfloor\sqrt n\rfloor}\right) $$

所以 Eratosthenes 筛法的时间复杂度为 $O\left(n\log\log n\right)$. 接下来证明 Mertens 第二定理的弱化版本，即

$$\sum_{k\leqslant\pi(n)}\dfrac1{p_k}=O\left(\log\log n\right) $$

根据 $\pi\left(n\right)=\Theta\left(\dfrac n{\log n}\right)$，可知第 $n$ 个素数的大小为 $\Theta\left(n\log n\right)$. 于是就有

$$O\left(\sum_{k=1}^{\pi(n)}\dfrac1{p_k}\right)=\displaystyle O\left(\sum_{k=2}^{\pi(n)}\dfrac1{k\log k}\right)=O\left(\displaystyle \int_2^{\pi(n)}\dfrac{\mathrm dx}{x\log x}\right)=O\left(\log\log\pi\left(n\right)\right)=O\left(\log\log n\right) $$$$\mathfrak{Q.E.D.}$$

由时间复杂度可知，上述做法仍然不够高效.

但是我们很容易考虑到，只需要筛选奇数就可以，因为除了 $2$ 以外的偶数都是合数，且除了 $3$ 以外的能被 $3$ 整除的数一定是合数，所以我们可以直接跳过它们，只用关心 $3$ 以外的奇数就好.

线性筛法

埃氏筛法仍有优化空间，它会将一个合数重复多次标记.

有没有什么办法省掉无意义的步骤呢？答案是肯定的.

如果能让每个合数都只被标记一次，那么时间复杂度就可以降到 $O\left(n\right)$ 了.

其代码如下：

int notPrime[N], prime[N], num = 0;
void Euler(int n) {
	for (int i = 2; i <= n; i ++) {
	    if (!notPrime[i]) {
		    prime[num ++] = i;
	    }
	    for (int j = 0; j < num; j ++) {
		    if (i * prime[j] > n) {
			    break;
  			}
		    notPrime[i * prime[j]] = 1;
		    if (i % prime[j] == 0) {
			    break;
  			}
        }
	}
    return;
}

上面的这种线性筛法也称为 Euler 筛法（欧拉筛法）.