题链:
题解:
FFT入门题。
(终于接触到迷一样的FFT了)
初学者在对复数和单位根有简单了解的基础上,可以直接看《再探快速傅里叶变换》(毛啸)。
(主要用于求两个序列的卷积)
代码:
递归版:
#include<bits/stdc++.h> #define MAXN 300000 using namespace std; const double Pi=acos(-1); struct Z{ double real,image; Z(double _real=0,double _image=0):real(_real),image(_image){} Z operator - ()const {return Z(-real,-image);} friend Z operator + (const Z &A,const Z &B){return Z(A.real+B.real,A.image+B.image);}; friend Z operator - (const Z &A,const Z &B){return A+(-B);} friend Z operator * (const Z &A,const Z &B){return Z(A.real*B.real-A.image*B.image,A.image*B.real+A.real*B.image);} }; void FFT(int n,Z *Y,int sn){ if(n==1) return; Z L[n>>1],R[n>>1]; for(int k=0;k<n;k+=2) L[k>>1]=Y[k],R[k>>1]=Y[k+1]; FFT(n/2,L,sn); FFT(n/2,R,sn); Z dw=Z(cos(2*Pi/n),sin(sn*2*Pi/n)),w=Z(1,0),tmp; for(int k=0;k<n/2;k++,w=w*dw) tmp=w*R[k],Y[k]=L[k]+tmp,Y[k+n/2]=L[k]-tmp; } int main(){ static Z A[MAXN],B[MAXN]; int n,m; scanf("%d%d",&n,&m); n++; m++; for(int i=0,x;i<n;i++) scanf("%d",&x),A[i]=Z(x,0); for(int i=0,x;i<m;i++) scanf("%d",&x),B[i]=Z(x,0); m=n+m-1; for(n=1;n<m;n<<=1); FFT(n,A,1); FFT(n,B,1); for(int i=0;i<n;i++) A[i]=A[i]*B[i]; FFT(n,A,-1); for(int i=0;i<m;i++) printf("%d ",(int)(A[i].real/n+0.5)); return 0; }
非递归版:
#include<bits/stdc++.h> #define MAXN 300000 using namespace std; const double Pi=acos(-1); struct Z{ double real,image; Z(double _real=0,double _image=0):real(_real),image(_image){} Z operator - ()const {return Z(-real,-image);} friend Z operator + (const Z &A,const Z &B){return Z(A.real+B.real,A.image+B.image);}; friend Z operator - (const Z &A,const Z &B){return A+(-B);} friend Z operator * (const Z &A,const Z &B){return Z(A.real*B.real-A.image*B.image,A.image*B.real+A.real*B.image);} }; int order[MAXN]; void FFT(int n,Z *Y,int sn){ for(int i=0;i<n;i++) if(i<order[i]) swap(Y[i],Y[order[i]]); for(int d=2;d<=n;d<<=1){ Z dw=Z(cos(2*Pi/d),sin(sn*2*Pi/d)),w,tmp; for(int i=0;w=Z(1,0),i<n;i+=d) for(int k=i;k<i+d/2;k++,w=w*dw) tmp=w*Y[k+d/2],Y[k+d/2]=Y[k]-tmp,Y[k]=Y[k]+tmp; } } int main(){ static Z A[MAXN],B[MAXN]; int n,m,len; scanf("%d%d",&n,&m); n++; m++; for(int i=0,x;i<n;i++) scanf("%d",&x),A[i]=Z(x,0); for(int i=0,x;i<m;i++) scanf("%d",&x),B[i]=Z(x,0); m=n+m-1; for(len=0,n=1;n<m;n<<=1) len++; for(int i=1;i<n;i++) order[i]=(order[i>>1]>>1)|((i&1)<<(len-1)); FFT(n,A,1); FFT(n,B,1); for(int i=0;i<n;i++) A[i]=A[i]*B[i]; FFT(n,A,-1); for(int i=0;i<m;i++) printf("%d ",(int)(A[i].real/n+0.5)); return 0; }