Note -「多项式」基础模板（FFT/NTT/多模 NTT）光速入门

  进阶篇戳这里。

目录

• $\omega$ 何为「多项式」$\omega$
  • $\omega$ 基本概念 $\omega$
  • $\omega$ 系数表示法 & 点值表示法 $\omega$
• $\omega$ 傅里叶（Fourier）变换 $\omega$
  • $\omega$ 概述 $\omega$
  • $\omega$ 前置知识 - 复数 $\omega$
  • $\omega$ 单位根 $\omega$
  • $\omega$ 快速傅里叶正变换（FFT）$\omega$
  • $\omega$ 快速傅里叶逆变换（IFFT）$\omega$
  • $\omega$ 迭代实现 $\omega$
  • $\omega$ 例题 $\omega$
    • $\omega$「洛谷 P3803」「模板」多项式乘法（FFT）$\omega$
      • $\omega$ 题意简述 $\omega$
      • $\omega$ 数据规模 $\omega$
      • $\omega$ $\text{Solution}$ $\omega$
      • $\omega$ $\text{Code}$ $\omega$
• $\omega$ 快速数论变换（NTT）$\omega$
  • $\omega$ 原根 $\omega$
  • $\omega$ 实现 $\omega$
  • $\omega$ NTT 模数 $\omega$
  • $\omega$ 奇怪的模数 - 任意模数 NTT $\omega$
    • $\omega$ 三模 NTT $\omega$
    • $\omega$ 拆系数 FFT（MTT）$\omega$
      • $\omega$ 七次转五次 $\omega$
      • $\omega$ 五次转四次 $\omega$
  • $\omega$ 例题 $\omega$
    • $\omega$ 「洛谷 P4245」「模板」任意模数 NTT
      • $\omega$ 题意简述 $\omega$
      • $\omega$ 数据规模 $\omega$
      • $\omega$ $\text{Solution}$ $\omega$
      • $\omega$ $\text{Code}$ $\omega$
    • $\omega$ 「CF 623E」Transforming Sequence $\omega$
      • $\omega$ 题意简述 $\omega$
      • $\omega$ 数据规模 $\omega$
      • $\omega$ $\text{Solution}$ $\omega$
• $\omega$ 参考资料 $\omega$

$\omega$ 何为「多项式」$\omega$

  多项式（polynomial）是指由变量（variable）、系数（coefficient）以及它们之间的加、减、乘、幂运算（非负整数次方）得到的表达式。多项式就是整式。

$\omega$ 基本概念 $\omega$

  在初中阶段就有所接触：整式可以包含加、减、乘、除、乘方五种对变量的运算，同时要求除数不含有变量。设 $x,y$ 为变量，其余字母为常量，则 $xy,a+by$ 是整式，$\sqrt{a}x+\pi y~(a\ge0)$ 也是整式，因为根号下的 $a$ 是常量；$\frac{x}y,\sqrt{xy}$ 则不是整式。

  在本文中，我们所说的多项式通常指一元多项式，即仅含一个变量 $x$ 的多项式。并简记一个多项式 $A(x)=a_0+a_1x+a_2x^2+\cdots+a_nx^n$，其中 $n$ 是 $A(x)$ 的次数。

$\omega$ 系数表示法 & 点值表示法 $\omega$

  在上文中 $A(x)=a_0+a_1x+a_2x^2+\cdots+a_nx^n$ 就是多项式的系数表示法。而可以证明，对于一个 $n$ 次多项式，可以用任意 $n+1$ 个其函数图像上不重合的点来唯一确定这个多项式。我们任取 $n+1$ 个横坐标的值 $x_0,x_1,\dots,x_n$，代入 $A(x)$，求得 $n+1$ 个点 $(x_0,y_0),(x_1,y_1),\dots,(x_n,y_n)$，就可以用这 $n+1$ 个点唯一确定 $A(x)$，这也被称为 $A(x)$ 的点值表示法。注意到 $x$ 是任意的，所以同一个多项式的点值表示法不唯一。

  接下来，抛出我们的问题：设有两个多项式 $A(x),B(x)$，如何计算它们的乘积（亦称作卷积） $C(x)=A(x)\cdot B(x)$ 呢？

  一个简单的想法，我们可以利用乘法分配律求到 $C(x)$ 每一项的系数。于是：

$$C(x)=\sum_{i=0}^{+\infty}\left(\sum_{j+k=i}a_jb_k\right)x^i$$

  注意求和中的 $+\infty$ 只是一种形式记法，事实上，当 $i$ 大于 $A(x)$ 的次数与 $B(x)$ 的次数之和时，$x^i$ 的系数 $c_i$ 显然为 $0$。

  不过呢……从算法设计的角度，这种求多项式乘法的方式的复杂度是 $\mathcal{O}(n^2)$ 的，好慢 qwq。

  联系到点值表示法，如果对于同一个横坐标 $x$，其在 $A(x)$ 上的坐标为 $(x,y_a)$，在 $B(x)$ 上的坐标为 $(x,y_b)$，那么代入最初的表达式，它在 $C(x)=A(x)\cdot B(x)$ 上的坐标就是 $(x,y_ay_b)$ 呐。所以，在已知 $A(x)$ 和 $B(x)$ 的点值表示法时，我们可以用 $\mathcal{O}(n)$ 的时间求出 $C(x)$ 的点值表示法！

