数论基础 — 模运算、CRT、二次剩余、Miller-Rabin

日期: 2026-11-05 方向: 密码学数学 阶段: Phase 4 - 密码学数学基础 (Day 181-194) 标签: #密码学 #数论 #CRT #Miller-Rabin #二次剩余

---

今日目标

---

一、动机

---

二、严格定义

2.1 模运算

定义 2.1：$a \equiv b \pmod n \Leftrightarrow n \mid (a - b)$。

性质：模运算保留 $+, -, \times$；除法仅当 $\gcd(b, n) = 1$ 时合法。

2.2 Fermat 与 Euler 定理

定理 2.2 (Fermat 小)：$p$ 素数，$\gcd(a,p)=1$ → $a^{p-1} \equiv 1 \pmod p$。

定理 2.3 (Euler)：$\gcd(a,n)=1$ → $a^{\varphi(n)} \equiv 1 \pmod n$。

Carmichael 函数 $\lambda(n)$：最小 $k$ 使 $a^k \equiv 1$ 对所有 $\gcd(a,n)=1$。$\lambda(n) \mid \varphi(n)$。

例：$\lambda(15) = \text{lcm}(\lambda(3), \lambda(5)) = \text{lcm}(2, 4) = 4$。

2.3 CRT (中国剩余定理)

定理 2.4 (CRT)：$n_1, \dots, n_k$ 两两互素，$N = \prod n_i$。则 $$\mathbb{Z}/N \cong \mathbb{Z}/n_1 \times \mathbb{Z}/n_2 \times \dots \times \mathbb{Z}/n_k$$

具体：方程组 $x \equiv a_i \pmod{n_i}$ 在 $\mathbb{Z}/N$ 中有唯一解。

构造解：$x = \sum a_i \cdot M_i \cdot M_i^{-1} \mod N$，其中 $M_i = N/n_i$，$M_i^{-1}$ 是 $M_i$ 模 $n_i$ 逆。

2.4 二次剩余 (QR)

定义 2.5：$a$ 是 $\bmod p$ 二次剩余 ($a \in QR_p$) 若存在 $x$ 使 $x^2 \equiv a \pmod p$。

Legendre 符号： $$\left(\frac{a}{p}\right) = \begin{cases} 1 & a \in QR_p, a \ne 0 \ -1 & a \notin QR_p \ 0 & p \mid a \end{cases}$$

Euler 准则：$\left(\frac{a}{p}\right) \equiv a^{(p-1)/2} \pmod p$。

2.5 素性测试

Fermat 测试：对若干 $a$ 测 $a^{n-1} \equiv 1 \pmod n$。素数必过；合数可能过（Carmichael 数 561, 1105, ...）。

Miller-Rabin：$n-1 = 2^s d$，测试 $a^d, a^{2d}, \dots, a^{2^{s-1}d}$ 是否经过 1 或 -1。Carmichael 数也无法躲。

复杂度：单次 $O(\log^3 n)$。$k$ 轮误判率 $\le 4^{-k}$。

---

三、完整推导

3.1 CRT 正确性证明

断言：$x = \sum a_i M_i (M_i^{-1} \mod n_i)$ 满足 $x \equiv a_j \pmod {n_j}$。

证明：

• 对 $i \ne j$：$M_i = N/n_i$ 含因子 $n_j$，故 $M_i \equiv 0 \pmod{n_j}$
• 对 $i = j$：$M_j \cdot (M_j^{-1} \mod n_j) \equiv 1 \pmod{n_j}$

代入：$x \equiv a_j \cdot 1 \equiv a_j \pmod{n_j}$。$\blacksquare$

唯一性：若 $x_1, x_2$ 都满足，则 $n_i \mid (x_1 - x_2)$ for all $i$，故 $N \mid (x_1 - x_2)$。

3.2 Euler 准则证明

断言：$a^{(p-1)/2} \equiv \pm 1 \pmod p$，且 $+1 \Leftrightarrow a \in QR_p$。

证明：

• $(a^{(p-1)/2})^2 = a^{p-1} \equiv 1$（FLT），所以 $a^{(p-1)/2} \equiv \pm 1$
• 若 $a = b^2$：$a^{(p-1)/2} = b^{p-1} \equiv 1$
• 反向：$\mathbb{F}_p^*$ 循环（生成元 $g$），$a = g^k$。$a^{(p-1)/2} = g^{k(p-1)/2}$。由 $\text{ord}(g) = p-1$，此 $\equiv 1 \Leftrightarrow (p-1) \mid k(p-1)/2 \Leftrightarrow k$ 偶数 $\Leftrightarrow a$ 是平方。$\blacksquare$

