FFT

多项式

定义

多项式是形如 $A(x) = \sum_{i = 0}^n a_i \times x_i$A(x)=∑i=0n​ai​×xi​ 的表达式，其中 $x$x 为不定元，这是多项式的系数表示法。

将多项式中不定元的最大次数称为多项式的次数，记做 $\partial(A(x))$∂(A(x))。

运算

$\begin{array}{c} A(x)=\sum_{i=0}^{n} a_{i} x^{i}, B(x)=\sum_{i=0}^{m} b_{i} x^{i} \\ \ \\ A(x) \pm B(x)=\sum_{i=0}^{\max \{n, m\}}\left(a_{i} \pm b_{i}\right) x^{i} \\ \ \\ A(x) \times B(x)=\sum_{i=0}^{n} \sum_{j=0}^{m}\left(a_{i} \times b_{j}\right) x^{i+j} \\ =\sum_{k=0}^{n+m} \sum_{i=0}^{n}\left(a_{i} \times b_{k-i}\right) x^{k} \end{array}$A(x)=∑i=0n​ai​xi,B(x)=∑i=0m​bi​xi A(x)±B(x)=∑i=0max{n,m}​(ai​±bi​)xi A(x)×B(x)=∑i=0n​∑j=0m​(ai​×bj​)xi+j=∑k=0n+m​∑i=0n​(ai​×bk−i​)xk​

多项式乘法满足交换律，结合律和分配率。

多项式的除法通常指带余除法。

$A(x) = q(x)B(x)+r(x) \quad \partial(r(x)) < \partial(B(x)) \\ \ \\ A(x) \equiv r(x) \quad(\bmod B(x))$A(x)=q(x)B(x)+r(x)∂(r(x))<∂(B(x)) A(x)≡r(x)(modB(x))

$q(x)$q(x) 称为 $B(x)$B(x) 除 $A(x)$A(x) 的商式，$r(x)$r(x) 称为 $B(x)$B(x) 除 $A(x)$A(x) 的余式。若 $r(x)=0$r(x)=0，则称 $B(x)$B(x) 整除 $A(x)$A(x)，记作 $B(x)|A(x)$B(x)∣A(x)。

点值表示法

$n$n 个点可以确定唯一的 $n-1$n−1 次多项式。如果多项式 $f(x)$f(x) 和 $g(x)$g(x) 的次数不超过 $n$n，且对于 $n+1$n+1 个不同的数都有相同的对应值，则 $f(x) = g(x)$f(x)=g(x)。