$\omega$ 傅里叶（Fourier）变换 $\omega$

  我们继续解决上文多项式乘法的问题——找到一个高效的算法，实现系数与点值表示法的转换。

  本节中，若非特别说明，$n=2^k,k\in\mathbb{N}$。

$\omega$ 概述 $\omega$

  离散傅里叶变换（Discrete Fourier Transform，DFT），是傅里叶变换在时域和频域上都呈离散的形式，将信号的时域采样变换为其 DTFT 的频域采样。

  FFT 是一种高效实现 DFT 的算法，称为快速傅立叶变换（Fast Fourier Transform，FFT）。它对傅里叶变换的理论并没有新的发现，但是对于在计算机系统或者说数字系统中应用离散傅立叶变换，可以说是进了一大步。

$\omega$ 前置知识 - 复数 $\omega$

  其实很简单：我们把形如 $z=a+bi$ 的数称为复数，其中 $a,b\in\mathbb R$； $a$ 为实部，$b$ 为虚部，$i$ 有 $i^2=-1$，为虚数单位。我们可以把复数与二维平面上的向量 $(a,b)$ 一一对应，此时 $x$ 轴称为实轴，$y$ 轴称为虚轴。并记其与 $x$ 轴正方向的夹角为辐角 $\theta$，向量的模长为 $|z|$。如图，对于 $z=3+2i$：

  其中 $l$ 为模长，$\theta$ 为辐角。接下来我们研究复数的一些基本运算：

• 加减法：$z_1\pm z_2=(a_1\pm a_2)+(b_1\pm b_2)i$。

• 乘法：$z_1z_2=(a_1a_2-b_1b_2)+(a_1b_2+a_2b_1)i$。

    进一步研究，发现：

  $$|z_1z_2|^2=(a_1^2a_2^2-2a_1a_2b_1b_2+b_1^2b_2^2)+(a_1^2b_2^2+2a_1a_2b_1b_2+a_2^2b_1^2)=(a_1^2+b_1^2)(a_2^2+b_2^2)$$

    最后一步用了简单的因式分解技巧。由于 $a_1^2+b_1^2=|z_1|^2,a_2^2+b_2^2=|z_2|^2$，所以有 $|z_1z_2|=|z_1||z_2|$。

    此外，运用 $\tan$ 的和角公式，还可以求到 $\theta_{z_1z_2}=\theta_{z_1}+\theta_{z_2}$。

    所以，复数的乘法有一个重要的几何意义：模长相乘，辐角相加。

    最后，我们定义 $\overline{z}$ 为 $z$ 关于 $x$ 轴对称的向量，称为 $z$ 的共轭复数。即若 $z=a+bi$，有 $\overline{z}=a-bi$。可以发现 $z\cdot\overline{z}=a^2+b^2\in\mathbb R$。 为后文方便叙述，记 $\operatorname{conj}(z)=\overline z$。

$\omega$ 单位根 $\omega$

  定义 $n$ 次单位根为所有满足 $z^n=1$ 的 $z$。根据乘法的几何意义，可知 $|z|=1$，$n\theta_z=2k\pi$，其中 $k\in\mathbb{Z}$。于是满足条件且小于 $2\pi$ 的 $\theta_z$ 集合为 $\{0,\frac{1\times2\pi}{n},\frac{2\times2\pi}{n},\dots,\frac{(n-1)\times2\pi}{n}\}$，共 $n$ 个。我们记 $\omega_n^k$ （$\omega$，/'əʊmɪɡə/）为第 $k$ 个 $n$ 次单位根（$k=0,1,\dots,n-1$）。如图，当 $n=8$ 时：

  运用简单的三角函数知识，得到：

$$\omega_n^k=(\cos\frac{2\pi k}{n},\sin\frac{2\pi k}{n})$$

  可以发现，单位根实际上体现了一个周期意义：$\omega_n^n=\omega_n^0$。所以此后单位根的上标（乘方或是编号，它们在运算意义上等价）默认对单位根次数取模。单位根有如下性质：

• $\omega_n^j\omega_n^k=\omega_n^{j+k}$。
• $\omega_{kn}^{kj}=\omega_n^j$。
• $\omega_{2n}^k=-\omega_{2n}^{k+n}$。

  可以结合乘法几何意义与公式理解。

$\omega$ 快速傅里叶正变换（FFT）$\omega$

  有必要重申我们的目标：实现系数与点值表示法的转换。

  现在，设有 $A(x)=a_0+a_1x+a_2x^2+\cdots+a_{n-1}x^{n-1}$，我们希望代入 $n$ 次单位根 $\omega_n^0,\omega_n^1,\dots,\omega_n^{n-1}$ 这 $n$ 个复数来取得 $A(x)$ 的点值表示法。

  分开奇数项与偶数项，记：

