FFT 学习笔记

FFT 学习笔记

要求多项式 \(f(x) \times g(x)\) 的乘积。

可以把 \(f(x)\) 和 \(g(x)\) 看做一个 \(n\) 次方程，并且过点 \((x_0, f(x_0)),(x_1, f(x_1)) \dots (x_n, f(x_n))\)。

现在已知 \(f, g\) 的系数表示法。

我们知道 \(f(x) \times g(x)\) 的点值表示法是 \((x_0, f(x_0) \times g(x_0)), (x_1, f(x_1) \times g(x_1)) \dots (x_n, f(x_n) \times g(x_n))\)。

所以 FFT 的核心思想就是将 \(f, g\) 转换成点值表示法，再把乘积的点值表示法转换成系数表示法。

单位根

考虑在复平面上的一个单位圆，进行 \(n\) 等分，得到 \(n\) 个复数。定义 \(\omega_n\) 表示幅角（\(n\) 等分的单位大小），有

\[\omega_n = \cos \frac{2\pi}{n} + i\sin \frac{2\pi}{n} \]

定义 \(\omega_n^k\) 表示从第一个点开始，逆时针的第 \(n\) 个点，那么

\[\omega_n^k = \cos \frac{2k\pi}{n} + i \sin \frac{2k\pi}{n} \]

单位根具有以下性质：

\[\omega_{rn}^{rk} = \omega_n^k \\ \omega_n^{k + \frac{n}{2}} = -\omega_n^k \\ (\omega_n^k) ^ 2 = \omega_n^{2k} \]

对于第一个式子，拆开即可证明。对于第二个式子，容易想到 \(\frac{n}{2}\) 是一个半圆，且 \(\omega_n^k\) 和 \(-\omega_n^k\) 关于原点对称，即可证明。

对于第三个式子的证明：

我们知道欧拉公式：

\[e^{i\theta} = \cos\theta + i\sin\theta \]

所以

\[\omega_n^k = e^{i\frac{2k\pi}{n}} = e^{i2k\pi/n} \\ (\omega_n^k) ^ 2 = e^{i4k\pi/n} \\ \omega_n^{2k} = e^{i4k\pi/n} = (\omega_n^k) ^ 2 \]

证毕。

快速傅里叶变换

对于多项式

\[f(x) = \sum_{i = 0}^{n-1} a_ix^i \in X_{n-1} \]

我们设 \(n = 2^s\)，一般 \(2^s \ge n\) 但是我们可以让高次项系数化为 \(0\)。所以令 \(n = 2^s\)

所以我们可以将 \(f(x)\) 按照系数的奇偶性分成两个部分。令

\[f_1(x) = a_0 + a_2x^1 + a_4x^2 + \dots + a_{n-2}x^{\frac{n}{2} - 1} \\ f_2(x) = a_1 + a_3x^1 + a_5x^2 + \dots + a_{n-1}x^{\frac{n}{2} - 1} \]

我们发现

\[f_1(x^2) = a_0 + a_2x^2 + a_4x^4 + \dots + a_{n-2}x^{n-2} \\ f_2(x^2) = a_1 + a_3x^2 + a_5x^4 + \dots + a_{n-1}x^{n-1} \]

所以

\[xf_2(x^2) = a_1x^1 + a_3x^3 + a_5x^5 + \dots a_{n-1}x^{n} \]

所以

\[f(x) = f_1(x^2) + xf_2(x^2) \]

令 \(i = \omega_n^k\)，所以

\[f(i) = f_1(i^2) + i\cdot f_2(i^2) \\ \because i = \omega_n^k \\ \therefore f(\omega_n^k) = f_1 (\omega_n^{2k}) + \omega_n^k \cdot f_2 (\omega_n^{2k}) \\ \therefore f(\omega_n^k) = f_1 (\omega_{\frac{n}{2}}^k) + \omega_n^k \cdot f_2 (\omega_{\frac{n}{2}}^k) \]

因为

\[w_n^{k + \frac{n}{2}} = -w_n^k \]

所以

\[f(w_n^{k + \frac{n}{2}}) = f_1(w_{\frac{n}{2}}^k) - w_n^k \cdot f_2(w_{\frac{n}{2}}^k) \]

因为平方的存在所以 \(f_1, f_2\) 括号中间的值不变。

前面说过，\(n = 2^s\)。根据上式，要求 \(f(\omega_n^k)\)，就要求 \(f_1(\omega_{\frac{n}{2}}^k), f_1(\omega_{\frac{n}{2}}^k)\)，对于子任务 \(f_1(\omega_{\frac{n}{2}}^k)\)，我们把 \(f_1\) 看做 \(f\)，继续做递归，直到 \(f_1(\omega_1^k) = f_2(\omega_1^k) = 1\)，最后回溯 + 更新答案。

然后我们知道了所有的 \(\omega_n^k\) 和 \(f(\omega_n^k)\)，就得到了 \(f\) 的点值表示。

快速傅里叶逆变换

啊这个东西看不懂啊……不管了，上代码！

#include <bits/stdc++.h>
using namespace std;
using ll = long long;
const int N = 1 << 22;
const double eps = 1e-6, pi = acos(-1.0);
complex<double> a[N], b[N];// 分别表示两个多项式的系数表达法 
int n, m;// 分别表示两个多项式的次数 

void FFT(complex<double> *f, int n, int inv){
// 分别表示要处理的 f，f 的长度，inv = 1 表示系数转点，inv = -1 表示点转系数 
	if(n == 1)	return;
	int mid = n >> 1;
	complex<double> f1[mid + 1], f2[mid + 1];
	for(int i = 0;i <= n;i += 2){// 拆分多项式 
		f1[i >> 1] = f[i];
		f2[i >> 1] = f[i + 1];
	}
	FFT(f1, mid, inv);// 分治处理 
	FFT(f2, mid, inv);
	complex<double> w0(1, 0), wn(cos(2 * pi / n), inv * sin(2 * pi / n));
	// 分别表示 w_n^k 和 w_n（单位根）
	for(int i = 0;i < mid;i++, w0 *= wn){
		f[i] = f1[i] + w0 * f2[i];
		f[i + n/2] = f1[i] - w0 * f2[i];
	}
}

signed main(){
//	ios::sync_with_stdio(false);
//	cin.tie(nullptr);
	scanf("%d%d", &n, &m);
	for(int i = 0;i <= n;i++){
		double x;	scanf("%lf", &x);
		a[i].real(x);// 表示将实部赋值为 x 
	}
	for(int i = 0;i <= m;i++){
		double x;	scanf("%lf", &x);
		b[i].real(x);
	}
	int len = 1 << max(int(ceil(log2(n + m))), 1);// FFT 的长度
	FFT(a, len, 1);// 将 a 转为点值表达
	FFT(b, len, 1);// 将 b 转为点值表达
	for(int i = 0;i <= len;i++)	a[i] = a[i] * b[i];// 点值乘积
	FFT(a, len, -1);// 将 a 转为系数表达 
	for(int i = 0;i <= n + m;i++){
		printf("%.0f ", a[i].real() / len + eps);
//		cout<<int(a[i].real() / len + eps)<<" ";
	}
	return 0;
}