9 素性测试

好,那我们走 Stage 9——换换口味,从”纯数学证明”转到”算法与计算”。这一块气质和前八个 Stage 很不一样:它是数论和计算机科学的交叉,是你密码学腹地的核心工具库,而且对你这个 USACO Platinum 选手来说,会有种”回到主场”的熟悉感。

Confidence: 0.9

---

Stage 9:计算数论与素性/分解算法

前八个 Stage 你在”用纸笔证明”,这一阶段你在”设计能跑的算法”。核心问题变成:面对一个几百位的大数,怎么高效地判断它是不是素数、怎么分解它、怎么求模平方根、怎么破解离散对数? 这些是 RSA、椭圆曲线密码学的底层引擎。

五块:① Miller-Rabin 素性测试 → ② Pollard’s rho 分解 → ③ Tonelli-Shanks 模平方根 → ④ Pohlig-Hellman 离散对数 → ⑤ Hensel 提升。每块:解决什么 → 核心思想 → 完整步骤/例子 → 代码 → 陷阱。

---

9.1 Miller-Rabin 素性测试 ★最重要

解决什么

快速判断一个大数(几百位)是不是素数。这是 RSA 生成密钥的第一步——你需要随机生成大素数,就得能快速验证”这个随机数是不是素数”。

为什么不能用试除法

判断 $n$ 是否素数,最朴素的方法是试除到 $\sqrt{n}$。但如果 $n$ 有 300 位,$\sqrt{n}$ 有 150 位,试除要 $10^{150}$ 次——宇宙毁灭都算不完。需要根本不同的方法。

核心思想:费马小定理的加强版

第一层(费马测试):费马小定理说,若 $p$ 素数,则 $a^{p-1}\equiv 1(\mathrm{mod}p)$。所以如果找到某个 $a$ 使 $a^{n-1}\not\equiv 1(\mathrm{mod}n)$,$n$ 铁定是合数。

问题:有些合数也”骗过”费马测试(叫卡迈克尔数,如 561),对所有 $a$ 都满足 $a^{n-1}\equiv 1$。所以费马测试会误判。

第二层(二次探测):Miller-Rabin 加了一道关键检查。把 $n-1$ 写成 $2^s\cdot d$($d$ 奇),然后利用一个事实:

若 $p$ 是素数,则方程 $x^2\equiv 1(\mathrm{mod}p)$ 只有两个解:$x\equiv \pm 1$。

(因为 $x^2-1=(x-1)(x+1)\equiv 0$,素数模下必有 $x\equiv 1$ 或 $x\equiv -1$。)

所以在计算 $a^{n-1}$ 的过程中,如果发现一个”非平凡的 1 的平方根”(既不是 $+1$ 也不是 $-1$,但平方后是 1),那 $n$ 铁定是合数——因为素数不允许这种平方根存在。

完整算法步骤

判断 $n$ 是否素数:

1. 把 $n-1=2^s\cdot d$($d$ 为奇数)
2. 随机选一个底数 $a$($2\le a\le n-2$)
3. 计算 $x=a^d\mathrm{mod}n$
4. 若 $x\equiv 1$ 或 $x\equiv n-1$,通过这一轮(可能是素数),换下一个 $a$
5. 否则,反复平方 $x$ 共 $s-1$ 次:
   • 每次 $x\leftarrow x^2\mathrm{mod}n$
   • 若某次 $x\equiv n-1$,通过这一轮
   • 若平方后变成 1 却没经过 $-1$ → 抓到非平凡平方根,$n$ 是合数!
6. 若 $s-1$ 次平方都没出现 $n-1$ → $n$ 是合数
7. 多换几个 $a$ 重复(测 $k$ 个底数),全通过则极大概率是素数

完整例子:测试 $n=561$(卡迈克尔数,费马测试会骗过)

$n=561$,$n-1=560=2^4\times 35$,所以 $s=4,d=35$。

