FFT快速傅里叶变换-递归版-带注释模板
题链:uoj#34多项式乘法
题意:给你两个多项式,请输出乘起来后的多项式。
//打法照抄hyc并复制了部分注释
真·存模板系列
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cmath>
#include<iostream>
#include<algorithm>
using namespace std;
#define N 301000
//!!! N>=2^18=262144>n+m>=200000
const double pi=acos(-1);
struct node
{
double x,y;
node(){x=y=0;}
node(double x,double y):x(x),y(y){}
}a[N],b[N];
node operator + (node x,node y) {return node(x.x+y.x,x.y+y.y);}
node operator - (node x,node y) {return node(x.x-y.x,x.y-y.y);}
node operator * (node x,node y) {return node(x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x);}
void fft(node *s,int n,int t)
//t就是标记它是不是逆的FFT
//(最后不是要从点值表达转回系数表达吗t=-1表示处理这个
//用ωn^(-1)替换ωn 并将计算结果的每个元素除以n。
{
if (n==1) return;
node a0[n>>1],a1[n>>1];
for (int i=0;i<=n;i+=2)
a0[i>>1]=s[i],a1[i>>1]=s[i+1];
fft(a0,n>>1,t);fft(a1,n>>1,t);
node wn(cos(2*pi/n),t*sin(2*pi/n)),w(1,0);//就是ωn和ω
//如果要算ωn^(-1)的话 根据负数的指数形式的定义e^(iu)=cos(u)+isin(u)
//ωn^(-1)=e^(-2πi/n) 设u=2πi/n e^(-2πi/n)=cos(-u)+isin(-u)=cos(u)-isin(u)
//所以如果反过来的话 让sin乘个t=-1就好了
for (int i=0;i<(n>>1);i++,w=w*wn)
s[i]=a0[i]+w*a1[i],s[i+(n>>1)]=a0[i]-w*a1[i];
//w^2=(w+n/2)^2 均匀分布在圆上面
//w[i^2,n]=w[i/2,n/2] 折半引理
//s[i]=a0’(i^2)+i*a1’(i^2)=a0(i)+i*a1(i)
//s[i+n/2]=a0’((i+n/2)^2)+i*a1’((i+n/2)^2)=a0’(i^2)-i*a1’(i^2)
//因为i=-(i+n/2) 折半引理
}
int main()
{
//freopen("a.in","r",stdin);
//freopen("a.out","w",stdout);
int n,m,fn,i;
scanf("%d%d",&n,&m);
for (i=0;i<=n;i++) scanf("%lf",&a[i].x);
for (i=0;i<=m;i++) scanf("%lf",&b[i].x);
fn=1;while (fn<=n+m) fn<<=1;
fft(a,fn,1);fft(b,fn,1);
for (i=0;i<=fn;i++) a[i]=a[i]*b[i];
fft(a,fn,-1);
for (i=0;i<=n+m;i++) printf("%d ",(int)(a[i].x/fn+0.5));
//强转int会自动向下取整 所以要加0.5让它四舍五入
printf("\n");
return 0;
}