概念
多项式乘法时,我们发现暴力乘十分缓慢,但是点值乘十分快速。考虑求 \(A\) 和 \(B\) 的卷积。
一个 \(n\) 次多项式可以被 \(n+1\) 个点确定。
设多项式 \(A(x)\) 的系数为 \((a_0,a_1,\cdots,a_n)\)
对其奇偶分类得 \(A(x)=\sum\limits a_{2i}*x^{2i}+\sum a_{2i+1}*x^{2i+1}\)
提取 \(n\) 得 \(A(x)=\sum\limits a_{2i}*x^{2i}+x*\sum a_{2i+1}*x^{2i}\)
化简得 \(A(x)=\sum\limits a_{2i}*(x^2)^i+x*\sum a_{2i+1}*(x^2)^i\)
设 \((\frac{n}{2}-1)\) 次多项式 \(A_0=\sum\limits a_{2i}*x^{2i},A_1=\sum\limits a_{2i+1}*x^{2i}\)
则 \(A(x)=A_0(x^2)+x*A_1(x^2)\)
我们发现可以递归
但是时间复杂度不对
然后是一些高深件
复数
定义 \(i^2=-1\)
则所有形如 \(z=a+b*i\) 的数 \(z\) 组成的集合为复数集,记为 \(C\)。而 \(a\) 称为 \(z\) 的实部, \(b\) 为虚部。
复平面
类似坐标系,纵虚横实。复数 \(z\) 对应从原点指向 \((a,b)\) 的向量。幅角为实轴横方向与向量的夹角,记为 \(θ\)。
加法:与向量加法一样
乘法:复数相乘,模长相乘,幅角相加。\((a+bi)*(c+di)=(ac-bd)+(bc+ad)i\)
单位根
在复平面上,以原点为圆心,\(1\) 为半径作圆,所得的圆叫单位圆。以圆点为起点,圆的 \(n\) 等分点为终点,做 \(n\) 个向量,设幅角为正且最小的向量对应的复数为 \(\omega_n\),称为 \(n\) 次单位根。剩余 \(n-1\) 个点则为 \(\omega_n^2,\omega_n^3,\omega_n^4,\cdots\)。
单位根的幅角为周角的 \(\dfrac{1}{n}\)。
欧拉公式有:\(w_n^k=\cos k\frac{2π}{n}+i\sin k\frac{2π}{n}\)
性质 \(1:\) \(\omega_{n}^{k}=\omega_{2n}^{2k}\)
性质 \(2:\) \(\omega _{n}^{k+\frac{n}{2}}=-\omega _{n}^{k}\)
性质 \(3:\) \(w_n^n=1\)
图上自证不难
FFT
承接高深件
高能推柿子
\(A(x)=A_0(x^2)+x*A_1(x^2)\)
带入 \(x=w_n^k(k<\frac{1}{2}n)\) 得
带入 \(x=w_n^{k+\frac{n}{2}}(k<\frac{1}{2}n)\) 得
不难发现,两式的值可以互相转换,且两室分别处理了一半区间。递归处理,复杂度 \(O(n\log n)\)。
IFFT
函数转点并相乘后,需要转化回系数。此时使用傅里叶逆变换。
拉差
拉差不行,还是要推柿子。
设向量 \(c_1,c_2,\cdots c_{n-1}\) 满足 \(c_k=\sum\limits^{n}_{j=1} y_j(w_n^{-k})^j\)
推柿子懒了,多项式早日消失
最终
\(c_k=n*a_k\)
\(a_k=\frac{c_k}{n}\)
代码
递归实现即可
使用STL:complex实现复数
但是但是,递归常数很大。
改成迭代就可以了。
尻尻板子:
#include<iostream>
#include<cstdio>
#include<complex>
#include<cmath>
#include<algorithm>
using namespace std;
double pi=acos(-1);
int tax[5000001];
void Rev(int n)
{
int d=n>>1;
int len=-1;
tax[++len]=0;
tax[++len]=d;
for(int i=2;i<=n;i<<=1)
{
d>>=1;
for(int p=0;p<i;p++)
{
tax[++len]=tax[p]|d;
}
}
}
void FFT(complex<double> *A,int N)
{
for(int i=1;i<N;i++)
{
if(tax[i]>i)
{
swap(A[i],A[tax[i]]);
}
}
for(int len=2,M=1;len<=N;M=len,len<<=1)
{
complex<double> W(cos(pi/M),sin(pi/M));
complex<double> w(1.0,0.0);
for(int l=0,r=len-1;r<=N;l+=len,r+=len)
{
complex<double> w0=w;
for(int p=l;p<l+M;p++)
{
complex<double> x=A[p]+w0*A[p+M];
complex<double> y=A[p]-w0*A[p+M];
A[p]=x;
A[p+M]=y;
w0*=W;
}
}
}
}
void IFFT(complex<double> *A,int N)
{
FFT(A,N);
reverse(A+1,A+N);
}
int n,m;
complex<double> F[5000001],G[5000001];
int main()
{
scanf("%d%d",&n,&m);
for(int i=0,_;i<=n;i++)
{
scanf("%d",&_);
F[i]=_;
}
for(int i=0,_;i<=m;i++)
{
scanf("%d",&_);
G[i]=_;
}
int p2=1;
while(p2<=m+n) p2<<=1;
Rev(p2);
FFT(F,p2);
FFT(G,p2);
for(int i=0;i<p2;i++)
{
F[i]*=G[i];
}
IFFT(F,p2);
for(int i=0;i<=n+m;i++)
{
printf("%d ",int(F[i].real()/p2+0.5));
}
}
NTT板子
#include<iostream>
#include<cstdio>
#include<complex>
#include<cmath>
#include<algorithm>
#define int long long
#define mod 998244353
using namespace std;
double pi=acos(-1);
int tax[5000001];
int n,m;
void Rev(int pn)
{
for(int i=0;i<pn;i++)
{
tax[i]=tax[i>>1]>>1|((i&1)<<(max((int)ceil(log2(n+m)),1ll)-1));
}
}
long long poww(int a,int b)
{
int re=1;
while(b)
{
if(b&1)
{
re*=a;
re%=mod;
}
a*=a;
a%=mod;
b>>=1;
}
return re;
}
void NTT(long long *a,int n,int type) //type:1正0反
{
for(int i=0;i<n;++i)
{
if(i<tax[i])
{
swap(a[i],a[tax[i]]);
}
}
for(int i=1;i<n;i<<=1)
{
long long gn=poww(3,(mod-1)/(i<<1));
for(int j=0;j<n;j+=(i<<1))
{
long long g0=1;
for(int k=0;k<i;++k,g0=g0*gn%mod)
{
long long x=a[j+k],y=g0*a[i+j+k]%mod;
a[j+k]=(x+y)%mod;
a[i+j+k]=(x-y+mod)%mod;
}
}
}
if(type!=1)
{
int nn=poww(n,mod-2);
reverse(a+1,a+n);
for(int i=0;i<n;i++)
{
a[i]=1ll*a[i]*nn%mod;
}
}
}
long long F[5000001],G[5000001];
signed main()
{
scanf("%lld%lld",&n,&m);
for(int i=0,_;i<=n;i++)
{
scanf("%lld",&_);
F[i]=_;
}
for(int i=0,_;i<=m;i++)
{
scanf("%lld",&_);
G[i]=_;
}
int p2=1;
while(p2<=m+n) p2<<=1;
Rev(p2);
NTT(F,p2,1);
NTT(G,p2,1);
for(int i=0;i<p2;i++)
{
F[i]*=G[i];
}
NTT(F,p2,0);
for(int i=0;i<=n+m;i++)
{
printf("%lld ",F[i]);
}
}