FFT(1)

FFT

Complex

struct complex{
	double re,im;
	complex(double r,double i){re=r,im=i;}
	complex(){re=0.0,im=0.0;}
	complex operator+(complex b){return complex(re+b.re,im+b.im);}
	complex operator-(complex b){return complex(re-b.re,im-b.im);}
	complex operator*(complex b){return complex(re*b.re-im*b.im,re*b.im+im*b.re);}
};

FFT

void fft(complex* a,int n,int* bitrev,int f){
	for(i=0;i<n;++i) if(i<bitrev[i]) swap(a[i],a[bitrev[i]]);
	for(i=1;i<n;i<<=1){
		complex wn(cos(pi/i),f*sin(pi/i));
		p=i<<1;
		for(j=0;j<n;j+=p){
			complex w(1,0);
			for(k=0;k<i;++k,w=w*wn){
				t=a[j+k],t2=w*a[j+k+i];
				a[j+k]=t+t2,a[j+k+i]=t-t2;
			}
		}
	}
	if(f<0) for(i=0;i<n;++i) a[i].re/=n,a[i].im/=n;
}

bitrev

for(N=1;N<pt;N<<=1,++L);
for(i=0;i<N;++i) bitrev[i]=(bitrev[i>>1]>>1)|((i&1)<<(L-1));

UR#34

这题我第一次写FFT.

#include <cstdio>
#include <cmath>
struct complex{
	double re,im;
	complex(double r,double i){re=r,im=i;}
	complex(){re=0.0,im=0.0;}
	complex operator+(complex b){return complex(re+b.re,im+b.im);}
	complex operator-(complex b){return complex(re-b.re,im-b.im);}
	complex operator*(complex b){return complex(re*b.re-im*b.im,re*b.im+im*b.re);}
};
#define pi 3.1415926535897932384626
complex t,t2;
int i,j,k,p;
void swap(complex& a,complex& b){
	t=a;
	a=b;
	b=t;
}
void fft(complex* a,int n,int* bitrev,int f){
	for(i=0;i<n;++i) if(i<bitrev[i]) swap(a[i],a[bitrev[i]]);
	for(i=1;i<n;i<<=1){
		complex wn(cos(pi/i),f*sin(pi/i));
		p=i<<1;
		for(j=0;j<n;j+=p){
			complex w(1,0);
			for(k=0;k<i;++k,w=w*wn){
				t=a[j+k],t2=w*a[j+k+i];
				a[j+k]=t+t2,a[j+k+i]=t-t2;
			}
		}
	}
	if(f<0) for(i=0;i<n;++i) a[i].re/=n,a[i].im/=n;
}
int n,m,bitrev[1000000];
complex a[1000000],b[1000000];
int N,L,pt;
int main(){
	scanf("%d%d",&n,&m);
	++n;++m;
	for(i=0;i<n;++i) scanf("%lf",&(a[i].re));
	for(i=0;i<m;++i) scanf("%lf",&(b[i].re));
	pt=n+m-1;
	for(N=1;N<pt;N<<=1,++L);
	for(i=0;i<N;++i) bitrev[i]=(bitrev[i>>1]>>1)|((i&1)<<(L-1));
	fft(a,N,bitrev,1);
	fft(b,N,bitrev,1);
	for(i=0;i<N;++i) a[i]=a[i]*b[i];
	fft(a,N,bitrev,-1);
	for(i=0;i<pt;++i) printf("%d ", (int)(a[i].re+0.5));
	return 0;
}

原文地址:https://www.cnblogs.com/tmzbot/p/4627360.html