【总结】FFT算法在信息竞赛中的应用

FFT算法本身就是一种优化，优化（类似）卷积运算的时间复杂度

（卷积：∑i,jai∗bj−i∑i,jai∗bj−i）。

FFT的本质，其实是利用复数的一些特殊性质，将一个多项式快速地在点值和系数两种表示方法间来回切换。再利用两个多项式点值表示法相乘的复杂度为O(n)，来达到降时间的目的。

FFT算法的前导概念

首先介绍关于复数的一些定义及性质

复数：一种形如a+bia+bi的数称为复数，其中aa为实部，bb为虚部，ii为虚部单位（i2=−1）（i2=−1）

对于每个a+bia+bi，a−bia−bi为它的共轭复数

其实可以将复数看作一个二项式，其运算法则与实数概念上的二项式几乎相同，只需要注意i即可。

另外，复数也可以用向量来表示：表示为在复平面上(a,b)的向量。向量的模长为a2+b2−−−−−−√a2+b2。

所以，同样地，复数还可以用三角函数来表示：即

a+bi=a2+b2−−−−−−√∗(cosθ+sinθi)a+bi=a2+b2∗(cosθ+sinθi)

在这种表示方法下，两个复数相乘，就可以表示为其角度相加，模长相乘

p∗(cosα+sinαi)∗q∗(cosβ+sinβi)p∗(cosα+sinαi)∗q∗(cosβ+sinβi)

=p∗q∗(sinα∗cosβi+cosα∗sinβi−sinα∗sinβ+cosα∗cosβ)=p∗q∗(sinα∗cosβi+cosα∗sinβi−sinα∗sinβ+cosα∗cosβ)

=p∗q∗(cos(α+β)+sin(α+β)i)=p∗q∗(cos(α+β)+sin(α+β)i)

然后是多项式的一些概念：

首先我们将A(x)=anxn+an−1xn−1+an−2xn−2...a1x+a0A(x)=anxn+an−1xn−1+an−2xn−2...a1x+a0称为一个多项式，同时这种形式也就是系数表示法。

点值表示法：

我们可以将n-1次（包括常数一共n项）的多项式，用位于A(x)上的n个点来表示。

很容易就会发现，点值表示时，两多项式相乘复杂度为O(n)

FFT算法介绍

首先还需要引入一个单位复根的概念:

单位复根就是ωn=1ωn=1的n个解

设ωkn=cosθ+sinθi,θ=2kπn（k=0,1,2...n−1）ωnk=cosθ+sinθi,θ=2kπn（k=0,1,2...n−1）

很容易发现：ωkn=ωk+nnωnk=ωnk+n即有周期性

离散傅里叶变换

下面开始玄学推理：

首先1n∑n−1k=0ωvkn=[vmodn=0]1n∑k=0n−1ωnvk=[vmodn=0]

将v换为p+q

[(p+q)modn=r]=[(p+q−r)=0]=1n∑n−1k=0ω(p+q−r)kn=1n∑n−1k=0ωpknωqknω−rkn[(p+q)modn=r]=[(p+q−r)=0]=1n∑k=0n−1ωn(p+q−r)k=1n∑k=0n−1ωnpkωnqkωn−rk

设crcr为多项式相乘后第i项的系数

有：cr=∑p,q[(p+q)mod n=r]apbq=1n∑n−1k=0ωpknωqknω−rknapbq=1n∑n−1k=0ω−rknωpknapωqknbqcr=∑p,q[(p+q)mod n=r]apbq=1n∑k=0n−1ωnpkωnqkωn−rkapbq=1n∑k=0n−1ωn−rkωnpkapωnqkbq

现在，令dm=∑n−1k=0ωmknakdm=∑k=0n−1ωnmkak

gm=∑n−1k=0ωmknbkgm=∑k=0n−1ωnmkbk

hm=dm×gmhm=dm×gm

则cm=1n∑n−1k=0ω−mkhkcm=1n∑k=0n−1ω−mkhk

将a,b转化为d,e以及将h转化为c,都是离散傅里叶变换，其复杂度均为O(n2)O(n2)。

快速傅里叶变换

现在我们要降低离散傅里叶变换的时间复杂度。

一个很简单的思路是这样的：

首先只针对于n为2的整次幂的情况

将A(x)按照奇偶次项分离，记为A0(x),A1(x)A0(x),A1(x)

对于任意m<n2m<n2

A(ωmn)=A0(ω2mn)+ωmnA1(ω2mn)=A0(ωmn/2)+ωmnA1(ω2mn)A(ωnm)=A0(ωn2m)+ωnmA1(ωn2m)=A0(ωn/2m)+ωnmA1(ωn2m)

A(ωm+n/2n)=A0(ω2mn)+ωm+n/2nA1(ω2mn)=A0(ωmn/2)−ωmnA1(ω2mn)A(ωnm+n/2)=A0(ωn2m)+ωnm+n/2A1(ωn2m)=A0(ωn/2m)−ωnmA1(ωn2m)

这样一来，我们就可以通过二分来计算离散傅里叶变换。这就是快速傅里叶变换

代码实现

首先是比较浅显易懂，与讲解结合较紧密的递归版本



缺点是常数过大。

因此，我们需要迭代版本：

观察递归时函数的参数：

(a0,a1,a2,a3,a4,a5,a6,a7)(a0,a1,a2,a3,a4,a5,a6,a7)

(a0,a2,a4,a6)(a1,a3,a5,a7)(a0,a2,a4,a6)(a1,a3,a5,a7)

(a0,a4)(a2,a6)(a1,a5)(a3,a7)(a0,a4)(a2,a6)(a1,a5)(a3,a7)

(a0)(a4)(a2)(a6)(a1)(a5)(a3)(a7)(a0)(a4)(a2)(a6)(a1)(a5)(a3)(a7)

在二进制下，可以表示为：

000，100，010，110，001，101，011，111

如果你从前向后看，会发现其实是一个倒序进位的二进制数！

这样，就能知道每一层的参数中系数的顺序了。

迭代版的模板：

void fft(cpx *a,int f){
int i,j,k;
for(i=1,j=0;i<N-1;i++){//N为大于两多项式长度之和的最小的2的整次幂
for(int d=N;j^=d>>=1,~j&d;);
if(i<j)
swap(a[i],a[j]);
}
for(i=1;i<N;i<<=1){
cpx wn(cos(Pi/i),f*sin(Pi/i));
for(j=0;j<N;j+=i<<1){
cpx w(1,0);
for(k=0;k<i;k++,w=w*wn){
cpx x=a[j+k],y=w*a[j+k+i];
a[j+k]=x+y;
a[j+k+i]=x-y;
}
}
}
if(f==-1)
for(i=0;i<N;i++)
a[i].r/=N;
}

看到这里，建议再往回看开头，也许能加深印象

提供一道模板题：UOJ34多项式乘法