在对称密码体制中，加密与解密共享同一把密钥；而公钥密码体制（public-key cryptography）的根基完全不同——它依赖于某些数学问题在”正向计算容易、逆向求解困难”这一方向上的不对称性。支撑这种不对称性的核心学科，正是数论（number theory）。

本文将从最基础的模运算出发，逐步引入群论的语言，探讨原根与离散对数问题，再延伸到中国剩余定理、二次剩余、大整数运算、素性测试与整数分解等关键议题。这些内容构成了理解 RSA、Diffie-Hellman、椭圆曲线密码等后续篇章不可或缺的数学工具箱。

先用一张图把本文的核心关系串起来。

graph TD
    A[模运算] --> B[欧几里得算法]
    B --> C[CRT]
    C --> D[素性测试]
    D --> E[公钥密码基础]

一、模运算基础

同余关系

给定正整数 n，若两个整数 a 与 b 满足 n | (a − b)，即 a − b 是 n 的倍数，则称 a 与 b 关于模 n 同余（congruent modulo n），记作 a ≡ b (mod n)。同余关系具有自反性、对称性与传递性，因此是一种等价关系，它将全体整数划分为 n 个等价类：{0, 1, 2, …, n − 1}，记作 Z/nZ 或简写为 Z_n。

模运算与日常的加减乘运算兼容。具体而言，若 a ≡ a’ (mod n) 且 b ≡ b’ (mod n)，则 a + b ≡ a’ + b’ (mod n)，a × b ≡ a’ × b’ (mod n)。这一性质保证了我们可以在中间步骤随时取模，而不影响最终结果——这对密码学中的大整数运算至关重要，因为它使我们能够在不让中间值爆炸式增长的前提下完成计算。

需要特别注意的是，模运算对除法并不直接兼容。a/b mod n 这一写法本身没有明确意义——必须先将除法转化为乘以逆元，即 a × b^(−1) mod n。这种区别在密码学协议的实现中经常容易引发错误，需要格外小心。

此外，模运算中的指数运算也有简化规则。由欧拉定理（后文详述），在计算 a^k mod n 时，可以先将指数 k 对 φ(n) 取模，从而显著减小指数的规模。这一技巧是快速模幂算法的理论基础之一。

计算示例

为了让读者更直观地体会模运算的规则，下面给出几个具体的计算示例。

示例 1：模乘法中的连续取模。计算 17 × 23 mod 7。方法一（直接计算）：17 × 23 = 391，391 = 55 × 7 + 6，因此 391 mod 7 = 6。方法二（先取模再相乘）：17 mod 7 = 3，23 mod 7 = 2，3 × 2 = 6，6 mod 7 = 6。两种方法结果一致。在密码学的大整数运算中，中间结果可能长达数千位，每一步都先取模可以将中间值控制在模数的范围内，极大减少计算开销。

