浅谈FFT(快速傅里叶变换)

本文主要简单写写自己在算法竞赛中学习FFT的经历以及一些自己的理解和想法。

FFT的介绍以及入门就不赘述了，网上有许多相关的资料，入门的话推荐这篇博客：FFT（最详细最通俗的入门手册），里面介绍得很详细。

为什么要学习FFT呢？因为FFT能将多项式乘法的时间复杂度由朴素的$O(n^2)$降到$O(nlogn)$，这相当于能将任意形如$f[k]=\sum\limits _{i+j=k}f[i]\cdot f[j]$的转移方程的计算在$O(nlogn)$的时间内完成。因此对于想要进阶dp的同学来说，FFT是必须掌握的技能之一。（虽然在赛场上可能没什么用武之地）

我学习FFT的过程也是比较曲折的，从接触到真正理解它的原理前前后后经历了半年的时间。（实际上我从去年接触了FFT之后就一直把它当做一个黑盒算法来用，研究的事就扔到一边了，只是偶尔简单推算过几次公式，直到这个月初才开始深入学习它的原理）

由于本人才疏学浅，所以自己的叙述若存在一些错误或者不足之处，敬请读者指正。

 首先FFT的作用是什么？可以将多项式的系数表达式转化成点值表达式（或者反过来，方法都是一样的）。FFT(a,n)的作用是将多项式a(系数表达式)从$w_{n}^{0}$到$w_{n}^{n-1}$的所有根对应的取值求出来。也就是说，设$f(x)=\sum\limits_{i=0}^{n-1}a[i]\cdot x^i$，经过FFT变换后，a[i]变成了$f(w_{n}^{i})$。

这个利用单位根来表示的点值表达式的一个好处是如果已知FFT(a0,n/2)以及FFT(a1,n/2)（a0为a的偶数次项所构成的多项式，a1为a的奇数次项所构成的多项式），则根据性质$\left\{\begin{matrix}\begin{aligned}&a[i]=a_0[i]+w_{n}^{i}\cdot a_1[i]\\&a[i+\frac{n}{2}]=a_0[i]-w_{n}^{i}\cdot a_1[i]\end{aligned} \end{matrix}\right.$可以在$O(n)$的时间内算出数组a的值。

 为什么要用单位根呢？因为对于任意的数组长度n，在FFT的过程中使用单位根都只需要计算n个不同变量的值，与数组长度是线性相关的，而且一定能保证取到n个不同的值。而假如取2,3,4这样的数的话，在对任意子数组进行FFT时仍需计算n个不同变量的值，这样的话总的复杂度仍为$O(n^2)$，没有丝毫降低。而假如取-1，1这样的数，虽然只需要计算常数个变量的值了，但无论如何只能取到一两个变量的值，也就是只能确定两点，无法确定一个具有n个维度的多项式。

接下来就是代码实现了。

首先我们做一下预处理：

1 typedef double db;
2 const db pi=acos(-1);

把double定义成db的作用，一是可以简化代码，二是需要调整精度的时候可以很方便地替换成其他变量类型，比如long double。

FFT的运算要用到复数，这就意味着我们必须找到一个能够代表复数的变量类型。图方便的话，C++库中内置的complex类就够用了。不过还是推荐自己写一个结构体，比C++自带的要快很多，而且也很好写。

由于复数是一个二元组，和二维平面上的点非常类似，因此可以直接套用二维几何中的点的结构体代码。加减数乘等操作都完全一样，只是多了个乘法。但这并不影响它的几何意义，因为在计算几何中两向量乘法我一般喜欢用dot(点积)和cross(叉积)两个函数来表示。此外，乘法运算符也可以表示坐标的旋转。

复数(点)的结构体代码如下：

1 struct P {
2     db x,y;
3     P operator+(const P& b) {
   return {x+b.x,y+b.y};}
4     P operator-(const P& b) {
   return {x-b.x,y-b.y};}
5     P operator*(const P& b) {
   return {x*b.x-y*b.y,x*b.y+y*b.x};}
6     P operator/(db b) {
   return {x/b,y/b};}
7 }

接下来就是FFT的实现了。有了FFT的基本概念和点的表示方法之后，我们不难写出这样的代码：（f为1代表正变换（取值），f为1代表逆变换（插值））

 1 void FFT(P* a,int n,int f) {
 2     if(n==1)return;
 3     static P b[N];
 4     for(int i=0; i<n; i+=2)b[i/2]=a[i],b[(i+n)/2]=a[i+1];
 5     for(int i=0; i<n; ++i)a[i]=b[i];
 6     FFT(a,n/2,f),FFT(a+n/2,n/2,f);
 7     P wn= {cos(2*pi/n),f*sin(2*pi/n)},w= {
   1,0};
 8     for(int i=0; i<n/2; ++i,w=w*wn) {
 9         P x=a[i],y=w*a[i+n/2];
10         a[i]=x+y,a[i+n/2]=x-y;
11     }