题解 P3803 【模板】多项式乘法（FFT）

感觉题解区不是写的太高深，就是写的太高深。所以给初中、小学和幼儿园的萌新准备一篇简单易懂的良心题解~

前置知识

一、多项式的系数表示法和点值表示法。$A(x)=\sum\limits_{i=0}^{n-1}a_i\cdot x^i$

系数：$(a_0,a_1,a_2...a_{n-2},a_{n-1})$。

点值：$(x_0,y_0),(x_1,y_1)...(x_{n-1},y_{n-1})$。易证，$n$ 个点确定一个多项式。

二、向量的加减乘。

三、圆的弧度制：$1^{\circ}=\frac{\pi}{180}rad,180^{\circ}=\pi rad$

四、复数：设 $a,b$ 为实数，$i^2=-1$，形如 $a+ib$ 的数叫复数。$x$ 为复平面的实数轴，$y$ 为虚数轴，从 $(0,0)$ 到 $(a,b)$ 的向量表示复数 $a+ib$，其中模长为 $\sqrt{a^2+b^2}$，幅角为以 $x$ 轴正半轴到已知向量的转角的有向角。

五、单位根。（下文默认 $n$ 为 $2$ 的正整数次幂）

在复平面上以 $(0,0)$ 为圆心，$1$ 为半径做圆。$n$ 等分后再做 $n$ 个向量。则幅角为正且最小的向量对应复数 $\omega_n$，称为 $n$ 次单位根。根据向量的乘法，其余 $n-1$ 个复数为 $\omega_n^2,\omega_n^3...\omega_n^n$。其中 $\omega_n^k=\cos{k\cdot\frac{2\pi}{n}}+\sin{k\cdot\frac{2\pi}{n}}$。

接下来证明两个之后会用到的公式。

1. 证明 $\omega_{2n}^{2k}=\omega_n^k$

$$\begin{aligned} \omega_{2n}^{2k} &=\cos{2k\cdot\frac{2\pi}{2n}}+\sin{2k\cdot\frac{2\pi}{2n}}\\ &=\cos{k\cdot\frac{2\pi}{n}}+\sin{k\cdot\frac{2\pi}{n}}\\ &=\omega_n^k \end{aligned}$$

2. 证明 $w_n^{k+\frac{n}{2}}=-\omega_n^k$

$$\begin{aligned} \omega_n^{k+\frac{n}{2}} &=\omega_n^k\cdot\omega_n^{\frac{n}{2}}\\ &=\omega_n^k\cdot(\cos{\frac{n}{2}\cdot\frac{2\pi}{n}}+\sin{\frac{n}{2}\cdot\frac{2\pi}{n}})\\ &=\omega_n^k\cdot(\cos\pi+\sin\pi)\\ &=-\omega_n^k \end{aligned}$$

特殊的，$w_n^0=w_n^n=1$。

正文

这里的 FFT 是要解决两个多项式相乘的问题。

如果用系数表示法直接暴力相乘，复杂度 $O(n^2)$，很难进一步优化，所以考虑点值表示法。

点值表示法的做法为，对于两个多项式分别找出足够的点，将其对应的函数值相乘，再通过某种方式转回系数表示法。

然而，如果我们随便选 $n$ 个点并求出其对应的值，单次复杂度线性，执行 $n$ 次 $O(n^2)$，还是不行。

因此我们所选的 $n$ 个点必须要有很好的性质，从而减少计算量。

而 $(\omega_n,\omega_n^2,\omega_n^3...\omega_n^n)$ 刚好有着很强的性质，可以作为所选的 $n$ 个值，通过分治可将复杂度降至 $O(n\log n)$。

接下来会简单介绍做法。

定义多项式 $A(x)=(a_0,a_1,a_2...a_{n-1})$。

按下标奇偶分类后得：

$$A(x)=(a_0+a_2\cdot x^2+a_4\cdot x^4+...+a_{n-2}\cdot x^{n-2})+(a_1\cdot x+a_3\cdot x^3+a_5\cdot x^5+...+a_{n-1}\cdot x^{n-1})$$

定义：

$$A_1(x)=a_0+a_2\cdot x+a_4\cdot x^2+a_6\cdot x^3+...+a_{n-2}\cdot x^{\frac{n}{2}-1}$$

$$A2(x)=a_1+a_3\cdot x+a_5\cdot x^2+a_7\cdot x^3+...+a_{n-1}\cdot x^{\frac{n}{2}-1}$$

那么可得：

$$A(x)=A_1(x^2)+x\cdot A_2(x^2)$$

带入 $\omega_n^k(k<\frac{n}{2})$ 得：

$$A(\omega_n^k)=A_1(\omega_n^{2k})+\omega_n^k\cdot A_2(\omega_n^{2k})$$

带入 $\omega_n^{k+\frac{n}{2}}$ 得：

$$\begin{aligned} A(\omega_n^{k+\frac{n}{2}}) &=A_1(\omega_n^{2k+n})+\omega_n^{k+\frac{n}{2}}\cdot A_2(\omega_n^{2k+n})\\ &=A_1(\omega_n^{2k}\cdot\omega_n^n)-\omega_n^k\cdot A_2(\omega_n^{2k}\cdot \omega_n^n)\\ &=A_1(\omega_n^{2k})-\omega_n^k\cdot A_2(\omega_n^{2k}) \end{aligned}$$

