快速傅里叶变换：从多项式乘法到信号处理的数学奇迹

引言：一个改变世界的算法

快速傅里叶变换（Fast Fourier Transform,）被誉为 20 世纪最重要的算法之一。从无线通信的 OFDM 调制，到 JPEG 图像的压缩编码；从 MRI 医疗成像，到引力波信号的频谱分析——FFT 无处不在。它之所以如此深远，根本在于它将一个看似不可能高效解决的问题（多项式乘法、信号处理、谱分析）从繁琐的 推演降低到了优雅。

本文将从多项式乘法的视角出发，逐步揭示 FFT 分治策略的数学本质，推导其递推结构，并探讨其在信号处理等现代技术中的应用。

一、从多项式乘法说起

1.1 多项式的两种表示

一个 n-1 次多项式 $A(x) = a_0 + a_1 x + a_2 x^2 + \cdots + a_{n-1} x^{n-1}$ 有两种等价表示方式：

• 系数表示法：直接用系数向量 $(a_0, a_1, \ldots, a_{n-1})$ 描述多项式。在这种表示下，多项式加法和求值操作的时间复杂度分别为 O(n) 和 O(n)（使用 Horner 法则），但乘法需要 O(n²) 的时间。
• 点值表示法：选取 n 个不同的点 $x_0, x_1, \ldots, x_{n-1}$，将多项式表示为 n 个点值对 $\{(x_0, A(x_0)), (x_1, A(x_1)), \ldots, (x_{n-1}, A(x_{n-1}))\}$。在这种表示下，多项式乘法只需 O(n)——只需将对应点值相乘即可！

1.2 朴素算法的瓶颈

如果我们能将系数表示法转换为点值表示法（求值），做乘法后，再将点值表示法转回系数表示法（插值），就能在 O(n log n) 时间内完成乘法！关键问题在于：求值和插值能否在 O(n log n) 时间内完成？直接计算需要 O(n²)，这正是 FFT 解决的问题。

二、FFT 的数学基础

2.1 单位根

FFT 选取的求值点是 n 次单位根（n 为 2 的幂，必要时补零）：$\omega_n^k = e^{2\pi i k / n}$，其中 $k = 0, 1, \ldots, n-1$。单位根具有关键的数学性质：

• 周期性：$\omega_n^{k+n} = \omega_n^k$
• 对称性：$\omega_n^{k+\frac{n}{2}} = -\omega_n^k$
• 消去引理：$\omega_{dn}^{dk} = \omega_n^k$
• 减半引理：$(\omega_n^{k+\frac{n}{2}})^2 = (\omega_n^k)^2$

特别是消去引理和减半引理，它们让 FFT 得以使用分治策略——将 n 个点的求值问题缩减为两个 n/2 个点的子问题。

2.2 分治求值（DFT 的 Cooley-Tukey 算法）

将多项式 $A(x)$ 按奇偶次幂拆分：$A(x) = A_{even}(x^2) + x \cdot A_{odd}(x^2)$，其中 $A_{even}$ 包含偶数下标系数，$A_{odd}$ 包含奇数下标系数。

那么在 $\omega_n^k$ 处求值可以写为：

$$A(\omega_n^k) = A_{even}(\omega_n^{2k}) + \omega_n^k \cdot A_{odd}(\omega_n^{2k})$$

由于 $\omega_n^{2k} = \omega_{n/2}^k$，我们只需要计算 $A_{even}$ 和 $A_{odd}$ 在 n/2 次单位根处的取值，然后合并即可。这产生了如下递归结构：

• 递归基：当 n=1 时，$A(\omega_1^0) = a_0$
• 递归步骤：分别计算 $A_{even}$ 和 $A_{odd}$ 在 n/2 次单位根处的取值，利用减半引理合并结果
• 蝴蝶操作（Butterfly Operation）：每一层的合并操作形如 $(a, b) \to (a + \omega b, a - \omega \omega b)$

2.3 复杂度推导

设 T(n) 为求值 n 个点的时间复杂度，由分治策略可得递推式：$T(n) = 2T(n/2) + O(n)$。根据主定理，$T(n) = O(n \log n)$。求值和插值（只需将 $\omega_n$ 替换为 $\omega_n^{-1}$）均为 $O(n \log n)$，因此多项式乘法的总复杂度为 $O(n \log n)$。

三、FFT 的迭代实现与位逆序

3.1 递归转迭代

虽然递归公式优雅，但在实际代码中，递归调用带来的函数调用开销不可忽视。Cooley-Tukey 算法可以通过迭代方式实现：从大小为 1 的子问题开始，逐层合并，直到覆盖整个多项式。这种方法不仅减少了调用开销，还能更好地利用 CPU 缓存。

3.2 位逆序置换（Bit-Reversal Permutation）

迭代版本的关键在于理解系数的排列规律。当我们将多项式按奇偶次幂反复拆分时，最终系数下标的排列恰好是原始下标的二进制位逆序。例如，当 n=8 时：