3.3 Miller-Rabin 正确性

关键事实：$n$ 素数 $\Rightarrow x^2 \equiv 1 \pmod n$ 仅有解 $x = \pm 1$。

合数 $n$ 可能有更多平方根 → Miller-Rabin 检测此异常。

算法：

n-1 = 2^s · d, d 奇
随机 a ∈ [2, n-2]
计算 x = a^d mod n
if x == 1 or x == n-1: return PROBABLY_PRIME
for i in 1..s-1:
  x = x² mod n
  if x == n-1: return PROBABLY_PRIME
return COMPOSITE

误判率：单次 ≤ 1/4（Rabin 1980）。常用 40 轮 → ≤ $2^{-80}$。

---

四、代码实现

"""
Day 188: 数论工具 — CRT, Miller-Rabin, Tonelli-Shanks
"""
from random import randrange
from math import gcd


# ============ 扩展欧几里得 + 模逆 ============
def egcd(a: int, b: int) -> tuple[int, int, int]:
    if b == 0:
        return a, 1, 0
    g, x1, y1 = egcd(b, a % b)
    return g, y1, x1 - (a // b) * y1


def modinv(a: int, n: int) -> int:
    g, x, _ = egcd(a % n, n)
    if g != 1:
        raise ValueError(f"{a} 与 {n} 不互素")
    return x % n


# ============ CRT ============
def crt(remainders: list[int], moduli: list[int]) -> int:
    """
    解 x ≡ a_i (mod n_i), n_i 两两互素
    返回 [0, N) 中唯一解, N = prod(n_i)
    """
    assert len(remainders) == len(moduli)
    N = 1
    for n in moduli: N *= n
    x = 0
    for a, n in zip(remainders, moduli):
        M = N // n
        x += a * M * modinv(M, n)
    return x % N


# ============ Miller-Rabin ============
def is_probable_prime(n: int, k: int = 40) -> bool:
    """Miller-Rabin: n 素数概率 ≥ 1 - 4^(-k)"""
    if n < 2: return False
    if n in (2, 3): return True
    if n % 2 == 0: return False

    # n - 1 = 2^s · d, d 奇
    d, s = n - 1, 0
    while d % 2 == 0:
        d //= 2
        s += 1

    for _ in range(k):
        a = randrange(2, n - 1)
        x = pow(a, d, n)
        if x == 1 or x == n - 1:
            continue
        for _ in range(s - 1):
            x = (x * x) % n
            if x == n - 1:
                break
        else:
            return False
    return True


def gen_prime(bits: int) -> int:
    """生成 bits 位概率素数"""
    while True:
        n = randrange(2 ** (bits - 1), 2 ** bits) | 1  # 奇
        if is_probable_prime(n):
            return n


# ============ Legendre 符号 + Tonelli-Shanks ============
def legendre(a: int, p: int) -> int:
    """Legendre 符号: 1 if QR, -1 if NQR, 0 if p|a"""
    if a % p == 0: return 0
    ls = pow(a, (p - 1) // 2, p)
    return -1 if ls == p - 1 else ls


def tonelli_shanks(n: int, p: int) -> int:
    """求 x² ≡ n (mod p), p 奇素数, n ∈ QR_p"""
    if legendre(n, p) != 1:
        raise ValueError("n 不是 QR")
    if p % 4 == 3:
        return pow(n, (p + 1) // 4, p)  # 简单情形

    # 一般情形
    q, s = p - 1, 0
    while q % 2 == 0:
        q //= 2
        s += 1

    # 找一个 NQR z
    z = 2
    while legendre(z, p) != -1:
        z += 1

    m = s
    c = pow(z, q, p)
    t = pow(n, q, p)
    r = pow(n, (q + 1) // 2, p)

    while True:
        if t == 1:
            return r
        # 找最小 i ∈ [0, m), t^(2^i) = 1
        t2i = t
        i = 0
        for j in range(1, m):
            t2i = (t2i * t2i) % p
            if t2i == 1:
                i = j
                break

        b = pow(c, 1 << (m - i - 1), p)
        m = i
        c = (b * b) % p
        t = (t * c) % p
        r = (r * b) % p


# ============ RSA + CRT 加速 ============
def rsa_decrypt_crt(c: int, p: int, q: int, d: int, n: int) -> int:
    """
    用 CRT 加速 RSA 解密 m = c^d mod n
    比直接 pow(c, d, n) 快 ~4×
    """
    dp = d % (p - 1)
    dq = d % (q - 1)
    qinv = modinv(q, p)

    mp = pow(c, dp, p)
    mq = pow(c, dq, q)

    h = (qinv * (mp - mq)) % p
    return (mq + h * q) % n


# ============ 演示 ============
if __name__ == "__main__":
    # 1. CRT
    rems, mods = [2, 3, 2], [3, 5, 7]
    x = crt(rems, mods)
    print(f"x ≡ 2 mod 3, ≡ 3 mod 5, ≡ 2 mod 7 → x = {x}")  # 23
    for r, m in zip(rems, mods):
        assert x % m == r

    # 2. Miller-Rabin
    print(f"\n997 素数? {is_probable_prime(997)}")
    print(f"561 (Carmichael 数) 素数? {is_probable_prime(561)}")  # False!
    print(f"1024-bit 素数: {gen_prime(1024).bit_length()} bit")

    # 3. Legendre + Tonelli
    p = 13
    for a in range(1, p):
        ls = legendre(a, p)
        if ls == 1:
            r = tonelli_shanks(a, p)
            assert (r * r) % p == a
            print(f"  √{a} mod {p} = {r}")

    # 4. RSA-CRT 加速演示
    p, q = gen_prime(512), gen_prime(512)
    n = p * q
    e = 65537
    phi = (p - 1) * (q - 1)
    d = modinv(e, phi)

    m = 12345
    c = pow(m, e, n)
    m_decoded = rsa_decrypt_crt(c, p, q, d, n)
    assert m == m_decoded
    print(f"\nRSA-CRT 解密验证通过, n = {n.bit_length()} bit")

输出：

x ≡ 2 mod 3, ≡ 3 mod 5, ≡ 2 mod 7 → x = 23
997 素数? True
561 (Carmichael 数) 素数? False
1024-bit 素数: 1024 bit
  √1 mod 13 = 1
  √3 mod 13 = 4
  √4 mod 13 = 2
  ...
RSA-CRT 解密验证通过, n = 1024 bit

---

五、密码学应用关联

RSA-CRT 加速分析

直接 $c^d \mod n$（$n$ 2048 bit, $d$ 2048 bit）：约 $3000$ 次模乘 in $\mathbb{Z}/n$。

CRT：分别计算 $c^{d_p} \mod p, c^{d_q} \mod q$（$d_p, d_q$ 1024 bit, $p, q$ 1024 bit）：每个约 $1500$ 模乘 in $\mathbb{Z}/p$。但 $\mathbb{Z}/p$ 模乘是 $\mathbb{Z}/n$ 的 1/4 成本。总：$2 \times 1500 \times 1/4 = 750$，4× 加速。

陷阱：CRT 实施漏算导致 Lenstra 1996 fault attack——故障注入让 $\mod p$ 算错，从 $\gcd(c^d - m, n) = p$ 直接因子分解。修复：CRT 后做完整验证 $m^e \stackrel{?}{=} c$。

---

六、真实参数

# RSA-2048 典型参数
N = 2048-bit composite (p × q, 各 1024 bit)
e = 65537 = 2^16 + 1
d = e^(-1) mod φ(N)
公钥 (N, e), 私钥 (p, q, d, dp, dq, qinv)

# EC 点压缩
secp256k1 压缩点: 0x02 || x  (y 偶) 或 0x03 || x  (y 奇)
解压: 算 x³+7 mod p, 用 Tonelli-Shanks 求 y

著名 Carmichael 数

• 561 = 3 × 11 × 17（最小）
• 1105, 1729 (Ramanujan), 2465, 2821, 6601

Carmichael 满足 $a^{n-1} \equiv 1 \pmod n$ 对所有 $\gcd(a,n)=1$，Fermat 测试无效，但 Miller-Rabin 仍能识别。

---

七、常见陷阱

1. "Fermat 测试足够"：不够。Carmichael 数躲过 Fermat。RSA 密钥生成必须 Miller-Rabin。
2. CRT-RSA fault attack：Lenstra 攻击。需做 CRT 重组后验证签名 / 解密结果。
3. Tonelli-Shanks 在 $p \equiv 3 \pmod 4$ 时简化：直接 $y = a^{(p+1)/4}$。secp256k1 满足此条件。$p \equiv 1 \pmod 4$ 需完整算法。
4. Legendre vs Jacobi：Legendre 仅对素 $p$；Jacobi 推广到合数 $n$，但 Jacobi = 1 不保证 QR。
5. 概率素数 vs 严格素数：$2^{-80}$ 误判率对密码学足够，但理论可证明素数（AKS, ECPP）成本高得多。
6. 小素数试除先做：Miller-Rabin 前用前 1000 素数试除可加速 10×（合数大概率含小因子）。

---

八、关键速查

公式

Python 函数

egcd(a, b) → (g, x, y)
modinv(a, n) → int
crt(rems, mods) → int
is_probable_prime(n, k=40) → bool
gen_prime(bits) → int
legendre(a, p) → int
tonelli_shanks(n, p) → int
rsa_decrypt_crt(c, p, q, d, n) → int

---

九、面试题

Q1: 解释 Carmichael 数为什么让 Fermat 测试失效。

答：Carmichael 数 $n$ 是合数但 $a^{n-1} \equiv 1 \pmod n$ 对所有 $\gcd(a,n)=1$。这是因为 $\lambda(n) \mid n-1$（Korselt 准则）。Fermat 测试从此错判为素。Miller-Rabin 进一步检查 $a^d, a^{2d}, \dots$ 序列中是否有非平凡 $\sqrt 1$（合数特征），即使 Carmichael 也会暴露。

Q2: RSA-CRT 加速比是多少？为什么？

答：约 4×。

• 直接 $c^d \mod n$ ($d, n$ 都 2048 bit)：模乘成本 $O(2048^2)$
• CRT 分两半在 $\mathbb{Z}/p, \mathbb{Z}/q$（1024 bit）算，模乘成本 $O(1024^2)$，是原来 1/4
• 算两次但每次只做一半指数（$d_p, d_q$ 都 1024 bit）
• 总计：$2 \times \text{half exp} \times \text{quarter cost} = 2 \times \frac{1}{2} \times \frac{1}{4} = \frac{1}{4}$ → 4×

Q3: 写出 Tonelli-Shanks 在 $p \equiv 3 \pmod 4$ 的简化形式。

答：直接 $y = a^{(p+1)/4} \mod p$。 证明：$y^2 = a^{(p+1)/2} = a^{(p-1)/2} \cdot a$。由 $a$ 是 QR，$a^{(p-1)/2} \equiv 1$，所以 $y^2 \equiv a$。$\blacksquare$

secp256k1 的 $p \equiv 3 \pmod 4$，所以点解压非常快。

Q4: CRT 在 Pohlig-Hellman 攻击中如何用？

答：群阶 $n = \prod p_i^{e_i}$。Pohlig-Hellman:

1. 在每个 $p_i^{e_i}$ 阶子群中解 $\log_{g} h \mod p_i^{e_i}$（用 BSGS / rho）
2. 得 $x \mod p_i^{e_i}$ for each $i$
3. 用 CRT 合成 $x \mod n$

复杂度 $\sum \sqrt{p_i}$ 远小于 $\sqrt n$。这就是为什么 secp256k1 必须选 $n$ 为大素数（无小因子可分）。

Q5: $\varphi(N) = (p-1)(q-1)$ 怎么用 Euler 推出？

答：$N = pq$，$p, q$ 互素。$\varphi$ 是积性函数（multiplicative），即 $\gcd(m,n)=1 \Rightarrow \varphi(mn) = \varphi(m)\varphi(n)$。$\varphi(p) = p-1$（除 $p$ 外都互素）。所以 $\varphi(N) = (p-1)(q-1)$。

RSA 安全性：知道 $\varphi(N)$ 等价于知道 $p+q$（联立 $pq=N, (p-1)(q-1)=\varphi$），等价于因子分解。

---

十、明日预告

Day 189: 多项式与 FFT — 多项式表示（系数/点值）、Lagrange 插值、FFT 复杂度推导（$O(n\log n)$ 蝴蝶变换）。这是 KZG / Plonk / FRI / STARK 的核心数据结构。

核心问题：怎么在 NTT 域里做「多项式承诺」让验证常数时间？

参考文献：

• Cormen et al, Introduction to Algorithms, Ch. 31 (Number Theory) + Ch. 30 (FFT).
• Galbraith Ch. 11-12.
• Rabin, "Probabilistic Algorithm for Testing Primality", J. Number Theory 1980.
• Lenstra, "Memo on RSA signature generation in the presence of faults", 1996.