P3338 [ZJOI2014] FFT

题意

传送门 P3338 [ZJOI2014]力

题解

对于左右两项，分别令 $k=j-i,k=i-j$
$$E_j=\sum _{k=1}^{j-1}\frac{q_{j-k}}{k^2}-\sum _{k=1}^{n-j}\frac{q_{j+k}}{k^2}$$ 令 $g_k=1/k^2,f_k=q_k$，上式转换为
$$\sum g_k{f}_i[k+i=j]+\sum g_k{f}_{k+j}$$ 左右两项都化为卷积的形式，FFT 求解即可。总时间复杂度$O(n\log n)$。

#include <bits/stdc++.h>
using namespace std;
using cp = complex<double>;
constexpr double PI = acos(-1);
vector<int> rev;
struct Poly : vector<cp>
{
    Poly() {}
    Poly(int n) : vector<cp>(n) {}
    Poly(const initializer_list<cp> &list) : vector<cp>(list) {}
    void fft(int n, bool inverse)
    {
        if ((int)rev.size() != n)
        {
            rev.resize(n);
            for (int i = 0; i < n; ++i)
                rev[i] = rev[i >> 1] >> 1 | (i & 1 ? n >> 1 : 0);
        }
        resize(n);
        for (int i = 0; i < n; ++i)
            if (i < rev[i])
                std::swap(at(i), at(rev[i]));

        for (int m = 1; m < n; m <<= 1)
        {
            int m2 = m << 1;
            cp _w = cp(cos(2 * PI / m2), sin(2 * PI / m2));
            if (inverse)
                _w = conj(_w);
            for (int i = 0; i < n; i += m2)
            {
                cp w = cp(1, 0);
                for (int j = 0; j < m; ++j, w *= _w)
                {
                    cp &x = at(i + j), &y = at(i + j + m);
                    cp t = w * y;
                    y = x - t;
                    x += t;
                }
            }
        }
    }
    void dft(int n) { fft(n, 0); };
    void idft(int n)
    {
        fft(n, 1);
        for (int i = 0; i < n; ++i)
            at(i) /= n;
    }
    Poly operator*(const Poly &p) const
    {
        auto a = *this, b = p;
        int k = 1, n = a.size() + b.size() - 1;
        while (k < n)
            k <<= 1;
        a.dft(k), b.dft(k);
        for (int i = 0; i < k; ++i)
            a[i] *= b[i];
        a.idft(k);
        a.resize(n);
        return a;
    }
};
constexpr int MAXN = 1E5 + 5;
int N;
double F[MAXN], G[MAXN], E[MAXN];

int main()
{
    ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
    cin >> N;
    for (int i = 1; i <= N; ++i)
        cin >> F[i];
    for (int i = 1; i <= N; ++i)
        G[i] = 1 / ((double)i * i);
    Poly g(N + 1), f(N + 1);
    for (int i = 0; i <= N; ++i)
        g[i].real(G[i]), f[i].real(F[i]);
    f = f * g;
    for (int i = 1; i <= N; ++i)
        E[i] += f[i].real();
    for (int i = 0; i <= N; ++i)
        g[i].real(G[i]), f[i].real(F[N - i]);
    f = f * g;
    for (int i = 1; i <= N; ++i)
        E[i] -= f[N - i].real();
    for (int i = 1; i <= N; ++i)
        cout << fixed << setprecision(3) << E[i] << '\n';
    return 0;
}