$$A_e(x)=\sum_{i\in[0,n)\land2|i}a_ix^{\frac{i}2}$$

$$A_o(x)=\sum_{i\in[0,n)\land2\not|i}a_ix^{\frac{i-1}2}$$

  利用 $A_e(x)$ 和 $A_o(x)$ 来表示 $A(x)$，有：

$$A(x)=A_e(x^2)+xA_o(x^2)$$

  设 $n=2m$，尝试代入一个 $\omega_n^k$：

$$A(\omega_n^k)=A_e(\omega_n^{2k})+\omega_n^kA_o(\omega_n^{2k})$$

  注意到 $\omega_n^{2k}=\omega_{2m}^{2k}=\omega_m^k$，所以：

$$A(\omega_n^k)=A_e(\omega_m^k)+\omega_n^kA_o(\omega_m^k)$$

  类似地，代入 $\omega_n^{k+m}=-\omega_n^k$，有：

$$A(\omega_n^{k+m})=A_e(\omega_m^k)-\omega_n^kA_o(\omega_m^k)$$

  所以，当已知 $A_e(\omega_m^k)$ 和 $A_o(\omega_m^k)$ 时，可以 $\mathcal{O}(1)$ 计算出 $A(\omega_n^k)$ 和 $A(\omega_n^{k+m})$。现在，就有了一个算法雏形：

Function name: FFT.

Input: a polynomial $A(x)=a_0+a_1x_1+\cdots$ which's length is $n~(n=2^k,k\in\mathbb N)$.

Output: a set $\mathcal A=\{A(\omega_n^0),A(\omega_n^1),\dots,A(\omega_n^{n-1})\}$.

if $n=1$ then

 return $\{(a_0,0)\}$

$m\leftarrow\frac{n}2$

$A_e(x)\leftarrow\sum_{i\in[0,n)\land2|i}a_ix^{\frac{i}2}$

$A_o(x)\leftarrow\sum_{i\in[0,n)\land2\not|i}a_ix^{\frac{i-1}2}$

$\mathcal A\leftarrow\varnothing$

$\mathcal{A}_e\leftarrow\operatorname{FFT}(A_e,m)$

$\mathcal{A}_o\leftarrow\operatorname{FFT}(A_o,m)$

$k\leftarrow0$

while $k<m$ do

 $\mathcal{A}_k\leftarrow{\mathcal{A}_e}_k+\omega_n^k{\mathcal{A}_o}_k$

​ $\mathcal{A}_{k+m}\leftarrow{\mathcal{A}_e}_k-\omega_n^k{\mathcal{A}_o}_k$

​ $k\leftarrow k+1$

return $\mathcal A$

  复杂度是 $\mathcal{T}(n)=2\mathcal{T}(\frac{n}2)+\mathcal{O}(n)=\mathcal{O}(n\log n)$。在后文中，“多项式 $A(x)$ 的点值表示”均指由上述 FFT 算法得到的 $\mathcal A$。

  别着急做题，这种实现方法 SPFA 了 qwq。

$\omega$ 快速傅里叶逆变换（IFFT）$\omega$

  利用 FFT，我们可以快速地将系数表示法转换成点值表示法，那怎么在同样优秀的复杂度内转回去呢？

  从本质入手，考虑到 $\mathcal A_j=\sum_{k=0}^{n-1}a_k(\omega_n^j)^k$，我们把它写成矩阵乘法的形式：

$$\begin{bmatrix} (\omega_n^0)^0&(\omega_n^0)^1&\cdots&(\omega_n^0)^{n-1}\\ (\omega_n^1)^0&(\omega_n^1)^0&\cdots&(\omega_n^1)^{n-1}\\ \vdots&\vdots&\ddots&\vdots\\ (\omega_n^{n-1})^0&(\omega_n^{n-1})^1&\cdots&(\omega_n^{n-1})^{n-1} \end{bmatrix} \begin{bmatrix} a_0\\a_1\\\vdots\\a_{n-1} \end{bmatrix}= \begin{bmatrix} \mathcal A_0\\ \mathcal A_1\\ \vdots\\ \mathcal A_{n-1} \end{bmatrix}$$

  令左侧的系数矩阵为 $U$，考虑系数矩阵 $V$：

$$V=\begin{bmatrix} (\omega_n^{-0})^0&(\omega_n^{-0})^1&\cdots&(\omega_n^{-0})^{n-1}\\ (\omega_n^{-1})^0&(\omega_n^{-1})^0&\cdots&(\omega_n^{-1})^{n-1}\\ \vdots&\vdots&\ddots&\vdots\\ (\omega_n^{-(n-1)})^0&(\omega_n^{-(n-1)})^1&\cdots&(\omega_n^{-(n-1)})^{n-1} \end{bmatrix}$$

  那么：

