Luogu3338 [ZJOI2014]力

https://www.luogu.com.cn/problem/P3338

\(FFT\)

\[E_{i}=\sum_{j=1}^{i-1} \frac{q_j}{(i-j)^2}-\sum_{j=i+1}^{n} \frac{q_j}{(i-j)^2}\\ 改成下标从0开始\\ E_{i}=\sum_{j=0}^{i-1} \frac{q_j}{(i-j)^2}-\sum_{j=i+1}^{n-1} \frac{q_j}{(i-j)^2}\\ 令f_i=q_i,g_i=\frac{1}{i^2}\\ E_{i}=\sum_{j=0}^{i-1} f_j*g_{i-j}-\sum_{j=i+1}^{n-1} f_j*g_{j-i}\\ \sum_{j=0}^{i-1} f_j*g_{i-j}直接FFT\\ \sum_{j=i+1}^{n-1} f_j*g_{j-i}\\ \sum_{j=i+1}^{n-1} f_{n-j-1}^{reverse} g_{j-i}\\ 又一个FFT,Over了 \]

\(C++ Code:\)

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cmath>
#define D double
#define N 400005
using namespace std;
int n,l,s;
D Pi=acos(-1.0);
struct virt
{
    D x,y;
    virt (D xx=0.0,D yy=0.0)
    {
        x=xx,y=yy;
    }
    virt operator + (virt b)
    {
        return virt(x+b.x,y+b.y);
    }
    virt operator - (virt b)
    {
        return virt(x-b.x,y-b.y);
    }
    virt operator * (virt b)
    {
        return virt(x*b.x-y*b.y,x*b.y+y*b.x);
    }
}ans,a[N],b[N],c[N];
int rev[N];
D q[N];
void FFT(virt *a,D t)
{
    for (int i=0;i<s;i++)
        if (i<rev[i])
            swap(a[i],a[rev[i]]);
    for (int mid=1;mid<s;mid <<=1)
    {
        virt wn=virt(cos(Pi/mid),t*sin(Pi/mid));
        for (int j=0;j<s;j+=(mid << 1))
        {
            virt w=virt(1.0,0.0);
            for (int k=0;k<mid;k++,w=w*wn)
            {
                virt x=a[j+k],y=w*a[j+k+mid];
                a[j+k]=x+y;
                a[j+k+mid]=x-y;
            }
        }
    }
}
int main()
{
    scanf("%d",&n);
    for (int i=0;i<n;i++)
        scanf("%lf",&q[i]);
    for (int i=0;i<n;i++)
        a[i]=virt(q[i],0.0);
    for (int i=1;i<n;i++)
        b[i]=virt(1.0/i/i,0.0);
    for (int i=0;i<n;i++)
        c[i]=virt(q[n-i-1],0.0);
    s=1,l=0;
    while (s<n+n)
    {
        s <<=1;
        l++;
    }
    for (int i=0;i<s;i++)
        rev[i]=(rev[i >> 1] >> 1) | ((i & 1) << (l - 1));
    FFT(a,1.0);
    FFT(b,1.0);
    FFT(c,1.0);
    for (int i=0;i<s;i++)
        a[i]=a[i]*b[i],c[i]=c[i]*b[i];
    FFT(a,-1.0);
    FFT(c,-1.0);
    for (int i=0;i<s;i++)
        a[i].x/=s,c[i].x/=s;
    for (int i=0;i<n;i++)
    {
        ans=a[i]-c[n-i-1];
        printf("%.3lf\n",ans.x);
    }
    return 0;
}