可以发现，当 $k$ 取遍 $[0,\frac{n}{2})$ 和 $[\frac{n}{2},n)$ 时，$A_1,A_2$ 转移项下标都取遍 $[\frac{n}{2},n)$。第二个式子与第一个式子只差一个符号，所以递归完后可以同时算，复杂度 $O(n\log n)$。

这步操作也被称为蝴蝶操作，即将 $A_1(\omega_n^i)$ 向 $A(\omega_n^i)$ 和 $A(\omega_n^{i+\frac{n}{2}})$ 连边，$A_2(\omega_n^i)$ 向 $A(\omega_n^i)$ 和 $A(\omega_n^{i+\frac{n}{2}})$ 连边，形成的图形像蝴蝶的双翅。（其实这里的连边是上标之间的连边，$A_2$ 接在 $A_1$ 右边）

可能会有人问：自变量为什么必须是 $\omega_n^k$，而不能是 $k$ 呢？

原因很简单：当 $k$ 取遍 $[0,\frac{n}{2})$ 和 $[\frac{n}{2},n)$ 时，两段对应的区间完全不同，导致搜索量无法减半，复杂度 $O(n^2)$。

这里想明白 FFT 基本就算学明白了。

现在我们已经得到新多项式的 $n$ 个点了，考虑怎么转化回系数，即 IFFT。

因为 $x_i=\omega_n^i,y_i=\sum\limits_{t=0}^{n-1}a_t\cdot x_i^t$。所以 $y_i=\sum\limits_{t=0}^{n-1}a_t\cdot(\omega_n^i)^t$。

定义 $c_k=\sum\limits_{i=0}^{n-1}y_i\cdot (w_n^{-k})^i$。

那么，

$$\begin{aligned} c_k &=\sum\limits_{i=0}^{n-1}\left[\left(\sum\limits_{j=0}^{n-1}a_j\cdot (\omega_n^i)^j\right)\cdot(w_n^{-k})^i\right]\\ &=\sum\limits_{i=0}^{n-1}\left(\sum\limits_{j=0}^{n-1}a_j\cdot(\omega_n^i)^j\cdot (\omega_n^{-k})^i\right)\\ &=\sum\limits_{j=0}^{n-1}a_j\cdot\sum\limits_{i=0}^{n-1}(w_n^j)^i\cdot(w_n^{-k})^i\\ &=\sum\limits_{j=0}^{n-1}a_j\cdot\sum\limits_{i=0}^{n-1}(w_n^{j-k})^i \end{aligned}$$

定义 $s(x)=\sum\limits_{i=0}^{n-1}x^i$。

两边同乘 $x$ 得：$x\cdot s(x)=\sum\limits_{i=1}^{n}x^i$。

两者相减得：$(x-1)s(x)=x^n-1$。

所以 $s(x)=\frac{x^n-1}{x-1}$。

带入 $\omega_n^k(k\neq n)$ 得：$s(x)=\frac{(w_n^k)^n-1}{w_n^k-1}=\frac{1-1}{w_n^k-1}=0$。

当 $j\neq k$ 时 $\sum\limits_{i=0}^{n-1}(w_n^{j-k})^i=0$。

当 $j=k$ 时 $\sum\limits_{i=0}^{n-1}(w_n^{j-k})^i=n$。

所以 $c_k=n\cdot a_k$，$a_k=\frac{c_k}{n}$。

FFT 算法到这里就结束了。

具体实现时，可对多项式 $A(x),B(x)$ 分别跑一遍 FFT，对应函数值相乘做一遍 IFFT 还原系数。

然后我们发现当递归到最底层时，序列各个元素编号二进制反转后是连续的。

所以可以先把最底层的弄出来，自底向上合并，常数很小。

code:

#include<bits/stdc++.h>
using namespace std;
typedef double db;
const int N=1e7+5;
const db pi=acos(-1.0);
int n,m,r[N],lim=1;
struct node{db x,y;node(db x=0.0,db y=0.0):x(x),y(y){}};
node operator+(node a,node b){return node(a.x+b.x,a.y+b.y);}
node operator-(node a,node b){return node(a.x-b.x,a.y-b.y);}
node operator*(node a,node b){return node(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}
void FFT(node *a,int op){
  for(int i=0;i<lim;++i)if(i<r[i])swap(a[i],a[r[i]]);
  for(int t=1;t<lim;t<<=1){
    node B=node(cos(pi/t),op*sin(pi/t)),mi,x,y;
    for(int len=t<<1,j=0;j<lim;j+=len){
      mi=node(1,0);
      for(int k=0;k<t;++k,mi=mi*B){
        x=a[j+k];y=mi*a[j+t+k];
        a[j+k]=x+y;a[j+t+k]=x-y;
      }
    }
  }
}
node a[N],b[N];
int main(){
  ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);
  cin>>n>>m;
  for(int i=0;i<=n;++i)cin>>a[i].x;
  for(int i=0;i<=m;++i)cin>>b[i].x;
  int l=0;
  while(lim<=n+m)lim<<=1,++l;
  for(int i=0;i<lim;++i)r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
  FFT(a,1);FFT(b,1);
  for(int i=0;i<lim;++i)a[i]=a[i]*b[i];
  FFT(a,-1);
  for(int i=0;i<=n+m;++i)cout<<(int)(a[i].x/lim+0.5)<<" ";
  cout<<endl;
  return 0;
}