bzoj3527: [Zjoi2014]力

这题算是FFT模板题。

一些小细节：

　　1.其中卷积的时候会发现一个多项式的下标会是负的。我们可以将这个多项式整体向后移n次，然后求出来的多项式的第n+i项就是原来第i项。

　　2.记得将重复用于FFT的数组先清零。。。

#include <iostream>
#include <cstdio>
#include <cmath>
using namespace std;
const int MAXN=300005;
const double pi=acos(-1);
struct complex
{
    double x,y;
    complex operator + (const complex & b)
    {
        return (complex){x+b.x,y+b.y};
    }
    complex operator - (const complex & b)
    {
        return (complex){x-b.x,y-b.y};
    }
    complex operator * (const complex & b)
    {
        return (complex){x*b.x-y*b.y,x*b.y+y*b.x};
    }
}a[MAXN],b[MAXN],c[MAXN],d[MAXN];
double q[MAXN];
int x,y,N,L,n,m;int rev[MAXN],Log[MAXN];
void DFT(complex A[],int N,int sign)
{
    static complex C[MAXN];
    for (int i=0;i<N;++i)
        if (i<rev[i])    swap(A[i],A[rev[i]]);
    int l;complex w,w0;
    for (int i=1;i<N;i<<=1)
    {
        l=i<<1;
        w=(complex){cos(pi/i),sin(pi/i)*sign};
        for (int j=0;j<N;j+=l)
        {
            w0=(complex){1,0};
            for (int k=0;k<i;++k)
            {
                C[j+k]=A[j+k]+w0*A[j+i+k];
                C[j+i+k]=A[j+k]-w0*A[j+i+k];
                w0=w0*w;
            }
        }
        for (int j=0;j<N;++j)    A[j]=C[j];
    }
    if (sign==-1)
        for (int i=0;i<N;++i)    A[i].x/=N;
}
int FFT(complex a[],complex b[],complex c[],int n,int m)
{
    int L,N;
    L=Log[n+m]+1;
    N=1<<L;
    for (int i=0;i<N;++i)
    {
        rev[i]=0;
        for (int j=0;j<L;++j)
            if (i&1<<j)    rev[i]+=1<<L-j-1;
    }
    DFT(a,N,1);DFT(b,N,1);
    for (int i=0;i<N;++i)    c[i]=a[i]*b[i];
    DFT(c,N,-1);
    return N;
}
int main()
{
    scanf("%d",&n);
    for (int i=1;i<=n;++i)    scanf("%lf",&q[i]);    
    for (int i=2;i<MAXN;++i)    Log[i]=Log[i>>1]+1;
    for (int i=1;i<=n;++i)    a[i].x=q[i];
    for (int i=1;i<=n;++i)    b[i].x=1.0/i/i;
    FFT(a,b,c,n,n);
    for (int i=0;i<MAXN;++i)
        a[i]=b[i]=(complex){0,0};
    for (int i=1;i<=n;++i)    a[i]=(complex){q[i],0};
    for (int i=n-1;i;--i)    b[i]=(complex){1.0/(n-i)/(n-i),0};
    FFT(a,b,d,n,n-1);
    for (int i=1;i<=n;++i)
        printf("%.6f
",c[i].x-d[i+n].x);
    return 0;
}