Skip to content
目录
(Updated: )•7 分钟阅读•0次浏览

数论变换（NTT）入门：数学公式、手算示例与 Python 实现

数学密码学#NTT#后量子密码学#数学#FFT
The butterfly data flow in FFT

**数论变换（Number Theoretic Transform，NTT）**是在有限域上进行的离散傅里叶变换。它将多项式系数变为单位根上的点值，通过逐点相乘和逆变换完成多项式乘法，所有结果都在选定模数下精确计算。

本文从模 17 的手算示例出发，解释正逆变换公式，并提供可直接运行的 Python 代码。理解基本的多项式乘法和模运算即可阅读。

DFT 和 FFT

DFT 是离散傅里叶变换，FFT 是快速计算 DFT 的算法。本文使用同类的蝶形算法计算 NTT，时间复杂度为 O(nlog⁡n)O(n \log n)。

FFT 在密码学中的问题

复数 FFT 的浮点实现需要控制舍入误差，才能恢复整数系数。NTT 使用有限域模运算，避免了这一步舍入，但需要选择合适的变换长度、素数模数和单位根。FFT 在误差界得到控制时也能用于精确的整数乘法；NTT 本身则不提供加密或安全保证。

NTT：基于整数的解决方案

数论变换（NTT）通过将复数的世界替换为模运算的世界来解决精度问题。NTT 不在无限的复数域上运算，而是在模素数 pp 的有限域整数上运算。

这样做有两个关键优势：

  1. 无浮点误差：所有计算都使用整数，确保完美精度。
  2. 高效性：计算速度极快，继承了 FFT 的 O(nlog⁡n)O(n \log n) 时间复杂度。

实现这一切的”魔法材料”是找到 FFT 中”单位原根”在整数域中的等价物。在 NTT 中，我们需要在有限域中找到一个特殊的整数 ω\omega（omega），它的行为方式与复数域中的对应物完全一致，从而可以执行正变换和逆变换。

卷积定理

长度为 nn 的 NTT 将循环卷积转换为逐点乘法：

NTT⁡(a∗nb)=NTT⁡(a)⊙NTT⁡(b).\operatorname{NTT}(a *_n b) = \operatorname{NTT}(a) \odot \operatorname{NTT}(b).

这里 ∗n*_n 表示模 xn−1x^n-1 的循环卷积，⊙\odot 表示逐点相乘，系数运算均模 pp。

逐步演算示例

我们计算 A(x)=2x+1A(x)=2x+1 和 B(x)=3x+4B(x)=3x+4 的乘积。

线性卷积与循环卷积

长度分别为 rr、ss 的系数向量，其普通乘积最多有 r+s−1r+s-1 个系数。本文的 radix-2 算法选择不小于该值的 2 的幂，并将输入补零，以免高次项在模 xn−1x^n-1 时回绕。这样得到线性卷积，但系数仍模 pp。若要恢复任意整数系数，需要足够大的模数或结合多个模数使用中国剩余定理；对于绝对值不超过 BB 的有符号系数，可使用大于 2B2B 的模数恢复。模 xn+1x^n+1 的负循环卷积需要不同构造，不由本文代码直接实现。

素数模数还须满足 n∣(p−1)n\mid(p-1)，单位根必须恰好具有阶 nn。对于 2 的幂 n>1n>1，需同时验证 ωn≡1\omega^n\equiv1 和 ωn/2≢1(modp)\omega^{n/2}\not\equiv1\pmod p，不能只验证前者。如果 gg 是乘法群的生成元，可以取 ω=g(p−1)/n mod p\omega=g^{(p-1)/n}\bmod p。

1. 参数设定

我们需要选择以下参数：

  • 模数（pp）：17。我们将在模 17 的整数域中运算。
  • 变换长度（nn）：4。两个 1 次多项式的乘积是 2 次多项式。为避免混叠，nn 必须是大于 2 的 2 的幂。我们将系数向量用零填充到长度 4。
  • 单位根（ω\omega）：13。我们需要一个元素 ω\omega 使得 ω4≡1(mod17)\omega^4 \equiv 1 \pmod{17}。131=1313^1 = 13，132≡16≡−113^2 \equiv 16 \equiv -1，133≡413^3 \equiv 4，134≡113^4 \equiv 1。
  • nn 的逆元（n−1n^{-1}）：13。因为 4×13=52=3×17+14 \times 13 = 52 = 3 \times 17 + 1，所以 4 模 17 的模逆元是 13。

