本文主要是介绍震惊！FFT竟然还能这么写？活到这么大没见过这么写FFT的！，希望对大家解决编程问题提供一定的参考价值，需要的开发者们随着小编来一起学习吧！
模板题
洛谷P3803 【模板】多项式乘法（FFT）
点值表达
大莉模拟的复杂度是 O ( n 2 ) O(n^2) O(n2)。
引入点值表达qwq（这部分摘抄自毒瘤czh的FFT讲稿）：
1. 1. 1.例子： A ( x ) = x 2 + 2 x − 1 A(x)=x^2+2x-1 A(x)=x2+2x−1可以被表达为 { ( 0 , − 1 ) , ( 1 , 2 ) , ( 2 , 7 ) } \{(0 , -1), (1 , 2),(2 , 7)\} {(0,−1),(1,2),(2,7)}
2. 2. 2.用处：在多项式乘法中，系数表达暴力模拟复杂度 O ( n 2 ) O(n^2) O(n2)，但是点值表达式可以 O ( n ) O(n) O(n)求出
更加具体地说，我们有两个点值表达
A ( x ) = { ( x 0 , A ( x 0 ) ) , ( x 1 , A ( x 1 ) ) , ( x 2 , A ( x 2 ) ) , . . . . , ( x n , A ( x n ) ) } A(x)=\{(x_0,A(x_0)),(x_1,A(x_1)),(x_2,A(x_2)),....,(x_n,A(x_n))\} A(x)={(x0,A(x0)),(x1,A(x1)),(x2,A(x2)),....,(xn,A(xn))}
B ( x ) = { ( x 0 , B ( x 0 ) ) , ( x 1 , B ( x 1 ) ) , ( x 2 , B ( x 2 ) ) , . . . . , ( x n , B ( x n ) ) } B(x)=\{(x_0,B(x_0)),(x_1,B(x_1)),(x_2,B(x_2)),....,(x_n,B(x_n))\} B(x)={(x0,B(x0)),(x1,B(x1)),(x2,B(x2)),....,(xn,B(xn))}
乘积 C ( x ) = A ( x ) × B ( x ) = { ( x 0 , A ( x 0 ) × B ( x 0 ) ) , ( x 1 , A ( x 1 ) × B ( x 1 ) ) , ( x 2 , A ( x 2 ) × B ( x 2 ) ) , . . . . , ( x n , A ( x n ) × B ( x n ) ) } C(x)=A(x)\times B(x)=\{(x_0,A(x_0)\times B(x_0)),(x_1,A(x_1)\times B(x_1)),(x_2,A(x_2)\times B(x_2)),....,(x_n,A(x_n)\times B(x_n))\} C(x)=A(x)×B(x)={(x0,A(x0)×B(x0)),(x1,A(x1)×B(x1)),(x2,A(x2)×B(x2)),....,(xn,A(xn)×B(xn))}
所以我们可以通过 O ( n ) O(n) O(n)复杂度，求出两个点值表达式乘积的结果（两个多项式的乘积又称卷积）
知道了点值表达式就珂以确定一个多项式了qwq（珂以列方程用高斯消元爆解一波qwq）
其中DFT算法是暴力 O ( n 2 ) O(n^2) O(n2)把系数表达转点值表达，FFT算法是 O ( n l o g n ) O(nlogn) O(nlogn)把系数表达转点值表达
DFT（离散傅里叶变换）
普通做法：随机取 n n n个值代入求值qwq。
DFT是考虑建立一个复平面，在上面怼出 画一个单位圆：
因为是单位圆，所以圆上的每一个点对应的复数模均为1。
把这个圆平均分成 n n n份，把每个点对应的复数代入。
具体地说，这个函数的点值表达是 { ( w n 0 , f ( w n 0 ) ) , ( w n 1 , f ( w n 1 ) ) , . . . , ( w n n − 1 , f ( w n n − 1 ) ) } \{(w_n^0,f(w_n^0)),(w_n^1,f(w_n^1)),...,(w_n^{n-1},f(w_n^{n-1}))\} {(wn0,f(wn0)),(wn1,f(wn1)),...,(wnn−1,f(wnn−1))}，
其中 w n k = c o s ( k n 2 π ) + i s i n ( k n 2 π ) w_n^k=cos(\frac{k}{n}2\pi)+isin(\frac{k}{n}2\pi) wnk=cos(nk2π)+isin(nk2π)
FFT（快速傅里叶变换）
w w w的性质
证明珂以看巨神mhy的博客qwq
这里只列出性质，就不证明了：
1. w n k = w 2 n 2 k 1.w_n^k=w_{2n}^{2k} 1.wnk=w2n2k
2. w n k + n 2 = − w n k 2.w_n^{k+\frac{n}{2}}=-w_n^k 2.wnk+2n=−wnk
3. w n 0 = w n n = 1 3.w_n^0=w_n^n=1 3.wn0=wnn=1
4. w n p ∗ w n q = w n p + q 4.w_n^p*w_n^q=w_n^{p+q} 4.wnp∗wnq=wnp+q
表达式的变形
首先把一个函数用系数表达式写出来（假设 n n n为偶数，如果不是就补一个 0 ∗ x n 0*x^n 0∗xn）：
f ( x ) = a 0 + a 1 x 1 + a 2 x 2 + a 3 x 3 + a 4 x 4 + . . . + a n − 1 x n − 1 f(x)=a_0+a_1x^1+a_2x^2+a_3x^3+a_4x^4+...+a_{n-1}x^{n-1} f(x)=a0+a1x1+a2x2+a3x3+a4x4+...+an−1xn−1
按照下标的奇偶性把 f f f分成两半：
f ( x ) = ( a 0 + a 2 x 2 + . . . + a n − 2 x n − 2 ) + x ( a 1 + a 3 x 2 + . . . + a n − 1 x n − 2 ) f(x)=(a_0+a_2x^2+...+a_{n-2}x^{n-2})+x(a_1+a_3x^2+...+a_{n-1}x^{n-2}) f(x)=(a0+a2x2+...+an−2xn−2)+x(a1+a3x2+...+an−1xn−2)
令 A 1 ( x ) = a 0 + a 2 x + a 4 x 2 + . . . + a n − 2 x n − 2 2 A_1(x)=a_0+a_2x+a_4x^2+...+a_{n-2}x^{\frac{n-2}{2}} A1(x)=a0+a2x+a4x2+...+an−2x2n−2， A 2 ( x ) = a 1 + a 3 x + a 5 x 2 + . . . + a n − 1 x n − 2 2 A_2(x)=a_1+a_3x+a_5x^2+...+a_{n-1}x^{\frac{n-2}{2}} A2(x)=a1+a3x+a5x2+...+an−1x2n−2
不难发现 f ( x ) = A 1 ( x 2 ) + x A 2 ( x 2 ) f(x)=A_1(x^2)+xA_2(x^2) f(x)=A1(x2)+xA2(x2)。
以下分别把 w n k w_n^k wnk和 w n k + n 2 w_n^{k+\frac{n}{2}} wnk+2n代入：
把 w n k w_n^k wnk代入，得到 f ( w n k ) = A 1 ( w n k ∗ w n k ) + x A 2 ( w n k ∗ w n k ) = A 1 ( w n 2 k ) + w n k A 2 ( w n 2 k ) f(w_n^k)=A_1(w_n^k*w_n^k)+xA_2(w_n^k*w_n^k)=A_1(w_n^{2k})+w_n^kA_2(w_n^{2k}) f(wnk)=A1(wnk∗wnk)+xA2(wnk∗wnk)=A1(wn2k)+wnkA2(wn2k)
由 w w w的性质1，珂以得到：
f ( w n k ) = A 1 ( w n 2 k ) + w n k A 2 ( w n 2 k ) f(w_n^k)=A_1(w_{\frac{n}{2}}^k)+w_n^kA_2(w_{\frac{n}{2}}^k) f(wnk)=A1(w2nk)+wnkA2(w2nk)
把 w n k + n 2 w_n^{k+\frac{n}{2}} wnk+2n代入，得到：
f ( w n k + n 2 ) = A 1 ( w n 2 k + n ) + w n k + n 2 A 2 ( w n 2 k + n ) = A 1 ( w n 2 k ∗ w n n ) − w n k A 2 ( w n 2 k ∗ w n n ) = A 1 ( w n 2 k ) − w n k A 2 ( w n 2 k ) f(w_n^{k+\frac{n}{2}})=A_1(w_n^{2k+n})+w_n^{k+\frac{n}{2}}A_2(w_n^{2k+n})=A_1(w_n^{2k}*w_n^n)-w_n^kA_2(w_n^{2k}*w_n^n)=A_1(w_n^{2k})-w_n^kA_2(w_n^{2k}) f(wnk+2n)=A1(wn2k+n)+wnk+2nA2(wn2k+n)=A1(wn2k∗wnn)−wnkA2(wn2k∗wnn)=A1(wn2k)−wnkA2(wn2k)
分治写法
认真观察……这 f ( w n k ) f(w_n^k) f(wnk)和 f ( w n k + n 2 ) f(w_n^{k+\frac{n}{2}}) f(wnk+2n)长得很像？
于是珂以递归（分治）解决，每次求解左半部分的答案，再推出右半部分的答案：
void FFT(complex<double> *ans,int n,int type) {//type=1表示FFT；type=-1表示IFFT if(n==1) return;//偶数存L，奇数存R complex<double> L[(n>>1)+1],R[(n>>1)+1];for(re i=0; i<n; i+=2) {L[i>>1]=ans[i];R[i>>1]=ans[i+1];}FFT(L,n>>1,type);FFT(R,n>>1,type);//玄学单位圆乱搞 complex<double> w(1,0);complex<double> k(cos(2.0*pi/n),type*sin(2.0*pi/n));int maxn=n>>1;for(re i=0; i<maxn; i++) {complex<double> tmp=w*R[i];ans[i]=L[i]+tmp;ans[i+maxn]=L[i]-tmp;w*=k;}
}
然后交到例题上，TLE？？
递归常数蜃大qwq，考虑使用迭代。
迭代写法
观察下面这张图：
发现原序列的每个数按对称轴翻转一下就是后序列的数qwq，称之为位逆序置换。
所以，我们可以预处理出系数位置递归到最底层的位置，自底向上进行递推。
处理出位逆序置换之后的方法（ R [ i ] R[i] R[i]表示 i i i进行位逆序置换之后得到的值， L L L表示二进制下的总位数）：
for(re i=0; i<MAXN; i++) { //MAXN表示要处理出来的数的总数//表示先把最后一位去掉，剩下的部分进行反转，然后把最后一位反转到最高位qwqR[i]=(R[i>>1]>>1)|((i&1)<<(L-1));
}
所以把递归强行用循环写出来：
第一维枚举区间长度，即递归版本中的 n n n。第二维枚举区间开始位置。第三维枚举当前数是区间中的第几个。
然后珂以大莉写一波qwq（我的版本中，区间长度实际上为 2 i 2i 2i）：
void FFT(Complex *ans,int n,int type) {//type=1表示FFT，type=-1表示IFFT //把序列交换成后序列的形式 for(re i=0; i<n; i++) if(i<R[i]) swap(ans[i],ans[R[i]]);for(re i=1; i<n; i<<=1) {//(2*pi)/(2*i)=pi/iComplex wn=(Complex){cos(pi/(double)i),type*sin(pi/(double)i)};for(re j=0; j<n; j+=(i<<1)) {Complex w=(Complex){1,0};for(re k=0; k<i; k++) {//类似递归版的处理const Complex x=ans[j+k],y=w*ans[i+j+k];ans[j+k]=x+y;ans[i+j+k]=x-y;w=w*wn;}}}
}
IFFT（快速傅里叶逆变换）
毒瘤czh：
美妙结论：一个多项式在分治的过程中乘上单位根的共轭复数，分治完的每一项除以n即为原多项式的每一项系数（具体看程序实现，证明略）
也就是在调用的时候把type设成-1就珂以了
例题完整代码
至此，就珂以写出例题了（注意预处理和边界条件qwq）：
#include<stdio.h>
#include<cstring>
#include<algorithm>
#include<math.h>
#define re register int
#define rl register ll
#define pi 3.14159265358979323846264338328
using namespace std;
int read() {re x=0,f=1;char ch=getchar();while(ch<'0' || ch>'9') {if(ch=='-') f=-1;ch=getchar();}while(ch>='0' && ch<='9') {x=10*x+ch-'0';ch=getchar();}return x*f;
}
inline void write(const int x) {if(x>9) write(x/10);putchar(x%10+'0');
}
const int Size=4000005;
//手写complex类，不过貌似也没快多少
struct Complex {double x,y;inline double real() {return x;}
};
inline Complex operator + (const Complex a,const Complex b) {return (Complex){a.x+b.x,a.y+b.y};
}
inline Complex operator - (const Complex a,const Complex b) {return (Complex){a.x-b.x,a.y-b.y};
}
inline Complex operator * (const Complex a,const Complex b) {return (Complex){a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x};
}
Complex A[Size],B[Size];
int R[Size];
void FFT(Complex *ans,const int n,const int type) {//type=1表示FFT，type=-1表示IFFTfor(re i=0; i<n; i++) if(i<R[i]) swap(ans[i],ans[R[i]]);for(re i=1; i<n; i<<=1) {//(2*pi)/(2*i)=pi/iconst Complex wn=(Complex){cos(pi/i),type*sin(pi/i)};for(re j=0; j<n; j+=i<<1) {Complex w=(Complex){1,0};for(re k=0; k<i; k++) {Complex x=ans[j+k],y=w*ans[i+j+k];ans[j+k]=x+y;ans[i+j+k]=x-y;w=w*wn;}}}
}
int ans[Size];
int main() {int n=read();int m=read();for(re i=0; i<=n; i++) {A[i]=(Complex){read(),0};}for(re i=0; i<=m; i++) {B[i]=(Complex){read(),0};}//now为大于n+m的第一个2的幂，因为总区间长度必须是2的幂qwqint now=1,L=0;while(now<=n+m) {now<<=1;L++;}for(re i=0; i<now; i++) {R[i]=(R[i>>1]>>1)|((i&1)<<(L-1));}FFT(A,now,1);FFT(B,now,1);for(re i=0; i<=now; i++) {A[i]=A[i]*B[i];}FFT(A,now,-1);for(re i=0; i<=n+m; i++) {write(int(A[i].real()/now+0.5)); //注意四舍五入putchar(' ');}return 0;
}
其他题目
洛谷P1919 A*B Problem升级版
即bzoj2179 （权限题）快速傅里叶变换
bzoj2194（权限题）快速傅里叶变换之二
bzoj3527 力
bzoj4827 礼物
……
…………
………………
题解？咕咕咕
还是有一篇的qwq
bzoj4332（权限题）
即洛谷P5075
题解传送门
这篇关于震惊！FFT竟然还能这么写？活到这么大没见过这么写FFT的！的文章就介绍到这儿，希望我们推荐的文章对编程师们有所帮助！
