深入理解 FFT

理论前置

知道啥是多项式（即 \(f(x)=\displaystyle\sum_{i=0}^{n-1}f_ix^i\) 这一类东西）。

知道啥是多项式的卷积（即 \((f\times g)(x)=h(x)\)，其中 \(h_i=\displaystyle\sum_{k=0}^if_kg_{i-k}\)）。

知道啥是复数（即 \(a+bi\) 这种，其中 \(i^2=-1\)）。

知道啥是单位根（即 \(x^n=1\) 的所有复数解，记作 \(\omega_n\)）。

引入

随着 OI 中毒瘤出题人和毒瘤数数题的变多，我们有时候会需要求两个多项式的卷积，根据定义暴力计算，我们有一个 \(\mathcal{O}(n^2)\) 的优秀做法。

显然毒瘤出题人不会止步于此，于是我们需要更快的算法。

这时就需要引入多项式的另外一种表示方法了。

点值表示法

对于 \(f(x)\)，我们不仅可以通过系数表示法（\(f(x)=\displaystyle\sum_{i=0}^{n-1}f_ix^i\)）来表示它，还可以通过点值表示法来表示它。

具体的，对于一个 \(n-1\) 次的多项式 \(f(x)\)，我们可以通过 \(n\) 个互不相同的点唯一确定这个多项式，即 \(f(x)=[(x_0,y_0),(x_1,y_1),\cdots,(x_{n-1},y_{n-1})]\)，其中 \(y_i=f(x_i)\)。

这种表示法对我们来说是很好的，因为对于两个多项式 \(f(x),g(x)\)，只要把点值表示中 \(y_i\) 两两相乘即可求出 \(h(x)\) 的点值表达（注意：由于 \(h(x)\) 是 \(n+m-2\) 次多项式，所以需要 \(n+m-1\) 个点才能唯一确定 \(h(x)\)）。这样，我们只需要 \(\mathcal{O}(n)\) 的时间复杂度即可求出两个多项式的卷积。

现在让我们整理一下算法流程：

1. 输入两个多项式 \(f(x),g(x)\)
2. 把 \(f(x),g(x)\) 各转成点值表达法（需要 \(n+m-1\) 个点）
3. 把点值表达法的 \(y_i\) 两两相乘
4. 把点值表达法转成 \(h(x)\) 的系数表达法
5. 输出多项式 \(h(x)\)。

现在为止一切都很好，但是我们会发现还存在一个问题：系数表达法和点值表达法的转换依旧是 \(O(n^2)\) 的，所以这个算法除了增大常数之外没有任何好处。

你先别急，先看完下面再说

系数表达法 -> 点值表达法

让我们先从简单的多项式入手，比如 \(f(x)=x^2\)，我们如何在较快的时间内求出它的点值表达呢？

很简单，我们可以代入一系列的正负值 \(\{\pm x_0,\pm x_1,\cdots,\pm x_{n/2-1}\}\)，因为我们有 \(f(x)=f(-x)\)，所以我们只需要求出 \(\dfrac{n}{2}\) 个点的点值表达了。

稍难一点，设 \(f(x)=x^3\)，由于 \(f(-x)=-f(x)\)，所以我们依旧可以参考上面的流程。

一般的，对于多项式 \(f(x)=\displaystyle\sum_{i=0}^{n-1}f_ix^i\)，我们可以把它分为两半，即 \(f_0(x)=\displaystyle\sum_{i=0}^{n/2-1}f_{2i}x^i\) 和 \(f_1(x)=\displaystyle\sum_{i=0}^{n/2-1}f_{2i+1}x^i\)，那么我们有

\[\begin{aligned} f(x_i)&=f_0(x_i^2)+x_if_1(x_i^2)\\ f(-x_i)&=f_0(x_i^2)-x_if_1(x_i^2)\\ \end{aligned} \]

可以发现对于 \(f_0\) 和 \(f_1\)，这等价于原问题：多点求值，因此我们可以分治下去求解。

看上去非常不错，但事实上这种做法暗藏着一个陷阱：经过一轮迭代后 \(x_i^2\) 均为正值，因此不能再套用这种做法。

怎么解决呢？

可以发现我们想要的是一系列值 \(x_0,x_1,\cdots,x_{n-1}\)，使得 \(x_{2i}=-x_{2i+1}\)，\(x_{4i}^2=-x_{4i+2}^2\)，\(x_{8i}^4=-x_{8i+4}^4\cdots\)

这时就需要引入复数了。

不妨设 \(x_0=1\)，则显然 \(x_1=-1\)，那么有 \(1=x_2^2\)，由于 \(x_i\) 互不相等，故 \(x_2,x_3=i,-i\)，以此类推。我们会发现 \(x_i\) 实际上就是 \(x^n=1\) 的全部复数解，或者说 \(\omega_n\)。

