洛谷P4721 分治FFT

给定序列 (g_{1,dots,n-1}) 求序列 (f_{0,dots,n-1})，其中

[f_i=sum_{j=1}^{i}f_{i-j}g_j\ f_0=1 ]

现在要求 (f_{[l,r]})，首先假设我们已经计算出 (f_{[l,mid]})，这一部分对 (f_{[mid+1,r]}) 的贡献为 (f_{[l,mid]}*g_{[1,r-l]})，只要先算出来再加到对应位置上就行了。时间复杂度 (O(nlog^2n))。

Code

#include <bits/stdc++.h>
using namespace std;

#define RG register int
#define LL long long

template<typename elemType>
inline void Read(elemType& T) {
    elemType X = 0, w = 0; char ch = 0;
    while (!isdigit(ch)) { w |= ch == '-';ch = getchar(); }
    while (isdigit(ch)) X = (X << 3) + (X << 1) + (ch ^ 48), ch = getchar();
    T = (w ? -X : X);
}

const int maxn = 2100000;

LL qpow(LL b, LL n, LL MOD) {
    if (MOD == 1) return 0;
    LL x = 1, Power = b % MOD;
    while (n) {
        if (n & 1) x = x * Power % MOD;
        Power = Power * Power % MOD;
        n >>= 1;
    }
    return x;
}

namespace Poly {
    int r[maxn];
    int L, limit;
    const LL P = 998244353, G = 3, Gi = 332748118;

    LL pinv(LL x) { return qpow(x, P - 2, P); }

    //快速数论变换 type=1:正变换 type=-1:逆变换
    void NTT(LL* A, int type) {
        for (int i = 0; i < limit; i++)
            if (i < r[i]) swap(A[i], A[r[i]]);
        for (int mid = 1; mid < limit; mid <<= 1) {
            LL Wn = qpow(type == 1 ? G : Gi, (P - 1) / (mid << 1), P);
            for (int j = 0; j < limit; j += (mid << 1)) {
                LL w = 1;
                for (int k = 0; k < mid; k++, w = (w * Wn) % P) {
                    int x = A[j + k], y = w * A[j + k + mid] % P;
                    A[j + k] = (x + y) % P;
                    A[j + k + mid] = (x - y + P) % P;
                }
            }
        }
        if (type == 1) return;
        LL inv_limit = pinv(limit);
        for (int i = 0; i < limit; ++i)
            A[i] = A[i] * inv_limit % P;
    }

    //多项式卷积 a(x): N-1次多项式 b(x): M-1次多项式
    void Conv(LL* a, int N, LL* b, LL M, LL* c) {
        L = 0; limit = 1;
        while (limit <= N + M) limit <<= 1, L++;
        fill(a + N, a + limit, 0);
        fill(b + M, b + limit, 0);
        fill(c, c + limit, 0);
        for (int i = 0; i < limit; i++) r[i] = (r[i >> 1] >> 1) | ((i & 1) << (L - 1));
        NTT(a, 1); NTT(b, 1);
        for (int i = 0; i < limit; i++) c[i] = a[i] * b[i] % P;
        NTT(c, -1);
    }

    LL bufA[maxn], bufB[maxn], bufC[maxn];

    void CDQ_NTT(LL* f, LL* g, int L, int R) {
        if (L == R) return;
        int mid = (L + R) >> 1;
        CDQ_NTT(f, g, L, mid);
        int p = 0; for (int i = L;i <= mid;++i) bufA[p++] = f[i];
        int q = 0; for (int i = 1;i <= R - L;++i) bufB[q++] = g[i];
        Conv(bufA, p, bufB, q, bufC);
        for (int i = 0;i + mid + 1 <= R;++i)
            f[i + mid + 1] = (f[i + mid + 1] + bufC[i + mid - L]) % P;
        CDQ_NTT(f, g, mid + 1, R);
    }
}

LL f[maxn], g[maxn];
int n;

int main() {
    Read(n);
    f[0] = 1;
    for (int i = 1;i <= n - 1;++i) Read(g[i]);
    Poly::CDQ_NTT(f, g, 0, n - 1);
    for (int i = 0;i < n;++i) {
        printf("%lld", f[i]);
        if (i + 1 < n) printf(" ");
    }
    printf("
");

    return 0;
}