$$(U\cdot V)_{ij}=\sum_{k=0}^{n-1}u_{ik}v_{kj}=\sum_{k=0}^{n-1}\omega_n^{k(j-i)}$$

  当 $i=j$ 时，显然 $(U\cdot V)_{ij}=n$。

  当 $i\not=j$ 时，发现是等比数列的形式。所以 $(U\cdot V)_{ij}=\frac{1-\omega_n^{n(j-i)}}{1-\omega_n^{j-i}}$。又因为 $\omega_n^n=\omega_n^0=1$，所以分子有 $1-1=0$，则原式 $=0$。

  综上，$(U\cdot V)_{ij}=[i=j]n$，即 $\frac{1}nV=U^{-1}$。

  以此，就有：

$$\begin{bmatrix} a_0\\a_1\\\vdots\\a_{n-1} \end{bmatrix}=\frac{1}{n}V\begin{bmatrix} \mathcal A_0\\ \mathcal A_1\\ \vdots\\ \mathcal A_{n-1} \end{bmatrix}= \frac{1}n\begin{bmatrix} (\omega_n^{-0})^0&(\omega_n^{-0})^1&\cdots&(\omega_n^{-0})^{n-1}\\ (\omega_n^{-1})^0&(\omega_n^{-1})^0&\cdots&(\omega_n^{-1})^{n-1}\\ \vdots&\vdots&\ddots&\vdots\\ (\omega_n^{-(n-1)})^0&(\omega_n^{-(n-1)})^1&\cdots&(\omega_n^{-(n-1)})^{n-1} \end{bmatrix} \begin{bmatrix} \mathcal A_0\\ \mathcal A_1\\ \vdots\\ \mathcal A_{n-1} \end{bmatrix}$$

  注意到 $\omega_n^{-k}=\omega_n^{n-k}$，所以依然能够用 FFT 的递归形式求解。

$\omega$ 迭代实现 $\omega$

  常数惊人呐……

  这里，我们将尝试把递归的 FFT Fast-Fast TLE 转换为仅用循环迭代的 FFT。

  观察递归时各系数的位置：

  观察第一层和最后一层：

原序列 dec	$0$	$1$	$2$	$3$	$4$	$5$	$6$	$7$
原序列 bin	$000$	$001$	$010$	$011$	$100$	$101$	$110$	$111$
新序列 dec	$0$	$4$	$2$	$6$	$1$	$5$	$3$	$7$
新序列 bin	$000$	$100$	$010$	$110$	$001$	$101$	$011$	$111$

  发现原序列下标的二进制翻转过来就是新序列下标。

  所以我们可以先安排好新序列下标，再从下往上迭代就可以了。

$\omega$ 例题 $\omega$

$\omega$「洛谷 P3803」「模板」多项式乘法（FFT）$\omega$

$\omega$ 题意简述 $\omega$

  link.

  给两个多项式的每项系数，求其相乘得到的多项式的每项系数。

$\omega$ 数据规模 $\omega$

  多项式长度 $n,m\le10^6$。

$\omega$ $\text{Solution}$ $\omega$

  直接给出完整代码，建议先模仿实现，以免常数惊人 owo。

  不过我不知道怎么折叠显示代码，影响阅读很抱歉 qwq。

$\omega$ $\text{Code}$ $\omega$

#include <cmath>
#include <cstdio>
#include <iostream>

inline int rint () {
	int x = 0, f = 1; char s = getchar ();
	for ( ; s < '0' || '9' < s; s = getchar () ) f = s == '-' ? -f : f;
	for ( ; '0' <= s && s <= '9'; s = getchar () ) x = x * 10 + ( s ^ '0' );
	return x * f;
}

template<typename Tp>
inline void wint ( Tp x ) {
	if ( x < 0 ) putchar ( '-' ), x = -x;
	if ( 9 < x ) wint ( x / 10 );
	putchar ( x % 10 ^ '0' );
}

const int MAXL = 1 << 21;
const double PI = acos ( -1 );
int n, m, rev[MAXL + 5];

struct Complex { // 复数结构体。
	double x, y;
	Complex () {} Complex ( const double tx, const double ty ): x ( tx ), y ( ty ) {}
	inline Complex operator + ( const Complex t ) const { return Complex ( x + t.x, y + t.y ); }
	inline Complex operator - ( const Complex t ) const { return Complex ( x - t.x, y - t.y ); }
	inline Complex operator * ( const Complex t ) const { return Complex ( x * t.x - y * t.y, x * t.y + y * t.x ); }
	inline Complex operator / ( const double t ) const { return Complex ( x / t, y / t ); }
} A[MAXL + 5], B[MAXL + 5];