如果选取 $n+1$n+1 个不同的数，$x_0 \sim x_n$x0​∼xn​ 对多项式求职，得到 $A(x_0), A(x_1), \cdots, A(x_n)$A(x0​),A(x1​),⋯,A(xn​)。那么称 $\{ (x_i,A(x_i) | 0 \leq i \leq n, i \in\mathbb{N}\}${(xi​,A(xi​)∣0≤i≤n,i∈N} 为多项式 $A(x)$A(x) 的点值表示。通过系数表示求点值表示的过程叫做求值，通过点值表示求系数表示的过程叫做插值。

设 $A(x)$A(x) 和 $B(x)$B(x) 是两个多项式，它们在 $x_0 \sim s_n$x0​∼sn​ 上的取值分别为 $c_0 \sim c_n$c0​∼cn​ 和 $D_0 \sim D_n$D0​∼Dn​。那么 $A(x) \times B(x)$A(x)×B(x) 在 $x_{0}, x_{1}, \ldots, x_{n}$x0​,x1​,…,xn​ 上的取值分别为 $c_{0} D_{0} \sim c_n D_n$c0​D0​∼cn​Dn​。

复数

• 定义

虚数单位 $i=\sqrt{-1}$i=−1​。

形如 $a + b_i (a, b \in \mathbb{R})$a+bi​(a,b∈R) 的数称为复数，$a$a 称为实部，$b$b 称为虚部。全体复数的集合记为 $\mathbb{C}$C。

• 运算

设 $z = a + bi, w = c+di$z=a+bi,w=c+di，有

$z \pm w=(a \pm c)+(b \pm d) i \\ \ \\ z \times w=a c+a d i+b c i+b d i^{2}=(a c-b d)+(a d+b c) i \\ \ \\ z \div w=\frac{a+b i}{c+d i}=\frac{(a+b i)(c-d i)}{(c+d i)(c-d i)}=\frac{(a c+b d)+(b c-a d) i}{c^{2}+d^{2}}$z±w=(a±c)+(b±d)i z×w=ac+adi+bci+bdi2=(ac−bd)+(ad+bc)i z÷w=c+dia+bi​=(c+di)(c−di)(a+bi)(c−di)​=c2+d2(ac+bd)+(bc−ad)i​

复数的相加满足平行四边形法则。

共轭复数：$\overline{z}=a-b i$z=a−bi，相当于对虚部取反。

• 几何意义

令横轴为实轴，纵轴为虚轴。复数 $a+bi$a+bi 就对应坐标为 $(a,b)$(a,b) 的点或向量。

复数对应的点到原点的距离称为复数的模长，记为 $r=|z|$r=∣z∣。 由勾股定理 $|z|=\sqrt{a^{2}+b^{2}}=\sqrt{z \overline{z}}$∣z∣=a2+b2​=zz​

复数与实轴正半轴形成的角称为复数的辐角，记为 $\varphi=\arg z$φ=argz。

复数也可以用模长和辐角表示，称为复数的极坐标形式。由极坐标与直角坐标的转化关系，$z=r(\cos \varphi+i \sin \varphi)$z=r(cosφ+isinφ)。

复数的乘法，按照模长相乘，辐角相加的运算规则。除法为乘法的逆运算。

复数的共轭，则为复数对应向量与实轴轴对称的结果。

单位根

欧拉公式

对于任意实数 $x$x，都有

$e^{i x}=\cos x+i \sin x$eix=cosx+isinx

证明略。这个公式可以说明当 $x$x 为实数时，函数 $f(x) = e^{ix}$f(x)=eix 可以再复数屏幕表述一个单位元，且 $x$x 为平面上的一个幅角。

定义

将方程 $z^n = 1$zn=1 的所有复数根称为 $n$n 次单位根，这样的单位根共有 $n$n 个。

$\large e^{\frac{2 \pi k i}{n}} \normalsize \quad(k = 0, 1, \cdots,n-1 )$en2πki​(k=0,1,⋯,n−1)

这些单位根分布在复数平面的单位圆上，构成了一个正 $n$n 边形，它们把单位圆等分成 $n$n 个部分。

简单的，定义 $\omega = e^{\frac{2 \pi k i}{n}}$ω=en2πki​，若有 $n$n 与 $k$k 互质，则 $\omega$ω 称为 $n$n 次本原单位根。$w^x\quad (x = 0, 1, \cdots, n-1)$wx(x=0,1,⋯,n−1) 构成了所有的 $n$n 次单位根。特殊的，当 $k=1$k=1 时，本原单位根记为 $\omega_{n}$ωn​。

性质

单位根的共轭复数为其倒数。

$\omega_{n}^{k}=\cos \frac{2 \pi k}{n} +i \sin \frac{2 \pi k}{n} \\ \ \\ \omega_{n}^{k}=\omega_{2 n}^{2 k} \\ \ \\ \omega_{n}^{\frac{n}{2}}=e^{i \pi}=-1 \\ \ \\ \omega_{n}^{k+\frac{n}{2}}=-\omega_{n}^{k}$ωnk​=cosn2πk​+isinn2πk​ ωnk​=ω2n2k​ ωn2n​​=eiπ=−1 ωnk+2n​​=−ωnk​

Cooley-Tukey 算法

离散傅里叶变换 DFT

多项式有两个表示方法——系数表示法和点值表示法。它们各有千秋，系数表示法适用范围广，但乘法操作为 $O(n^2)$O(n2)；点值表示法适用范围窄，但乘法操作为线性。那么能不能在两个表示法相互转换呢？

假设现在有一个 $n-1$n−1 次多项式，$A(x) = \sum_{i = 0}^n a_i \times x_i$A(x)=∑i=0n​ai​×xi​。为了方便，设 $n = 2^m, m \in \mathbb{N}$n=2m,m∈N，不足在高位补 $0$0。

将 $n$n 个 $n$n 次单位根 $\omega_{n}^{0}, \omega_{n}^{1}, \ldots, \omega_{n}^{n-1}$ωn0​,ωn1​,…,ωnn−1​ 代入 $A(x)$A(x) 将其转换成点值表达 $A\left(\omega_{n}^{k}\right)=\sum_{i=0}^{n-1} a_{i} \omega^{k i}$A(ωnk​)=∑i=0n−1​ai​ωki。 点值向量 $\vec{y}=\left(A\left(\omega_{n}^{0}\right), A\left(\omega_{n}^{1}\right), \ldots, A\left(\omega_{n}^{n-1}\right)\right)$y​=(A(ωn0​),A(ωn1​),…,A(ωnn−1​)) 称为系数向量 $\vec{a}=\left(a_{0}, a_{1}, \ldots, a_{n-1}\right)$a=(a0​,a1​,…,an−1​) 的离散傅里叶变换，记为 $\vec{y}= \operatorname{DFT}_{n}(\vec{a})$y​=DFTn​(a)。

可以发现离散傅里叶变换的复杂度为 $O(n^2)$O(n2) 的。Cooley-Tukey 算法对 DFT 进行了优化，它将每一项进行了奇偶分类。当 $k < \frac{n}{2}$k<2n​ 时

根据性质 $\omega_{n}^{k}=\omega_{2 n}^{2 k}$ωnk​=ω2n2k​

$A\left(\omega_{n}^{k}\right) =\sum_{i=0}^{n-1} a_{i} \omega_{n}^{k i} \\ =\sum_{i=0}^{n-1} a_{2 i} \omega_{n}^{2 k i}+\omega_{n}^{k} \sum_{i=0}^{\frac{n}{2}-1} a_{2 i+1} \omega_{n}^{2 k i} \\ = \sum_{i=0}^{\frac{n}{2}-1} a_{2 i} \omega_{\frac{n}{2}}^{k i}+\omega_{n}^{k} \sum_{i=0}^{\frac{n}{2}-1} a_{2 i+1} \omega_{\frac{n}{2}}^{k i}$A(ωnk​)=i=0∑n−1​ai​ωnki​=i=0∑n−1​a2i​ωn2ki​+ωnk​i=0∑2n​−1​a2i+1​ωn2ki​=i=0∑2n​−1​a2i​ω2n​ki​+ωnk​i=0∑2n​−1​a2i+1​ω2n​ki​

根据性质 $\omega_{n}^{k+\frac{n}{2}}=-\omega_{n}^{k}$ωnk+2n​​=−ωnk​

$\begin{aligned} A\left(\omega_{n}^{k+\frac{n}{2}}\right) &=\sum_{i=0}^{\frac{n}{2}-1} a_{2 i} \omega_{\frac{n}{2}}^{k i}+\omega_{n}^{k+\frac{n}{2}} \sum_{i=0}^{\frac{n}{2}-1} a_{2 i+1} \omega_{\frac{n}{2}}^{k i} \\ &=\sum_{i=0}^{\frac{n}{2}-1} a_{2 i} \omega_{\frac{n}{2}}^{k i}-\omega_{n}^{k} \sum_{i=0}^{\frac{n}{2}-1} a_{2 i+1} \omega_{\frac{n}{2}}^{k i} \end{aligned}$A(ωnk+2n​​)​=i=0∑2n​−1​a2i​ω2n​ki​+ωnk+2n​​i=0∑2n​−1​a2i+1​ω2n​ki​=i=0∑2n​−1​a2i​ω2n​ki​−ωnk​i=0∑2n​−1​a2i+1​ω2n​ki​​

可以发现左右两个式子相似。令左式的系数表达式为 $A_1$A1​，右式的系数表达式为 $A_2$A2​。

$A\left(\omega_{n}^{k}\right) = A_1(\omega_{\frac{n}{2}}^k) + A_2(\omega_{\frac{n}{2}}^k) \\ \ \\ A\left(\omega_{n}^{k+\frac{n}{2}}\right) = A_1(\omega_{\frac{n}{2}}^k) - A_2(\omega_{\frac{n}{2}}^k)$A(ωnk​)=A1​(ω2n​k​)+A2​(ω2n​k​) A(ωnk+2n​​)=A1​(ω2n​k​)−A2​(ω2n​k​)

这意味着，如果能求出 $A_1(\omega_{\frac{n}{2}}^k)$A1​(ω2n​k​) 和 $A_2(\omega_{\frac{n}{2}}^k)$A2​(ω2n​k​)，那么可以直接求出 $A\left(\omega_{n}^{k}\right)$A(ωnk​) 和 $A\left(\omega_{n}^{k+\frac{n}{2}}\right)$A(ωnk+2n​​)

于是，我们只要算出 $a_{0}, a_{2}, \ldots, a_{2 n-2}$a0​,a2​,…,a2n−2​ 和 $a_{1}, a_{3}, \ldots, a_{2 n-1}$a1​,a3​,…,a2n−1​ 的 DFT 就可以算出 $a_{0}, a_{1}, \ldots a_{2 n-1}$a0​,a1​,…a2n−1​ 的 DFT，递归下去即可。时间复杂度 $O\left(n \log _{2} n\right)$O(nlog2​n)。

傅里叶逆变换 IDFT

刚刚计算的是 $\vec{y}=\operatorname{DFT}_{n}(\vec{a}),$y​=DFTn​(a), 可以将多项式转化成点值表示，现在为了将点值表示转化成系数表示，需要计算 IDFT，它是 DFT 的逆运算，也就是插值的过程。提起插值，可能先想起拉格朗日插值，不过 FFT 不需要它。

考虑解一个线性同余方程组

$\left\{\begin{array}{ccccc} a_{0}\left(\omega_{n}^{0}\right)^{0} & +\cdots & +a_{n-2}\left(\omega_{n}^{0}\right)^{n-2} & +a_{n-1}\left(\omega_{n}^{0}\right)^{n-1} & =A\left(\omega_{n}^{0}\right) \\ a_{0}\left(\omega_{n}^{1}\right)^{0} & +\cdots & +a_{n-2}\left(\omega_{n}^{1}\right)^{n-2} & +a_{n-1}\left(\omega_{n}^{1}\right)^{n-1} & =A\left(\omega_{n}^{1}\right) \\ \vdots & \vdots & \vdots & \vdots \\ a_{0}\left(\omega_{n}^{n-1}\right)^{0} & +\cdots & +a_{n-2}\left(\omega_{n}^{n-1}\right)^{n-2} & +a_{n-1}\left(\omega_{n}^{n-1}\right)^{n-1} & =A\left(\omega_{n}^{n-1}\right) \end{array}\right.$⎩⎪⎪⎪⎪⎨⎪⎪⎪⎪⎧​a0​(ωn0​)0a0​(ωn1​)0⋮a0​(ωnn−1​)0​+⋯+⋯⋮+⋯​+an−2​(ωn0​)n−2+an−2​(ωn1​)n−2⋮+an−2​(ωnn−1​)n−2​+an−1​(ωn0​)n−1+an−1​(ωn1​)n−1⋮+an−1​(ωnn−1​)n−1​=A(ωn0​)=A(ωn1​)=A(ωnn−1​)​

写成矩阵的形式

$\left[\begin{array}{cccc} \left(\omega_{n}^{0}\right)^{0} & \left(\omega_{n}^{0}\right)^{1} & \cdots & \left(\omega_{n}^{0}\right)^{n-1} \\ \left(\omega_{n}^{1}\right)^{0} & \left(\omega_{n}^{1}\right)^{1} & \cdots & \left(\omega_{n}^{1}\right)^{n-1} \\ \vdots & \vdots & \ddots & \vdots \\ \left(\omega_{n}^{n-1}\right)^{0} & \left(\omega_{n}^{n-1}\right)^{1} & \cdots & \left(\omega_{n}^{n-1}\right)^{n-1} \end{array}\right]\left[\begin{array}{c} a_{0} \\ a_{1} \\ \vdots \\ a_{n-1} \end{array}\right]=\left[\begin{array}{c} A\left(\omega_{n}^{0}\right) \\ A\left(\omega_{n}^{1}\right) \\ \vdots \\ A\left(\omega_{n}^{n-1}\right) \end{array}\right]$⎣⎢⎢⎢⎢⎡​(ωn0​)0(ωn1​)0⋮(ωnn−1​)0​(ωn0​)1(ωn1​)1⋮(ωnn−1​)1​⋯⋯⋱⋯​(ωn0​)n−1(ωn1​)n−1⋮(ωnn−1​)n−1​⎦⎥⎥⎥⎥⎤​⎣⎢⎢⎢⎡​a0​a1​⋮an−1​​⎦⎥⎥⎥⎤​=⎣⎢⎢⎢⎡​A(ωn0​)A(ωn1​)⋮A(ωnn−1​)​⎦⎥⎥⎥⎤​

左矩阵即为系数矩阵 $V$V，要求它的逆矩阵 $V^{-1}$V−1。可以发现，这个矩阵可以使用范德蒙德矩阵的逆矩阵公式，令 $D$D 为 $V$V 的伴随矩阵。

${V}^{-1}=\frac{1}{n} {D}=\frac{1}{n}\left[\begin{array}{cccc} \left(\omega_{n}^{-0}\right)^{0} & \left(\omega_{n}^{-0}\right)^{1} & \cdots & \left(\omega_{n}^{-0}\right)^{n-1} \\ \left(\omega_{n}^{-1}\right)^{0} & \left(\omega_{n}^{-1}\right)^{1} & \cdots & \left(\omega_{n}^{-1}\right)^{n-1} \\ \vdots & \vdots & \ddots & \vdots \\ \left(\omega_{n}^{-(n-1)}\right)^{0} & \left(\omega_{n}^{-(n-1)}\right)^{1} & \cdots & \left(\omega_{n}^{-(n-1)}\right)^{n-1} \end{array}\right]$V−1=n1​D=n1​⎣⎢⎢⎢⎢⎢⎡​(ωn−0​)0(ωn−1​)0⋮(ωn−(n−1)​)0​(ωn−0​)1(ωn−1​)1⋮(ωn−(n−1)​)1​⋯⋯⋱⋯​(ωn−0​)n−1(ωn−1​)n−1⋮(ωn−(n−1)​)n−1​⎦⎥⎥⎥⎥⎥⎤​

证明较为复杂，我们简单验证一下，设 $E = V \times D$E=V×D，则有

$E_{i, j}=\sum_{k=0}^{n-1} D_{i k} V_{k j}=\sum_{k=0}^{n-1} \omega_{n}^{-i k} \omega_{n}^{k j}=\sum_{k=0}^{n-1} \omega_{n}^{k(j-i)}$Ei,j​=k=0∑n−1​Dik​Vkj​=k=0∑n−1​ωn−ik​ωnkj​=k=0∑n−1​ωnk(j−i)​

当 $i=j$i=j 时，显然 $E_{i,j}=n$Ei,j​=n。当 $i \neq j$i​=j 时，有

$E_{i, j}=\sum_{k=0}^{n-1}\left(\omega_{n}^{j-i}\right)^{k}=\frac{1-\left(\omega_{n}^{j-i}\right)^{n}}{1-\omega_{n}^{j-i}}=0$Ei,j​=k=0∑n−1​(ωnj−i​)k=1−ωnj−i​1−(ωnj−i​)n​=0

因此 $D$D 的主对角线都是 $n$n，其它位置都是 $0$0，除以 $n$n 后即为单位矩阵，证毕。

综上所述，可以得到

$\left[\begin{array}{c} a_{0} \\ a_{1} \\ \vdots \\ a_{n-1} \end{array}\right]=\frac{1}{n}\left[\begin{array}{cccc} \left(\omega_{n}^{-0}\right)^{0} & \left(\omega_{n}^{-0}\right)^{1} & \cdots & \left(\omega_{n}^{-0}\right)^{n-1} \\ \left(\omega_{n}^{-1}\right)^{0} & \left(\omega_{n}^{-1}\right)^{1} & \cdots & \left(\omega_{n}^{-1}\right)^{n-1} \\ \vdots & \vdots & \ddots & \vdots \\ -(n-1) & \left(\omega_{n}^{-(n-1)}\right)^{1} & \cdots & \left(\omega_{n}^{-(n-1)}\right)^{n-1} \end{array}\right]\left[\begin{array}{c} A\left(\omega_{n}^{0}\right) \\ A\left(\omega_{n}^{1}\right) \\ \vdots \\ A\left(\omega_{n}^{n-1}\right) \end{array}\right]$⎣⎢⎢⎢⎡​a0​a1​⋮an−1​​⎦⎥⎥⎥⎤​=n1​⎣⎢⎢⎢⎢⎢⎡​(ωn−0​)0(ωn−1​)0⋮−(n−1)​(ωn−0​)1(ωn−1​)1⋮(ωn−(n−1)​)1​⋯⋯⋱⋯​(ωn−0​)n−1(ωn−1​)n−1⋮(ωn−(n−1)​)n−1​⎦⎥⎥⎥⎥⎥⎤​⎣⎢⎢⎢⎡​A(ωn0​)A(ωn1​)⋮A(ωnn−1​)​⎦⎥⎥⎥⎤​

IDFT 就相当于把 DFT 过程中的 $w_n^i$wni​ 换成 $w_n^{-i}$wn−i​，然后做一次 DFT，将结果除以 $n$n 即可。

递归实现

迭代实现

NTT

原根