输入向量为：

  • a=[1,2,0,0]a = [1, 2, 0, 0]（表示 1+2x1 + 2x）
  • b=[4,3,0,0]b = [4, 3, 0, 0]（表示 4+3x4 + 3x）

2. 正向 NTT

我们用公式 Ak=∑j=0n−1ajωjk(mod17)A_k = \sum_{j=0}^{n-1} a_j \omega^{jk} \pmod{17} 将向量 aa 变换为点值形式 AA。

  • A0=1+2(13)0+0+0=3A_0 = 1 + 2(13)^0 + 0 + 0 = 3

  • A1=1+2(13)1=1+26=27≡10A_1 = 1 + 2(13)^1 = 1 + 26 = 27 \equiv \mathbf{10}

  • A2=1+2(13)2=1+2(16)=33≡16A_2 = 1 + 2(13)^2 = 1 + 2(16) = 33 \equiv \mathbf{16}

  • A3=1+2(13)3=1+2(4)=9≡9A_3 = 1 + 2(13)^3 = 1 + 2(4) = 9 \equiv \mathbf{9}

    A=[3,10,16,9]A = [3, 10, 16, 9]

同理，变换向量 bb：

  • B0=4+3(13)0=7B_0 = 4 + 3(13)^0 = 7

  • B1=4+3(13)1=4+39=43≡9B_1 = 4 + 3(13)^1 = 4 + 39 = 43 \equiv \mathbf{9}

  • B2=4+3(13)2=4+3(16)=52≡1B_2 = 4 + 3(13)^2 = 4 + 3(16) = 52 \equiv \mathbf{1}

  • B3=4+3(13)3=4+3(4)=16≡16B_3 = 4 + 3(13)^3 = 4 + 3(4) = 16 \equiv \mathbf{16}

    B=[7,9,1,16]B = [7, 9, 1, 16]

3. 逐点相乘

我们将变换后的向量逐元素相乘，模 17，得到 CC。

  • C0=3×7=21≡4C_0 = 3 \times 7 = 21 \equiv \mathbf{4}

  • C1=10×9=90≡5C_1 = 10 \times 9 = 90 \equiv \mathbf{5}

  • C2=16×1=16C_2 = 16 \times 1 = \mathbf{16}

  • C3=9×16=144≡8C_3 = 9 \times 16 = 144 \equiv \mathbf{8}

    C=[4,5,16,8]C = [4, 5, 16, 8]

4. 逆 NTT（INTT）

现在我们用逆变换公式将 CC 转回系数。我们使用 ω−1=4\omega^{-1} = 4（因为 13×4=52≡113 \times 4 = 52 \equiv 1），并将最终结果乘以 n−1=13n^{-1} = 13。

首先，计算求和 Sj=∑k=0n−1Ck(ω−1)jk(mod17)S_j = \sum_{k=0}^{n-1} C_k (\omega^{-1})^{jk} \pmod{17}：

  • S0=4(1)+5(1)+16(1)+8(1)=33≡16S_0 = 4(1) + 5(1) + 16(1) + 8(1) = 33 \equiv 16
  • S1=4(1)+5(4)+16(16)+8(13)≡4+3+1+2=10S_1 = 4(1) + 5(4) + 16(16) + 8(13) \equiv 4 + 3 + 1 + 2 = 10
  • S2=4(1)+5(16)+16(256)+8(169)≡4+12+16+9=41≡7S_2 = 4(1) + 5(16) + 16(256) + 8(169) \equiv 4 + 12 + 16 + 9 = 41 \equiv 7
  • S3=4(1)+5(13)+16(169)+8(2197)≡4+14+1+15=34≡0S_3 = 4(1) + 5(13) + 16(169) + 8(2197) \equiv 4 + 14 + 1 + 15 = 34 \equiv 0