inline void FFT ( const int n, Complex* A, const int tp ) {
	// n为多项式长度，保证是2的整数幂；A为每项系数；tp=1表示正变换，tp=-1表示逆变换。
	for ( int i = 0; i < n; ++ i ) if ( i < rev[i] ) std :: swap ( A[i], A[rev[i]] ); // 安排系数位置。
	for ( int i = 2, stp = 1; i <= n; i <<= 1, stp <<= 1 ) {
		Complex w ( cos ( PI / stp ), tp * sin ( PI / stp ) ); // 当前正在处理$\omega_i$。
		for ( int j = 0; j < n; j += i ) {
			Complex r ( 1, 0 );
			for ( int k = 0; k < stp; ++ k, r = r * w ) {
				// $r=\omega_i^k$。
				Complex ev ( A[j + k] ), ov ( r * A[j + k + stp] );
				A[j + k] = ev + ov, A[j + k + stp] = ev - ov;
			}
		}
	}
	if ( ! ~ tp ) for ( int i = 0; i < n; ++ i ) A[i] = A[i] / n; // 逆变换，最后每项/n。
}

inline void ployConv ( const int n, const int m, const Complex* A, const Complex* B, Complex* C ) {
	// 计算A(x)与B(x)的卷积，答案为C(x)。
	int lg = 0, len = 1;
	static Complex tmpA[MAXL + 5], tmpB[MAXL + 5];
	for ( ; len < n + m - 1; len <<= 1, ++ lg ); // C(x)长度为n+m-1，找到适合的长度len。
	for ( int i = 0; i < len; ++ i ) {
		tmpA[i] = A[i], tmpB[i] = B[i];
		rev[i] = ( rev[i >> 1] >> 1 ) | ( ( i & 1 ) << lg >> 1 ); // rev[i]为i翻转后的位置，自行理解。
	}
	FFT ( len, tmpA, 1 ), FFT ( len, tmpB, 1 ); // 正变换求到点值表示。
	for ( int i = 0; i < len; ++ i ) C[i] = tmpA[i] * tmpB[i]; // 求到C(x)的点值表示。
	FFT ( len, C, -1 ); // 逆变换求到C(x)。
}

int main () {
	n = rint () + 1, m = rint () + 1; // 注意这里要+1。
	for ( int i = 0; i < n; ++ i ) A[i].x = rint ();
	for ( int i = 0; i < m; ++ i ) B[i].x = rint ();
	ployConv ( n, m, A, B, A );
	for ( int i = 0; i < n + m - 1; ++ i ) {
		wint ( ( int ) round ( A[i].x ) );
		putchar ( i == n + m - 2 ? '\n' : ' ' );
	}
	return 0;
}

$\omega$ 快速数论变换（NTT）$\omega$

  事实上，由于 FFT 存在大量浮点数运算，其精度和常数都不尽人意。那么能否在整数域下找到单位根的替代品呢？

  我们来总结一下单位根应具有的性质：

• $\omega_n^j\omega_n^k=\omega_n^{j+k}$，注意这里的上标就意义上来说不是乘方，而是编号（当然运算上就是乘方 www）。
• $\omega_{kn}^{kj}=\omega_n^j$，我们需要这个性质进行递归处理。
• $\omega_{2n}^k=-\omega_{2n}^{k+n}$，这一性质体现在了合并计算的步骤。
• $(\forall j\not=k)(\omega_n^j\not=\omega_n^k)$，它保证了点值表示法的坐标不重合。
• $\omega_n^0=1$，它保证了逆变换顺利进行。

  考虑在模意义下……

$\omega$ 原根 $\omega$

  OI-Wiki 讲得好高深 qwq……

  对于素数 $p$，定义其原根 $g$ 为满足 $(\forall j,k\in[0,p-1),j\not=k)\left(g^j\not\equiv g^k\pmod p\right)$ 的数。模仿 $\omega$ 的定义，我们令 $\alpha_n^k=(g^\frac{p-1}n)^k$（这里 $\alpha$ 是随便取的名字 www），考虑其是否满足上述几条性质：

• $\alpha_n^j\alpha_n^k\equiv\alpha_n^{jk}\pmod p$，乘方性质，显然。
• $\alpha_{kn}^{kj}=(g^\frac{p-1}{kn})^{kj}=(g^\frac{p-1}n)^j=\alpha_n^j$，成立。
• $-\alpha_{2n}^{k+n}=-g^\frac{(p-1)(k+n)}{2n}=(-g^\frac{p-1}2)g^\frac{(p-1)k}{2n}\equiv g^\frac{(p-1)k}{2n}=\alpha_{2n}^k\pmod p$。注意因为 $(g^\frac{p-1}{2})^2\equiv1\pmod p$ 且根据 $g$ 的定义， $g^\frac{p-1}2=\alpha_{2n}^n\not=\alpha_{2n}^0\land\alpha_{2n}^0\equiv1\pmod p$，所以 $g^\frac{p-1}2\equiv-1\pmod p$，等式成立。
• 由定义，成立。
• 显然 $g\not=0$，成立。

  于是，把 FFT 中的所有 $\omega_n^k$ 替换为 $\alpha_n^k$ 就变成 NTT 啦！记得取模 owo。

$\omega$ 实现 $\omega$