取 $a=2$:

• $x=2^{35}\mathrm{mod}561$。算出 $2^{35}\equiv 263(\mathrm{mod}561)$
• $263\ne 1$ 且 $263\ne 560$,进入平方循环:
  • $263^2\mathrm{mod}561=166$(不是 560)
  • $166^2\mathrm{mod}561=67$(不是 560)
  • $67^2\mathrm{mod}561=1$ ← 平方变成了 1,但前一步是 67 不是 560!
• 抓到非平凡平方根($67^2\equiv 1$ 但 $67\ne \pm 1$),561 是合数! ✓

(费马测试在这里会被骗:$2^{560}\equiv 1(\mathrm{mod}561)$,但 Miller-Rabin 通过二次探测抓住了它。)

代码

import random
 
def miller_rabin(n, k=20):
    if n < 2: return False
    # 小素数快速排除
    for p in [2,3,5,7,11,13,17,19,23,29,31,37]:
        if n % p == 0:
            return n == p
    # 写 n-1 = 2^s * d
    d, s = n - 1, 0
    while d % 2 == 0:
        d //= 2
        s += 1
    # 测 k 个随机底数
    for _ in range(k):
        a = random.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   # 没找到 n-1,且抓到非平凡平方根 → 合数
    return True            # 所有底数都通过 → 极大概率素数
 
print(miller_rabin(561))            # False (卡迈克尔数被抓)
print(miller_rabin(10**18 + 9))     # True (大素数瞬间判定)
print(miller_rabin(2**61 - 1))      # True (梅森素数)

准确性

Miller-Rabin 是概率算法:测 $k$ 个底数,每个合数被误判为素数的概率 $<4^{-k}$。测 20 个底数,误判概率 $<10^{-12}$,实用上足够。

确定性版本:对 $n<3.3\times 10^{24}$,只需测固定的几个底数(2,3,5,7,11,13,17,19,23,29,31,37)就能确定判素,无需随机。

陷阱:⚠️ 关键是步骤 5 那个 for-else——只有”循环正常结束(没 break)“才返回合数。这个 Python 的 for-else 语法精确对应”没找到 $n-1$“的逻辑,容易写错。

---

9.2 Pollard’s rho 整数分解

解决什么

Miller-Rabin 判断”是不是素数”,但如果是合数,怎么把它分解?Pollard’s rho 负责这个——尤其擅长找出大数的中等大小因子。

核心思想:生日悖论 + 伪随机游走

生日悖论:23 个人里就有较大概率两人同生日——因为”碰撞”比直觉快得多($\sqrt{N}$ 量级)。

应用到分解:设 $n$ 有一个未知质因子 $p$。我们造一串伪随机数 $x_0,x_1,x_2,\dots$(用 $x_{i+1}=x_i^2+c\mathrm{mod}n$)。这串数模 $p$ 会比模 $n$ 更早出现重复(因为 $p<n$,模 $p$ 的”抽屉”更少,约 $\sqrt{p}$ 步就碰撞)。

关键:一旦 $x_i\equiv x_j(\mathrm{mod}p)$ 但 $x_i\ne x_j(\mathrm{mod}n)$,则 $p\mid (x_i-x_j)$ 但 $n\nmid (x_i-x_j)$,所以: $\gcd (x_i-x_j,n)=p\text{ (或 p 的倍数)}$ 捞出了因子 $p$!

名字来历

那串数最终会进入循环(因为模 $n$ 只有有限个值),画出来像希腊字母 ρ (rho) 的形状——一条尾巴接一个环。

Floyd 判圈(龟兔赛跑)

为了检测碰撞又不存储所有值,用 Floyd 判圈:一只”乌龟”走一步、一只”兔子”走两步,它们的差 $\gcd (|x_{\text{兔}}-x_{\text{龟}}|,n)$ 一旦 $>1$ 就找到因子。

完整例子:分解 $n=8051$

