快速傅里叶变换 (FFT) — 从零精通算法与数据结构——Google 面试系统备战 第9篇

第9章:快速傅里叶变换 (FFT)

本章目标

读完本章你会:

  • 解释为什么多项式乘法能从 O(n²) 降到 O(n log n)
  • 理解单位复根 ω_n 的对称性及其在 FFT 中的作用
  • 手写 Cooley-Tukey FFT 和逆 FFT
  • 应用 FFT 解决大整数乘法、信号卷积等问题
  • 将 FFT 模块添加到 algo_toolkit

知识讲解

从一个生活例子开始

两个 n 次多项式相乘:

A(x) = a₀ + a₁x + a₂x² + ... + a_{n-1}x^{n-1}
B(x) = b₀ + b₁x + b₂x² + ... + b_{n-1}x^{n-1}

C(x) = A(x) × B(x) = c₀ + c₁x + ... + c_{2n-2}x^{2n-2}
其中 c_k = Σ_{i=0}^{k} a_i · b_{k-i}   (卷积)

暴力法:每个 c_k 需要约 k 次乘法,总共 O(n²)。

核心洞察: 换个角度看多项式!

一个 n-1 次多项式由两种方式唯一确定:

  • 系数表示: n 个系数 (a₀, a₁, ..., a_{n-1}) —— 做乘法慢
  • 点值表示: n 个点 (x_i, y_i) —— 做乘法只需 O(n)!
两个多项式在相同 x 处的点值:
A 在 x 处: A(x) = y_A
B 在 x 处: B(x) = y_B
→ C(x) = y_A · y_B   ← 一次乘法!

对所有 2n 个点这么做 → O(n)

FFT 的思路就是在这两种表示之间快速转换:

系数表示 ──FFT  O(n log n)──→ 点值表示
                                  │  O(n) 逐点相乘
系数表示 ←─逆FFT O(n log n)─  点值表示

总复杂度:O(n log n)。这就是分治法的力量——与归并排序完全相同的复杂度结构。

工作原理

9.1 单位复根:FFT 的"魔法"来源

n 次单位复根 ω_n 满足 ω_n^n = 1。在复平面上,它们均匀分布在单位圆上:

ω_n^k = e^{2πi·k/n} = cos(2πk/n) + i·sin(2πk/n)

性质:
- ω_n^n = 1              (n 次后回到 1)
- ω_n^{n/2} = -1          (折半性质——关键!)
- ω_n^{k + n/2} = -ω_n^k  (对称性——关键!)
- ω_{dn}^{dk} = ω_n^k     (消去引理)

9.2 FFT 的分治结构

手工示例——n=4 的 FFT 完整计算过程:

对 A(x) = 1 + 2x + 3x² + 4x³,4 次单位复根是 ω₄⁰=1, ω₄¹=i, ω₄²=-1, ω₄³=-i。

Step 1: 抽出偶/奇数项
  A_even(x) = 1 + 3x      (系数来自 a₀, a₂)
  A_odd(x)  = 2 + 4x      (系数来自 a₁, a₃)

Step 2: 递归求 A_even 和 A_odd 在 ω₂ 处的值
  ω₂⁰=1, ω₂¹=-1
  A_even(1) = 1+3·1 = 4   A_even(-1) = 1+3·(-1) = -2
  A_odd(1)  = 2+4·1 = 6   A_odd(-1)  = 2+4·(-1) = -2

Step 3: 蝴蝶操作合成 n=4 的结果
  A(ω₄⁰) = A_even(ω₂⁰) + ω₄⁰ · A_odd(ω₂⁰) = 4 + 1·6 = 10
  A(ω₄¹) = A_even(ω₂¹) + ω₄¹ · A_odd(ω₂¹) = -2 + i·(-2) = -2-2i
  A(ω₄²) = A_even(ω₂⁰) - ω₄⁰ · A_odd(ω₂⁰) = 4 - 1·6 = -2
  A(ω₄³) = A_even(ω₂¹) - ω₄¹ · A_odd(ω₂¹) = -2 - i·(-2) = -2+2i

