bzoj3527 [Zjoi2014]力

http://www.elijahqi.win/2018/03/15/bzoj3527/
Description
给出n个数qi，给出Fj的定义如下：

令Ei=Fi/qi，求Ei.

Input
第一行一个整数n。
接下来n行每行输入一个数，第i行表示qi。
n≤100000，0 < qi < 1000000000

Output
n行，第i行输出Ei。与标准答案误差不超过1e-2即可。

Sample Input
5
4006373.885184
15375036.435759
1717456.469144
8514941.004912
1410681.345880
Sample Output
-16838672.693
3439.793
7509018.566
4595686.886
10903040.872
首先将公式里的q[i]消去 然后直接写就可以发现前半个sigma已经是卷积的形式了 那么针对后面的sigma我们直接采取换元的形式发现如果把式子倒着来 然后再* 那个平方项就又是一个卷积形式了 做两次fft即可
$Ei=\sum_{i>j}^{}{qj \over (i-j)^2}-\sum_{i<j}^{}{qj \over (j-i)^2}$
$Ei=\sum_{j=1}^{i-1}{qj \over (i-j)^2}-\sum_{j=i+1}^{n}{qj \over (j-i)^2}$
$Ei=\sum_{j=1}^{i-1}{a[j]*b[i-j]}-\sum_{j=i+1}^{n}{a[j]*b[j-i]}$
$a[i]=qi,b[i]={1 \over i^2}$
后面一块换元
$\sum_{j=0}^{n-1-i}{a[n-1-j]*b[j-i]}$

#include<cmath>
#include<cstdio>
#include<complex>
#include<algorithm>
#define N 110000
#define pi acos(-1)
using namespace std;
struct C{
    double a,b;
    inline C operator +(const C &y){return (C){a+y.a,b+y.b};}
    inline C operator -(const C &y){return (C){a-y.a,b-y.b};}
    inline C operator *(const C &y){return (C){a*y.a-b*y.b,a*y.b+b*y.a};}
    inline void operator *=(const C &y){*this=*this*y;}
}a[N<<2],a1[N<<2],b[N<<2];
double q[N];int n,R[N<<2];
inline void fft(C *x,int f){
    for (int i=0;i<n;++i) if (i<R[i]) swap(x[i],x[R[i]]);
    for (int i=1;i<n;i<<=1){
        C wn={cos(pi/i),f*sin(pi/i)};
        for (int j=0;j<n;j+=(i<<1)){
            C w={1,0},t1,t2;
            for (int k=0;k<i;++k,w*=wn) {
                t1=x[j+k],t2=x[i+j+k]*w;
                x[j+k]=t1+t2,x[i+j+k]=t1-t2;
            }
        }
    }if(f==-1) for (int i=0;i<n;++i) x[i].a/=n;
}
int main(){ 
    freopen("bzoj3527.in","r",stdin);
    scanf("%d",&n);for (int i=0;i<n;++i) scanf("%lf",&q[i]),a[i]=(C){q[i],0};
    for (int i=1;i<n;++i) b[i]=(C){double(1)/i/i,0};
    for (int i=0;i<n;++i) a1[i]=(C){q[n-i-1],0};int m=n;
    for (n=1;n<=m;n<<=1);n<<=1;int l=log2(n);
    for (int i=0;i<n;++i)   R[i]=(R[i>>1]>>1)|(i&1)<<l-1;
    fft(a,1);fft(b,1);fft(a1,1);
    for (int i=0;i<n;++i) a[i]*=b[i];fft(a,-1);
    for (int i=0;i<n;++i) a1[i]*=b[i];fft(a1,-1);//a1-> from the back
    for (int i=0;i<m;++i) printf("%f\n",a[i].a-a1[m-i-1].a);
    return 0;
}