FFT
多项式
定义：多项式是指由变量、系数以及它们之间的加、减、乘、幂运算（非负整数次方）得到的表达式。
- 项：多项式中的每个单项式叫做多项式的项。
- 多项式的次数：这些单项式中的最高项次数，就是这个多项式的次数。
一个 \(n\) 次多项式可以表示为 \(\sum\limits_{i=0}^na_ix^i\),其中 \(a_i\) 为系数。
复数
定义：我们把形如 \(z=a+bi\)（\(a,b\) 均为实数）的数称为 复数，其中 \(a\) 称为实部，\(b\) 称为虚部，\(i\) 称为虚数单位，（\(i^2=-1\)）
设 \(z_1=a+bi\),\(z_2=c+di\)
单位根
定义：记方程 \(\omega^n=1\left(n\in N^{+}\right)\) 的 \(n\) 个复数解为 \(n\) 次单位根。
可以解出 \(\omega=\mathrm{e}^{\frac{2k\pi i}{n}}\)
定义 主次单位根：\(\omega_n=\mathrm{e}^{\frac{2\pi i}{n}}\)
对于 所有单位根，记作：\(\omega_n^k=\mathrm{e}^{\frac{2k\pi i}{n}}\)
多项式乘法卷积
对于 \(n\) 次多项式 \(A,B\)，\(A=\sum\limits_{i=0}^na_ix^i\),\(B=\sum\limits_{i=0}^nb_ix^i\)
若 \(A\) 与 \(B\) 的乘法卷积为 \(C\)，\(C=\sum\limits_{i=0}^{n}\sum\limits_{j=0}^{n}a_ib_jx^{i+j}\)
多项式的点值表示
众所周知，一个 \(n\) 次函数可以由 \(n+1\) 个点对 \(\left(x_i,y_i\right)\) 唯一确定（其中 \(x_i\) 互不相同）
因此可以用 \(n+1\) 个点对表示出一个 \(n\) 次的多项式。
若 \(A=\left(\left(x_1,y_1\right),\left(x_2,y_2\right),\dots,\left(x_n,y_n\right)\right)\)，\(B=\left(\left(x_1,y'_1\right),\left(x_2,y_2'\right),\dots,\left(x_n,y_n'\right)\right)\)
则 \(A\times B=\left(\left(x_1,y_1y_1'\right),\left(x_2,y_2y_2'\right),\dots,\left(x_n,y_ny_n'\right)\right)\)
所以如果多项式的点值已知，那么多项式的乘法是 \(O\left(n\right)\) 的。
若把点对的数量扩大至 \(2n+1\)，则可以唯一确定 \(C\)。
DFT & IDFT
DFT(离散傅里叶变换) 是系数数列 \(a_0,a_1,\cdots,a_n\) 转为点值数列的过程，IDFT(逆离散傅里叶变换) 是 DFT 的逆运算。
注意到单位根存在以下性质
考虑用单位根的性质优化 DFT。
对于多项式 \(A\left(x\right)=\sum_{i=0}^{2n-1}a_ix^i\)，其中 \(2n=2^k\left(k\in N^{+}\right)\)，我们将 \(\omega_{2n}^0,\omega_{2n}^1,\cdots,\omega_{2n}^{{2n}-1}\) 代入公式计算点值。
现在重定义两个多项式
显然
将单位复根代入上式
发现对于 \(k\in\left[0,1,\cdots,n-1\right]\) \(A\left(\omega_{2n}^k\right)\) 和 \(A\left(\omega_{2n}^{n+k}\right)\) 是可以在一起计算的。
于是有以下伪代码
function FFT(complex A[], int Length)
if Length == 1 return
complex A0[Length / 2], A1[Length / 2]
for int i = 0 to Length / 2 - 1 with i += 1
A0[i] = A[i * 2]
A1[i] = A[i * 2 + 1]
FFT(A0, Length / 2)
FFT(A1, Length / 2)
complex wn = (cos(2 * Pi / Length), sin(2 * Pi / Length))
complex w = (1, 0)
for int i = 0 to Length / 2 - 1 with i += 1, w *= wn
A[i] = A0[i] + A1[i] * w
A[i + Length / 2] = A0[i] - A1[i] * w
IDFT 是 DFT 的逆运算，令 \(DFT\left(\left\{a_i\right\}\right)=\left\{b_i\right\}\)，已知
存在结论
证明：将前式代入后式
考虑 \(\sum_{i=0}^{n-1}\omega_n^{i\left(j-k\right)}\)
当 \(j=k\)，
当 \(j\neq k\)，
所以，
发现 DFT 和 IDFT 的公式形式一样，所以可以用一个函数实现。
蝴蝶操作
考虑 FFT 的分治过程，以 \(n=16\) 为例
其二进制下表示，
观察发现，若干次蝴蝶操作（奇偶分治）后，其数位二进制翻转后是连续的。
对于二进制翻转，可用递推计算，rev[i] = rev[i >> 1] >> 1 | (i & 1 ? n >> 1 : 0)。
考虑用蝴蝶定理使 FFT 的过程避免递归，可以先将 \(\left\{a_i\right\}\) 按 \(rev\) 序重新排列，在分治树上从下往上，
function FFT(complex A[], int Length, int on)
for int i = 0 to Length - 1 with i += 1
if i < Rev[i]
swap(A[i], A[Rev[i]])
for int k = 1 to n - 1 with k *= 2
complex wn = (cos(Pi / k), sin(on * Pi / k))
for int i = 0 to n with i += k * 2
complex w = (1, 0)
for int j = i to i + k - 1 with j += 1
complex x = a[j]
complex y = a[j + k]
a[j] = x + y
a[j + k] = x - y
if on == -1
for int i = 0 to n - 1 with i += 1
a[i] /= n
然后枚举当前合并的长度 \(2k\)，枚举合并区间开始位置 \(i\)，枚举区间中的元素 \(a_j\)。
struct node{
double x,y;
}A[MAXN],B[MAXN];
node operator+ (node p,node q){return (node){p.x+q.x,p.y+q.y};}
node operator- (node p,node q){return (node){p.x-q.x,p.y-q.y};}
node operator* (node p,node q){return (node){p.x*q.x-p.y*q.y,p.x*q.y+p.y*q.x};}
void FFT(node *X,int flag){
for (int i=0;i<len;i++) if (i<rev[i]) swap(X[i],X[rev[i]]);
for (int l=1;l<len;l<<=1){
node T=(node){cos(Pi/l),flag*sin(Pi/l)};
for (int s=0;s<len;s+=(l<<1)){
node t=(node){1,0};
for (int k=0;k<l;k++,t=t*T){
node Nx=X[s+k],Ny=t*X[s+k+l];
X[s+k]=Nx+Ny,X[s+k+l]=Nx-Ny;
}
}
}
for (register int i=0;i<len;i++) A[i].x/=len;
return;
}
int main(){
n=qr(),m=qr();
for (int i=0;i<=n;i++) A[i].x=qr();
for (int i=0;i<=m;i++) B[i].x=qr();
while (len<=n+m) len<<=1,++L;
for (int i=0;i<len;i++) R[i]=(R[i>>1]>>1)|((i&1)<<(L-1));
FFT(A,1);FFT(B,1);
for (int i=0;i<=len;i++) A[i]=A[i]*B[i];
FFT(A,-1);
for (int i=0;i<=n+m;i++) printf("%d ",(int)(A[i].x+0.5));
return 0;
}