inline void NTT ( const int n, int* A, const int tp ) {
	for ( int i = 0; i < n; ++ i ) if ( i < rev[i] ) A[i] ^= A[rev[i]] ^= A[i] ^= A[rev[i]];
	for ( int i = 2, stp = 1, w; i <= n; i <<= 1, stp <<= 1 ) {
		w = qkpow ( G, ( MOD - 1 ) / i );
              // G为原根，qkpow(a,b)为a的b次方模MOD的结果。
		if ( ! ~ tp ) w = qkpow ( w, MOD - 2 );
		for ( int j = 0; j < n; j += i ) {
			for ( int r = 1, k = j; k < j + stp; ++ k, r = 1ll * r * w % MOD ) {
				int ev = A[k], ov = 1ll * r * A[k + stp] % MOD;
				A[k] = ( ev + ov ) % MOD, A[k + stp] = ( ev - ov + MOD ) % MOD;
			}
		}
	}
	if ( ! ~ tp ) {
		int invn = qkpow ( n, MOD - 2 ); // invn为n的逆元。
		for ( int i = 0; i < n; ++ i ) A[i] = 1ll * A[i] * invn % MOD;
	}
}

$\omega$ NTT 模数 $\omega$

  注意到，$n$ 的大小是随多项式长度变化的，而我们必须保证 $n|(p-1)$，所以 $p-1$ 所含 $2$ 的因子个数决定了 NTT 能够接受的最大多项式长度。对于素数 $p$ 的选取自然也有了考究。

  记住一点，$998244353-1=2^{23}\times119$，适合做 NTT 模数，原根为 $3$；$10^9+7-1=2\times500000003$，不适合做 NTT 模数。

  更多 NTT 模数可见这篇博客。

  到此，多项式最最基础的知识就结束啦。

$\omega$ 奇怪的模数 - 任意模数 NTT $\omega$

  你 NTT 对模数有要求对吧？出题人就笑了……/xyx

  好吧，在不接受 FFT 精度误差的情况下，我们来想想办法应对模数不是 NTT 模数，甚至不是素数的毒瘤题吧！

$\omega$ 三模 NTT $\omega$

  我偏要用 NTT 模数 www。

  顾名思义，我们任取三个 NTT 模数 $p_1,p_2,p_3$，分别求出多项式乘积的第 $k$ 项模 $p_1,p_2,p_3$ 的结果，然后用中国剩余定理（CRT）解出该项模题目给定模数 $M$ 的结果就行了。不熟悉 CRT 的读者可以参考这篇博客。

  我个人认为不太优美懒，因此没有具体实现。在抱歉之余推荐这篇博客的模板 owo。

$\omega$ 拆系数 FFT（MTT）$\omega$

  我偏要用 FFT。

  缩写释义：MTT - MaoXiao Theory(?) Transformation。快速毛论变换（雾。

  假设多项式 $A(x)=a_0+a_1x+a_2x^2+\cdots,B(x)=b_0+b_1x+b_2x^2+\cdots$，长度均为 $n$，给定的模数 $m$，我们要求它们的卷积 $R(x)=A(x)\cdot B(x)$。如果直接 FFT，最大的系数可达到 $nm^2$ 的级别。我们假定 $n\le10^5,m\le10^9$，那么这样的精度误差是不被接受的。我们取 $M=\lceil\sqrt m\rceil$，并定义以下几个多项式：

$$\begin{aligned} C(x)&=\sum_{k=0}^{n-1}(a_k\bmod M)x^k\\ D(x)&=\sum_{k=0}^{n-1}\frac{a_k-(a_k\bmod M)}{M}x^k\\ E(x)&=\sum_{k=0}^{n-1}(b_k\bmod M)x^k\\ F(x)&=\sum_{k=0}^{n-1}\frac{a_k-(a_k\bmod M)}{M}x^k \end{aligned}$$

  可以发现，以上每个多项式的系数都是不超过 $M$ 的，那么任意两个多项式卷积的最大系数不超过 $nM^2\le10^{14}$，在 double 的承受范围以内。（保险起见，建议开 long double，不过对于某些 g++ 版本貌似并没有作用 owo。）我们用这四个多项式表示 $R(x)$：

$$R(x)=M^2\cdot D(x)F(x)+M\cdot\left(C(x)F(x)+D(x)E(x)\right)+C(x)E(x)$$

  这样，在 FFT 过程中，系数大小能被接受；作系数拼凑时步步取模亦能保证不溢出，我们就用 FFT 解决了任意模数 NTT 的问题。一般来说，取 $M=32768=2^{15}$，可以用位运算节省常数；还应利用 $\omega_n^k$ 的通项预处理出所有单位根以保障精度。

$\omega$ 七次转五次 $\omega$

  数一数，我们对 $C(x),D(x),E(x),F(x)$ 分别进行了一次 FFT，对 $D(x)F(x),$ $C(x)F(x)+D(x)E(x),$ $C(x)E(x)$ 分别进行了一次 IFFT，一共有七次变换操作，常数有点点大 qwq。

  接下来，我们假设有实系数多项式 $A(x)=a_0+a_1x+a_2x^2+\cdots,B(x)=b_0+b_1x+b_2x^2+\cdots$，长度均为 $n~(n=2^k,k\in\mathbb N)$，我们希望分别求到它们的点值表示 $\mathcal A$ 和 $\mathcal B$。于是定义：

