(各种)FFT模板

先来一个标准的归并版FFT 2881ms

#include<math.h>
#include<stdio.h>
#include<string.h>
#include<algorithm>
#define D double
using namespace std;
struct Z{ D x,y; } a[280010],b[280010],z[280010];
inline Z operator+ (Z x,Z y){ return (Z){x.x+y.x,x.y+y.y}; }
inline Z operator- (Z x,Z y){ return (Z){x.x-y.x,x.y-y.y}; }
inline Z operator* (Z x,Z y){ return (Z){x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x}; }
void FFT(Z* a,Z* b,int l,int r,int g){
	if(l==r) return;
	int m=l+r>>1,n=r-l+1,i,j=l,k=m+1;
	Z w=(Z){cos(2*M_PI/n),g*sin(2*M_PI/n)},z=(Z){1,0};
	for(i=l;i<=r;++i) b[~i&1?j++:k++]=a[i];
	FFT(b,a,l,m,g); FFT(b,a,m+1,r,g);
	for(i=l,j=m-l+1;i<=m;++i,z=z*w){
		a[i]=b[i]+z*b[i+j];
		a[i+j]=b[i]-z*b[i+j];
	}
}
int n,m,M=1,s[280010];
int main(){
	scanf("%d%d",&n,&m);
	for(int i=0;i<=n;++i) scanf("%lf",&a[i].x);
	for(int i=0;i<=m;++i) scanf("%lf",&b[i].x);
	while(M<=max(n,m)) M<<=1; M<<=1;
	FFT(a,z,0,M-1,1);
	FFT(b,z,0,M-1,1);
	for(int i=0;i<M;++i) a[i]=a[i]*b[i];
	FFT(a,z,0,M-1,-1);
	for(int i=0;i<M;++i) s[i]+=int(a[i].x/M+0.5);
	for(int i=0;i<=n+m;++i) printf("%d ",s[i]);
}
让后就可以写一个蝴蝶变换版本的 1097ms

#include<math.h>
#include<stdio.h>
#include<string.h>
#include<algorithm>
#define D double
using namespace std;
struct Z{ D x,y; } a[280010],b[280010],z[280010];
inline Z operator+ (Z x,Z y){ return (Z){x.x+y.x,x.y+y.y}; }
inline Z operator- (Z x,Z y){ return (Z){x.x-y.x,x.y-y.y}; }
inline Z operator* (Z x,Z y){ return (Z){x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x}; }
inline Z rt(double X){ return (Z){cos(X),sin(X)}; }
void FFT(Z* a,Z* b,int n,int g){
	for(int i=0,j,k,t;i<n;++i){
		for(j=0,k=i,t=n-1;t;t>>=1,k>>=1) j=(j<<1)|(k&1);
		b[j]=a[i];
	}
	for(int m=2;m<=n;m<<=1){
		Z w=rt(g*M_PI/m),z=(Z){1,0},u,v;
		for(int i=0,k=m>>1;i<k;++i){
			for(int j=i;j<n;j+=m){
				u=b[j]; v=b[j+k]*z;
				b[j]=u+v; b[j+k]=u-v;
			}
			z=z*w;
		}
	}
	for(int i=0;i<n;++i) a[i]=b[i];
}
int n,m,M=1,s[280010];
int main(){
	scanf("%d%d",&n,&m);
	for(int i=0;i<=n;++i) scanf("%lf",&a[i].x);
	for(int i=0;i<=m;++i) scanf("%lf",&b[i].x);
	while(M<=max(n,m)) M<<=1; M<<=1;
	FFT(a,z,M,2);
	FFT(b,z,M,2);
	for(int i=0;i<M;++i) a[i]=a[i]*b[i];
	FFT(a,z,M,-2);
	for(int i=0;i<M;++i) s[i]+=int(a[i].x/M+0.5);
	for(int i=0;i<=n+m;++i) printf("%d ",s[i]);
}
加速蝴蝶变换版本的 620ms

#include<math.h>
#include<stdio.h>
#include<string.h>
#include<algorithm>
#define D double
using namespace std;
int n,m,M=1,L,s[280010],r[280010];
struct Z{ D x,y; } a[280010],b[280010],z[280010];
inline Z rt(double X){ return (Z){cos(X),sin(X)}; }
inline Z operator+ (Z x,Z y){ return (Z){x.x+y.x,x.y+y.y}; }
inline Z operator- (Z x,Z y){ return (Z){x.x-y.x,x.y-y.y}; }
inline Z operator* (Z x,Z y){ return (Z){x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x}; }
void FFT(Z* a,Z* b,int n,int g){
	for(int i=0;i<n;++i) b[r[i]]=a[i];
	for(int m=2;m<=n;m<<=1){
		Z w=rt(g*M_PI/m),u,v,z;
		for(int i,k=m>>1,j=0;j<n;j+=m){
			for(z=(Z){1,i=0};i<k;++i,z=z*w){
				u=b[i+j]; v=z*b[i+j+k];
				b[i+j]=u+v; b[i+j+k]=u-v;
			}
		}
	}
	memcpy(a,b,n*sizeof(Z));
}
int main(){
	scanf("%d%d",&n,&m);
	for(int i=0;i<=n;++i) scanf("%lf",&a[i].x);
	for(int i=0;i<=m;++i) scanf("%lf",&b[i].x);
	for(;M<=max(n,m);++L) M<<=1; M<<=1;
	for(int i=0;i<M;++i) r[i]=(r[i>>1]>>1)|((i&1)<<L);
	FFT(a,z,M,2);
	FFT(b,z,M,2);
	for(int i=0;i<M;++i) a[i]=a[i]*b[i];
	FFT(a,z,M,-2);
	for(int i=0;i<M;++i) s[i]+=int(a[i].x/M+0.5);
	for(int i=0;i<=n+m;++i) printf("%d ",s[i]);
}

