【[ZJOI2014]力】

题目

好神仙啊

\[F_{j}=\sum_{i<j}\frac{q_iq_j}{(i-j)^2}-\sum_{j<i}\frac{q_iq_j}{(i-j)^2}\]

求\(\frac{F_j}{q_j}\)

显然

\[\frac{F_j}{q_j}=\sum_{i<j}\frac{q_i}{(i-j)^2}-\sum_{j<i}\frac{q_i}{(i-j)^2}\]

先来看前面的那个柿子如何去搞

设\(x=j-i\)

那么

\[\sum_{i<j}\frac{q_i}{(i-j)^2}=\sum_{x=1}^{j-1}\frac{q_{j-x}}{x^2}\]

我们搞出来两个多项式，\(G(x)=\frac{1}{x^2},H(x)=q_x\)

那么就会发现

\[G\times H(j)=\sum_{i=1}^{j-1}G(i)H(j-i)=\sum_{i=1}^{j-1}\frac{q_{j-i}}{i^2}\]

哎这不就是了我们要求的东西了吗

我们发现还有后面那个东西，我们只需要把\(q\)反置再来一遍\(FFT\)就好了

代码

#include<algorithm>
#include<iostream>
#include<cstring>
#include<cstdio>
#include<cmath>
#define maxn 500005
#define re register
#define LL long long
#define max(a,b) ((a)>(b)?(a):(b))
#define min(a,b) ((a)<(b)?(a):(b))
#define double long double
const double Pi=acos(-1);
inline int read()
{
    char c=getchar();int x=0;while(c<'0'||c>'9') c=getchar();
    while(c>='0'&&c<='9') x=(x<<3)+(x<<1)+c-48,c=getchar();
    return x;
}
struct complex
{
    double r,c;
    complex (double a=0,double b=0) {r=a,c=b;}
}f[maxn],g[maxn],og,og1,t;
complex operator +(complex a,complex b) {return complex(a.r+b.r,a.c+b.c);}
complex operator -(complex a,complex b) {return complex(a.r-b.r,a.c-b.c);}
complex operator *(complex a,complex b) {return complex(a.r*b.r-a.c*b.c,a.r*b.c+a.c*b.r);}
int rev[maxn],len,n;
double q[maxn],ans[maxn];
inline void FFT(complex *f,int v)
{
    for(re int i=0;i<=len;i++) if(i<rev[i]) std::swap(f[i],f[rev[i]]);
    for(re int i=2;i<=len;i<<=1)
    {
        int ln=i>>1;
        og1=complex(cos(Pi/ln),v*sin(Pi/ln));
        for(re int l=0;l<len;l+=i)
        {
            og=complex(1,0);
            for(re int x=l;x<l+ln;x++)
            {
                t=og*f[x+ln];
                f[x+ln]=f[x]-t;
                f[x]=f[x]+t;
                og=og*og1;
            }
        }
    }
}
int main()
{
    scanf("%d",&n);
    for(re int i=1;i<=n;i++) scanf("%Lf",&q[i]);
    for(re int i=1;i<=n;i++) f[i].r=q[i],f[i].c=0;
    for(re int i=1;i<=n;i++) g[i].r=1.0/i/i,g[i].c=0;
    len=1;while(len<2*n+2) len<<=1;
    for(re int i=0;i<=len;i++) rev[i]=(rev[i>>1]>>1)|((i&1)?(len>>1):0);
    FFT(f,1),FFT(g,1);
    for(re int i=0;i<len;i++) f[i]=f[i]*g[i];
    FFT(f,-1);
    for(re int i=0;i<len;i++) ans[i]=f[i].r/len;
    memset(f,0,sizeof(f)),memset(g,0,sizeof(g));
    for(re int i=1;i<=n;i++) f[i].r=q[n-i+1],f[i].c=0;
    for(re int i=1;i<=n;i++) g[i].r=1.0/i/i,g[i].c=0;
    FFT(f,1),FFT(g,1);
    for(re int i=0;i<len;i++) f[i]=f[i]*g[i];
    FFT(f,-1);
    for(re int i=1;i<=n;i++) printf("%.3Lf\n",ans[i]-f[n-i+1].r/len);
    return 0;
}

转载于:https://www.cnblogs.com/asuldb/p/10277333.html