这是一个悲伤的故事。。
打算从今天开始重学FFT,然后顺便学NTT
因为很多基本操作都不记得了。。
一开始是想在学校教了复数之后再学的。。
但是似乎并没有这个时间了。。
于是就去自己看了高中的复数内容。。
然后发现和高中的微积分类容一样,都是没什么用的东西。。
和我想要的根本不是一个东西,好吧。。
然后去苦战了一波算法导论
发现我看不懂。。
投降。。
于是被下几个结论,就和当时第二次学的时候一样草草收尾了。。
这是第三次了,都没有学会TAT
于是根据结论,我还是只会递归版本
Update4.18
递归版的实现虽然复杂度是一样的,但是常数实在是太大了,太大了,太大了。实在受不了,于是搞了一些迭代的版本
一般来说,速度可能是3,4,5,6倍等等
迭代版的
NTT
void ntt(LL *a,LL n,LL op)
{
for (LL u=0;u<n;u++) bin[u]=(bin[u>>1]>>1)|((u&1)*(n>>1));
for (LL u=0;u<n;u++) if (u<bin[u]) swap(a[u],a[bin[u]]);
for (LL u=1;u<n;u<<=1)
{
LL wn=pow(op==1?g:gi,(MOD-1)/(u<<1)),w,t;
for (LL i=0;i<n;i=i+(u<<1))
{
w=1;
for (LL k=0;k<u;k++)
{
t=w*a[u+i+k]%MOD;
a[u+i+k]=(a[i+k]-t+MOD)%MOD;
a[i+k]=(a[i+k]+t)%MOD;
w=w*wn%MOD;
}
}
}
if(op==-1)
{
LL Inv=pow(n,MOD-2);
for(LL i=0;i<n;i++) a[i]=a[i]*Inv%MOD;
}
}
FFT
void fft (qq *a,int n,int o)
{
for (int u=0;u<n;u++) bin[u]=(bin[u>>1]>>1)|((u&1)*(n>>1));
for (int u=1;u<n;u++)
if (u<bin[u])
swap(a[u],a[bin[u]]);
for (int u=1;u<n;u<<=1)
{
qq wn=qq(cos(2*pi/(u<<1)),o*sin(2*pi/(u<<1)));
for (int i=0;i<n;i+=(u<<1))
{
qq w=qq(1,0);
for (int j=0;j<u;j++)
{
qq t=w*a[u+i+j];
a[u+i+j]=a[i+j]-t;
a[i+j]=a[i+j]+t;
w=w*wn;
}
}
}
}
递归版的
FFT
#include<cstdio>
#include<iostream>
#include<cstring>
#include<cmath>
#include<complex>
using namespace std;
const double pi=M_PI;
const int N=1000005;
int n,m;
complex<double> a[N],b[N];
void fft (complex<double> *a,int n,int o)
{
if (n==1) return ;
int k=(n>>1);
complex<double> w=1,wn(cos(2*pi/n),o*sin(2*pi/n)),a0[k],a1[k];
for (int u=0;u<k;u++)
{
int i=u*2;
a0[u]=a[i];
a1[u]=a[i+1];
}
fft(a0,k,o);fft(a1,k,o);
for (int u=0;u<k;u++)
{
a[u]=a0[u]+w*a1[u];
a[u+k]=a0[u]-w*a1[u];
w=w*wn;
}
}
int main()
{
scanf("%d%d",&n,&m);
for (int u=0;u<=n;u++) scanf("%lf",&a[u]);
for (int u=0;u<=m;u++) scanf("%lf",&b[u]);
m=m+n;n=1;while (n<=m) n<<=1;
fft(a,n,1);fft(b,n,1);
for (int u=0;u<=n;u++) a[u]*=b[u];
fft(a,n,-1);
for (int u=0;u<=m;u++)
printf("%d ",(int)(a[u].real()/n+0.5));
return 0;
}
NTT
就是用原根代替单位复数根
使用范围一点都不广。。
什么?你问我原根怎么找?
暴力找啊!一般都不大的
基本上和上面这个代码是一样的
CODE:
#include<cstdio>
#include<iostream>
#include<cstring>
#include<cmath>
#include<complex>
using namespace std;
typedef long long LL;
const LL MOD=998244353,g=3,gi=332748118;
const LL N=1000005;
LL n,m;
LL a[N],b[N];
LL pow (LL x,LL y)
{
if (y==1) return x;
LL lalal=pow(x,y>>1);
lalal=lalal*lalal%MOD;
if (y&1) lalal=lalal*x%MOD;
return lalal;
}
void ntt (LL *a,LL n,LL o)
{
if (n==1) return ;
LL k=(n>>1);
LL w=1,wn=pow(o==1?g:gi,(MOD-1)/n),a0[k],a1[k];
for (LL u=0;u<k;u++)
{
LL i=u*2;
a0[u]=a[i];
a1[u]=a[i+1];
}
ntt(a0,k,o);ntt(a1,k,o);
for (LL u=0;u<k;u++)
{
a[u]=a0[u]+w*a1[u]%MOD;
a[u]=(a[u]%MOD+MOD)%MOD;
a[u+k]=a0[u]-w*a1[u];
a[u+k]=(a[u+k]%MOD+MOD)%MOD;
w=w*wn%MOD;
}
}
int main()
{
scanf("%I64d%I64d",&n,&m);
for (LL u=0;u<=n;u++) scanf("%I64d",&a[u]);
for (LL u=0;u<=m;u++) scanf("%I64d",&b[u]);
m=m+n;n=1;while (n<=m) n<<=1;
ntt(a,n,1);ntt(b,n,1);
for (LL u=0;u<=n;u++) a[u]*=b[u];
ntt(a,n,-1);
LL inv=pow(n,MOD-2);
for (LL u=0;u<=m;u++) printf("%I64d ",a[u]*inv%MOD);
return 0;
}