最后，乘以 n−1=13n^{-1} = 13：

  • c0=16×13=208≡4c_0 = 16 \times 13 = 208 \equiv \mathbf{4}
  • c1=10×13=130≡11c_1 = 10 \times 13 = 130 \equiv \mathbf{11}
  • c2=7×13=91≡6c_2 = 7 \times 13 = 91 \equiv \mathbf{6}
  • c3=0×13=0c_3 = 0 \times 13 = \mathbf{0}

最终结果为 c=[4,11,6,0]c = [4, 11, 6, 0]，对应多项式 6x2+11x+46x^2 + 11x + 4。

与直接展开一致：(2x+1)(3x+4)=6x2+11x+4(2x+1)(3x+4) = 6x^2 + 11x + 4。

NTT 与 INTT 算法

如果你想进一步了解，让我们来看看正变换和逆变换的细节。NTT 将多项式的系数向量 aa 变换为点值表示 AA，INTT 则执行逆过程。

正向 NTT 定义为： Ak=∑j=0n−1ajωjk(modN)A_k = \sum_{j=0}^{n-1} a_j \omega^{jk} \pmod{N}

逆 NTT 定义为： aj=(n−1)∑k=0n−1Akω−jk(modN)a_j = (n^{-1}) \sum_{k=0}^{n-1} A_k \omega^{-jk} \pmod{N}

其中：

  • aa 是输入系数向量。
  • AA 是变换后的输出向量。
  • NN 是素数模数。
  • nn 是变换长度（本文的 radix-2 实现要求为 2 的幂）。
  • ω\omega 是模 NN 的 nn 次单位原根。
  • n−1n^{-1} 是 nn 模 NN 的模乘法逆元。

为什么逆变换要乘以 n−1n^{-1}？

单位原根满足：∑k=0n−1ω(j−ℓ)k\sum_{k=0}^{n-1}\omega^{(j-\ell)k} 在 j=ℓj=\ell 时为 nn，其余情况为零。使用逆单位根求和后，每个原始系数都被放大了 nn 倍，因此还须乘以 nn 的模逆元。示例中 4⋅13≡1(mod17)4\cdot13\equiv1\pmod{17}，所以乘以 13。

以下是一个 NTT 的 Python 参考实现，其结构与 Cooley-Tukey FFT 算法类似。

