前置知识

  • FFT
  • 阶与原根

介绍

数论变换（Number-Theoretic Transform, NTT）是离散傅里叶变换（DFT）在数论基础上的实现。它利用模算术替代了复数运算，从而避免了 DFT 的精度误差问题，并且在整数运算环境下通常速度更快。

实际上，NTT 指的是数论变换，而 FNTT 指的是快速数论变换（Fast Number-Theoretic Transform, FNTT）。但竞赛圈中的 NTT 一般指的是 FNTT。注意本文会严格区分。

在算法竞赛中，当题目要求在特定模数下计算多项式卷积时，我们就需要使用 NTT。

FNTT 的算法结构与 FFT 完全一致，但其数学基础从复数域转移到了有限域（模算术）。

NTT 的使用受限于模数的选择。一个常用的 NTT 模数是 998244353，我们在文末会解释其特殊性。

阶与原根

在讲解 NTT 之前，需要理解原根的性质，因为它将扮演 DFT 中单位根的角色。

阶

若正整数 \(a,p\) 互素且 \(p > 1\)，则满足 \(a^n \equiv 1 \pmod p\) 的最小正整数 \(n\)，称为 \(a\) 模 \(p\) 的阶，记作 \(\text{ord}_p(a)\)。

性质：对于 \(i \in [0, \text{ord}_p(a) - 1]\)，所有的 \(a^i \pmod p\) 的结果互不相同。

证明

假设存在两个整数 \(j, k\) 满足 \(0 \le j < k < \text{ord}_p(a)\) 且 \(a^j \equiv a^k \pmod p\)。 因为 \(a, p\) 互素，所以 \(a\) 存在模 \(p\) 的逆元。两边同乘以 \((a^{-1})^j\) 可得 \(a^{k-j} \equiv 1 \pmod p\)。 因为 \(0 < k-j < \text{ord}_p(a)\)，这与阶的最小性定义矛盾。

原根

设 \(p\) 是一个素数，整数 \(g\) 是模 \(p\) 的原根，当且仅当 \(g\) 模 \(p\) 的阶为 \(\varphi(p) = p-1\)。

性质：若 \(g\) 是模 \(p\) 的原根，则它的幂次 \(g^0, g^1, \dots, g^{p-2}\) 在模 \(p\) 意义下构成了模 \(p\) 的既约剩余系，即遍历了 \(1, 2, \dots, p-1\) 的所有值。

数论变换 (NTT)

DFT 的核心是利用单位根 \(\omega_n\) 的性质。在 NTT 中，我们将寻找一个在模意义下具有类似性质的替代品——这就是原根的作用。

设 \(g\) 是模 \(p\) 的一个原根。我们定义 NTT 中的 主 \(n\) 次单位根 为： \[ g_n = g^{(p-1)/n} \pmod p \]

这个 \(g_n\) 具有和 \(\omega_n\) 非常相似的性质：

  1. 周期性 (消去引理) \[(g_n)^n = g^{(p-1)/n \cdot n} = g^{p-1} \equiv 1 \pmod p \quad \text{(根据费马小定理)}\]

  2. 折半引理 \[g_{2n}^{2k} = g^{2k(p-1)/2n} = g^{k(p-1)/n} = g_n^k \pmod p\]

  3. 对称性 \[(g_n)^{n/2} = g^{(p-1)/2} \equiv -1 \pmod p \quad\]

    因为 \((g^{(p-1)/2})^2 \equiv 1 \pmod p\ \text{且} \ g^{(p-1)/2} \not\equiv 1 \pmod p\)

    因此，可以推导出最重要的蝴蝶变换性质： \[ \begin{align*} (g_n)^{k+n/2} &= (g_n)^k \cdot (g_n)^{n/2} \\ &\equiv - (g_n)^k \pmod p \end{align*} \]

我们发现，\(g_n\) 完美地复刻了单位根 \(\omega_n\) 在 DFT 中所需的全部性质。因此，我们可以将 DFT 算法中的 \(\omega_n\) 替换为 \(g_n\)，复数运算替换为模 \(p\) 意义下的整数运算，从而得到 NTT。

