快速傅里叶变换(FFT)相关内容汇总
作者：HocRiser。最终核对并完稿于2018年1月10日。支持转载。
目录：
概述
1. 前置技能：数学基础
1.1 多项式概念与运算。
1.2 单元微积分初步与泰勒展开
1.3 普通型生成函数与指数型生成函数
1.4 线性代数相关（矩阵，行列式与特征多项式）
1.5 组合数与伯努利数
1.6 常系数齐次线性递推
1.7 初等数论与初等群论
1.8 卷积概念与O(n^2)求法
1.9 拉格朗日插值法
1.10 复数及其运算
2. FFT：快速傅里叶变换算法总述
2.1 概述
2.2 算法简介
2.2.1 多项式求值与插值
2.2.2 单位复数根及其性质
2.2.3 DFT&IDFT：离散与逆离散傅里叶变换
2.2.4 蝴蝶操作与二进制优化
3. NTT&FWT：快速数论变换与快速沃尔什变换
3.1 原根与模意义下复数运算
3.2 快速数论变换模板
3.3 快速沃尔什变换简介
4. 基础应用
4.1 大整数乘法与多项式乘法
4.2 循环卷积与简单例题
4.3 伯努利数与自然数幂和
5. 形式幂级数与多项式运算
5.1 多项式逆元
5.2 多项式开方
5.3 牛顿迭代法与拉格朗日插值
5.4 多项式ln
5.5 多项式exp
5.6 多项式求幂
5.7 多项式除法
6. 生成函数概述
6.1 普通型生成函数
6.2 指数型生成函数
6.3 生成函数的应用与FFT优化
7. 线性递推，分治FFT与任意模数FFT
7.1 常系数齐次线性递推
7.2 FFT与CDQ分治
7.3 任意模数FFT
8. 拓展与总结
8.1 WC与CTSC以上级别算法与知识内容
8.2 补充说明
8.3 结语
概述
快速傅里叶变换（Fast Fourier Transform），是信息通讯领域的一个重要内容，是傅里叶级数的应用之一。近年来，OI系列赛事中，FFT出现的频率越来越高，应用范围也越来越广。本文致力于介绍FFT在OI中的运用，但不会花大量篇幅介绍算法，只提供笔者认为较好的学习资料与习题，希望能给正在学习相关算法的OIers带来帮助。
傅里叶级数：与本文无关，有兴趣可以学习《高等数学 同济》工科版下册。
傅里叶分析：与本文无关，提供一个介绍网址https://zhuanlan.zhihu.com/p/19763358
傅里叶变换：求两个多项式的卷积的有力工具。
第一部分 前置技能：数学基础
1.1 多项式概念与运算
多项式的概念应该比较熟悉。
但在学习FFT之前，请学习：
(1)多项式运算（加减乘除与开方）
(2)多项式的系数表达与点值表达
(3)多项式求值与插值
(4)高精度运算与程序实现
内容较简单，不再赘述
1.2 微积分初步与泰勒展开
微积分是现代数学的基础，应用很广。
主要运用了“极限”，“化曲为直”，“无穷”的思想，以直代曲。
微积分在物理方程，几何切线，曲率与曲边形面积等方面均有重要运用。
国内FFT题目一般不要求掌握多元微积分，常微分，偏微分和微分方程解法，需掌握的有（后四项只需大致了解）：
(1) 导数与极限概念
(2) 常用函数求导与微分
(3) 微分中值定理：罗尔定理，柯西中值定理与拉格朗日中值定理。
(4) 洛必达公式
(5) 泰勒展开与麦克劳林展开式
(6) 积分与广义积分
(7) 简单曲线积分与曲边形面积求法
(8) 辛普森积分
(9) 微积分的简单物理与几何运用
(10) 曲率，曲率圆与曲率半径
需要注意的是，多项式的求导运算其实也是抽象线性空间的线性变换。（符合线性变换的两条公理）
1.3 普通型生成函数与指数型生成函数
生成函数，又称母函数，能用化离散为连续的方法化简很多求值问题，这在后面会有详细介绍。需掌握的有：两种生成函数的概念与简单运用。
1.4 线性代数相关
线性代数最基础的三个知识点：矩阵，行列式与线性方程。
要掌握：
(1) 矩阵运算（矩阵加减，矩阵乘法，矩阵求幂）
(2) 行列式计算（余子式，代数余子式）
(3) 范得蒙行列式与克莱姆法则
(4) 向量空间与几何意义
(5) 矩阵初等变换与线性方程组
(6) 高斯消元法及其程序实现
(7) 线性递推与矩阵的关系
1.5 组合数与伯努利数
阶乘求法，模意义下阶乘逆元，杨辉三角（帕斯卡三角）
排列组合求法，伯努利数
内容较简单，不再赘述
1.6 常系数齐次线性递推
主要是特征系列问题。
矩阵特征值，特征向量，特征多项式，特征方程，特征根与相关问题。
可以查阅相关资料，提供一个网址https://www.zhihu.com/question/51662733
其余不再赘述
1.7 初等数论与初等群论
群论方面只需掌握概念（封闭性，结合性，幺元与逆元）即可，这方面的应用主要在Burnside与Polya法则上，与FFT关系不大
初等数论需要掌握的有：
(1) 最大公约数与最小公倍数
(2) 中国剩余定理与扩展中国剩余定理（不互质情况）
(3) 欧拉函数及其性质与欧拉函数公式
(4) 欧几里得算法与扩展欧几里得算法
(5) Lucas定理与扩展Lucas定理
(6) BSGS算法与扩展BSGS算法
(7) 费马小定理，威尔逊定理与欧拉定理
(8) 快速幂算法（二进制初步）
(9) 模意义与逆元求法
(10) 更相减损术
1.8 卷积概念与O(n^2)求法
信号处理方面较常用，与FFT算法有密切联系，可自行查找资料，略
1.9 拉格朗日插值法 略
1.10 复数及其运算
复数是数系的扩充，需要了解的内容有：
(1) 复数的概念
(2) 复数的运算：加减乘除
(3) 复数运算的几何意义：向量加减与“模长相乘，幅角相加”
(4) 高次方程的复数根
(5) 共轭复数及其性质
(6) 单位复数根及其性质
其余内容不再赘述。
第二部分 FFT：快速傅里叶变换算法总述
2.1 概述
快速傅里叶变换，是求两个多项式卷积的算法，其时间复杂度为O(nlogn)，优于普通卷积求法，且根据有关证明，快速傅里叶变换是基于变换求卷积的理论最快算法。
关于FFT的介绍，最详细易懂的是《算法导论》上的内容。
其大致介绍与代码在这里：https://www.geek-share.com/detail/2713625660.html.
2.2 算法简介
2.2.1 多项式求值与插值
现在我们已经有了多项式系数表达与点值表达的基础，不难发现，乘法运算在点值表达下只需要O(n)的复杂度即可完成，这样我们会自然的想到，可以先将多项式系数表达转变为点值表达，再在点值表达下做运算，最后转换回系数表达。
按照这个流程，我们得到了快速傅里叶变换的基本思想。其中第一步与第三步的复杂度为O(nlogn)，第二步复杂度为O(n)。在插值的过程中，虽然可以用著名的拉格朗日插值法，但其复杂度仍未O(n^2)，所以我们必须使用更有的方法。
其中，求值的过程叫作离散傅里叶变换(Discrete Fourier Transform)，插值的过程叫作逆离散副叶变换(IDFT)。
2.2.2 单位复数根及其性质
我们知道复数将数系的范围从数轴扩展到整个二维平面（复平面），而复数单位根实际上本质是单位圆的半径，它有三个重要定理：消去引理，折半引理与求和引理。它们的几何意义可以通过圆的无穷对称性理解。
单位复数根的求法：根据cos，sin函数和exp函数的泰勒展开式（麦克劳林）和复数运算的几何意义可以轻松得出。
2.2.3 DFT&IDFT：离散与逆离散傅里叶变换
首先要找到一个保证不会被超过的次数界，为了保证圆的对称性，实际次数界n通常为大于最大运算所需次数界的最小的2的次幂。
DFT：通过递归将系数表达化点值表达的过程利用分治得到，分支过程中次数界每次减半。其中递归时主要用到了消去引理，回溯时主要用到了折半引理，而求和引理将在优化中用到。显然易知DFT分治的时间复杂度为O(nlogn).
IDFT：DFT的逆运算，插值过程。根据多项式与单位复数根的性质，IDFT与DFT的唯一不同在于复数根枚举是反向的以及最终向量的每一个单位都要除以n。
这样，我们先将运算的两个参数（向量a，b）通过DFT化为点值表达，再通过O(n)的乘法运算化为a*b的点值表达，最后用IDFT将a*b化回系数表达。
2.2.4 蝴蝶操作与二进制优化
前面的算法是递归形式的，常数较大。
由于前面的过程是蝴蝶操作，可以根据二进制的性质优化常数并改为非递归形式，具体方法见上面的网址链接。
至此，我们运用一系列的数学知识，将两个多项式的卷积运算成功优化到了O(nlogn)的复杂度。
最后给出完整FFT代码：UOJ#34 多项式乘法
#include<cstdio>
#include<complex>
#include<algorithm>
#define rep(i,l,r) for (int i=l; i<=r; i++)
using namespace std;
typedef complex<double> C;
const int N=262145;
const double pi=acos(-1.);
int n,m,L,x,rev
;
C a
,b
;
void DFT(C a[],int f){
rep(i,0,n-1) if (i<rev[i]) swap(a[i],a[rev[i]]);
for (int i=1; i<n; i<<=1){
C wn(cos(pi/i),f*sin(pi/i));
for (int p=i<<1,j=0; j<n; j+=p){
C w(1,0);
for (int k=0; k<i; k++,w*=wn){
C x=a[j+k],y=w*a[i+j+k]; a[j+k]=x+y; a[i+j+k]=x-y;
}
}
}
}
int main(){
scanf("%d%d",&n,&m);
rep(i,0,n) scanf("%d",&x),a[i]=x;
rep(i,0,m) scanf("%d",&x),b[i]=x;
m=n+m; for (n=1; n<=m; n<<=1) L++;
rep(i,0,n-1) rev[i]=(rev[i>>1]>>1)|((i&1)<<(L-1));
DFT(a,1); DFT(b,1);
rep(i,0,n-1) a[i]=a[i]*b[i];
DFT(a,-1);
rep(i,0,m) printf("%d ",(int)(a[i].real()/n+0.5));
return 0;
}
第三部分 NTT&FWT：快速数论变换与快速沃尔什变换
3.1 原根与模意义下复数运算
原根定义与性质及其求法https://www.geek-share.com/detail/2578640841.html
在NTT中有重要地位，但需要掌握的内容并不多，不再赘述。
3.2 快速数论变换模板
与FFT几乎没有区别，只是将复数运算换成了整数模意义下运算。
#include<cstdio>
#include<complex>
#include<algorithm>
#define rep(i,l,r) for (int i=l; i<=r; i++)
using namespace std;
const int N=262145,md=(119<<23)+1,G=3;//998244353
int n,m,L,x,rev
,a
,b
;
int ksm(int a,int b){
int res;
for (res=1; b; a=(1ll*a*a)%md,b>>=1)
if (b & 1) res=(1ll*res*a)%md;
return res;
}
void NTT(int a[],int f){
rep(i,0,n-1) if (i<rev[i]) swap(a[i],a[rev[i]]);
for (int i=1; i<n; i<<=1){
int wn=ksm(G,(f==1) ? (md-1)/(i<<1) : (md-1)-(md-1)/(i<<1));
for (int p=i<<1,j=0; j<n; j+=p){
int w=1;
for (int k=0; k<i; k++,w=1ll*w*wn%md){
int x=a[j+k],y=1ll*w*a[i+j+k]%md;
a[j+k]=(x+y)%md; a[i+j+k]=(x-y+md)%md;
}
}
}
if (f==1) return;
int inv=ksm(n,md-2);
rep(i,0,n-1) a[i]=1ll*a[i]*inv%md;
}
int main(){
scanf("%d%d",&n,&m);
rep(i,0,n) scanf("%d",&x),a[i]=x;
rep(i,0,m) scanf("%d",&x),b[i]=x;
m=n+m; for (n=1; n<=m; n<<=1) L++;
rep(i,0,n-1) rev[i]=(rev[i>>1]>>1)|((i&1)<<(L-1));
NTT(a,1); NTT(b,1);
rep(i,0,n-1) a[i]=1ll*a[i]*b[i]%md;
NTT(a,-1);
rep(i,0,m) printf("%d ",a[i]);
return 0;
}
3.3 快速沃尔什变换简介
快速沃尔什变换(Fast Walsh-Hadamard Transform)是基于位运算卷积的多项式算法。现在可以运用的运算有：与运算and，或运算or，异或运算xor。
适用范围并不是非常广，UOJ上有部分相关应用且难度较大，有兴趣可自行查阅资料。
例题：BZOJ4589，BZOJ4911
第四部分 基础应用
4.1 大整数乘法与多项式乘法
显然，大整数乘法的本质就是多项式乘法，所以高精度算法显然可以使用FFT加速。例题BZOJ2179，裸题。
对于其扩展，有BZOJ2194 快速傅里叶之二，可以考察对快速傅里叶算法的理解是否深刻。
4.2 循环卷积与简单例题
循环卷积就是将多项式（向量）看作一个可循环的环，在这上面的运算，只需要将其中一个多项式复制一份接在自己后面即可。
例题：BZOJ4927 HNOI2017 礼物
4.3 伯努利数与自然数幂和
伯努利数分为两种：B-（第一类伯努利数）与B+（第二类伯努利数），区别就在于前者B1=-1/2，后者B1=1/2，其余项完全一样。
但虽然只有一项不同，二者的递推式与生成函数定义是不同的，应用在自然数幂和等方面的公式推导也有不同，一般来说前者是国际公认的标准。
伯努利数在求自然数幂和方面的应用主要是加速幂和过程，能将通过预处理将O(k^2*T)的复杂度优化为O(k^2+k*T)。而通过生成函数与FFT可以在对数时间复杂度内求值。
一篇很好的介绍文章https://www.geek-share.com/detail/2620536568.html
关于其对数时间复杂度的求法，在文章的结尾有介绍，原网页：http://picks.logdown.com/posts/189620-the-inverse-element-of-polynomial
第五部分 形式幂级数与多项式运算
5.1 多项式逆元
多项式逆元，是指在多项式运算中，在模x的某个次幂的意义下，乘上多项式逆元的效果与除以此多项式是一样的
https://www.geek-share.com/detail/2703764480.html
https://www.geek-share.com/detail/2653934120.html
http://blog.miskcoo.com/2015/05/polynomial-inverse
很好理解，代码也很好实现
可见2015国集论文《生成函数的运算和组合计数问题》
void inverse(int a[],int b[],int l){
int n=l<<1;
if (l==1) { b[0]=ksm(a[0],p-2); return; }
inverse(a,b,l>>1);
rep(i,0,l) t[i]=a[i],t[i+l]=0;
DFT(t,n,1); DFT(b,n,1);
rep(i,0,n) b[i]=1ll*b[i]*(2-1ll*t[i]*b[i]%p+p)%p;
DFT(b,n,-1);
rep(i,l,n) b[i]=0;
}
5.2 多项式开方
顾名思义，对多项式在模x的某个次幂的意义下作开方操作
http://picks.logdown.com/posts/189620-inverse-element-of-polynomial
http://picks.logdown.com/posts/202388-square-root-of-polynomial
代码类似
void sqrt(int a[],int b[],int l){
int n=l<<1;
if (l==1) { b[0]=1; return; }
sqrt(a,b,l>>1);
rep(i,0,l) t[i]=a[i],t[i+l]=0,ib[i]=ib[i+l]=0;
inverse(b,ib,l); DFT(t,n,1); DFT(b,n,1); DFT(ib,n,1);
rep(i,0,n) b[i]=1ll*inv[2]*(b[i]+1ll*t[i]*ib[i]%p)%p;
DFT(b,n,-1);
rep(i,l,n) b[i]=0;
}
5.3 牛顿迭代法与拉格朗日插值
需要认识到的是，泰勒展开本身就是牛顿迭代法的拓展。
http://blog.miskcoo.com/2015/06/polynomial-with-newton-method
http://picks.logdown.com/posts/197262-polynomial-division
5.4 多项式ln
前面的网址已经介绍了ln的求法（导数与积分的妙用）
void ln(int a[],int b[],int l){
int n=l<<1;
rep(i,0,n) da[i]=ia[i]=0;
rep(i,0,l) da[i]=1ll*(i+1)*a[i+1]%p;
inverse(a,ia,l); DFT(da,n,1); DFT(ia,n,1);
rep(i,0,n) b[i]=1ll*da[i]*ia[i]%p;
DFT(b,n,-1); b[0]=0;
for (int i=l; i; i--) b[i]=1ll*inv[i]*b[i-1]%p;
rep(i,l,n) b[i]=0;
}
5.5 多项式exp
同上，不再赘述
void exp(int a[],int b[],int l){
int n=l<<1;
if (l==1) { b[0]=1; return; }
exp(a,b,l>>1);
rep(i,0,n) t[i]=0; ln(b,t,l);
rep(i,0,l) t[i]=(a[i]-t[i]+p)%p;
t[0]=(t[0]+1)%p; DFT(b,n,1) DFT(t,n,1);
rep(i,0,n) b[i]=ll*b[i]*t[i]%p;
DFT(b,n,-1) rep(i,l,n) b[i]=0;
}
5.6 多项式求幂
现在我们已经有了ln与exp的求法，求出幂就很简单了。
void power(int a[],int k,int l){
int n=l<<1;
rep(i,0,n) t[i]=0;
ln(a,t,l);
rep(i,0,l) t[i]=1ll*k*t[i]%p;
rep(i,0,n) a[i]=0;
exp(t,a,l);
}
5.7 多项式除法
并不是很常用但了解一下并不费时间。
http://picks.logdown.com/posts/197262-polynomial-division
http://blog.miskcoo.com/2015/05/polynomial-division
#include<cstdio>
#include<algorithm>
#define rep(i,l,r) for (int i=l; i<=r; i++)
using namespace std;
const int N=300100,mod=998244353,g=3;
int n,m,p,p2,a
,b
,c
,d
,e
,lg[N<<1],fac
,fin
,tmp
,tmp2
,rev
;
int ksm(int a,int b){
int res;
for (res=1; b; a=1ll*a*a%mod,b>>=1)
if (b&1) res=1ll*res*a%mod;
return res;
}
void DFT(int a[],int n,int f){
for (int i=0; i<n; i++) if (i<rev[i]) swap(a[i],a[rev[i]]);
for (int i=1; i<n; i<<=1){
int wn=ksm(g,(f==1) ? (mod-1)/(i<<1) : (mod-1)-(mod-1)/(i<<1));
for (int p=i<<1,j=0; j<n; j+=p){
int w=1;
for (int k=0; k<i; k++,w=1ll*w*wn%mod){
int x=a[j+k],y=1ll*w*a[i+j+k]%mod; a[j+k]=(x+y)%mod; a[i+j+k]=(x-y+mod)%mod;
}
}
}
if (f==1) return;
int inv=ksm(n,mod-2);
for (int i=0; i<n; i++) a[i]=1ll*a[i]*inv%mod;
}
void get(int a[],int b[],int l){
if (l==1){ b[0]=ksm(a[0],mod-2); return ;}
get(a,b,l>>1); int n=l<<1;
for (int i=0; i<l; i++) tmp[i]=a[i],tmp[i+l]=0;
for (int i=0; i<n; i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(lg
-1));
DFT(tmp,n,1); DFT(b,n,1);
for (int i=0; i<n; i++) tmp[i]=1ll*b[i]*(2-1ll*tmp[i]*b[i]%mod+mod)%mod;
DFT(tmp,n,-1);
for (int i=0; i<l; i++) b[i]=tmp[i],b[i+l]=0;
}
void mul(int a[],int b[],int c[],int n){
for (int i=0; i<n; i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(lg
-1));
for (int i=0; i<n; i++) tmp[i]=a[i],tmp2[i]=b[i];
DFT(tmp,n,1); DFT(tmp2,n,1);
for (int i=0; i<n; i++) tmp[i]=1ll*tmp[i]*tmp2[i]%mod;
DFT(tmp,n,-1);
for (int i=0; i<n; i++) c[i]=tmp[i];
}
int main(){
freopen("div.in","r",stdin);
freopen("div.out","w",stdout);
scanf("%d%d",&n,&m);
lg[1]=0; rep(i,2,n<<3) lg[i]=lg[i>>1]+1;
rep(i,0,n) scanf("%d",&a[i]); reverse(a,a+n+1);
rep(i,0,m) scanf("%d",&b[i]); reverse(b,b+m+1);
p=1; while (p<=n) p<<=1;
p2=p<<1; get(b,c,p);
mul(c,a,c,p2); reverse(c,c+(n-m)+1);
rep(i,n-m+1,p2) c[i]=0;
reverse(a,a+n+1); reverse(b,b+m+1);
mul(c,b,d,p2); rep(i,0,m) e[i]=(a[i]-d[i]+mod)%mod;
rep(i,0,n-m) printf("%d ",c[i]); puts("");
rep(i,0,m-1) printf("%d ",e[i]); puts("");
return 0;
}
多项式除法
另外，多点求值与快速插值，可以参考
http://blog.miskcoo.com/2015/04/polynomial-multiplication-and-fast-fourier-transform
例题：BZOJ3625，BZOJ4555
第六部分 生成函数概述
6.1 普通型生成函数(Ordinary Generating Function, OGF)
很常用的生成函数，大量例题在这里https://www.cnblogs.com/candy99/category/946506.html，这里不再赘述。
6.2 指数型生成函数(Exponential Generating Function, EGF)
同理，但难度大一些，FFT中用的比重少一些。
6.3 生成函数的应用与FFT优化
生成函数是FFT的一大用途，在子集计数，容斥原理等方面，大量的生成函数运算可以用FFT实现。生成函数将复杂度从指数复杂度变为多项式复杂度，而FFT将其优化为对数复杂度。
普通型生成函数是FFT的主要用途之一，按照笔者的个人理解，它是一种以连续（或者说是形式幂级数）的函数方法解决离散计数问题的思想。事实上，快速傅里叶变换也是运用了这种思想，本来是两个数列（或者说两个离散函数）的卷积，通过一个没有具体意义的"x"，构造出连续的函数，从而使能应用于连续函数的那一套理论成功应用于此类问题。
容斥是生成函数的一个部分，有关容斥的计数问题应该很快想到FFT，这类问题要特别注意是否容斥完全了（有没有重和漏）。另外有个小技巧，在多项式连续进行几个操作时，可以先求值，然后在点值表达下进行多个运算（修改y坐标即可），最后再一起插值，大大减小了常数。
例题：BZOJ3509，BZOJ4503
第七部分 线性递推，分治FFT与任意模数FFT
7.1 常系数齐次线性递推
https://www.geek-share.com/detail/2704711402.html
系数为定值的一次递推式构成的数列就是常系数齐次线性递推数列。
暴力求值的复杂度一般是O(nk)的，当n非常大而k较小时可以使用矩阵快速幂优化（比较基础的应用，此处不再展开），复杂度为O(k^3*log n)。但对于某些题目，k的三次方的复杂度还有更多的优化空间，能被优化为O(k^2*log n)。这里涉及到一定的矩阵对角化的知识。利用相似矩阵的性质，我们先根据Hamilton-Carley定理列出矩阵的特征方程，然后取k+1个值代入方程，解出特征根，再通过特征根找到特征向量并拼成一个矩阵，最后就可以直接将矩阵快速幂化成k个特征值的快速幂了。
当然对于一般的矩阵，我们可以通过拉格朗日插值法求出特征根，但这样也是有局限性的，因为拉格朗日插值法本身的复杂度就是O(k^4)。
7.2 FFT与CDQ分治
https://www.geek-share.com/detail/2701881800.html
CDQ分治是NOI选手陈丹琦在国家队论文中提出的分治算法，在处理偏序问题等方面有着重要运用，而FFT也可以使用这种方法。
很多情况下我们会得到一个数列的递推式（当然这不是常系数齐次线性的），而这个递推式可以化成两个已知的数列的卷积，这个时候就可以用CDQ分治解决。
设我们先在要处理[l,r]的数列值，我们先递归处理[l,mid]的值，然后计算[l,mid]的值对(mid,r]的值的影响（贡献），这通常也是一个卷积的形式，这个时候只要直接构造出两个卷积数列，然后套用FFT即可。
复杂度为O(nlog^2n)且常数较大，但与暴力相比还是有很大优势的。
例题：BZOJ4555
7.3 任意模数FFT
有时会要求在模意义下做FFT，而模数不是NTT的专用模数，用FFT又会卡精度和爆long long，这是就需要任意模数FFT（matthew99在2016论文中讲过三模数FFT，但笔者还没有研究）
http://www.cnblogs.com/zhoushuyu/p/8763215.html#4015017
任意模数FFT的思想很简单，选M=sqrt(mod)，然后将两个卷积式子分别拆开，然后利用卷积的分配率合并，一共7次DFT。
#include<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define rep(i,l,r) for (int i=l; i<=r; i++)
typedef long long ll;
typedef long double ld;
using namespace std;
const int N=400010;
const ld pi=acos(-1.);
int n,m,mod,a1
,a2
,rev
,res
;
struct Com{
ld x,y;
Com(ld _x=0,ld _y=0):x(_x),y(_y){};
}A
,B
,C
,D
,s1
,s2
,s3
;
Com operator +(Com a,Com b){ return Com(a.x+b.x,a.y+b.y); }
Com operator -(Com a,Com b){ return Com(a.x-b.x,a.y-b.y); }
Com operator *(Com a,Com b){ return Com(a.x*b.x-a.y*b.y,a.y*b.x+a.x*b.y); }
void FFT(Com a[],int n,int f){
for (int i=0; i<n; i++) if (i<rev[i]) swap(a[i],a[rev[i]]);
for (int i=1; i<n; i<<=1){
Com wn(cos(pi/i),f*sin(pi/i));
for (int p=i<<1,j=0; j<n; j+=p){
Com w(1,0);
for (int k=0; k<i; k++,w=w*wn){
Com x=a[j+k],y=w*a[i+j+k]; a[j+k]=x+y; a[i+j+k]=x-y;
}
}
}
if (f==-1) for (int i=0; i<n; i++) a[i].x/=n;
}
void solve(){
int s=n+m; int nn=1,L=0;
for (; nn<s; nn<<=1,L++);
for (int i=0; i<nn; i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(L-1));
for (int i=0; i<nn; i++) A[i]=a1[i]>>15,B[i]=a1[i]&0x7fff;
for (int i=0; i<nn; i++) C[i]=a2[i]>>15,D[i]=a2[i]&0x7fff;
FFT(A,nn,1); FFT(B,nn,1); FFT(C,nn,1); FFT(D,nn,1);
for (int i=0; i<nn; i++) s1[i]=A[i]*C[i],s2[i]=A[i]*D[i]+B[i]*C[i],s3[i]=B[i]*D[i];
FFT(s1,nn,-1); FFT(s2,nn,-1); FFT(s3,nn,-1);
for (int i=0; i<nn; i++)
res[i]=((((ll)round(s1[i].x)%mod)<<30)%mod+(((ll)round(s2[i].x)%mod)<<15)%mod+(ll)round(s3[i].x)%mod)%mod;
}
int main(){
freopen("mtt.in","r",stdin);
freopen("mtt.out","w",stdout);
scanf("%d%d%d",&n,&m,&mod);
for (int i=0; i<=n; i++) scanf("%d",&a1[i]);
for (int i=0; i<=m; i++) scanf("%d",&a2[i]);
solve();
for (int i=0; i<=n+m; i++) printf("%d ",res[i]);
return 0;
}
第八部分 拓展与总结
8.1 WC与CTSC以上级别算法与知识内容
(1) 二项式反演与子集反演
(2) 圆反演与拉格朗日反演
(3) 快速莫比乌斯变换
(4) 子集计数与容斥原理
还有少量新兴算法也与FFT与FWT相关，也有部分题目解题推导中会用到FWT但最终代码里并没有体现，可见快速傅里叶与沃尔什算法有时也会成为解题的纯数学工具。
8.2 补充说明
FFT的模板应该背熟并能在5min内正确默写完。还有几道FFT相关题目的代码可能可以在本博客的“快速傅里叶变换”标签中找到（包括上面例题中的几道），关于FFT的一些其它用法这篇文章不再涉及，有兴趣可以自行查阅资料。
FFT的常数较大，除了二进制优化之外，还有不少效果较好的优化方法，如matthew99在集训队论文中提出的方法，但由于实现较复杂，所以考场上的用处并不大。一个比较常用的技巧是，在作复数运算时自己写一个复数类而不是使用STL库中自带的complex类。除此之外，FFT的精度损失也是比较严重的，所以看到类似998244353,1004535809（这两个素数的原根均为3）之类的费马素数时，应尽量使用NTT。
FFT可能会用到的素数及其原根：http://blog.miskcoo.com/2014/07/fft-prime-table
8.3 结语
快速傅里叶变换算法本身并没有很大的扩展，但由于非常多的数学定理与应用都会用到多项式卷积运算，这就使FFT在现在比赛中的运用占比十分大。对于这方面的理解，要通过做题与自己演算才能加深。最后希望各位正在学习FFT的OIers能熟练运用并在比赛中的相关题目上获得高分。
鸣谢：给我讲解过FFT的同学以及上面引用的所有网页的作者，如有侵权请及时与笔者联系。
由于时间匆忙，可能有个别笔误与错误的地方，希望各位发现后在评论区给予指正，感激不尽。
- Flash Player 10 相关内容汇总
- appledoc相关内容汇总（非原创）
- Android开发环境——Eclipse ADT相关内容汇总
- IOS相关培训机构课程内容汇总
- appledoc相关内容汇总
- CRF相关内容汇总（转自其他博客）
- Android开发环境——SDK相关内容汇总
- Android开发环境——SDK相关内容汇总
- Android开发环境——SDK相关内容汇总
- 2010年会计从业资格考试相关内容汇总（全国各地）
- Android开发环境——Eclipse ADT相关内容汇总
- Android开发环境——连接驱动ADB相关内容汇总
- Android百度地图相关内容汇总
- Android开发环境——连接驱动ADB相关内容汇总
- Android百度地图相关内容汇总
- Android开发环境——调试器 DDMS相关内容汇总
- Android百度地图相关内容汇总
- Android开发环境——模拟器AVD相关内容汇总
- Android百度地图相关内容汇总
- Android百度地图相关内容汇总