用 $f(x)=x^2+1\mathrm{mod}8051$,从 $x=2$ 开始:

步	龟 $x$	兔 $y$	$\gcd (\Vert x-y\Vert ,8051)$
1	5	26	$\gcd (21,8051)=1$
2	26	7474	$\gcd (7448,8051)=1$
3	7474	871	$\gcd (6603,8051)=97$

找到因子 97!$8051=97\times 83$。✓

代码

from math import gcd
 
def pollard_rho(n):
    if n % 2 == 0: return 2
    x = y = 2
    c = 1
    d = 1
    f = lambda v: (v*v + c) % n
    while d == 1:
        x = f(x)           # 乌龟走一步
        y = f(f(y))        # 兔子走两步
        d = gcd(abs(x - y), n)
    return d if d != n else None   # 返回找到的因子
 
print(pollard_rho(8051))      # 97 (或 83)
print(pollard_rho(10403))     # 101 (10403=101×103)

陷阱:⚠️ Pollard’s rho 找到的是某个因子(不一定是素因子),需要递归分解;且偶尔失败(返回 $n$ 本身),要换参数 $c$ 重试。它对有小因子的数很快,但对两个大素数的积(如 RSA 模数)依然慢——这正是 RSA 安全的基础。

---

9.3 Tonelli-Shanks 模平方根

解决什么

求 $x^2\equiv a(\mathrm{mod}p)$ 的解(模平方根)。你 Stage 4 学了 $p\equiv 3(\mathrm{mod}4)$ 的捷径($x\equiv a^{(p+1)/4}$),但 $p\equiv 1(\mathrm{mod}4)$ 时捷径失效,需要 Tonelli-Shanks。

核心思想(通俗版)

把 $p-1=2^s\cdot q$($q$ 奇)。算法分两步:

1. 先得到一个”接近”的候选 $x\approx a^{(q+1)/2}$
2. 反复修正:这个候选的平方和 $a$ 差一个”2 的幂次方根因子”,用一个二次非剩余构造的”修正因子”逐步调整,像调音一样把候选逼近真正的平方根

直觉:在循环群里,误差总是 $2^k$ 阶的元素,每次修正消掉一层 2 的幂,有限步内归零。

完整例子:求 $x^2\equiv 2(\mathrm{mod}17)$

$17\equiv 1(\mathrm{mod}4)$(捷径失效!),$p-1=16=2^4\times 1$,$s=4,q=1$。

(这里 $q=1$ 比较特殊,我给你看小例子的结果:)

• 找一个非剩余,比如 3($3^8\equiv -1(\mathrm{mod}17)$)
• 经过 Tonelli-Shanks 的迭代修正……
• 得到 $x\equiv 6(\mathrm{mod}17)$

验证:$6^2=36\equiv 36-34=2(\mathrm{mod}17)$ ✓。另一根 $-6\equiv 11$。

代码

