【BZOJ3527】[ZJOI2014]力
题面
题解
易得
[E_i=sum_{j<i}frac{q_j}{(i-j)^2}-sum_{j>i}frac{q_j}{(i-j)^2}
]
设(f_i=q_i)，(g_i=i^2)
[E_i=sum_{j<i}f_jg_{i-j}-sum_{j>i}f_jg_{i-j}
]
将(f)翻转得到(h)
[E_i=sum_{j<i}f_jg_{i-j}-sum_{j<i}h_jg_{i-j}
]
这™就是两个卷积啊。。。
就分别(FFT)求出两个卷积然后一减即可
代码
#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <complex>
using namespace std;
#define sqr(x) (1.0 * (x) * (x))
typedef complex<double> Complex;
const double PI = acos(-1.0);
const int MAX_N = 3e5 + 5;
int n, N, M, P, r[MAX_N];
double q[MAX_N], ans[MAX_N];
Complex a[MAX_N], b[MAX_N];
void FFT(Complex *p, int op) {
for (int i = 0; i < N; i++) if (i < r[i]) swap(p[i], p[r[i]]);
for (int i = 1; i < N; i <<= 1) {
Complex rot(cos(PI / i), op * sin(PI / i));
for (int j = 0; j < N; j += (i << 1)) {
Complex w(1, 0);
for (int k = 0; k < i; k++, w *= rot) {
Complex x = p[j + k], y = w * p[i + j + k];
p[j + k] = x + y, p[i + j + k] = x - y;
}
}
}
}
int main () {
scanf("%d", &n);
for (int i = 1; i <= n; i++) scanf("%lf", &q[i]);
N = M = n - 1;
for (int i = 0; i <= N; i++) a[i] = q[i + 1], b[i] = 1.0 / sqr(i + 1);
M += N;
for (N = 1; N <= M; N <<= 1, ++P) ;
for (int i = 0; i < N; i++) r[i] = (r[i >> 1] >> 1) | ((i & 1) << (P - 1));
FFT(a, 1); FFT(b, 1);
for (int i = 0; i <= N; i++) a[i] = a[i] * b[i];
FFT(a, -1);
for (int i = 2; i <= n; i++) ans[i] += (double)(a[i - 2].real() / N);
for (int i = 0; i <= N; i++) a[i].real() = b[i].real() = a[i].imag() = b[i].imag() = 0;
for (int i = 0; i < n; i++) a[i] = q[n - i], b[i] = 1.0 / sqr(i + 1);
FFT(a, 1); FFT(b, 1);
for (int i = 0; i <= N; i++) a[i] = a[i] * b[i];
FFT(a, -1);
for (int i = n - 1; i; i--) ans[i] -= (double)(a[n - i - 1].real() / N);
for (int i = 1; i <= n; i++) printf("%0.3lf
", ans[i]);
return 0;
}