验证:A(1) = 1+2+3+4 = 10 ✓,这就是 ω₄⁰ 处的值。

FFT 的核心是逐次抽取偶数项和奇数项

A(x) = a₀ + a₁x + a₂x² + a₃x³ + ... + a_{n-1}x^{n-1}

A_even(x) = a₀ + a₂x + a₄x² + ...           (偶数项)
A_odd(x)  = a₁ + a₃x + a₅x² + ...           (奇数项)

则: A(x) = A_even(x²) + x · A_odd(x²)

现在在单位复根上求值:

对于 k = 0, 1, ..., n/2 - 1:

A(ω_n^k) = A_even(ω_{n/2}^k) + ω_n^k · A_odd(ω_{n/2}^k)
A(ω_n^{k+n/2}) = A_even(ω_{n/2}^k) - ω_n^k · A_odd(ω_{n/2}^k)
                                       ↑ 因为 ω_n^{n/2} = -1

美妙之处: 两个点值用了相同的子问题计算结果,只差一个符号!这就把规模 n 的问题变成了两个规模 n/2 的问题,合并代价 O(n)。

递推式: T(n) = 2T(n/2) + O(n) → Θ(n log n) —— 主定理情况 2。

9.3 逆 FFT

点值回到系数也需要 O(n log n),使用几乎相同的算法:

逆 FFT 的魔法:将 ω_n 替换为 ω_n^{-1},最后除以 n

这意味着同一份代码就能做 FFT 和逆 FFT,只需传一个参数控制方向。代码复用之美。

9.4 FFT 的应用

大整数乘法: 两个 n 位数字相乘,视为多项式在 x=10 处的值。用 FFT 做 O(n log n)。

信号卷积: 信号处理中滤波器本质上就是卷积。FFT → 逐点相乘 → IFFT,加速卷积从 O(n²) 到 O(n log n)。

快速多项式乘法: 是许多计算机代数系统(如 SymPy、Mathematica)的底层。


代码实战

include/algo/fft.h

#ifndef ALGO_FFT_H_
#define ALGO_FFT_H_

#include <complex>
#include <vector>

namespace algo {

using Complex = std::complex<double>;

// FFT:coeff → point-value。invert=true 时做逆 FFT
void Fft(std::vector<Complex>& a, bool invert);

// 多项式乘法:卷积
std::vector<Complex> MultiplyPolynomials(
    const std::vector<Complex>& a,
    const std::vector<Complex>& b);

// 大整数乘法(字符串输入输出)
std::string MultiplyBigIntegers(const std::string& num1,
                                const std::string& num2);

}  // namespace algo

#endif  // ALGO_FFT_H_

include/algo/fft_impl.h

#ifndef ALGO_FFT_IMPL_H_
#define ALGO_FFT_IMPL_H_

#include <algorithm>
#include <cmath>
#include <string>

