ZJOI2014力

题意

给出 nn 个数 q1,q2,…qnq_1,q_2, \dots q_n ，定义

Fj=∑i=1j−1qi×qj(i−j)2−∑i=j+1nqi×qj(i−j)2F_j=\sum_{i=1}^{j-1} \frac{q_i\times q_j}{(i-j)^2} - \sum_{i=j+1}^n \frac{q_i\times q_j}{(i-j)^2}

Ei=FiqiE_i=\frac{F_i}{q_i}

对于 1≤i≤n1\le i\le n ，求 EiE_i 的值。1≤n≤105,0<qi<1091\le n \le 10^5 ,0\lt q_i \lt 10^9 .

题解

Ej=∑i=1j−1qi(i−j)2−∑i=j+1nqi(i−j)2Ej=∑i=0jqi(i−j)2−∑i=jnqi(i−j)2E_j=\sum_{i=1}^{j-1}\frac{q_i}{(i-j)^2}-\sum_{i=j+1}^{n}\frac{q_i}{(i-j)^2}\\ E_j=\sum_{i=0}^{j}\frac{q_i}{(i-j)^2}-\sum_{i=j}^{n}\frac{q_i}{(i-j)^2}\\

设 gi=1i2g_i=\frac{1}{i^2}

Ej=∑i=0jqi⋅gj−i−∑i=jnqi⋅gi−jE_j=\sum_{i=0}^j q_i\cdot g_{j-i} -\sum_{i=j}^nq_i\cdot g_{i-j}\\

左边已经是卷积形式了，考虑化简右边：

∑i=jnqi⋅gi−j=∑i=0n−jqi+j⋅giqi′=qn−i=∑i=0n−jqn−i−j′⋅gi\begin{aligned} &\sum_{i=j}^n q_i\cdot g_{i-j}\\ =&\sum_{i=0}^{n-j} q_{i+j}\cdot g_{i}\\ &q'_i=q_{n-i}\\ =&\sum_{i=0}^{n-j} q'_{n-i-j}\cdot g_i \end{aligned}

也化成了卷积的形式，那么

Ej=∑i=0jqi⋅gj−i−∑i=0n−jqn−i−j′⋅giE_j=\sum_{i=0}^j q_i\cdot g_{j-i}-\sum_{i=0}^{n-j} q'_{n-i-j}\cdot g_i\\

A=q∗g,  B=q′∗gEj=Aj−Bn−jA=q\ast g,\;B=q'\ast g\\ E_j=A_j-B_{n-j}

上 FFT 即可.

code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
#include <bits/stdc++.h>
#define _max(x,y) ((x>y)?x:y)
#define _min(x,y) ((x<y)?x:y)
using namespace std;
typedef long long ll;
template<typename _Tp>
inline void red(_Tp &x) {
    x=0;bool fg=0;char ch=getchar();
    while(ch<'0'||ch>'9') {if(ch=='-') fg^=1;ch=getchar();}
    while(ch>='0'&&ch<='9') x=(x<<1)+(x<<3)+(_Tp)(ch^48),ch=getchar();
    if(fg) x=-x;
}
template<typename _Tp> bool umax(_Tp &x,_Tp y) {return (x<y)?(x=y,true):(false);}
template<typename _Tp> bool umin(_Tp &x,_Tp y) {return (x>y)?(x=y,true):(false);}
const int N = 262144;
typedef double db;
const db Pi = acos(-1);
struct com {
    db r,i;
    com():r(0),i(0) {}
    com(db _r,db _i):r(_r),i(_i){}
    com(int n):r(cos(Pi*2/n)),i(sin(Pi*2/n)){} // w_n^1
    com operator+(const com &b)const{return com(r+b.r,i+b.i);} 
    com operator-(const com &b)const{return com(r-b.r,i-b.i);} 
    com operator*(const com &b)const{return com(r*b.r-i*b.i,i*b.r+r*b.i);}
}A[N],B[N],g[N];
int rev[N],n,m;
void FFT(com f[],int n,int fg) {
    for(int i=0;i<n;++i) rev[i]=(rev[i>>1]>>1)|((i&1)?(n>>1):0);
    for(int i=0;i<n;++i) (i<rev[i])&&(swap(f[i],f[rev[i]]),1);
    for(int h=2;h<=n;h<<=1) {
        com Dt(h),w,tt; Dt.i*=fg;
        int len=(h>>1);
        for(int j=0;j<n;j+=h) {
            w=com(1,0);
            for(int k=j;k<j+len;++k) {
                tt=f[k+len]*w;
                f[k+len]=f[k]-tt;
                f[k]=f[k]+tt;
                w=w*Dt;
            } 
        }
    }
    if(fg==-1) for(int i=0;i<n;++i) f[i].r/=n;
}
int main() {
    red(n);
    for(int i=1;i<=n;++i){
        scanf("%lf",&A[i].r);
        B[n-i]=A[i];
    } 
    for(int i=1;i<=n;++i) g[i].r=(db)(1.0/i/i);
    m=ceil(log2(n*2)); m=1<<m;
    FFT(A,m,1);FFT(B,m,1);FFT(g,m,1);
    for(int i=0;i<m;++i) A[i]=A[i]*g[i], B[i]=B[i]*g[i];
    FFT(A,m,-1);FFT(B,m,-1);
    for(int i=1;i<=n;++i){
        printf("%.3lf\n",A[i].r-B[n-i].r);
    }
    return 0;
}