OI-wiki 数论变换(NTT / FNTT)完全指南:原理、素数模数选型与 C++ 实现
OI-wiki 数论变换NTT / FNTT完全指南原理、素数模数选型与 C 实现【免费下载链接】OI-wiki:star2: Wiki of OI / ICPC for everyone. 某大型游戏线上攻略内含炫酷算术魔法项目地址: https://gitcode.com/GitHub_Trending/oi/OI-wiki数论变换Number-Theoretic Transform, NTT是离散傅里叶变换DFT在数论基础上的实现快速数论变换Fast NTT, FNTT则是快速傅里叶变换FFT在数论基础上的实现。本文以 OI-wiki 的 ntt.md 为主体结合仓库内完整可运行的 NTT 模板代码 与官方 测试样例系统讲解 NTT 的数学原理、常用素数模数的选择依据、DFT/FFT/NTT/FNTT 四者之间的关系以及一个可直接提交评测的迭代实现。读完本文你将掌握为什么竞赛中最常见的模数是998244353并具备独立写出、验证和迁移 NTT 模板的能力。从 FFT 的痛点说起为什么需要 NTT快速傅里叶变换可以在 $O(n\log n)$ 时间内完成两个 $n$ 次多项式的乘法比朴素的 $O(n^2)$ 逐项卷积快得多。但 FFT 在实现上存在一些固有的缺点数据向量必须乘以复系数矩阵加以处理每个复系数的实部和虚部是一个正弦及余弦函数因此大部分系数都是浮点数必须进行复数且是浮点数的运算计算量较大浮点数运算会产生不可避免的舍入误差当结果需要精确的整数时尤其危险。数论变换正是为解决这些问题而生NTT 解决的是多项式乘法带模数的情况所有运算都在模 $p$ 意义下以整数进行结果精确无误差。当然它也有一些限制——受模数限制且涉及的数通常比较大。目前最常见的模数是998244353后文会解释这一选择背后的数学原因。前置知识学习数论变换需要以下前置知识相关知识可以在 OI-wiki 对应页面中学习离散傅里叶变换DFTNTT 是 DFT 在有限域上的对应物参见 快速傅里叶变换生成子群NTT 依赖的「原根」本质上是模 $p$ 既约剩余系乘法群循环群的生成元原根primitive rootNTT 用模 $p$ 意义下的本原单位根替换复数单位根参见 原根离散对数与「阶」的概念密切相关是理解原根性质的辅助工具。多项式与生成函数的整体框架可参见 多项式与生成函数。数论变换的定义与核心公式定义有限域上的 DFT在数学中NTT 是关于任意环上的离散傅里叶变换。在有限域的情况下通常称为数论变换NTT。数论变换通过将离散傅里叶变换化为 $F\mathbb{Z}/p$——即整数模质数 $p$——来构造。这是一个有限域只要 $n$ 可除 $p-1$就存在本原 $n$ 次方根因此我们有 $p\xi n1$其中 $\xi$ 为正整数。具体来说对于质数 $pqn1$$n2^m$原根 $g$ 满足$$ g^{qn} \equiv 1 \pmod p $$将 $g_n g^q \pmod p$ 看作 $\omega_n$复数本原 $n$ 次单位根的等价物则它满足相似的性质$$ g_n^n \equiv 1 \pmod p,\qquad g_n^{n/2} \equiv -1 \pmod p $$这正是 FFT 中单位根最核心的两条性质乘方循环与半角取负也是蝶形运算能够成立的基础。这里有一个记号上的说明因为涉及到数论变换这里的 $n$为了区分 FFT 中的 $n$下文记作 $N$可以比 FFT 中的 $n$ 大但只需要把 $\dfrac{qN}{n}$ 看作这里的 $q$ 即可这样能够避免大小问题。常用素数模数速查表OI-wiki 的 ntt.md 给出了五个经过验证的常用素数模数$pqn1$ 形式其中 $2$ 的幂次部分决定了可支持的最大变换长度模数 $p$分解形式原根 $g$可支持的最大 $2$ 幂阶167772161$5 \times 2^{25}1$3$2^{25}$469762049$7 \times 2^{26}1$3$2^{26}$754974721$3^2 \times 5 \times 2^{24}1$11$2^{24}$998244353$7 \times 17 \times 2^{23}1$3$2^{23}$1004535809$479 \times 2^{21}1$3$2^{21}$这些模数都满足「$p-1$ 含有足够大的 $2$ 的幂因子」这一条件因为 NTT 的变换长度必须是 $2$ 的幂。其中998244353支持到 $2^{23}8388608$对于绝大多数算法竞赛题目的多项式乘法场景已经足够。以表格第二行为例$4697620497\times 2^{26}1$即 $q7$$n2^{26}$原根 $g3$那么 $g_n g^q 3^7 \pmod p$ 就是阶为 $2^{26}$ 的本原单位根。迭代过程中的单位根计算在迭代蝶形合并到长度 $l$ 时需要使用的单位根为$$ g_l g^{\frac{p-1}{l}} $$或者等价地写作$$ \omega_n g_l g_N^{\frac{N}{l}} g_N^{\frac{p-1}{l}} $$也就是说$N$ 阶单位根 $g_N g^{\frac{p-1}{N}}$而长度为 $l$ 的变换所需单位根可以直接由 $g_N$ 的幂得到这避免了每一层都重新做一次大指数快速幂。不过常见的简洁模板如仓库中的 ntt_1.cpp直接对每层计算 $g^{(p-1)/l}$实现上更加直观。快速数论变换分治加速FNTT 与 NTT 的关系快速数论变换FNTT是数论变换NTT增加分治操作之后的快速算法。快速数论变换使用的分治办法与快速傅里叶变换使用的分治办法完全一致——这意味着只需在快速傅里叶变换的代码基础上进行简单修改即可得到快速数论变换的代码。在算法竞赛中常提到的「NTT」一词往往实际指的是快速数论变换一般默认「数论变换」是指「快速数论变换」。这样简写是有逻辑依据的类似地「快速傅里叶变换」FFT实际指「快速离散傅里叶变换」FDFT但由于「快速」只能作用于离散情形甚至是本原单位根阶数为 $2$ 的幂的特殊情形不能作用于连续情形因此「离散」一词被省略FDFT 变为 FFT——FFT 永远指特殊的离散情形数论变换或快速数论变换是在取模意义下进行的操作不存在连续的情形永远是离散的自然也无需提到「离散」一词在算法领域不进行提速的操作是无意义的。FFT 文档中专门介绍 DFT是因为 DFT 在信号处理、图像处理领域有独立应用且 DFT 是 FFT 的原理与前置知识。在不引起混淆的情形下常用 NTT 代指 FNTT。为了严谨下文将 NTT 与 FNTT 两个词进行分离。DFT、FFT、NTT、FNTT 的关系四个概念之间的转换关系可以概括为两条路径增加分治操作在 DFT 与 NTT 的基础上增加分治操作分别得到 FFT 与 FNTT。分治操作的办法与原理奇偶分组、位逆序置换、蝶形运算可以参见 快速傅里叶变换替换运算与单位根在 DFT 与 FFT 的基础上将复数加法与复数乘法替换为模 $p$ 意义下的加法和乘法数值大小一般限制在 $0$ 到 $p-1$ 之间将本原单位根替换为模 $p$ 意义下相同阶数阶数为 $2$ 的幂的本原单位根即可得到 NTT 与 FNTT。由于替换的运算只涉及加法和乘法因此 DFT、FFT、NTT、FNTT 拥有相同的原理它们均在满足加法与乘法的环上进行无需域上满足除法运算的更严格条件。环与域的一般定义可参见 抽象代数基本概念。适用范围原根存在即可行事实上只要拥有原根——即群论中的生成元——该模数下的 NTT 或 FNTT 即可进行。考虑到模数为 $1$、$2$ 和 $4$ 的情形太小不具有实际意义对于奇素数 $p$ 和正整数 $\alpha$只要给出模数为 $p^\alpha$ 和 $2p^\alpha$ 的原根 $g$采用同样的办法NTT 或 FNTT 仍然可以进行。这与 原根 一文中的原根存在定理模 $m$ 原根存在当且仅当 $m1,2,4,p^e,2p^e$完全吻合。模板代码逐行解析OI-wiki 的 ntt.md 给出了 Library Checker - Convolution 对应的完整实现仓库中对应的源文件为 docs/math/code/poly/ntt/ntt_1.cpp。这个实现是非递归迭代版的蝶形运算整体流程与 fft.md 中「倍增法实现」完全一致。完整代码#include algorithm #include cstdio using namespace std; constexpr int N 1 20, P 998244353; int qpow(int x, int y) { int res 1; while (y) { if (y 1) res 1ll * res * x % P; x 1ll * x * x % P; y 1; } return res; } int r[N]; void ntt(int *x, int lim, int opt) { for (int i 0; i lim; i) if (r[i] i) swap(x[i], x[r[i]]); for (int m 2; m lim; m 1) { int k m 1; int gn qpow(3, (P - 1) / m); for (int i 0; i lim; i m) { int g 1; for (int j 0; j k; j, g 1ll * g * gn % P) { int tmp 1ll * x[i j k] * g % P; x[i j k] (x[i j] - tmp P) % P; x[i j] (x[i j] tmp) % P; } } } if (opt -1) { reverse(x 1, x lim); int inv qpow(lim, P - 2); for (int i 0; i lim; i) x[i] 1ll * x[i] * inv % P; } } int A[N], B[N], C[N]; int main() { int n, m; scanf(%d %d, n, m); for (int i 0; i n; i) scanf(%d, A[i]); for (int i 0; i m; i) scanf(%d, B[i]); int N max(n, m), lim 1; while (lim (N 1)) lim 1; for (int i 0; i lim; i) r[i] (i 1) * (lim 1) (r[i 1] 1); ntt(A, lim, 1); ntt(B, lim, 1); for (int i 0; i lim; i) C[i] 1ll * A[i] * B[i] % P; ntt(C, lim, -1); for (int i 0; i n m - 1; i) printf(%d%c, C[i], \n[i n m - 2]); return 0; }关键组成部分说明常量与快速幂。P 998244353即前文速查表中的首选模数N 1 20是数组上限变换长度lim不能超过它实际长度按需取 2 的幂。qpow(x, y)是标准的模意义快速幂用于计算单位根 $g^{(P-1)/m}$ 与逆元 $x^{-1} \bmod P$由费马小定理$x^{-1}x^{P-2}$。位逆序置换。非递归实现的第一个步骤是 bit-reversal permutation把元素按照二进制位反转后的下标重排从而将递归分治改写成自底向上的迭代。第 50 行的递推式r[i] (i 1) * (lim 1) (r[i 1] 1);正是 fft.md 中推导出的 $O(n)$ 递推公式 $R(x)\left\lfloor R(\lfloor x/2\rfloor)/2 \right\rfloor (x\bmod 2)\times \frac{len}{2}$ 的直接代码化(i 1) * (lim 1)处理最低位翻转到最高位r[i 1] 1处理其余位的翻转。蝶形合并。外层循环m从 2 倍增到lim代表当前合并的段长k m 1是半段长。gn qpow(3, (P - 1) / m)正是前文公式 $g_l g^{(p-1)/l}$此处 $lm$。内层两行x[i j k] (x[i j] - tmp P) % P; // f(ω^{kn/2}) G - ω^k·H x[i j] (x[i j] tmp) % P; // f(ω^k) G ω^k·H对应蝶形运算的两个式子 $f(\omega_n^k)G(\omega_{n/2}^k)\omega_n^k H(\omega_{n/2}^k)$ 与 $f(\omega_n^{kn/2})G(\omega_{n/2}^k)-\omega_n^k H(\omega_{n/2}^k)$。(x - tmp P) % P中的 P是为了避免模运算出现负数。逆变换。与 FFT 逆变换完全相同的手法reverse(x 1, x lim)将所有非零下标反转等价于把单位根换成其逆 $g^{-1}$对应 fft.md 方法二再对全体元素乘以 $\lim^{-1} \bmod P$ 完成归一化。因此同一个ntt函数通过opt参数同时承担正变换与逆变换两种角色opt 1时为正变换opt -1时为逆变换。主流程。标准的「点值相乘」三步对 $A$、$B$ 分别做正变换在频域逐点相乘C[i] A[i] * B[i] % P再做逆变换得到卷积结果。变换长度取lim ≥ n m - 1的最小的 2 的幂第 48–49 行这是因为循环卷积长度必须超过结果多项式次数高阶项不足需要补零。用仓库样例验证模板仓库为这份模板配套了官方的输入输出样例可以直接验证实现的正确性输入 ntt_1.in4 5 1 2 3 4 5 6 7 8 9期望输出 ntt_1.ans5 16 34 60 70 70 59 36手动验证$(12x3x^24x^3)\times(56x7x^28x^39x^4)$ 展开后系数依次为 $1\times55$$1\times62\times516$$1\times72\times63\times534$$1\times82\times73\times64\times560$后续为 $70,70,59,36$与期望输出完全一致。在本地或在线评测系统如 Library Checker 的convolution_mod题中编译运行该模板即可验证 NTT 的正确性。NTT 在多项式运算体系中的角色NTT 是 OI-wiki 多项式算法体系的地基。在 多项式与生成函数 中可以看到以快速傅里叶变换为基石的多项式算法赋予了选手直接操纵生成函数的能力而在模意义下这份「基石」的角色正是由 NTT 承担的。仓库中大量多项式高级操作的代码都以 NTT 为底层算子多项式求逆elementary-func.md 中的polyinv使用倍增法迭代过程中对每个长度做两次 DFT 与一次 IDFT其模数同样为998244353多项式开方sqrt_1.cpp 中实现了一个风格略有差异的 NTT使用change位逆序函数 qpow(3, (mod - 1) / q)逐层求单位根并在开方迭代中反复调用NTT、inv其他如 初等函数、牛顿迭代 等页面描述的多项式指数、对数、除法等操作复杂度均建立在 $O(n\log n)$ 的 NTT 乘法之上。这说明了掌握 NTT 模板的实用性它不仅是解决「裸卷积」题目的工具更是后续一切模意义多项式算法的通用加速器。从源码结构看仓库中 docs/math/code/poly/ntt/ 目录维护着这份最基础的模板而 docs/math/code/poly/ 下其他目录如sqrt/、comp-rev/、czt/的代码均直接使用998244353及其原根 3印证了该模数在算法竞赛实践中的主流地位。常见问题与工程建议为什么首选 998244353它满足 $p-17\times17\times2^{23}$含 $2^{23}$ 因子能支持到长度 8388608 的变换原根很小$g3$计算单位根时快速幂代价低且 $p 2^{30}$在 32 位有符号整数范围内一次int乘法用 64 位long long承接即可无需高精度或__int128。这是它成为竞赛事实标准的原因。什么时候需要其他模数当答案需要对某个特定素数取模或需要 CRT中国剩余定理合并多个 NTT 结果来突破单模数的上限时可以使用速查表中的其他模数如 469762049、1004535809 等它们拥有互素的结构便于三模数 NTT CRT 恢复出更大范围内的精确值。精度与误差。与 FFT 的浮点运算不同NTT 全程是整数模运算结果精确无舍入误差这是它在「答案必须精确取模」的题目中不可替代的原因。长度限制。每个模数能支持的变换长度受限于 $p-1$ 中 $2$ 的幂因子变换长度lim必须是 2 的幂且不超过该上限否则 $g^{(p-1)/m}$ 不再是本原 $m$ 次单位根结果将出错。参考资料与拓展阅读OI-wiki 的 ntt.md 文末推荐了以下扩展资料供想深入了解 NTT/FFT 推导细节的读者参考FWT快速沃尔什变换零基础详解 qaqACM/OIFFT快速傅里叶变换0 基础详解附 NTTACM/OINumber-theoretic transform (NTT) — WikipediaTutorial on FFT/NTT — The tough made simple (Part 1)NTT 模板 — CSDN 博客相关页面快速傅里叶变换的完整推导与两种非递归实现见 fft.md多项式求逆、开方等基于 NTT 的高级操作见 elementary-func.md原根与阶的数论基础见 primitive-root.md。【免费下载链接】OI-wiki:star2: Wiki of OI / ICPC for everyone. 某大型游戏线上攻略内含炫酷算术魔法项目地址: https://gitcode.com/GitHub_Trending/oi/OI-wiki创作声明:本文部分内容由AI辅助生成(AIGC),仅供参考