快速傅立叶变换（FFT）及其应用

 1 #include<complex>
 2 #include<cmath>
 3 #include<cstdio>
 4 #define pi acos(-1)
 5 using namespace std;
 6 const int maxn=131072;
 7 int n,m,L;
 8 int R[maxn];
 9 complex<double> a[maxn],b[maxn];
10 int read()
11 {
12     int x=0,f=1;char ch=getchar();
13     while(ch<‘0‘||ch>‘9‘){if(ch==‘-‘)f=-1;ch=getchar();}
14     while(ch>=‘0‘&&ch<=‘9‘){x=x*10+ch-‘0‘;ch=getchar();}
15     return x*f;
16 }
17 void fft(complex<double> *a,int f)
18 {
19     for(int i=0;i<n;i++) if(i<R[i]) swap(a[i],a[R[i]]);
20     for(int i=1;i<n;i<<=1)
21     {
22         complex<double> wn(cos(pi/i),f*sin(pi/i));
23         for(int j=0;j<n;j+=(i<<1))
24         {
25             complex<double> w(1,0);
26             for(int k=0;k<i;k++,w*=wn)
27             {
28                 complex<double> x=a[j+k],y=w*a[j+k+i];
29                 a[j+k]=x+y;a[j+k+i]=x-y;
30             }
31         }
32     }
33     if(f==-1) for(int i=0;i<n;i++) a[i]/=n;
34 }
35 int main()
36 {
37     scanf("%d",&n);
38     n--;
39     for(int i=0;i<=n;i++)
40     {
41         a[i]=read();
42         b[n-i]=read();
43     }
44     m=2*n;
45     for(n=1;n<=m;n<<=1) L++;
46     for(int i=0;i<n;i++) R[i]=(R[i>>1]>>1)|((i&1)<<(L-1));
47     fft(a,1);fft(b,1);
48     for(int i=0;i<=n;i++) a[i]*=b[i];
49     fft(a,-1);
50     for(int i=m/2;i<=m;i++) printf("%d\n",(int)(a[i].real()+0.1));
51     return 0;
52 }