• 下标 1（二进制 001）→ 位逆序后为 4（二进制 100）
• 下标 3（二进制 011）→ 位逆序后为 6（二进制 110）

因此迭代 FFT 通常以一次位逆序置换开始，然后从大小为 2 的子问题开始逐层合并。

3.3 经典迭代代码结构

迭代 FFT 的核心代码段如下（C 语言伪代码）：

// 1. 位逆序重排
for (int i = 0; i < n; i++)
    if (i < rev[i]) swap(a[i], a[rev[i]]);

// 2. 逐层合并，蝴蝶操作
for (int len = 2; len <= n; len <<= 1) {
    double ang = 2 * PI / len * (invert ? -1 : 1);
    wlen = complex(cos(ang), sin(ang));
    for (int i = 0; i < n; i += len) {
        w = complex(1);
        for (int j = 0; j < len / 2; j++) {
            u = a[i+j];
            v = a[i+j+len/2] * w;
            a[i+j] = u + v;
            a[i+j+len/2] = u - v;
            w *= wlen;
        }
    }
}

四、逆 FFT：从点值回到系数

求值完成后，我们需要通过逆 FFT（IFFT）将点值表示转换回系数表示。数学上，DFT 矩阵的逆是 DFT 矩阵的共轭转置除以 n。因此 IFFT 的实现只需：

1. 对输入序列进行 FFT，但单位根取 $\omega_n^{-k}$
2. 将结果除以 n

这意味着 FFT 和 IFFT 共享绝大部分代码。在多项式乘法的场景中，IFFT 之后需要将各系数四舍五入到最近整数（因为浮点运算引入了微小误差）。

五、工程挑战与优化

5.1 复数运算的浮点精度

经典的迭代 FFT 使用 double 精度浮点数，但对于 n > 10⁶ 的大规模应用，浮点误差可能使结果失效。常见的解决方案包括：

• 使用 NTT（数论变换）：在模素数域上进行类似 FFT 的卷积运算，避免了浮点误差。要求 p 是费马素数（如 $p = 998244353 = 119 \cdot 2^{23} + 1$），确保单位根存在。
• 高精度拆分：将大数拆分为高低两部分分别卷积，合并时处理进位。

5.2 缓存友好性

FFT 的内存访问模式并非完全连续。当 n 较大时，蝴蝶操作在距离较远的地址间进行数据交换，可能导致缓存未命中。实际的高性能 FFT 库（如 FFTW）采用分块策略，将大问题分解为能完全放入 L1/L2 缓存的子问题，并使用代码生成（code generation）和运行时调优（runtime adaptation）来自动选择最优的分解策略。

5.3 SIMD 加速

现代 CPU 的向量指令（AVX/NEON）可以并行处理多个蝴蝶操作。实数 FFT 可以进一步压缩：通过将两个实数序列打包为复数序列进行一次 FFT，利用对称性分离结果，实现两倍的加速。

六、信号处理中的应用

6.1 频谱分析与滤波

在信号处理中，FFT 将时域信号转换到频域，使得滤波操作变得极为简单：在频域中乘以滤波器的频率响应，再进行 IFFT 即可得到滤波后的信号。这比时域卷积更高效，且能直观地理解信号的频率特性。

6.2 卷积定理

FFT 最重要的应用之一是快速卷积：时域中的卷积等价于频域中的逐点乘法。对于长度为 n 的序列和长度为 m 的序列，直接卷积为 O(nm)，而通过 FFT 实现只需 O((n+m) log(n+m))。在图像处理中，这意味着大尺寸高斯模糊能在与图像大小成线性对数的时间内完成。

6.3 大整数乘法

FFT 还可用于超大整数的乘法（如 Python 大整数和 GMP 库）。将整数视为多项式（基为 10 或 2³²），用 FFT 计算卷积，然后统一处理进位。对于百万位级别的整数，FFT 乘法比传统的 Karatsuba 和 Toom-Cook 算法更快。

七、从 FFT 到多项式技术的延伸

FFT 是多项式计算的基础构件。在此基础上，可以构建：

• NTT-friendly 算法：竞赛编程中的模意义下多项式乘法、生成函数求值
• Schönhage-Strassen 算法：利用 FFT 在 $\mathbb{Z}[\omega]$ 环上实现超大规模整数乘法
• FWT（快速沃尔什变换）：用于处理异或卷积和子集动态规划
• Chirp Z-Transform：将 FFT 推广到任意螺旋路径上的采样

总结

FFT 的精妙之处在于：它通过选取单位根作为求值点，利用数学上的对称性和周期性，将 O(n²) 的问题规约为 O(n log n)。这个平凡而深刻的分治策略，在诞生五十年后的今天，仍然是信号处理、数据压缩和高效计算的核心基石。理解 FFT，不仅是掌握一个算法，更是体会数学之美与工程智慧的完美融合。