【HEOI / TJOI 2016】求和（CDQ 分治 + FFT）

题目链接：【HEOI / TJOI 2016】求和

转载自：algocode 算法博客

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn = 1 << 18 | 5;
const ll P = 998244353, g = 3;
int n, r[maxn];
ll inv[maxn], fact[maxn], finv[maxn];
ll f[maxn], a[maxn], b[maxn];
void prework(int n) {
    inv[1] = 1;
    for (int i = 2; i <= n; i++) {
        inv[i] = -(P / i) * inv[P % i] % P;
    }
    fact[0] = finv[0] = 1;
    for (int i = 1; i <= n; i++) {
        fact[i] = fact[i - 1] * i % P;
        finv[i] = finv[i - 1] * inv[i] % P;
    }
}
ll mpow(ll a, ll b) {
    ll res = 1;
    b = (b % (P - 1) + P - 1) % (P - 1);
    for (; b; b >>= 1, a = a * a % P) {
        if (b & 1) {
            res = res * a % P;
        }
    }
    return res;
}
void ntt(ll *a, int n, int type) {
    for (int i = 0; i < n; i++) {
        if (i < r[i]) {
            swap(a[i], a[r[i]]);
        }
    }
    for (int k = 1; k < n; k <<= 1) {
        ll w = mpow(g, type * (P - 1) / (k << 1));
        for (int i = 0; i < n; i += k << 1) {
            ll v = 1;
            for (int j = i; j < i + k; j++, v = v * w % P) {
                ll x = a[j], y = v * a[j + k] % P;
                a[j] = (x + y) % P, a[j + k] = (x - y) % P;
            }
        }
    }
    if (type == -1) {
        ll x = mpow(n, P - 2);
        for (int i = 0; i < n; i++) {
            a[i] = a[i] * x % P;
        }
    }
}
void cdq(int lb, int rb) {
    if (lb == rb) {
        return;
    }
    int mid = (lb + rb) / 2;
    cdq(lb, mid);
    int len = rb - lb + 1, n, k = 0;
    for (n = 1; n < len; n <<= 1) k++;
    for (int i = 0; i < n; i++) {
        r[i] = (r[i >> 1] >> 1) | ((i & 1) << (k - 1));
    }
    for (int i = 0; i < n; i++) {
        a[i] = b[i] = 0;
    }
    for (int i = lb; i <= mid; i++) {
        a[i - lb] = f[i];
    }
    for (int i = 0; i <= rb - lb; i++) {
        b[i] = finv[i];
    }
    ntt(a, n, 1), ntt(b, n, 1);
    for (int i = 0; i < n; i++) {
        a[i] = a[i] * b[i] % P;
    }
    ntt(a, n, -1);
    for (int i = mid + 1; i <= rb; i++) {
        f[i] = (f[i] + 2 * a[i - lb]) % P;
    }
    cdq(mid + 1, rb);
}
int main () {
    scanf("%d", &n);
    prework(n);
    f[0] = 1;
    cdq(0, n);
    ll ans = 0;
    for (int i = 0; i <= n; i++) {
        ans = (ans + f[i] * fact[i]) % P;
    }
    printf("%lld\n", (ans % P + P) % P);
    return 0;
}