def ntt(a, N, w, n):
    """Radix-2 NTT；调用者须保证 N 为素数。"""
    if n < 1 or n & (n - 1) or len(a) != n:
        raise ValueError("Use a power-of-two length matching the input")
    if N < 2 or (N - 1) % n:
        raise ValueError("The prime modulus must satisfy n | (N - 1)")
    if pow(w, n, N) != 1 or (n > 1 and pow(w, n // 2, N) == 1):
        raise ValueError("The root must have exact order n")
    # 1. 比特反转置换
    A = [value % N for value in a]
    for i in range(n):
        rev_i = int(f'{i:0{n.bit_length()-1}b}'[::-1], 2)
        if i < rev_i:
            A[i], A[rev_i] = A[rev_i], A[i]
 
    # 2. 蝶形运算
    s = 1
    while (m := 2**s) <= n:
        w_m = pow(w, n // m, N)
        for k in range(0, n, m):
            w_pow = 1
            for j in range(m // 2):
                u = A[k + j]
                t = (A[k + j + m // 2] * w_pow) % N
                A[k + j] = (u + t) % N
                A[k + j + m // 2] = (u - t + N) % N
                w_pow = (w_pow * w_m) % N
        s += 1
    return A

以下是逆 NTT 的 Python 代码。注意它几乎完全相同，只是使用了逆根并对最终结果进行了缩放。

def intt(A, N, w_inv, n):
    """逆数论变换"""
    # INTT 算法与 NTT 相同，只需使用逆根
    a = ntt(A, N, w_inv, n)
 
    # 用 n_inv 缩放结果
    n_inv = pow(n, N - 2, N)
    for i in range(n):
        a[i] = (a[i] * n_inv) % N
    return a

比特反转步骤是一个预处理置换，它将输入数组重新排列为后续迭代蝶形运算所需的正确顺序。正是这种重排使得算法如此高效。

运行完整 Python 示例

使用 Python 3.8 或更新版本，先复制上面的两个函数，再运行以下代码。系数按常数项到高次项排列。这是教学实现，没有实现密码学用途所需的常数时间运算。

p, n, w = 17, 4, 13
a = [1, 2, 0, 0]
b = [4, 3, 0, 0]
A = ntt(a, p, w, n)
B = ntt(b, p, w, n)
C = [x * y % p for x, y in zip(A, B)]
c = intt(C, p, pow(w, -1, p), n)
 
assert intt(A, p, pow(w, -1, p), n) == a
assert c == [4, 11, 6, 0]
print("NTT(a):", A)
print("NTT(b):", B)
print("Product coefficients:", c)
NTT(a): [3, 10, 16, 9]
NTT(b): [7, 9, 1, 16]
Product coefficients: [4, 11, 6, 0]

标准 NTT 之外：变体与近亲

虽然标准数论变换是现代密码学的主力，但它属于一个更大的变换家族，每个成员都针对特定的约束或数学结构进行了优化。

  • 离散傅里叶变换（DFT）：所有变换的鼻祖。NTT 在有限域的整数上运算，而 DFT 在复数域上运算。快速傅里叶变换（FFT）只是计算 DFT 的一种高效算法。
  • 离散余弦变换（DCT）：DFT 的近亲，只使用实数和余弦函数。它以信号处理中的”能量压缩”特性闻名（支撑了 JPEG 和 MP3），但其数论变体也存在于专门的整数算术应用中。
  • 离散伽罗瓦变换（DGT）：一种在有限域的伽罗瓦扩张（如模 pp 的高斯整数）上运算的专用变体。它有效地将更多信息打包到每个元素中，有可能将给定多项式次数所需的变换长度减半。这在优化 Kyber 等后量子方案时尤为重要。
  • 离散加权变换（DWT）：一种对输入和输出向量施加加权因子的变体。常用于执行”负循环”卷积（在格密码中很常见），或者在不需要标准循环卷积那种完整零填充开销的情况下计算线性卷积。
  • 沃尔什-哈达玛变换（WHT）：一种仅使用加法和减法（无需乘法）的简化变换。虽然对于一般卷积不如 NTT/FFT 强大，但其速度使它在量子算法、纠错码和特定的密码布尔函数中很有价值。
  • 无理基离散加权变换（IBDWT）：DWT 的一种变体，用于涉及超大数的计算（如互联网梅森素数大搜索）。它允许使用浮点硬件进行高效乘法，同时通过精确的误差分析保持精确的整数结果。

快速对比

变换域主要用途核心优势
DFT/FFT复数信号处理、音频/图像成熟的硬件支持，直观的频率分析
NTT有限域（整数）密码学、大整数算术精确计算（无舍入误差），高效的模运算
DCT实数压缩（JPEG/MP3）能量压缩，实值运算（无需复数存储）
DGT伽罗瓦扩张高级后量子密码（如 Kyber）缩短变换长度，更高的信息密度
DWT加权环格密码、超大数乘法高效的负循环卷积，减少填充开销
WHT二值/整数量子算法、纠错码极速（无需乘法运算）

为什么 NTT 是密码学的颠覆者？

NTT 在密码学中的主要应用是执行快速多项式乘法。许多现代密码方案，特别是基于格的方案，其构建依赖于涉及大型多项式的运算。

如果用朴素方法进行这些乘法运算，速度慢到不切实际。NTT 提供了巨大的加速，使这些先进的密码方案变得可行。

NIST 的 ML-KEM 标准（FIPS 203）源自 CRYSTALS-Kyber，规定了密钥封装机制；ML-DSA 标准（FIPS 204）源自 CRYSTALS-Dilithium，规定了数字签名方案。它们使用各自特定的 NTT 构造与参数，本文的教学代码不能替代标准实现。

相关代数基础可参阅密码学公式指南。

结语

数论变换是现代密码学的一块基石，它提供了构建未来安全系统所需的速度和精度。虽然它在幕后默默运作，但这个优雅的数学工具使得创建不仅安全于当今计算机、也能抵御未来量子威胁的密码系统成为可能。这是一个绝佳的例子，展示了抽象数学如何在保障我们数字世界的安全中找到了强大而实用的应用。