$$P(x)=A(x)+iB(x)$$

$$Q(x)=A(x)-iB(x)$$

  其中 $i$ 是虚数单位。设 $P$ 与 $Q$ 的点值表示分别为 $\mathcal P,\mathcal Q$，并记 $\omega_n^k$ 的辐角 $\frac{2\pi k}n$ 为 $\theta_n^k$。那么：

$$\begin{aligned} \mathcal P_k&=A(\omega_n^k)+iB(\omega_n^k)\\ &=\sum_{j=0}^{n-1}a_j\omega_n^k+ib_j\omega_n^k\\ &=\sum_{j=0}^{n-1}(a_j+ib_j)(\cos\theta_n^k+i\sin\theta_n^k)\\\\ \mathcal Q_k&=A(\omega_n^k)-iB(\omega_n^k)\\ &=\sum_{j=0}^{n-1}a_j\omega_n^{jk}-ib_j\omega_n^{jk}\\ &=\sum_{j=0}^{n-1}(a_j-ib_j)(\cos\theta_n^{jk}+i\sin\theta_n^{jk})\\ &=\sum_{j=0}^{n-1}(a_j\cos\theta_n^{jk}+b_j\sin\theta_n^{jk})+i(a_j\sin\theta_n^{jk}-b_j\cos\theta_n^{jk})\\ &=\operatorname{conj}\left(\sum_{j=0}^{n-1}(a_j\cos\theta_n^{jk}+b_j\sin\theta_n^{jk})-i(a_j\sin\theta_n^{jk}-b_j\cos\theta_n^{jk})\right)\\ &=\operatorname{conj}\left(\sum_{j=0}^{n-1}\left(a_j\cos(-\theta_n^{jk})-b_j\sin(-\theta_n^{jk})\right)+i\left(a_j\sin(-\theta_n^{jk})-b_j\cos(-\theta_n^{jk})\right)\right)\\ &=\operatorname{conj}\left(\sum_{j=0}^{n-1}(a_j+ib_j)\left(\cos(-\theta_n^{jk})+i\sin(-\theta_n^{jk})\right)\right)\\ &=\operatorname{conj}\left(\sum_{j=0}^{n-1}(a_j+ib_j)\left(\cos\theta_n^{-jk}+i\sin\theta_n^{-jk}\right)\right)\\ &=\operatorname{conj}\left(\sum_{j=0}^{n-1}(a_j+ib_j)\left(\cos\theta_n^{(n-k)j}+i\sin\theta_n^{(n-k)j}\right)\right)\\ &=\operatorname{conj}(\mathcal P_{n-k}) \end{aligned}$$

  推导时，时刻留意 $\theta_n$ 的上标在模 $n$ 意义下嗷！

  最终，发现 $\mathcal Q$ 可以在求得 $\mathcal P$ 的前提下 $\mathcal O(1)$ 求到。只需要用 $\mathcal P$ 和 $\mathcal Q$ 表示出 $\mathcal A$ 和 $\mathcal B$：

$$\mathcal A_k=\frac{\mathcal P_k+\mathcal Q_k}2\\ \mathcal B_k=i\frac{\mathcal Q_k-\mathcal P_k}2$$

  就可以用一次 FFT 求出两个多项式的点值表示啦！（我认为毛啸的论文在 $\mathcal B_k$ 的表达式处有误，右侧分子应为本文的 $\mathcal Q_k-\mathcal P_k$ 而非 $\mathcal P_k-\mathcal Q_k$，欢迎指正 qwq。）

  于是，我们可以把对 $C(x),D(x),E(x),F(x)$ 的 FFT 合并为两组，就只需要五次 FFT 啦！

$\omega$ 五次转四次 $\omega$

  正变换已经足够优秀了，我们接下来考虑逆变换的合并。按照上文的推导，我们只需要求到 $P(x)$ 或者 $Q(x)$ 就可以了——因为原多项式是实系数的，所以对于 $P(x)$ 或 $Q(x)$ 的每一项，我们能够区分出其实部就是 $A(x)$ 的系数而虚部就是 $B(x)$ 的系数。既然已知 $\mathcal A$ 和 $\mathcal B$，直接利用定义式，有：

$$\mathcal P_k=\mathcal A_k+i\mathcal B_k$$

  逆变换 $\mathcal P$ 得到 $P(x)$，就可以把两次逆变换合并为一次了。最终，我们把四次正变换合并为两次，三次逆变换合并为两次，仅用四次 FFT 就解决了这一问题。

$\omega$ 例题 $\omega$

$\omega$ 「洛谷 P4245」「模板」任意模数 NTT

$\omega$ 题意简述 $\omega$

  link.

  给两个多项式，求其在模 $p$ 意义下的卷积。