最终版本:490ms

#include<cctype>
#include<math.h>
#include<stdio.h>
#include<string.h>
#include<algorithm>
#define D double
#define Rg register
#define N 280010
using namespace std;
int n,m,M=1,L,s[N],r[N];
struct Z{ D x,y; } a[N],b[N],z[N];
inline Z rt(D X){ return (Z){cos(X),sin(X)}; }
inline Z operator+ (Z x,Z y){ return (Z){x.x+y.x,x.y+y.y}; }
inline Z operator- (Z x,Z y){ return (Z){x.x-y.x,x.y-y.y}; }
inline Z operator* (Z x,Z y){ return (Z){x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x}; }
inline int rd(){  
	int X=0,w=0; char ch=0;
    while(!isdigit(ch)) {w|=ch=='-';ch=getchar();}
    while(isdigit(ch)) X=(X<<3)+(X<<1)+(ch^48),ch=getchar();
    return w?-X:X; 
} 
void FFT(Z* a,Z* b,int n,int g){
	for(int i=0;i<n;++i) b[r[i]]=a[i];
	for(int m=2;m<=n;m<<=1){
		Rg Z w=rt(g*M_PI/m),u,v,z;
		for(Rg int i,k=m>>1,j=0;j<n;j+=m){
			for(z=(Z){1,i=0};i<k;++i,z=z*w){
				u=b[i+j]; v=z*b[i+j+k];
				b[i+j]=u+v; b[i+j+k]=u-v;
			}
		}
	}
	memcpy(a,b,n*sizeof(Z));
}
int main(){
	n=rd(); m=rd();
	for(int i=0;i<=n;++i) a[i].x=rd();
	for(int i=0;i<=m;++i) b[i].x=rd();
	for(;M<=max(n,m);++L) M<<=1; M<<=1;
	for(int i=0;i<M;++i) r[i]=(r[i>>1]>>1)|((i&1)<<L);
	FFT(a,z,M,2);
	FFT(b,z,M,2);
	for(int i=0;i<M;++i) a[i]=a[i]*b[i];
	FFT(a,z,M,-2);
	for(int i=0;i<M;++i) s[i]+=int(a[i].x/M+0.5);
	for(int i=0;i<=n+m;++i) printf("%d ",s[i]);
}

数论版(NTT) 649ms

#include<cctype>
#include<math.h>
#include<stdio.h>
#include<string.h>
#include<algorithm>
#define N 280010
#define Rg register
#define LL long long
#define M 1004535809 
using namespace std;
int n,m,T,L,r[N],a[N],b[N],z[N]; LL W[N],inv;
inline LL pow(LL x,LL k,LL& s){
	for(s=1;k;x=x*x%M,k>>=1) k&1?s=s*x%M:0;
}
inline int rd(){  
	int X=0; char ch=0;
    while(!isdigit(ch)) ch=getchar();
    while(isdigit(ch)) X=(X<<3)+(X<<1)+(ch^48),ch=getchar();
    return X; 
}
void NTT(int* a,int* b,int n,int g){
	for(int i=0;i<n;++i) b[r[i]]=a[i];
	for(int m=2;m<=n;m<<=1){
		LL w=W[(g?(n/m):(n-n/m))],u,v,z;
		for(int i,k=m>>1,j=0;j<n;j+=m)
			for(z=1,i=0;i<k;++i,z=z*w%M){
				u=b[i+j]; v=z*b[i+j+k]%M;
				b[i+j]=(u+v)%M; b[i+j+k]=(u-v+M)%M;
			}
	}
	memcpy(a,b,n<<2);
}
int main(){
	n=rd(); m=rd();
	for(int i=0;i<=n;++i) a[i]=rd();
	for(int i=0;i<=m;++i) b[i]=rd();
	T=n+m; m=max(n,m); n=1;
	for(;n<=m;++L) n<<=1; n<<=1;
	for(int i=0;i<n;++i) r[i]=(r[i>>1]>>1)|((i&1)<<L);
	pow(3,(M-1)/n,W[*W=1]);
	for(int i=2;i<=n;++i) W[i]=W[i-1]*W[1]%M;
	NTT(a,z,n,1);
	NTT(b,z,n,1);
	for(int i=0;i<n;++i) a[i]=(LL)a[i]*b[i]%M;
	NTT(a,z,n,0);
	pow(n,M-2,inv);
	for(int i=0;i<=T;++i) printf("%lld ",a[i]*inv%M);
}


原文地址:https://www.cnblogs.com/Extended-Ash/p/9477140.html