示例 2：模幂运算的指数约化。计算 3^100 mod 7。由后文将要介绍的费马小定理，当 p = 7 时有 3^6 ≡ 1 (mod 7)。因此 100 = 6 × 16 + 4，3^100 = (3⁶⁾16 × 3^4 ≡ 1^16 × 81 (mod 7)。而 81 = 11 × 7 + 4，故 3^100 mod 7 = 4。这一技巧将一个天文数字级的幂运算化简为了小整数运算——在 RSA 中，模数和指数都是 2048 位的大整数，如果不借助指数约化和快速幂算法，解密一条消息所需的时间将超过宇宙的年龄。

示例 3：模逆元的手动计算。求 5^(−1) mod 17，即找 x 使得 5x ≡ 1 (mod 17)。使用扩展欧几里得算法：17 = 3 × 5 + 2，5 = 2 × 2 + 1，回代得 1 = 5 − 2 × 2 = 5 − 2 × (17 − 3 × 5) = 7 × 5 − 2 × 17，因此 5 × 7 ≡ 1 (mod 17)，即 5^(−1) mod 17 = 7。验证：5 × 7 = 35 = 2 × 17 + 1 ≡ 1 (mod 17)，正确。在 RSA 密钥生成中，求解 d ≡ e^(−1) (mod φ(n)) 本质上就是这个过程的大整数版本。

模逆元

在实数域中，若 a ≠ 0 则存在 a 的乘法逆元 1/a。类比地，在 Z_n 中，若存在 x 使得 a × x ≡ 1 (mod n)，则称 x 为 a 模 n 的乘法逆元（modular multiplicative inverse），记作 a^(−1) mod n。模逆元存在的充要条件是 gcd(a, n) = 1，即 a 与 n 互素（coprime）。

扩展欧几里得算法

计算模逆元最经典的方法是扩展欧几里得算法（Extended Euclidean Algorithm）。普通的欧几里得算法通过反复取余来求 gcd(a, b)；扩展版本则在此过程中同时维护系数 x, y，使得 a × x + b × y = gcd(a, b)。当 gcd(a, n) = 1 时，等式变为 a × x + n × y = 1，对其两边取模 n，得到 a × x ≡ 1 (mod n)，故 x 即为所求的模逆元。

下面给出 Python 实现：

def extended_gcd(a: int, b: int) -> tuple[int, int, int]:
    """返回 (g, x, y) 使得 a*x + b*y = g = gcd(a, b)"""
    if b == 0:
        return a, 1, 0
    g, x1, y1 = extended_gcd(b, a % b)
    # 由 b*x1 + (a % b)*y1 = g
    # 即 b*x1 + (a - (a//b)*b)*y1 = g
    # 即 a*y1 + b*(x1 - (a//b)*y1) = g
    return g, y1, x1 - (a // b) * y1


def mod_inverse(a: int, n: int) -> int:
    """计算 a 模 n 的乘法逆元，不存在则抛出异常"""
    g, x, _ = extended_gcd(a % n, n)
    if g != 1:
        raise ValueError(f"模逆元不存在: gcd({a}, {n}) = {g}")
    return x % n

正确性证明：算法的不变量是在每一层递归中，返回的 (g, x, y) 满足 a × x + b × y = g。基础情况 b = 0 时，gcd(a, 0) = a，此时 a × 1 + 0 × 0 = a 成立。归纳步中，假设下层递归已正确返回 g, x1, y1 使得 b × x1 + (a mod b) × y1 = g，将 a mod b 展开为 a − ⌊a/b⌋ × b，经代数整理即得上层的系数关系，不变量得以维持。又因为欧几里得算法在有限步后终止（每次递归 b 严格减小且始终非负），算法正确性得证。

该算法的时间复杂度为 O(log(min(a, b)))，在密码学实践中效率极高。

值得一提的是，扩展欧几里得算法还可以用于求解一般形式的线性同余方程 a × x ≡ c (mod n)。当 gcd(a, n) | c 时，方程有 gcd(a, n) 个不同的解；否则无解。在密码协议的密钥生成阶段，求解这类方程是非常常见的操作。

二、群论入门

群的定义与公理

群（group）是一个集合 G 配上一个二元运算 ·，满足以下四条公理：

1. 封闭性（closure）：对任意 a, b ∈ G，有 a · b ∈ G。
2. 结合律（associativity）：对任意 a, b, c ∈ G，有 (a · b) · c = a · (b · c)。
3. 单位元（identity element）：存在 e ∈ G，使得对任意 a ∈ G，有 e · a = a · e = a。
4. 逆元（inverse）：对任意 a ∈ G，存在 a^(−1) ∈ G，使得 a · a^(−1) = a^(−1) · a = e。

若运算还满足交换律 a · b = b · a，则称该群为阿贝尔群（abelian group）。

元素的阶与循环群

群元素 a 的阶（order）是使得 a^k = e 成立的最小正整数 k，记作 ord(a)。若群 G 中存在某个元素 g，使得 G 中每个元素都可以写成 g 的幂次，即 G = {g^0, g^1, g^2, …}，则称 G 为循环群（cyclic group），g 为 G 的生成元（generator）。

循环群在密码学中地位极为重要。Diffie-Hellman 密钥交换、ElGamal 加密、数字签名算法（DSA）等方案都需要在一个大的循环群中工作。

笔者认为，理解「模运算背后的群结构」是学习公钥密码学最关键的一次思维跳跃。初学者往往只把 mod 当作一种取余数的算术操作，却没有意识到这些余数在乘法下构成了一个具有严格代数结构的群。一旦跨越了这道门槛，Diffie-Hellman、ElGamal、椭圆曲线等看似截然不同的方案就会在同一个框架下统一起来——它们本质上都是在某个群中利用「正向操作容易、逆向求解困难」的不对称性。群论不是公钥密码学的装饰，而是它的骨架。理解了这一点，你在阅读后续关于 Z_p* 的原根结构、椭圆曲线点群、乃至格密码中模格结构的讨论时，都会感受到同一条思想主线贯穿始终。

子群与拉格朗日定理

群 G 的子集 H，若在同一运算下也构成群，则称 H 为 G 的子群（subgroup）。判定子群的一个常用充分条件是：H 非空，且对任意 a, b ∈ H，有 a · b^(−1) ∈ H。

拉格朗日定理（Lagrange’s theorem）指出：有限群 G 的任意子群 H 的阶（即元素个数）必整除 G 的阶。证明思路是利用陪集（coset）划分：对任意 a ∈ G，左陪集 aH = {a · h : h ∈ H} 的大小与 H 相同，且不同的左陪集互不相交，它们的并恰好覆盖整个 G，因此 |G| = [G : H] × |H|，其中 [G : H] 为陪集的个数，称为 H 在 G 中的指数（index）。

作为直接推论，群中任意元素的阶也必须整除群的阶。这意味着在一个阶为 n 的有限群中，对任意元素 a，必有 a^n = e。

这一定理在密码学中有着深远的影响。例如，它直接导出了费马小定理与欧拉定理，也是理解椭圆曲线群结构的关键。在实际密码协议中，利用拉格朗日定理可以验证某个元素是否属于预期的子群——这是防范小子群攻击（small subgroup attack）的重要手段。

三、Z_p* 的结构

有了群论的基本语言，现在我们可以回到最具体的数论对象——模素数的乘法群——来观察抽象理论如何变成密码学的实用工具。前一节介绍的循环群、元素的阶、拉格朗日定理等概念，将在这里获得完全具体的实例化：Z_p* 不仅是一个有限阿贝尔群，更是一个循环群，它的生成元（原根）和离散对数结构直接支撑了 Diffie-Hellman 与 ElGamal 方案的安全性。而我们对群阶的分析，也将自然地延伸到第五节中国剩余定理的”分而治之”思想——CRT 本质上就是利用互素模数下群结构的直积分解。

模素数的乘法群

取素数 p，集合 Z_p* = {1, 2, …, p − 1} 在模 p 乘法下构成一个群。由于 p 为素数，1 到 p − 1 中的每个整数都与 p 互素，因此都有模逆元。该群的阶为 p − 1，且它是阿贝尔群。

更为重要的是，Z_p* 是循环群。这是数论中一个深刻的结论：对任何素数 p，Z_p* 中必定存在生成元，即原根（primitive root）。

原根

若 g ∈ Z_p* 满足 ord(g) = p − 1，即 g 的幂次遍历 Z_p* 的所有元素，则 g 为模 p 的原根。例如，取 p = 7，元素 3 是一个原根，因为 3^1 = 3, 3^2 = 2, 3^3 = 6, 3^4 = 4, 3^5 = 5, 3^6 = 1 (mod 7)，恰好生成了 {1, 2, 3, 4, 5, 6} 的一个排列。

让我们通过一个更完整的例子来展示原根的寻找与验证过程。取 p = 13，则 Z_13* = {1, 2, …, 12}，群的阶为 12 = 2^2 × 3。按照上述验证方法，需要检验 g^(12/2) = g^6 ≢ 1 (mod 13) 且 g^(12/3) = g^4 ≢ 1 (mod 13)。

先试 g = 2：2^4 = 16 ≡ 3 (mod 13)，不等于 1，通过第二个检验；2^6 = 64 ≡ 64 − 4 × 13 = 12 ≡ −1 (mod 13)，也不等于 1，通过第一个检验。因此 2 是模 13 的原根。验证：2^1 = 2, 2^2 = 4, 2^3 = 8, 2^4 = 3, 2^5 = 6, 2^6 = 12, 2^7 = 11, 2^8 = 9, 2^9 = 5, 2^10 = 10, 2^11 = 7, 2^12 = 1 (mod 13)，确实遍历了 Z_13* 的所有 12 个元素。

再试 g = 3：3^4 = 81 ≡ 81 − 6 × 13 = 3 (mod 13)，通过第二个检验；3^6 = 3^4 × 3^2 ≡ 3 × 9 = 27 ≡ 1 (mod 13)，未通过第一个检验——g^6 ≡ 1 意味着 3 的阶整除 6 而非 12，因此 3 不是原根。事实上 ord(3) = 3，它只能生成子群 {1, 3, 9}。这个例子清楚地展示了拉格朗日定理的力量：3 的阶 3 确实整除群的阶 12，而 2 的阶恰好等于 12，所以只有 2 才能生成整个群。

Z_p* 中原根的个数恰为 φ(p − 1)，其中 φ 是欧拉函数。当 p − 1 含有大素因子时，原根的个数通常相当可观。在密码学协议中选取生成元时，常需验证所选元素确实是原根或者属于某个大素数阶子群的生成元。

验证 g 是否为模 p 的原根的标准方法是：对 p − 1 的每个素因子 q_i，检验 g^((p−1)/q_i) ≢ 1 (mod p)。如果所有检验都通过，则 g 是原根。这比逐一计算 g 的所有幂次要高效得多，前提是已知 p − 1 的分解。

原根的存在使得 Z_p* 中的每个元素 a 都可以唯一地表示为 a = g^x mod p（0 ≤ x ≤ p − 2），此处的 x 称为 a 关于基 g 的离散对数（discrete logarithm），记作 x = log_g(a) mod p。计算离散对数被认为是一个困难问题——这正是 Diffie-Hellman 密钥交换与 ElGamal 加密方案安全性的基石。

欧拉函数

欧拉函数 φ(n)（Euler’s totient function）定义为小于等于 n 的正整数中与 n 互素的个数。对于素数 p，φ(p) = p − 1；对于两个不同素数的乘积 n = p × q，φ(n) = (p − 1)(q − 1)。更一般地，对 n 的标准分解 n = p_1^(a_1) × p_2^(a_2) × … × p_k^(a_k)，有 φ(n) = n × ∏(1 − 1/p_i)。

φ(n) 在 RSA 密码体制中扮演核心角色：RSA 的公钥与私钥之间的数学关系正是通过 φ(n)（或更精确地说，λ(n) = lcm(p − 1, q − 1)）来建立的。

四、欧拉定理与费马小定理

欧拉定理

欧拉定理（Euler’s theorem）声称：若 gcd(a, n) = 1，则 a^φ(n) ≡ 1 (mod n)。

该定理可由群论优雅地导出：a 是群 Z_n* 的元素，由拉格朗日定理，a 的阶整除群的阶 φ(n)，因此 a^φ(n) = (a^(ord(a)))(φ(n)/ord(a)) = e^(φ(n)/ord(a)) = 1。

费马小定理

当 n = p 为素数时，φ(p) = p − 1，欧拉定理退化为费马小定理（Fermat’s little theorem）：若 p 不整除 a，则 a^(p−1) ≡ 1 (mod p)。等价形式为：对任意整数 a，a^p ≡ a (mod p)。

在计算模逆元中的应用

由欧拉定理，a^φ(n) ≡ 1 (mod n)，即 a × a^(φ(n)−1) ≡ 1 (mod n)，因此 a^(φ(n)−1) mod n 就是 a 的模逆元。特别地，当模数为素数 p 时，a^(−1) ≡ a^(p−2) (mod p)。在实际密码库中，当 φ(n) 已知时，这种基于快速幂（square-and-multiply）的求逆方法有时比扩展欧几里得算法更便于实现常数时间执行，从而抵抗侧信道攻击（side-channel attack）。

快速幂算法的思想是将指数写成二进制表示，然后通过反复平方与乘法来计算幂次。例如计算 a^13 mod n：13 的二进制为 1101，从最高位到最低位扫描，遇到 1 则平方后乘以 a，遇到 0 则仅平方。整个过程只需 O(log k) 次模乘，其中 k 为指数。这一算法是密码学中几乎所有公钥操作的基本构件。

在工程实现中，快速幂与 Montgomery 乘法（后文详述）的结合是标准做法。在每一步平方或乘法操作中，使用 Montgomery 乘法来避免除法运算，使得整个模幂过程获得显著的性能提升。

五、中国剩余定理

定理陈述

中国剩余定理（Chinese Remainder Theorem, CRT）指出：若正整数 n_1, n_2, …, n_k 两两互素，则对任意整数 a_1, a_2, …, a_k，同余方程组

x ≡ a_1 (mod n_1)
x ≡ a_2 (mod n_2)
…
x ≡ a_k (mod n_k)

在模 N = n_1 × n_2 × … × n_k 下有唯一解。

构造性证明

令 N_i = N / n_i，由于 gcd(N_i, n_i) = 1，N_i 在模 n_i 下有逆元 M_i ≡ N_i^(−1) (mod n_i)。那么解为 x ≡ Σ(a_i × N_i × M_i) (mod N)。

验证：对第 j 个方程，当 i ≠ j 时 N_i 含因子 n_j，故 a_i × N_i × M_i ≡ 0 (mod n_j)；当 i = j 时，a_j × N_j × M_j ≡ a_j × 1 (mod n_j)。各项求和后 x ≡ a_j (mod n_j)，唯一性由模 N 的等价关系保证。

从代数角度看，CRT 实际上给出了一个环同构：Z_N ≅ Z_(n_1) × Z_(n_2) × … × Z_(n_k)。这一同构不仅保持加法结构，也保持乘法结构。在此同构下，模 N 的运算可以被分解为若干独立的、模数更小的运算，分别执行后再组合回来。这种”分而治之”的思想贯穿了整个密码学的高效实现。

在 RSA-CRT 优化中的应用

RSA 解密时需计算 m = c^d mod n，其中 n = p × q。直接在模 n 下做幂运算代价很高。利用 CRT，可以分别在模 p 和模 q 下计算：

m_p = c^(d mod (p−1)) mod p
m_q = c^(d mod (q−1)) mod q

然后用 CRT 合并得到 m mod n。由于模 p 和模 q 的运算各自处理的数字大小约为模 n 的一半位长，而模幂运算的代价与位长的立方成正比，因此 CRT 版本的 RSA 解密速度约为朴素版本的四倍。这也是实际 RSA 实现中几乎都采用 CRT 优化的原因。

此外，模 p 与模 q 的计算完全独立，可以在硬件中并行执行，进一步降低延迟。

需要注意的是，CRT 优化引入了额外的安全风险。Bellcore 攻击表明，如果 CRT-RSA 的计算过程中出现故障（fault），攻击者可以利用错误的签名值通过简单的 gcd 运算直接提取出 p 或 q，从而完全攻破密钥。因此，实际实现中必须对 CRT 计算的结果进行验证——通常在最终输出前检查签名是否正确，这仅需一次模幂验证即可完成。

CRT 的设计模式视角

从工程实践来看，CRT 远不止是一个数论定理——笔者更愿意把它看作密码学中一种反复出现的「设计模式」。RSA-CRT 优化只是最经典的应用；在 Shamir 秘密共享方案中，CRT 的思想被用来将一个秘密拆分为多个份额并在需要时重建；在某些多方计算（MPC）协议中，参与方各自在不同模数下执行运算，最后通过 CRT 合并结果；在全同态加密（FHE）的 RNS（Residue Number System）表示中，CRT 使得大整数运算可以被分解为多路并行的小整数运算，显著提升了吞吐量。一个经常被忽视的观点是，CRT 本质上揭示了一种「分而治之」的结构性分解——只要你的问题可以在互素的模数下独立处理，CRT 就能保证合并的正确性与唯一性。认识到这一点，你会在密码学的各个角落不断与 CRT 重逢。

下面通过一个具体例子来展示 CRT 的构造过程。求解同余方程组：x ≡ 2 (mod 3)，x ≡ 3 (mod 5)，x ≡ 2 (mod 7)。此处 N = 3 × 5 × 7 = 105。计算 N_1 = 105/3 = 35，N_2 = 105/5 = 21，N_3 = 105/7 = 15。求各模逆：M_1 = 35^(−1) mod 3 ≡ 2^(−1) mod 3 = 2（因为 2 × 2 = 4 ≡ 1 (mod 3)）；M_2 = 21^(−1) mod 5 ≡ 1^(−1) mod 5 = 1；M_3 = 15^(−1) mod 7 ≡ 1^(−1) mod 7 = 1。故 x ≡ 2 × 35 × 2 + 3 × 21 × 1 + 2 × 15 × 1 = 140 + 63 + 30 = 233 ≡ 233 mod 105 = 23。验证：23 mod 3 = 2，23 mod 5 = 3，23 mod 7 = 2，全部正确。

六、二次剩余

定义与勒让德符号

设 p 为奇素数，整数 a 满足 gcd(a, p) = 1。若存在 x 使得 x^2 ≡ a (mod p)，则 a 为模 p 的二次剩余（quadratic residue），否则为二次非剩余（quadratic non-residue）。

勒让德符号（Legendre symbol）(a/p) 定义为：若 a 为模 p 的二次剩余则取值 1，若为二次非剩余则取值 −1，若 p | a 则取值 0。

欧拉准则

欧拉准则（Euler’s criterion）给出了一种高效判断二次剩余性的方法：(a/p) ≡ a^((p−1)/2) (mod p)。

其证明直接利用费马小定理：a^(p−1) ≡ 1 (mod p)，即 (a^((p−1)/2))2 ≡ 1 (mod p)，故 a^((p−1)/2) 只能为 1 或 −1。通过原根的幂次表示可以进一步证明，恰好在 a 为二次剩余时取值 1，否则取值 −1。

雅可比符号

当模数 n 为合数时，勒让德符号推广为雅可比符号（Jacobi symbol）(a/n)，定义为 n 的素因子分解中各勒让德符号的乘积。需注意的是，雅可比符号为 1 并不意味着 a 是模 n 的二次剩余——这一微妙差别在密码学方案的安全性分析中至关重要。

与密码学的关联

二次剩余假设（Quadratic Residuosity Assumption）认为在不知道 n 的分解的情况下，判断某个数是否为模 n 的二次剩余是困难的。基于此假设，Goldwasser 与 Micali 在 1982 年构造了第一个被证明具有语义安全性的概率加密方案。Rabin 密码体制的解密过程也直接涉及模合数的平方根计算——而这一计算等价于整数分解。

在椭圆曲线密码中，判断一个元素是否为二次剩余同样是点解压缩（point decompression）操作的关键步骤。椭圆曲线上的点可以只存储 x 坐标和 y 的符号位（即一个比特），恢复完整坐标时需要计算 y^2 = x^3 + ax + b mod p 的平方根——而模素数的平方根计算本身就依赖于二次剩余的判定。当 p ≡ 3 (mod 4) 时，平方根可以直接通过 y = a^((p+1)/4) mod p 计算；一般情况下则需要 Tonelli-Shanks 算法。

七、大整数运算

密码学中的典型操作数长度为 2048 位甚至更长。普通 CPU 的字长仅为 64 位，因此需要专门的多精度算术（multi-precision arithmetic）库来处理这些大整数。

乘法算法

最基本的多精度乘法是教科书算法（schoolbook multiplication），其时间复杂度为 O(n^2)，其中 n 为字数。当操作数较大时，Karatsuba 算法通过分治策略将复杂度降至 O(n^1.585)。对更大的操作数，基于快速傅里叶变换（FFT）的 Schonhage-Strassen 算法可达到 O(n log n log log n) 的渐近复杂度。在实际密码学库（如 OpenSSL、GMP）中，通常根据操作数大小动态选择算法：小数使用教科书法，中等数使用 Karatsuba，极大数使用 FFT 方法。

Montgomery 乘法

在模幂运算中，每一步乘法后都需要取模。传统的取模运算需要一次除法，代价很高。Montgomery 于 1985 年提出的 Montgomery 乘法（Montgomery multiplication）巧妙地避开了除法运算。

核心思想是引入 Montgomery 表示：将普通整数 a 映射为 aR mod n，其中 R = 2^k 是一个略大于 n 的 2 的幂。在这一表示下，两个 Montgomery 形式的数相乘后，通过一次”Montgomery 约化”（仅涉及乘法、加法与右移操作，不需要除法）即可得到乘积的 Montgomery 形式。

def montgomery_reduce(t: int, n: int, n_inv: int, r_bits: int) -> int:
    """
    Montgomery 约化：将 t 从 [0, n*R) 约化到 [0, n) 的 Montgomery 域
    n_inv 满足 n * n_inv ≡ -1 (mod R)，R = 2^r_bits
    """
    r = 1 << r_bits
    mask = r - 1
    # 关键步骤：m = (t * n_inv) mod R，仅取低位
    m = (t * n_inv) & mask
    # t = (t + m * n) / R，右移即除以 R
    result = (t + m * n) >> r_bits
    if result >= n:
        result -= n
    return result


def montgomery_multiply(a_mont: int, b_mont: int, n: int,
                        n_inv: int, r_bits: int) -> int:
    """Montgomery 域中的乘法：输入输出均为 Montgomery 形式"""
    t = a_mont * b_mont
    return montgomery_reduce(t, n, n_inv, r_bits)

在整个模幂运算过程中，所有中间值都保持在 Montgomery 形式下，仅在最终结果时进行一次转换回普通表示。这使得整个运算过程完全避免了代价高昂的除法指令，显著提升了性能。

Barrett 约化

Barrett 约化（Barrett reduction）是另一种避免除法的取模方法。其思路是预计算 n 的倒数的近似值 ⌊2^(2k) / n⌋，然后通过乘法与移位来估算商。具体而言，对于要约化的值 a，先计算 q = ⌊a × ⌊2^(2k) / n⌋ / 2^(2k)⌋ 作为 ⌊a / n⌋ 的近似值，然后计算 r = a − q × n。由于近似的误差至多为 1，最终可能只需一次减法修正。Barrett 方法适用于模数固定的场景，常与 Montgomery 方法互为补充。

在实际密码库中，Montgomery 乘法更常用于模幂运算的内层循环（因为连续乘法可以保持在 Montgomery 域中），而 Barrett 约化则更适用于单次取模操作或与其他算法的接口处。两者的选择取决于具体的应用场景与硬件平台特性。

八、素性测试

在密码学中，生成大素数是关键操作。RSA 需要两个大素数 p 和 q；Diffie-Hellman 需要在一个大素数阶群中工作。因此高效的素性测试（primality test）至关重要。

试除法

最朴素的素性测试是试除法（trial division）：依次用 2, 3, 5, …, √n 去除 n。这对小素数筛选足够，但对密码学级别的千位数则完全不可行。实际中通常先用试除法排除小因子（如前几千个素数），再进入更高级的测试。

费马测试

费马素性测试（Fermat primality test）基于费马小定理的逆否命题：若 a^(n−1) ≢ 1 (mod n)，则 n 一定不是素数。随机选取若干基 a 进行测试，如果所有测试都通过，则 n “很可能”是素数。

但费马测试有一个致命弱点：卡迈克尔数（Carmichael number）。这类合数对所有与之互素的基 a 都满足费马条件，因此无论选多少个基都无法用费马测试识别出来。最小的卡迈克尔数是 561 = 3 × 11 × 17。虽然卡迈克尔数相对稀少，但它们有无穷多个——Alford、Granville 与 Pomerance 在 1994 年证明了这一点。因此，费马测试在理论上不够可靠，实际中不推荐单独使用。

Miller-Rabin 测试

Miller-Rabin 素性测试（Miller-Rabin primality test）在费马测试的基础上增加了对”非平凡平方根”的检测，有效解决了卡迈克尔数的问题。

算法思路如下：将 n − 1 写成 2^s × d 的形式（d 为奇数）。对随机选取的基 a，先计算 x = a^d mod n，然后反复平方 s 次。在这个过程中，如果某一步出现 x ≡ −1 (mod n)，则可以提前判定”可能为素数”；如果最终 x ≠ 1，则 n 一定是合数。

每一轮 Miller-Rabin 测试将合数误判为素数的概率不超过 1/4。进行 t 轮独立测试后，误判概率降至 4^(−t)。在密码学实践中，通常取 t 在 20 到 64 之间，此时出错的概率远低于硬件故障的概率。

import random

def miller_rabin(n: int, rounds: int = 40) -> bool:
    """
    Miller-Rabin 素性测试
    返回 True 表示 n 极大概率为素数，False 表示 n 确定为合数
    """
    if n < 2:
        return False
    if n == 2 or n == 3:
        return True
    if n % 2 == 0:
        return False

    # 将 n-1 写成 2^s * d，d 为奇数
    s, d = 0, n - 1
    while d % 2 == 0:
        s += 1
        d //= 2

    for _ in range(rounds):
        a = random.randrange(2, n - 1)
        x = pow(a, d, n)  # a^d mod n

        if x == 1 or x == n - 1:
            continue

        for _ in range(s - 1):
            x = pow(x, 2, n)
            if x == n - 1:
                break
        else:
            # 内层循环正常结束，未遇到 n-1
            return False

    return True

算法中利用了 Python 内置的三参数 pow(a, d, n) 来高效计算模幂，其底层实现使用了快速幂算法。

值得注意的是，当 n 为待测合数时，Miller-Rabin 测试中使 n 被判定为”可能为素数”的基 a 称为强伪证（strong liar）。对于任何奇合数 n，强伪证的个数不超过 (n − 1)/4——这就是每轮误判概率不超过 1/4 的理论保证。而当 n 为卡迈克尔数时，费马测试中的伪证可以多达 n − 2 个，两者的差异极为显著。

在某些确定性应用中，可以选取固定的基集合。例如，对于 64 位以内的整数，使用基集 {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37} 就足以保证完全正确——不存在通过所有这些基的检验的合数。这一结论是通过穷举计算验证的，使得 Miller-Rabin 在特定范围内可以用作确定性测试。

AKS 确定性测试

2002 年，Agrawal、Kayal 与 Saxena 发表了 AKS 素性测试，首次证明素性判定可以在多项式时间内确定性地完成。这一结果在理论上具有里程碑意义，将 PRIMES 问题归入复杂度类 P。然而，AKS 的实际运行时间远慢于 Miller-Rabin，因此在密码学工程中几乎不使用。它的价值更多在理论层面——证明了素性测试不需要依赖随机性假设。

安全素数与强素数

在密码学协议中，通常对素数的选择有额外要求。安全素数（safe prime）是指形如 p = 2q + 1 的素数，其中 q 也是素数。使用安全素数可以确保 Z_p* 中存在一个阶为 q 的大素数阶子群，这对 Diffie-Hellman 类协议的安全性非常重要——它排除了 Pohlig-Hellman 攻击利用小因子分解群阶的可能性。

强素数（strong prime）的定义则更为严格，要求 p − 1 含有大素因子（抵御 Pollard p − 1 算法）、p + 1 含有大素因子（抵御 Williams p + 1 算法）、且 p − 1 的大素因子 r 本身的 r − 1 也含有大素因子。不过，随着 RSA 密钥长度增加到 2048 位以上，使用强素数的额外安全收益在实际中已不那么显著。

关于「随机素数」的常见误解

笔者认为，初学者关于 RSA 密钥生成最常见的误解，是把「随机选取大素数」理解为某种均匀随机的过程。实际上，密码学中的素数生成是一个「随机采样加确定性验证」的迭代过程：先用密码学安全的随机数生成器（CSPRNG）产生一个奇数候选值，用小素数做快速筛选（排除能被前几千个素数整除的候选），再用 Miller-Rabin 多轮测试判断其是否（极大概率）为素数，若不是则递增到下一个候选值继续测试。由素数定理，n 附近的素数密度约为 1/ln(n)，对于 1024 位的数，大约每 710 个奇数中就有一个素数，因此平均只需几百次测试即可找到。

真正的安全关键不在于素数本身的「随机性」，而在于初始候选值的熵源质量——如果随机数生成器有偏差或可预测性，即使最终找到的素数本身完全合法，攻击者也可以大幅缩小搜索空间。2012 年 Lenstra 等人的大规模研究发现，互联网上约 0.2% 的 RSA 公钥共享了一个素因子，根本原因正是密钥生成时的随机数源质量不足（特别是在嵌入式设备和虚拟机首次启动时熵池尚未充分填充的情况下）。这一案例深刻地说明了：密码系统的安全瓶颈往往不在数学算法本身，而在于工程实现中的随机性保障。从工程实践来看，确保高质量的熵源（如硬件随机数生成器 HRNG、操作系统的 /dev/urandom 或 getrandom 系统调用），比选择哪种素性测试算法要重要得多。

九、整数分解

公钥密码（尤其是 RSA）的安全性直接依赖于大整数分解（integer factorization）的困难性。攻击者若能将 RSA 模数 n 分解为 p × q，就可以直接计算私钥 d，从而完全攻破密码系统。

试除法

最朴素的分解方法是试除法：依次用小素数去除 n。其复杂度为 O(√n)，对于 2048 位的 n 而言完全不可行（√n 仍有约 1024 位），但试除法可以快速发现小因子。

Pollard rho 算法

Pollard 的 rho 算法（Pollard’s rho algorithm）利用伪随机序列中的碰撞来找因子。它的期望复杂度为 O(n^(1/4))，对于具有中等大小因子的合数非常高效。算法的核心思想是利用生日悖论：在模 p（n 的某个未知因子）下，一个伪随机序列在大约 O(√p) 步后就会出现碰撞，通过 Floyd 的龟兔算法可以高效检测碰撞。

二次筛法

二次筛法（Quadratic Sieve, QS）由 Pomerance 在 1981 年提出。它的核心思想是找到一组满足 x^2 ≡ y^2 (mod n) 的关系，然后利用 gcd(x − y, n) 来提取因子。通过系统地筛选光滑数（smooth number）并用高斯消元法在 GF(2) 上组合出完全平方，二次筛法达到了亚指数复杂度 L_n(1/2, 1)。在 1990 年代，它是分解百位十进制数的主力算法。

数域筛法

数域筛法（Number Field Sieve, NFS）是目前已知最快的通用整数分解算法，其启发式复杂度为 L_n(1/3, (64/9)^(1/3))，约为 exp(1.923 × (ln n)^(1/3) × (ln ln n)^(2/3))。NFS 的思想是在两个不同的数域中同时寻找光滑关系，然后通过代数关系将它们组合在一起。

NFS 分为多项式选择、筛选、线性代数和平方根四个主要阶段。多项式选择的质量直接影响后续步骤的效率；筛选阶段是最耗时的部分，但高度可并行化；线性代数阶段使用块 Lanczos 或块 Wiedemann 算法处理 GF(2) 上的超大稀疏矩阵。

当前分解记录

截至目前，使用 NFS 分解的最大 RSA 类型整数（即两个大致相同大小的素数之积）为 RSA-250，这是一个 829 位（250 位十进制数）的整数，于 2020 年被分解。这一记录耗费了约 2700 个 CPU 核年的计算量。

对比之下，当前推荐的 RSA 密钥长度为 2048 位（约 617 位十进制数）。从 250 位到 617 位的跨越，在 NFS 的亚指数复杂度下意味着计算量增长数十个数量级——这就是 RSA-2048 在可预见的未来被认为是安全的根本原因。

然而，需要注意的是，Shor 量子算法可以在多项式时间内完成整数分解。一旦大规模容错量子计算机实现，RSA 将不再安全。这也是后量子密码学（post-quantum cryptography）成为当前研究热点的原因。

分解困难性与 RSA 安全

严格来说，RSA 的安全性依赖于 RSA 问题（给定 n 和 e，从 c = m^e mod n 求 m）的困难性，而非直接等价于整数分解。目前尚未证明 RSA 问题与分解问题等价，但也没有发现不分解 n 就能攻破 RSA 的高效算法。在实践中，分解 n 是对 RSA 最有效的已知攻击途径。

本文建立了公钥密码学所需的数论工具箱。从最基本的模运算与扩展欧几里得算法出发，经由群论的抽象语言理解了循环群、原根与离散对数的概念，再深入到中国剩余定理、二次剩余、大整数运算、素性测试与整数分解等核心议题。这些知识为后续理解 RSA、Diffie-Hellman 与椭圆曲线密码学奠定了坚实的数学基础。下一篇将聚焦于公钥密码的开山之作——RSA 密码体制。

密码学百科系列 · 第 14 篇

← 上一篇：密钥派生函数 | 系列目录 | 下一篇：RSA →