快速数论变换 (FNTT)

FNTT (Fast Number-Theoretic Transform) 就是 NTT 的快速实现，其算法结构与 FFT 完全一致，同样是基于分治思想和蝴蝶变换，时间复杂度为 \(O(n \log n)\)。

快速逆数论变换 (IFNTT)

IFNTT (Inverse Fast Number-Theoretic Transform) 是 FNTT 的逆变换，用于将点值表示法还原为系数表示法。

与 IFFT 类似，IFNTT 的实现只需要在 FNTT 的基础上做两处修改：

  1. 将变换过程中使用的 “单位根” \(g_n\) 替换为其在模 \(p\) 意义下的逆元 \((g_n)^{-1}\)。
  2. 变换结束后，将结果数组的每个元素都乘以 \(n\) 在模 \(p\) 意义下的逆元 \(n^{-1}\)。

模数的性质

因为 FNTT 的实现依赖分治算法，要求多项式的长度必须为 \(2\) 的幂次。而 主 \(n\) 次单位根 \(g_n = g^{(p-1)/n} \pmod p\) 要求指数一定要是整数，即 \(n | P-1\)，所以 \(P\) 必须满足：

\[ P = c \cdot 2^k + 1 \]

其中 \(k\) 要足够大，至少要大于题目中最大多项式长度的幂次。

而 \(998244353 = 113 \cdot 2^{23} + 1\)，符合条件。

最后有一些常见的 NTT 模数：

  • \(167772161 = 5 \cdot 2^{25} + 1, \ g = 3\)

  • \(469762049 = 7 \cdot 2^{26} + 1, \ g = 3\)

  • \(754974721 = 45 \cdot 2^{24} + 1, \ g = 11\)

  • \(1004535809 = 479 \cdot 2^{21} + 1, \ g = 3\)

实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
#include <bits/stdc++.h>
using namespace std;
using ll = long long;
constexpr int N = 3e6 + 5;
constexpr int P = 998244353, G = 3, IG = 332748118;
int n, m;
int r[N];
ll a[N], b[N];
int limit = 1, l = 0;

ll fp(ll x, ll p = P - 2) {
    int res = 1;
    while (p) {
        if (p & 1) res = res * x % P;
        x = x * x % P;
        p >>= 1;
    }
    return res;
}

ll add(ll x, ll y) { return x + y < P ? x + y : x + y - P; }
ll sub(ll x, ll y) { return x - y >= 0 ? x - y : x - y + P; }

void FNTT(ll *a, int type) {
    for (int i = 1; i < limit; i++)
        if (i < r[i]) swap(a[i], a[r[i]]);
    for (int mid = 1; mid < limit; mid <<= 1) {
        int gn = fp((type == 1) ? G : IG, (P - 1) / (mid << 1));
        for (int i = 0; i < limit; i += mid << 1) {
            ll g = 1;
            for (int j = 0; j < mid; j++, g = g * gn % P) {
                ll x = a[i + j], y = a[i + j + mid] * g % P;
                a[i + j] = add(x, y), a[i + j + mid] = sub(x, y);
            }
        }
    }
    if (type == -1) {
        int inv = fp(limit, P - 2);
        for (int i = 0; i < limit; i++) a[i] = a[i] * inv % P;
    }
}

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr), cout.tie(nullptr);
    cin >> n >> m;
    for (int i = 0; i <= n; i++) cin >> a[i];
    for (int i = 0; i <= m; i++) cin >> b[i];
    while (limit <= n + m) limit <<= 1, l++;
    for (int i = 1; i < limit; i++)
        r[i] = (r[i >> 1] >> 1) | ((i & 1) << (l - 1));
    FNTT(a, 1), FNTT(b, 1);
    for (int i = 0; i < limit; i++) a[i] = a[i] * b[i] % P;
    FNTT(a, -1);
    for (int i = 0; i <= n + m; i++) cout << a[i] << ' ';
    cout << '\n';
    return 0;
}