这题算是FFT模板题。
一些小细节：
1.其中卷积的时候会发现一个多项式的下标会是负的。我们可以将这个多项式整体向后移n次，然后求出来的多项式的第n+i项就是原来第i项。
2.记得将重复用于FFT的数组先清零。。。
1 #include <iostream> 2 #include <cstdio> 3 #include <cmath> 4 using namespace std; 5 const int MAXN=300005; 6 const double pi=acos(-1); 7 struct complex 8 { 9 double x,y; 10 complex operator + (const complex & b) 11 { 12 return (complex){x+b.x,y+b.y}; 13 } 14 complex operator - (const complex & b) 15 { 16 return (complex){x-b.x,y-b.y}; 17 } 18 complex operator * (const complex & b) 19 { 20 return (complex){x*b.x-y*b.y,x*b.y+y*b.x}; 21 } 22 }a[MAXN],b[MAXN],c[MAXN],d[MAXN]; 23 double q[MAXN]; 24 int x,y,N,L,n,m;int rev[MAXN],Log[MAXN]; 25 void DFT(complex A[],int N,int sign) 26 { 27 static complex C[MAXN]; 28 for (int i=0;i<N;++i) 29 if (i<rev[i]) swap(A[i],A[rev[i]]); 30 int l;complex w,w0; 31 for (int i=1;i<N;i<<=1) 32 { 33 l=i<<1; 34 w=(complex){cos(pi/i),sin(pi/i)*sign}; 35 for (int j=0;j<N;j+=l) 36 { 37 w0=(complex){1,0}; 38 for (int k=0;k<i;++k) 39 { 40 C[j+k]=A[j+k]+w0*A[j+i+k]; 41 C[j+i+k]=A[j+k]-w0*A[j+i+k]; 42 w0=w0*w; 43 } 44 } 45 for (int j=0;j<N;++j) A[j]=C[j]; 46 } 47 if (sign==-1) 48 for (int i=0;i<N;++i) A[i].x/=N; 49 } 50 int FFT(complex a[],complex b[],complex c[],int n,int m) 51 { 52 int L,N; 53 L=Log[n+m]+1; 54 N=1<<L; 55 for (int i=0;i<N;++i) 56 { 57 rev[i]=0; 58 for (int j=0;j<L;++j) 59 if (i&1<<j) rev[i]+=1<<L-j-1; 60 } 61 DFT(a,N,1);DFT(b,N,1); 62 for (int i=0;i<N;++i) c[i]=a[i]*b[i]; 63 DFT(c,N,-1); 64 return N; 65 } 66 int main() 67 { 68 scanf("%d",&n); 69 for (int i=1;i<=n;++i) scanf("%lf",&q[i]); 70 for (int i=2;i<MAXN;++i) Log[i]=Log[i>>1]+1; 71 for (int i=1;i<=n;++i) a[i].x=q[i]; 72 for (int i=1;i<=n;++i) b[i].x=1.0/i/i; 73 FFT(a,b,c,n,n); 74 for (int i=0;i<MAXN;++i) 75 a[i]=b[i]=(complex){0,0}; 76 for (int i=1;i<=n;++i) a[i]=(complex){q[i],0}; 77 for (int i=n-1;i;--i) b[i]=(complex){1.0/(n-i)/(n-i),0}; 78 FFT(a,b,d,n,n-1); 79 for (int i=1;i<=n;++i) 80 printf("%.6f ",c[i].x-d[i+n].x); 81 return 0; 82 }