[Luogu P3338] [BZOJ 3527] [ZJOI2014]力

洛谷传送门

BZOJ传送门

###题目描述

给出$n$个数$q_i$，给出$F_j$的定义如下：

$F_j={\sum }_{i<j}\frac{q_i{q}_j}{(i-j{)}^2}-{\sum }_{i>j}\frac{q_i{q}_j}{(i-j{)}^2}$

令$E_i=F_i/q_i$，求$E_i$.

###输入输出格式

####输入格式：

第一行一个整数$n$。

接下来$n$行每行输入一个数，第i行表示qi。

####输出格式：

$n$行，第$i$行输出$E_i$。

与标准答案误差不超过$1e-2$即可。

###输入输出样例

####输入样例#1：

5
4006373.885184
15375036.435759
1717456.469144
8514941.004912
1410681.345880

####输出样例#1：

-16838672.693
3439.793
7509018.566
4595686.886
10903040.872

###说明

对于30%的数据，$n\le 1000$。

对于50%的数据，$n\le 60000$。

对于100%的数据，$n\le 100000$，$0<qi<1000000000$。

解题分析

设$in{v}_i=\frac{1}{i^2}$, 则所求变成了：
$$E_j=\sum _{i=1}^{j-1}in{v}_{j-i}\times q_i+\sum _{i=j+1}^{n}in{v}_{i-j}\times q_i$$
前面可以卷积卷一卷， 后面看起来不是很爽， 考虑将$q_i$中元素倒过来得到$q_i^{\prime }$， 那么可以这样写：
$$E_j=\sum _{i=1}^{j-1}in{v}_{j-i}\times q_i+\sum _{i=1}^{n-j}in{v}_i\times q_{n-j-i+1}^{\prime }$$

这下可以用FFT两次卷积搞出来了…

代码如下：

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cctype>
#include <cmath>
#include <algorithm>
#define R register
#define IN inline
#define gc getchar()
#define MX 400050
#define W while
#define db double
const db PI = std::acos((db)-1);
struct Complex {db re, im;} a[MX], b[MX];
IN Complex operator + (const Complex &x, const Complex &y) {return {x.re + y.re, x.im + y.im};}
IN Complex operator - (const Complex &x, const Complex &y) {return {x.re - y.re, x.im - y.im};}
IN Complex operator * (const Complex &x, const Complex &y) {return {x.re * y.re - x.im * y.im, x.re * y.im + x.im * y.re};}
db inv[MX], q[MX], ans1[MX], ans2[MX], re[MX];
int rev[MX], len, tot = 1, lg;
IN void FFT(Complex *dat, const int &typ)
{
    for (R int i = 0; i < tot; ++i) if(rev[i] < i) std::swap(dat[i], dat[rev[i]]);
    R int seg, bd, step, now, cur;
    Complex tar1, tar2, base, deal;
    for (seg = 1; seg < tot; seg <<= 1)
    {
        step = seg << 1; base = {std::cos(PI / seg), typ * std::sin(PI / seg)};
        for (now = 0; now < tot; now += step)
        {
            bd = now + seg; deal = {1, 0};
            for (cur = now; cur < bd; ++cur, deal = deal * base)
            {
                tar1 = dat[cur], tar2 = dat[cur + seg] * deal;
                dat[cur] = tar1 + tar2, dat[cur + seg] = tar1 - tar2;
            }
        }
    }
}
IN void FFT_init(db *dat1, db *dat2, db *ans)
{
    for (R int i = 0; i < tot; ++i)
    a[i] = {dat1[i], 0}, b[i] = {dat2[i], 0};
    FFT(a, 1), FFT(b, 1);
    for (R int i = 0; i < tot; ++i) a[i] = a[i] * b[i];
    FFT(a, -1);
    for (R int i = 1; i <= len; ++i) ans[i] = a[i].re / tot;
}
int main(void)
{
    scanf("%d", &len);
    for (R int i = 1; i <= len; ++i) scanf("%lf", &q[i]), inv[i] = 1.0 / i / i, re[i] = q[i];
    std::reverse(re + 1, re + 1 + len);
    W (tot <= len) lg++, tot <<= 1; tot <<= 1, ++lg;
    for (R int i = 0; i < tot; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << lg - 1);
    FFT_init(q, inv, ans1);
    FFT_init(re, inv, ans2);
    for (R int i = 1; i <= len; ++i) printf("%.3lf\n", ans1[i] - ans2[len + 1 - i]);
}