代码

#include <bits/stdc++.h>
using namespace std;
using i64 = int64_t;
using u64 = uint64_t;


// constexpr
constexpr int MOD = 998244353, g = 3, g_inv = 332748118;
constexpr int MAXB = 21, MAXN = 1 << MAXB;

void addAssign(int &x, int y) { x += y; if (x >= MOD) x -= MOD; }
int modAdd(int x, int y) { x += y; return x < MOD ? x : x - MOD; }
void subAssign(int &x, int y) { x -= y; if (x < 0) x += MOD; }
int modSub(int x, int y) { x -= y; return x < 0 ? x + MOD : x; }
void mulAssign(int &x, int y) { x = i64(x) * y % MOD; }
int modMul(int x, int y) { return i64(x) * y % MOD; }
int neg(int x) { return x ? MOD - x : 0; }

int power(i64 a, i64 b) {
    a %= MOD; i64 res = 1;
    while (b) {
        if (b & 1) res = modMul(res, a);
        a = modMul(a, a);
        b >>= 1;
    }
    return res;
}

i64 inv(i64 x) { return power(x, MOD - 2); }
void divAssign(int &x, int y) { x = modMul(x, inv(y)); }
int modDiv(int x, int y) { return modMul(x, inv(y)); }

using Z = int;

vector<Z> bfg = []() {
    vector<Z> res(MAXN + 1); res[0] = 1;
    int t = MAXB; res[1 << t] = 31;
    for (int i = t; i; i--) res[1 << (i - 1)] = modMul(res[1 << i], res[1 << i]);
    for (int i = 1; i < MAXN; i++) res[i] = modMul(res[i & -i], res[i & (i - 1)]);
    return res;
}();

void ntt(vector<Z>& a, int n, int type = 1) {
    static array<u64, MAXN> ca{};
    for (int i = 0; i < n; i++) ca[i] = a[i];
    if (type == 1) {
        for (int m = n >> 1; m; m >>= 1) {
            for (int j = 0, *g = &bfg[0]; j < n; j += (m << 1), g++) {
                for (int k = j; k < j + m; k++) {
                    Z o = ca[k + m] * (*g) % MOD;
                    ca[k + m] = ca[k] + MOD - o, ca[k] += o;
                }
            }
            if (m == 1 << 10) for (int i = 0; i < n; i++) ca[i] %= MOD;
        }
        for (int i = 0; i < n; i++) a[i] = ca[i] % MOD;
    } else {
        for (int m = 1; m < n; m <<= 1) {
            for (int j = 0, *g = &bfg[0]; j < n; j += (m << 1), g++) {
                for (int k = j; k < j + m; k++) {
                    int o = ca[k + m] % MOD;
                    ca[k + m] = (ca[k] + MOD - o) * (*g) % MOD, ca[k] += o;
                }
            }
            if (m == 1 << 10) for (int i = 0; i < n; i++) ca[i] %= MOD;
        }
        Z invn = inv(n);
        for (int i = 0; i < n; i++) a[i] = modMul(ca[i] % MOD, invn);
        reverse(a.begin() + 1, a.begin() + n);
    }
}

vector<Z> convolution(vector<Z> a, vector<Z> b, int size = INT_MAX) {
    int n = a.size(), m = b.size();
    int k = 1; while (k < n + m - 1) k <<= 1;
    if (min(n, m) <= 66) {
        if (n < m) swap(n, m), swap(a, b);
        vector<Z> ans(n + m - 1);
        for (int i = 0; i < n; i++) for (int j = 0; j < m; j++)
            addAssign(ans[i + j], modMul(a[i], b[j]));
        return ans;
    }
    a.resize(k), b.resize(k);
    ntt(a, k), ntt(b, k);
    for (int i = 0; i < k; i++) mulAssign(a[i], b[i]);
    ntt(a, k, -1);
    a.resize(n + m - 1);
    return a;
}

void poly_multiply(unsigned *a, int n, unsigned *b, int m, unsigned *c) {
    vector<Z> va(a, a + n + 1), vb(b, b + m + 1);
    auto vc = convolution(va, vb);
    copy(vc.begin(), vc.end(), (Z*)c);
}

评测结果