怎么计算 \(\omega_n\) 呢？可以发现单位根其实就是把单位圆进行了 \(n\) 等分，即 \(\omega_n^k=e^{i(2k\pi/n)}=\cos(2k\pi/n)+i\sin(2k\pi/n)\)。

至此，我们就可以通过分治将此算法做到 \(\mathcal{O}(n\log n)\) 了。

点值表达法 -> 系数表达法

考虑系数表示法变成点值表示法的另外一种方式：线性代数。

我们有

\[\begin{cases} f_0+x_0f_1+x_0^2f_2+\cdots+x_0^{n-1}f_{n-1}=f(x_0)\\ f_0+x_1f_1+x_1^2f_2+\cdots+x_1^{n-1}f_{n-1}=f(x_1)\\ \vdots\\ f_0+x_{n-1}f_1+x_{n-1}^2f_2+\cdots+x_{n-1}^{n-1}f_{n-1}=f(x_{n-1})\\ \end{cases} \]

使用矩阵表示即为

\[\begin{bmatrix} f(x_0)\\ f(x_1)\\ \vdots\\ f(x_{n-1}) \end{bmatrix} = \begin{bmatrix} 1 & x_0 & x_0^2 & x_0^3 & \cdots & x_0^{n-1}\\ 1 & x_1 & x_1^2 & x_1^3 & \cdots & x_1^{n-1}\\ 1 & x_2 & x_2^2 & x_2^3 & \cdots & x_2^{n-1}\\ \vdots&\vdots&\vdots&\vdots&\ddots&\vdots\\ 1 & x_{n-1} & x_{n-1}^2 & x_{n-1}^3 & \cdots & x_{n-1}^{n-1}\\ \end{bmatrix} \begin{bmatrix} f_0\\ f_1\\ \vdots\\ f_{n-1} \end{bmatrix} \]

代入 \(x_k=\omega_n^k\) 可得

\[\begin{bmatrix} f(x_0)\\ f(x_1)\\ \vdots\\ f(x_{n-1}) \end{bmatrix} = \begin{bmatrix} 1 & 1 & 1 & 1 & \cdots & 1\\ 1 & \omega_n & \omega_n^2 & \omega_n^3 & \cdots & \omega_n^{n-1}\\ 1 & \omega_n^2 & \omega_n^4 & \omega_n^6 & \cdots & \omega_n^{2(n-1)}\\ \vdots&\vdots&\vdots&\vdots&\ddots&\vdots\\ 1 & \omega_n^{n-1} & \omega_n^{2(n-1)} & \omega_n^{3(n-1)} & \cdots & \omega_n^{(n-1)(n-1)}\\ \end{bmatrix} \begin{bmatrix} f_0\\ f_1\\ \vdots\\ f_{n-1} \end{bmatrix} \]

那么如果我们想要进行 IFFT，本质上就是对上面这个大矩阵求逆。

由于这是一个有一些很好的性质的矩阵，所以它的逆也非常的漂亮：

\[\begin{bmatrix} 1 & 1 & 1 & 1 & \cdots & 1\\ 1 & \omega_n & \omega_n^2 & \omega_n^3 & \cdots & \omega_n^{n-1}\\ 1 & \omega_n^2 & \omega_n^4 & \omega_n^6 & \cdots & \omega_n^{2(n-1)}\\ \vdots&\vdots&\vdots&\vdots&\ddots&\vdots\\ 1 & \omega_n^{n-1} & \omega_n^{2(n-1)} & \omega_n^{3(n-1)} & \cdots & \omega_n^{(n-1)(n-1)}\\ \end{bmatrix}^{-1}=\frac{1}{n}\begin{bmatrix} 1 & 1 & 1 & 1 & \cdots & 1\\ 1 & \omega_n^{-1} & \omega_n^{-2} & \omega_n^{-3} & \cdots & \omega_n^{-(n-1)}\\ 1 & \omega_n^{-2} & \omega_n^{-4} & \omega_n^{-6} & \cdots & \omega_n^{-2(n-1)}\\ \vdots&\vdots&\vdots&\vdots&\ddots&\vdots\\ 1 & \omega_n^{-(n-1)} & \omega_n^{-2(n-1)} & \omega_n^{-3(n-1)} & \cdots & \omega_n^{-(n-1)(n-1)}\\ \end{bmatrix} \]

不难发现逆矩阵和原矩阵几乎一模一样，也可以使用 FFT 的方式计算。

至此，我们可以在 \(\mathcal{O}(n\log n)\) 的时间内计算两个多项式的乘积。

参考视频