$\omega$ 数据规模 $\omega$

  多项式长度 $n\le10^5$，系数 $a_k,b_k\le10^9$，$p\le10^9+9$。

$\omega$ $\text{Solution}$ $\omega$

  这里仅给出四次 FFT 的 MTT 做法。

$\omega$ $\text{Code}$ $\omega$

#include <cmath>
#include <cstdio>
#include <iostream>

typedef long long LL;

inline int rint () { /* 快读 */ }

inline void wint ( const int x ) { /* 快输 */ }

const int MAXN = 1 << 18;
const double PI = acos ( -1 );
int n, m, p, A[MAXN + 5], B[MAXN + 5], rev[MAXN + 5];

struct Complex {
	/* 复数结构体 */
} omega[MAXN + 5], P[MAXN + 5], Q[MAXN + 5], C[MAXN + 5], D[MAXN + 5], E[MAXN + 5], F[MAXN + 5];

inline void FFT ( const int n, Complex* A, const int tp ) { /* FFT，应使用预处理好的单位根omega[]。 */ }

int main () {
	n = rint () + 1, m = rint () + 1, p = rint ();
	for ( int i = 0; i < n; ++ i ) A[i] = rint () % p, P[i] = Complex ( A[i] & 0x7fff, A[i] >> 15 );
	for ( int i = 0; i < m; ++ i ) B[i] = rint () % p, Q[i] = Complex ( B[i] & 0x7fff, B[i] >> 15 );
	int lg = 0, len = 1;
	for ( ; len < n + m - 1; len <<= 1, ++ lg );
	for ( int i = 0; i < len; ++ i ) rev[i] = ( rev[i >> 1] >> 1 ) | ( ( i & 1 ) << lg >> 1 );
	for ( int i = 1; i < len; i <<= 1 ) { // 预处理单位根。
		for ( int k = 0; k < i; ++ k ) {
			omega[len / i * k] = Complex ( cos ( PI * k / i ), sin ( PI * k / i ) );
		}
	}
	FFT ( len, P, 1 ), FFT ( len, Q, 1 );
	for ( int i = 0; i < len; ++ i ) {
		Complex t ( P[( len - i ) % len].x, -P[( len - i ) % len].y );
		C[i] = ( P[i] + t ) / 2, D[i] = Complex ( 0, 1 ) * ( t - P[i] ) / 2;
	}
	for ( int i = 0; i < len; ++ i ) {
		Complex t ( Q[( len - i ) % len].x, -Q[( len - i ) % len].y );
		E[i] = ( Q[i] + t ) / 2, F[i] = Complex ( 0, 1 ) * ( t - Q[i] ) / 2;
	}
	for ( int i = 0; i < len; ++ i ) {
		Complex c ( C[i] ), d ( D[i] ), e ( E[i] ), f ( F[i] );
		C[i] = c * e, D[i] = c * f + d * e, E[i] = d * f;
		P[i] = C[i] + Complex ( 0, 1 ) * E[i];
	}
	FFT ( len, D, -1 ), FFT ( len, P, -1 );
	for ( int i = 0; i < n + m - 1; ++ i ) {
		int c = ( LL ( P[i].x + 0.5 ) % p + p ) % p, d = ( LL ( D[i].x + 0.5 ) % p + p ) % p, e = ( LL ( P[i].y + 0.5 ) % p + p ) % p;
		wint ( ( c + ( 1ll << 15 ) % p * d % p + ( 1ll << 30 ) % p * e % p ) % p );
		putchar ( i + 1 == n + m - 1 ? '\n' : ' ' );
	}
	return 0;
}

$\omega$ 「CF 623E」Transforming Sequence $\omega$

$\omega$ 题意简述 $\omega$

  link.

  求有多少个长度为 $n$ 的序列 $\{a_n\}$，满足 $a_k\in(0,2^k)\land a_k\not\subseteq\bigcup_{j=1}^{k-1}a_j$。这里的集合运算把 $a_k$ 作为了二进制集合。对 $10^9+7$ 取模。

$\omega$ 数据规模 $\omega$

  $n\le10^{18},k\in[1,3\times10^4]$。

$\omega$ $\text{Solution}$ $\omega$

  详见我的题解。题意描述略有不同，请自行理解。为什么我的 MTT 有八次变换啊……

  专门提起这道题的目的，除了让大家实战运用 MTT 之外，还想让大家感受倍增这种几乎能融入所有 OI 知识点的思想。

  本来还有多项式全家桶，但因为太懒所以烂尾咯 qwq。

$\omega$ 参考资料 $\omega$

  按在文中第一次借鉴或引用的位置顺序排列：

• 多项式 @百度百科.
• 整式 @百度百科.
• 快速傅里叶变换 @OI-Wiki.
• 多项式乘法 @Tiw_Air_OAO.
• FFT 用到的各种素数 @Miskcoo's Space.
• 中国剩余定理 CRT @yyb.
• 再探快速傅里叶变换 @毛啸.（网上资源）