def tonelli_shanks(a, p):
    # 求 x^2 ≡ a (mod p), 假设解存在
    if pow(a, (p-1)//2, p) != 1:
        return None   # a 不是二次剩余,无解
    if p % 4 == 3:
        return pow(a, (p+1)//4, p)   # 捷径
    # 一般情况
    q, s = p - 1, 0
    while q % 2 == 0:
        q //= 2; s += 1
    # 找一个二次非剩余 z
    z = 2
    while pow(z, (p-1)//2, p) != p - 1:
        z += 1
    m, c, t, r = s, pow(z, q, p), pow(a, q, p), pow(a, (q+1)//2, p)
    while t != 1:
        # 找最小 i 使 t^(2^i)=1
        i, temp = 0, t
        while temp != 1:
            temp = temp*temp % p; i += 1
        b = pow(c, 1 << (m-i-1), p)
        m, c, t, r = i, b*b % p, t*b*b % p, r*b % p
    return r
 
print(tonelli_shanks(2, 17))   # 6 (或 11)
x = tonelli_shanks(2, 17)
print(f"验证: {x}^2 mod 17 = {x*x % 17}")   # 2

用途:椭圆曲线上”给定 x 求 y”(因为 $y^2=x^3+ax+b$,要开平方),所以是 ECC 的核心子程序。

陷阱:⚠️ 用前要先确认 $a$ 是二次剩余(用欧拉判别 $a^{(p-1)/2}\equiv 1$),否则无解。

---

9.4 Pohlig-Hellman 离散对数

解决什么

求离散对数 $g^x\equiv h(\mathrm{mod}p)$ 的 $x$。你 Stage 4 学了 BSGS($O(\sqrt{n})$),但当群的阶 $n$ 能分解成小质因子时,Pohlig-Hellman 快得多。

核心思想:CRT 用在离散对数上

这是 CRT 思想的绝妙应用。如果群的阶 $n=p_1^{e_1}\cdots p_k^{e_k}$:

1. 把”求模 $n$ 的离散对数”这个大问题,拆成在每个 $p_i^{e_i}$ 子群里求小的离散对数
2. 每个小问题规模只有 $p_i^{e_i}$,用 BSGS 快速解
3. 用 CRT 把各个 $x\mathrm{mod}{p}_i^{e_i}$ 拼回完整的 $x$

直觉:大的离散对数难,但如果能”分而治之”成几个小的离散对数,再用 CRT 合并,总难度由最大的质因子决定,而不是 $n$ 本身。

对密码学的致命启示 ★

这个算法揭示了一个密码学的核心安全原则:

如果群的阶有很多小质因子,离散对数就不安全!

所以 Diffie-Hellman、ECC 选参数时,必须让群的阶含一个大质因子(否则 Pohlig-Hellman 直接攻破)。这就是为什么密码学里要用”安全素数”$p=2q+1$($q$ 也是素数)——保证群的阶 $p-1=2q$ 有大质因子 $q$。

完整例子(骨架)

求 $g^x\equiv h$ 在阶 $n=p_1^{e_1}\cdots$ 的群里:

• 对每个 $p_i$:把方程”投影”到 $p_i$ 阶子群(两边取 $n/p_i^{e_i}$ 次幂),解出 $x\mathrm{mod}{p}_i^{e_i}$
• CRT 合并所有 $x\mathrm{mod}{p}_i^{e_i}$ 得 $x\mathrm{mod}n$

代码(简化版)

from sympy.ntheory import discrete_log
from sympy import factorint
 
# sympy 的 discrete_log 内部就用了 Pohlig-Hellman + BSGS
# 求 g^x ≡ h (mod p)
p, g, h = 1019, 2, 5
x = discrete_log(p, h, g)
print(f"离散对数 x = {x}")
print(f"验证: {g}^{x} mod {p} = {pow(g, x, p)}")   # 应等于 h=5
 
# 展示群阶的因子分解(决定 Pohlig-Hellman 的效率)
print(f"群阶 p-1 = {p-1} = {factorint(p-1)}")

陷阱:⚠️ Pohlig-Hellman 只在”群阶光滑(smooth,即只有小质因子)“时快;若群阶有大质因子,它退化成 BSGS 的速度。密码学正是靠”故意让群阶有大质因子”来防御它。

---

9.5 Hensel 提升 (Hensel’s Lemma)

解决什么

解高次同余 $f(x)\equiv 0(\mathrm{mod}{p}^k)$。你解过线性的 $ax\equiv b$,但 $x^2\equiv 2(\mathrm{mod}{7}^3)$ 这种”高次 + 质数幂模”需要 Hensel。

核心思想:牛顿迭代的数论版

“提升”:先在简单的模 $p$ 下找到解,然后像牛顿迭代一样,把这个解一步步”提升”到模 $p^2$、$p^3$……每次提升用导数信息修正。

公式:若 $f(r)\equiv 0(\mathrm{mod}{p}^k)$ 且 $f^{\prime }(r)\not\equiv 0(\mathrm{mod}p)$,则模 $p^{k+1}$ 的解是: $r^{\prime }=r-f(r)\cdot [f^{\prime }(r){]}^{-1}(\mathrm{mod}{p}^{k+1})$ (和牛顿法 $x_{n+1}=x_n-\frac{f(x_n)}{f^{\prime }(x_n)}$ 一模一样,只是在模 $p$ 意义下)

完整例子:解 $x^2\equiv 2(\mathrm{mod}{7}^2=49)$

第一步(模 7):$x^2\equiv 2(\mathrm{mod}7)$,试得 $x\equiv 3$($3^2=9\equiv 2$)或 $x\equiv 4$。取 $r=3$。

第二步(提升到模 49):

• $f(x)=x^2-2$,$f^{\prime }(x)=2x$,$f^{\prime }(3)=6$
• $f(3)=9-2=7$
• $6^{-1}(\mathrm{mod}7)$?这里要模 49 算:$r^{\prime }=3-7\cdot 6^{-1}$。$6^{-1}\equiv 6(\mathrm{mod}7)$(因 $6\times 6=36\equiv 1$)
• $r^{\prime }=3-7\times 6=3-42=-39\equiv 10(\mathrm{mod}49)$

验证:$10^2=100\equiv 100-98=2(\mathrm{mod}49)$ ✓

所以 $x\equiv 10(\mathrm{mod}49)$(另一根 $-10\equiv 39$)。

代码

def hensel_lift(f, df, p, k, r):
    # 把模 p 的解 r 提升到模 p^k
    mod = p
    for _ in range(k - 1):
        mod *= p
        fr = f(r) % mod
        dfr_inv = pow(df(r), -1, mod)   # f'(r) 的逆
        r = (r - fr * dfr_inv) % mod
    return r
 
f  = lambda x: x*x - 2
df = lambda x: 2*x
# 从模 7 的解 3 提升到模 7^2
root = hensel_lift(f, df, 7, 2, 3)
print(f"x^2≡2 mod 49 的解: {root}")           # 10
print(f"验证: {root}^2 mod 49 = {root*root % 49}")   # 2

用途:p-adic 数的基础、解高次同余、某些密码学构造。

陷阱:⚠️ Hensel 提升要求 $f^{\prime }(r)\not\equiv 0(\mathrm{mod}p)$(导数非零),否则提升可能失败或有多个/无解——这对应牛顿法在导数为零处失效。

---

Stage 9 全景速查表

算法	解决什么	核心思想	复杂度
9.1 Miller-Rabin	判素数	费马 + 二次探测	$O(k{\log }^{3}n)$
9.2 Pollard’s rho	分解	生日悖论 + 判圈	$O(n^{1/4})$
9.3 Tonelli-Shanks	模平方根	逐层消 2 的幂	$O({\log }^{2}p)$
9.4 Pohlig-Hellman	离散对数	CRT 分而治之	由最大质因子定
9.5 Hensel	高次同余	牛顿迭代提升	$O(k)$ 次提升

---

我的观点:Stage 9 对你而言是一个特别”对味”的阶段,我想说清它为什么和你高度契合,以及它在你大局中的关键位置。

第一,这一阶段是你的”主场”。 前八个 Stage 是纯数学证明,而 Stage 9 是算法——对你这个 USACO Platinum 选手来说,“设计高效算法”正是你的看家本领。你会发现这里的思维方式很熟悉:Pollard’s rho 用生日悖论(概率 + 复杂度分析)、Floyd 判圈(你 USACO 一定写过快慢指针)、Miller-Rabin 的位运算优化——这些都是竞赛编程的核心套路。前八个 Stage 你是在”学数学家怎么想”,这一阶段你是在用”你本来就擅长的算法思维”解决数论问题。 这种”回到主场”的感觉,会让你学得格外顺。

第二,也是最重要的——这一阶段是你密码学梦想的”工具锻造车间”。 我之前说 Stage 9 严格属于”计算数论”而非初等数论,所以从”主线”角度跳过了它。但从你的密码学兴趣角度,Stage 9 才是真正的核心装备库:Miller-Rabin 是 RSA 生成大素数的引擎、Pollard’s rho 是分解攻击的武器、Tonelli-Shanks 是椭圆曲线”求点”的必备子程序、Pohlig-Hellman 直接揭示了”为什么密码学参数要那样选”。你研究的 ECDLP、你向往的椭圆曲线密码,底层全是这些算法在支撑。 所以你选择走 Stage 9,其实是走向了比 Stage 10(二次互反)更贴近你真实目标的方向——这个选择很对。

第三,Stage 9 让你看到数论的”另一张面孔”。 你前八个 Stage 看到的数论是”优美的、证明的、纸笔的”;Stage 9 让你看到它”实用的、计算的、能跑的”一面。同一个费马小定理——在 Stage 3 它是个要证明的定理,在 Stage 9 (Miller-Rabin) 它变成了判素数的算法引擎;同一个 CRT——在 Stage 3 它是解同余方程组,在 Stage 9 (Pohlig-Hellman) 它变成了攻破离散对数的策略。这种”同一个数学,既能证明又能计算”的双重性,正是计算数论的魅力,也是连接”纯数学”和”密码学应用”的桥梁。你站在这座桥上,一边是你走过的八个 Stage,一边是你向往的椭圆曲线密码。

我特别想让你记住 9.1 和 9.4 这两个:Miller-Rabin 是”最该掌握的算法”(密码学、CP 都用),而 Pohlig-Hellman 揭示的”群阶要有大质因子”原则,是理解一切离散对数密码学安全性的钥匙——你将来学椭圆曲线密码,会反复回到这个原则。

关键提醒 (Caveats):

• 我为每个算法给了可运行的代码,但有些是教学简化版(尤其 Tonelli-Shanks、Pohlig-Hellman 的完整工程实现更复杂),理解核心思想比记代码重要
• Miller-Rabin 是概率算法(测 $k$ 底数误判率 $<4^{-k}$);但对小范围有确定性版本(固定底数),两者别混淆
• Pollard’s rho 找到的是某个因子(可能是合数),需递归分解到素因子;且偶尔失败需换参数 $c$
• 9.4 的代码我用了 sympy 的 discrete_log(内部含 Pohlig-Hellman),没手写完整算法——完整实现较长,理解”CRT 分而治之”的思想是关键
• Tonelli-Shanks 我给的 $x^2\equiv 2(\mathrm{mod}17)$ 例子,中间迭代细节略过了(直接给结果),完整迭代可以对照代码理解
• 这些算法的”安全启示”(如 Pohlig-Hellman → 群阶要有大质因子)对你理解密码学至关重要,比算法本身更值得记
• 本块代码用到 sympy,安装:pip install sympy -i https://pypi.tuna.tsinghua.edu.cn/simple

Stage 9 讲完了——你看到了数论的”算法面孔”:Miller-Rabin 判素数、Pollard’s rho 分解、Tonelli-Shanks 求模平方根、Pohlig-Hellman 攻离散对数、Hensel 提升解高次同余。这是你密码学腹地的核心工具库,也是你算法思维的主场。

下一步你挑:

1. 出一套 Stage 9 练习题(算法实现 + 手算追踪,会偏编程一些,延续 PDF 格式)
2. 进 Stage 10(二次互反完整理论,初等数论主线的最后收尾)
3. 正式开启椭圆曲线/密码学新旅程(你的应许之地——Stage 9 刚锻造的这些工具,正是为它准备的;Miller-Rabin、Tonelli-Shanks 会在那里直接派上用场)