namespace algo {

namespace {

constexpr double kPi = 3.14159265358979323846;

// 位反转——Cooley-Tukey 的迭代实现需要
int BitReverse(int x, int log_n) {
  int result = 0;
  for (int i = 0; i < log_n; ++i) {
    result = (result << 1) | (x & 1);
    x >>= 1;
  }
  return result;
}

}  // namespace

inline void Fft(std::vector<Complex>& a, bool invert) {
  int n = static_cast<int>(a.size());
  if (n <= 1) return;

  // 1. 确保 n 是 2 的幂(补零)
  int log_n = 0;
  while ((1 << log_n) < n) ++log_n;
  if ((1 << log_n) != n) {
    n = 1 << log_n;
    a.resize(n);
  }

  // 2. 位反转置换——将偶数/奇数项分组的效果
  for (int i = 0; i < n; ++i) {
    int j = BitReverse(i, log_n);
    if (i < j) std::swap(a[i], a[j]);
  }

  // 3. 迭代 Cooley-Tukey —— 自底向上合并
  for (int len = 2; len <= n; len <<= 1) {
    double angle = 2.0 * kPi / len * (invert ? -1.0 : 1.0);
    Complex wlen(std::cos(angle), std::sin(angle));

    for (int i = 0; i < n; i += len) {
      Complex w(1.0);
      for (int j = 0; j < len / 2; ++j) {
        Complex u = a[i + j];
        Complex v = a[i + j + len / 2] * w;

        // 蝴蝶操作:两个结果源自相同的两个输入
        a[i + j] = u + v;
        a[i + j + len / 2] = u - v;  // 仅此处的符号!
        w *= wlen;
      }
    }
  }

  // 4. 逆 FFT 时除以 n
  if (invert) {
    for (auto& x : a) x /= static_cast<double>(n);
  }
}

inline std::vector<Complex> MultiplyPolynomials(
    const std::vector<Complex>& a,
    const std::vector<Complex>& b) {
  int n = 1;
  int result_size = static_cast<int>(a.size() + b.size()) - 1;
  while (n < result_size) n <<= 1;

  std::vector<Complex> fa(a.begin(), a.end());
  std::vector<Complex> fb(b.begin(), b.end());
  fa.resize(n);
  fb.resize(n);

  Fft(fa, false);
  Fft(fb, false);

  // O(n) 逐点相乘
  for (int i = 0; i < n; ++i) fa[i] *= fb[i];

  Fft(fa, true);
  fa.resize(result_size);
  return fa;
}

inline std::string MultiplyBigIntegers(const std::string& num1,
                                       const std::string& num2) {
  if (num1 == "0" || num2 == "0") return "0";

  std::vector<Complex> a, b;
  for (auto it = num1.rbegin(); it != num1.rend(); ++it)
    a.emplace_back(*it - '0');
  for (auto it = num2.rbegin(); it != num2.rend(); ++it)
    b.emplace_back(*it - '0');

  auto c = MultiplyPolynomials(a, b);

  // 处理进位
  std::string result;
  int carry = 0;
  for (std::size_t i = 0; i < c.size(); ++i) {
    int digit = static_cast<int>(std::round(c[i].real())) + carry;
    carry = digit / 10;
    result.push_back(static_cast<char>('0' + (digit % 10)));
  }
  while (carry > 0) {
    result.push_back(static_cast<char>('0' + (carry % 10)));
    carry /= 10;
  }

  // 反转(因为我们从低位存到高位)
  std::reverse(result.begin(), result.end());
  return result;
}

}  // namespace algo

#endif  // ALGO_FFT_IMPL_H_

代码关键点:

  • 位反转置换(bit-reversal permutation)将系数按偶/奇分组重新排列——这是迭代 Cooley-Tukey 的精髓
  • 蝴蝶操作:u+vu-v 是 FFT 的效率核心——两个结果从相同的两次读操作导出
  • 逆 FFT 只改两个地方:wlen 的方向(angle 取负)和最后的 ÷ n

本章小结

  1. FFT 的核心洞察:多项式乘法的瓶颈在系数表示,换到点值表示就能 O(n) 乘
  2. FFT 在 O(n log n) 时间内完成系数 ↔ 点值的转换
  3. 单位复根的性质(ω_n^{n/2} = -1)使得蝴蝶操作把 n 规模变成两个 n/2 子问题
  4. 逆 FFT 只比 FFT 多一个符号变化和除以 n——代码完全复用
  5. 应用:多项式乘法、大整数乘法、信号卷积——都是 O(n log n)

关键术语

术语 释义
单位复根 ω_n = e^{2πi/n},满足 ω_n^n = 1 的复数
蝴蝶操作 FFT 的原子操作:从 (u, v) 计算 (u + v·ω, u - v·ω)
位反转 将索引的二进制位倒序,用于 Cooley-Tukey 迭代版
NTT 数论变换——用模素数下的原根替代复数域上的单位复根
posted @ 2026-06-22 01:09  Yobeeo  阅读(18)  评论(0)    收藏  举报