BZOJ3527: [Zjoi2014]力

n<=100000个实数$q_i$求数组f：

$f(i)=sum_{j<i}frac{q_j}{(j-i)^2}-sum _{j>i}frac{q_j}{(j-i)^2},iin [0,n)$

题目的i和j是反过来的。。难受，不管啦。

两部分，先来第一个：

$sum_{j<i}frac{q_j}{(j-i)^2}= sum_{j=0}^{i-1}q_jfrac{1}{(i-j)^2}$，一个卷积，简单。

再来第二个：

等等你会<号不会>？？把数列反过来不就完了？？等等我们要优雅点，假装是自己“推”出来的：

$sum _{j>i}frac{q_j}{(j-i)^2}=sum_{j=i+1}^{n-1}q_jfrac{1}{(j-i)^2}$

闪一句：令k=j-i。

$=sum_{k=1}^{n-i-1}q_{k+i}frac{1}{k^2}=sum_{k=1}^{n-i-1}q'_{n-i-1-k}frac{1}{k^2}$

很好又一卷积。

写的时候注意下，一个数组乘完就废了。。别再变回去，精度伤不起

#include<cstdio>
#include<cstdlib>
#include<algorithm>
#include<cstdlib>
#include<cstring>
#include<complex>
#include<math.h>
//#include<iostream>
using namespace std;

int n,m,wei;
#define maxn 262222
const double pi=acos(-1);
typedef complex<double> cp;
cp a[maxn],b[maxn],c[maxn],d[maxn],e[maxn],f[maxn]; int rev[maxn];

void dft(cp *a,int n,int type)
{
    for (int i=0;i<n;i++) if (i<rev[i]) {cp t=a[i]; a[i]=a[rev[i]]; a[rev[i]]=t;}
    for (int i=1;i<n;i<<=1)
    {
        cp base=cp(cos(pi/i),type*sin(pi/i));
        for (int j=0,p=(i<<1);j<n;j+=p)
        {
            cp t=cp(1,0);
            for (int k=0;k<i;k++,t*=base)
            {
                cp tmp=t*a[j+k+i];
                a[j+k+i]=a[j+k]-tmp;
                a[j+k]+=tmp;
            }
        }
    }
}

void mul(cp *a,cp *b,cp *c)
{
    if (!rev[1]) for (int i=1;i<n;i++)
    {
        rev[i]=0;
        for (int j=0;j<wei;j++) rev[i]|=((i>>j)&1)<<(wei-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);
    for (int i=0;i<n;i++) c[i]/=n;
}

int main()
{
    scanf("%d",&n); double x;
    for (int i=0;i<n;i++) scanf("%lf",&x),a[i]=d[n-1-i]=x; a[n]=d[n]=0;
    b[0]=e[0]=0; for (int i=1;i<=n;i++) b[i]=e[i]=1.0/i/i;
    m=n+n; for (n=1,wei=0;n<=m;n<<=1,wei++); mul(a,b,c); m>>=1;
    m<<=1; mul(d,e,f); m>>=1;
    for (int i=0;i<m;i++) printf("%.3lf
",c[i].real()-f[m-i-1].real());
    return 0;
}