说在前面
鸣谢
h
s
y
\mathcal{h{\color {red}sy}}
hsy神仙的word文档。
虽然您说您被zxh虐,但是您还是把我调着打/kel…
Solution
直接上式子了。。。
E
j
=
F
j
q
j
=
∑
i
=
0
j
−
1
q
i
(
i
−
j
)
2
−
∑
i
=
j
+
1
n
−
1
q
i
(
i
−
j
)
2
E_j=\frac{F_j}{q_j}=\sum_{i=0}^{j-1}\frac{q_i}{(i-j)^2}-\sum_{i=j+1}^{n-1}\frac{q_i}{(i-j)^2}
Ej=qjFj=i=0∑j−1(i−j)2qi−i=j+1∑n−1(i−j)2qi
设
f
[
i
]
=
s
g
n
(
i
)
i
2
(
i
!
=
0
)
,
f
[
i
]
=
0
(
i
=
0
)
f[i]=\frac{sgn(i)}{i^2}(i!= 0),f[i]=0(i=0)
f[i]=i2sgn(i)(i!=0),f[i]=0(i=0)
所以
E
j
=
∑
i
=
0
n
−
1
q
[
i
]
f
[
j
−
i
]
E_j=\sum_{i=0}^{n-1}q[i]f[j-i]
Ej=i=0∑n−1q[i]f[j−i]
这就是一个卷积了,所以赶紧
F
a
F
a
T
a
\mathcal{FaFaTa}
FaFaTa
Warning
FFT的空间在一般情况下要开4倍,但是本体因为 j − i j-i j−i可能小于 0 0 0所以右开大了一倍
Code
#include <cstdio>
#include <algorithm>
#include <cmath>
#define N 100010
typedef long double LD;
const LD pi = std :: acos(-1);
struct cp{
LD a, b;
cp(LD a_ = 0.0, LD b_ = 0.0) {
a = a_; b = b_;
}
inline cp operator + (const cp &nxt) const {return cp(a + nxt.a, b + nxt.b); }
inline cp operator - (const cp &nxt) const {return cp(a - nxt.a, b - nxt.b); }
inline cp operator * (const cp &nxt) const {return cp(a * nxt.a - b * nxt.b, a * nxt.b + b * nxt.a); }
}q[N << 3], f[N << 3];
int re[N << 3], len, l;
inline void fft(cp *a, int fg) {
for (int i = 0; i < len; ++i) {
if (i < re[i]) std :: swap(a[i], a[re[i]]);
}
for (int l = 1; l < len; l <<= 1) {
cp bs = cp(std :: cos(pi / l), std :: sin(pi / l) * fg);
for (int j = 0; j < len; j += (l << 1)) {
cp now = cp(1.0, 0.0);
for (int k = 0; k < l; ++k) {
cp nx = a[j + k], ny = a[j + k + l] * now;
a[j + k] = nx + ny;
a[j + k + l] = nx - ny;
now = now * bs;
}
}
}
}
int main() {
int n;
scanf("%d", &n);
for (int i = 0; i < n; ++i) {
scanf("%Lf", &q[i].a);
}
for (int i = - n + 1; i < n; ++i) {
if (!i) f[n - 1].a = 0;
else f[i + n - 1].a = 1.0 * abs(i) / i / i / i;
}
for (l = 0, len = 1; len <= 3 * n - 2; len <<= 1, ++l);
for (int i = 0; i < len; ++i)
re[i] = (re[i >> 1] >> 1 | ((i & 1) << (l - 1)));
fft(q, 1);
fft(f, 1);
for (int i = 0; i < len; ++i) q[i] = q[i] * f[i];
fft(q, -1);
for (int i = 0; i < n; ++i) {
printf("%.3Lf\n", q[i + n - 1].a / len);
}
return 0;
}