const int SZ = 1 << 19, P = 998244353; int pw(int a, int m) { if (m < 0) m += P - 1; int res = 1; while (m) m & 1 ? res = 1ll * res * a % P : 0, a = 1ll * a * a % P, m >>= 1; return res; } int tr[SZ]; int init(int n) { // n指多项式的长度,非度数 int len = 1; while (len < n) len <<= 1; FOR(i, 0, len - 1) tr[i] = tr[i >> 1] >> 1 | ((i & 1) * (len / 2)); return len; } int dft_w[2][19][SZ]; void init_G() { FOR(j2, 0, 18) { int j = 1 << j2, wn = pw(3, (P - 1) / (j << 1)), *W = dft_w[0][j2]; W[0] = 1; FOR(i, 1, j - 1) W[i] = W[i - 1] * 1ll * wn % P; wn = pw(3, -(P - 1) / (j << 1)), W = dft_w[1][j2]; W[0] = 1; FOR(i, 1, j - 1) W[i] = W[i - 1] * 1ll * wn % P; } } void dft(int *f, int len, int tag) { FOR(i, 0, len - 1) if (i < tr[i]) swap(f[i], f[tr[i]]); FOR(i, 0, len - 1) f[i] = (f[i] + P) % P; for (int j = 1, j2 = 0; j < len; j <<= 1, ++j2) for (int i = 0, *W = dft_w[tag == -1][j2] /*pw(3,(P-1)/(j<<1)*tag)*/; i < len; i += j << 1) for (int k = i, u, v, w = 1; k < i + j; k++, w = W[k - i] /*1ll*w*wn%P*/) u = f[k], v = f[k + j] * 1ll * w % P, f[k] = u + v, f[k] -= f[k] >= P ? P : 0, f[k + j] = u - v, f[k + j] += f[k + j] < 0 ? P : 0; if (tag == -1) for (int i = 0, ilen = pw(len, P - 2); i < len; i++) f[i] = 1ll * f[i] * ilen % P; } int inv_h[SZ]; void inv(const int *f, int *g, int n) { // f*g=1 mod x^n if (n == 1) return g[0] = pw(f[0], P - 2), void(); inv(f, g, (n + 1) / 2); int len = init(2 * n); FOR(i, 0, len - 1) inv_h[i] = 0; FOR(i, (n + 1) / 2, len - 1) g[i] = 0; FOR(i, 0, n - 1) inv_h[i] = f[i]; dft(inv_h, len, 1), dft(g, len, 1); FOR(i, 0, len - 1) inv_h[i] = 1ll * inv_h[i] * g[i] % P * g[i] % P; FOR(i, 0, len - 1) g[i] = (g[i] * 2ll - inv_h[i] + P) % P; dft(g, len, -1); // FOR(i,n,len-1)g[i]=0; } int ln_h[SZ]; void ln(const int *f, int *g, int n) { // ln f=g mod x^n ,f[0]=1 inv(f, g, n); ln_h[n - 1] = 0; FOR(i, 1, n - 1) ln_h[i - 1] = f[i] * 1ll * i % P; int len = init(n * 2); FOR(i, n, len - 1) ln_h[i] = g[i] = 0; dft(ln_h, len, 1), dft(g, len, 1); FOR(i, 0, len - 1) g[i] = 1ll * g[i] * ln_h[i] % P; dft(g, len, -1); ROF(i, n - 1, 1) g[i] = g[i - 1] * 1ll * pw(i, P - 2) % P; g[0] = 0; // f[0]=1 } int exp_t[SZ], exp_h[SZ]; void exp(const int *f, int *g, int n) { // f[0]==0 if (n == 1) return g[0] = 1, void(); exp(f, g, (n + 1) / 2); FOR(i, (n + 1) / 2, n - 1) g[i] = 0; ln(g, exp_h, n); FOR(i, 0, n - 1) exp_t[i] = (f[i] - exp_h[i] + P) % P; exp_t[0] = (exp_t[0] + 1) % P; int len = init(n * 2); FOR(i, n, len) g[i] = exp_t[i] = 0; dft(g, len, 1), dft(exp_t, len, 1); FOR(i, 0, len - 1) g[i] = g[i] * 1ll * exp_t[i] % P; dft(g, len, -1); } int dm_f_r[SZ], dm_g_r[SZ], dm_h[SZ]; void div_mod(const int *f, int n, const int *g, int m, int *q, int &lq, int *r, int &lr) { // f = g*q + r if (n < m) return q[0] = 0, lq = 1, memcpy(r, g, sizeof(int) * m), lr = m, void(); memcpy(dm_f_r, f, sizeof(int) * n), memcpy(dm_g_r, g, sizeof(int) * m); reverse(dm_f_r, dm_f_r + n), reverse(dm_g_r, dm_g_r + m); int MOD = n - m + 1; inv(dm_g_r, dm_h, MOD); int len = init(MOD + MOD); FOR(i, MOD, len - 1) dm_f_r[i] = dm_h[i] = 0; dft(dm_f_r, len, 1), dft(dm_h, len, 1); FOR(i, 0, len - 1) dm_h[i] = dm_h[i] * 1ll * dm_f_r[i] % P; dft(dm_h, len, -1); FOR(i, 0, MOD - 1) q[i] = dm_h[i]; reverse(q, q + MOD); lq = MOD; memcpy(dm_f_r, g, sizeof(int) * m), memcpy(dm_g_r, q, sizeof(int) * lq); len = init(lq + m); FOR(i, m, len - 1) dm_f_r[i] = 0; FOR(i, lq, len - 1) dm_g_r[i] = 0; dft(dm_f_r, len, 1), dft(dm_g_r, len, 1); FOR(i, 0, len - 1) dm_g_r[i] = 1ll * dm_f_r[i] * dm_g_r[i] % P; dft(dm_g_r, len, -1); FOR(i, 0, n - 1) r[i] = (f[i] - dm_g_r[i] + P) % P; lr = m - 1; } const int PLSZ = 1 << 27, SZ2 = SZ << 2; int POOL[PLSZ], LP; int *newArray(int len) { int *p = POOL + LP; memset(p, 0, sizeof(int) * len); LP += len; // assert(LP<PLSZ); return p; } int *mtic_F[SZ2], mtic_L[SZ2]; int mtic_t1[SZ], mtic_t2[SZ]; void mtic_dfs1(int u, int l, int r, const int *a) { if (l == r) return mtic_F[u] = newArray(mtic_L[u] = 2), mtic_F[u][0] = P - a[l], mtic_F[u][1] = 1, void(); int mid = (l + r) >> 1; mtic_dfs1(u << 1, l, mid, a), mtic_dfs1(u << 1 | 1, mid + 1, r, a); memcpy(mtic_t1, mtic_F[u << 1], sizeof(int) * mtic_L[u << 1]); memcpy(mtic_t2, mtic_F[u << 1 | 1], sizeof(int) * mtic_L[u << 1 | 1]); int len = init(mtic_L[u << 1] + mtic_L[u << 1 | 1]); FOR(i, mtic_L[u << 1], len - 1) mtic_t1[i] = 0; FOR(i, mtic_L[u << 1 | 1], len - 1) mtic_t2[i] = 0; dft(mtic_t1, len, 1), dft(mtic_t2, len, 1); FOR(i, 0, len - 1) mtic_t2[i] = mtic_t2[i] * 1ll * mtic_t1[i] % P; dft(mtic_t2, len, -1); mtic_F[u] = newArray(mtic_L[u] = mtic_L[u << 1] + mtic_L[u << 1 | 1] - 1); memcpy(mtic_F[u], mtic_t2, sizeof(int) * mtic_L[u]); } void mtic_dfs2(int u, int l, int r, int *f, int n, const int *a, int *b) { if (l == r) { b[l] = 0; ROF(i, n - 1, 0) b[l] = (1ll * b[l] * a[l] + f[i]) % P; return; } int mid = (l + r) >> 1; int lp, lq, *g; div_mod(f, n, mtic_F[u << 1], mtic_L[u << 1], mtic_t1, lp, mtic_t2, lq); g = newArray(lq); FOR(i, 0, lq - 1) g[i] = mtic_t2[i]; mtic_dfs2(u << 1, l, mid, g, lq, a, b); div_mod(f, n, mtic_F[u << 1 | 1], mtic_L[u << 1 | 1], mtic_t1, lp, mtic_t2, lq); g = newArray(lq); FOR(i, 0, lq - 1) g[i] = mtic_t2[i]; mtic_dfs2(u << 1 | 1, mid + 1, r, g, lq, a, b); } /** * 多点求值 * * @param f 多项式 * @param n 多项式长度(度数 + 1) * @param a 点值序列,下标从 1 开始 * @param m 点值序列长度 * @param b 记录答案的序列 */ void multi_index_calc(int *f, int n, const int *a, int m, int *b) { LP = 0; mtic_dfs1(1, 1, m, a); mtic_dfs2(1, 1, m, f, n, a, b); } int *mtii_M[SZ2], mtii_LM[SZ2]; int mtii_t1[SZ], mtii_t2[SZ]; void mtii_dfs1(int u, int l, int r, const int *X) { if (l == r) return mtii_M[u] = newArray(mtii_LM[u] = 2), mtii_M[u][0] = P - X[l], mtii_M[u][1] = 1, void(); int mid = (l + r) >> 1; mtii_dfs1(u << 1, l, mid, X), mtii_dfs1(u << 1 | 1, mid + 1, r, X); memcpy(mtii_t1, mtii_M[u << 1], sizeof(int) * mtii_LM[u << 1]); memcpy(mtii_t2, mtii_M[u << 1 | 1], sizeof(int) * mtii_LM[u << 1 | 1]); int len = init(mtii_LM[u << 1] + mtii_LM[u << 1 | 1]); FOR(i, mtii_LM[u << 1], len - 1) mtii_t1[i] = 0; FOR(i, mtii_LM[u << 1 | 1], len - 1) mtii_t2[i] = 0; dft(mtii_t1, len, 1), dft(mtii_t2, len, 1); FOR(i, 0, len - 1) mtii_t2[i] = mtii_t2[i] * 1ll * mtii_t1[i] % P; dft(mtii_t2, len, -1); mtii_M[u] = newArray(mtii_LM[u] = mtii_LM[u << 1] + mtii_LM[u << 1 | 1] - 1); memcpy(mtii_M[u], mtii_t2, sizeof(int) * mtii_LM[u]); } int mtii_M1[SZ], mtii_L1, mtii_V[SZ]; int *mtii_F[SZ2], mtii_LF[SZ2]; void mtii_dfs2(int u, int l, int r, const int *Y) { if (l == r) return mtii_F[u] = newArray(mtii_LF[u] = 1), mtii_F[u][0] = Y[l] * 1ll * pw(mtii_V[l], P - 2) % P, void(); int mid = (l + r) >> 1; mtii_dfs2(u << 1, l, mid, Y), mtii_dfs2(u << 1 | 1, mid + 1, r, Y); memcpy(mtii_t1, mtii_M[u << 1], sizeof(int) * mtii_LM[u << 1]); memcpy(mtii_t2, mtii_F[u << 1 | 1], sizeof(int) * mtii_LF[u << 1 | 1]); int l1 = mtii_LM[u << 1] + mtii_LF[u << 1 | 1] - 1, l2 = mtii_LF[u << 1] + mtii_LM[u << 1 | 1] - 1; mtii_F[u] = newArray(mtii_LF[u] = max(l1, l2)); int len = init(mtii_LM[u << 1] + mtii_LF[u << 1 | 1]); FOR(i, mtii_LM[u << 1], len - 1) mtii_t1[i] = 0; FOR(i, mtii_LF[u << 1 | 1], len - 1) mtii_t2[i] = 0; dft(mtii_t1, len, 1), dft(mtii_t2, len, 1); FOR(i, 0, len - 1) mtii_t2[i] = mtii_t2[i] * 1ll * mtii_t1[i] % P; dft(mtii_t2, len, -1); FOR(i, 0, l1 - 1) mtii_F[u][i] = mtii_t2[i]; memcpy(mtii_t1, mtii_F[u << 1], sizeof(int) * mtii_LF[u << 1]); memcpy(mtii_t2, mtii_M[u << 1 | 1], sizeof(int) * mtii_LM[u << 1 | 1]); len = init(mtii_LF[u << 1] + mtii_LM[u << 1 | 1]); FOR(i, mtii_LF[u << 1], len - 1) mtii_t1[i] = 0; FOR(i, mtii_LM[u << 1 | 1], len - 1) mtii_t2[i] = 0; dft(mtii_t1, len, 1), dft(mtii_t2, len, 1); FOR(i, 0, len - 1) mtii_t2[i] = mtii_t2[i] * 1ll * mtii_t1[i] % P; dft(mtii_t2, len, -1); FOR(i, 0, l2 - 1) mtii_F[u][i] = (mtii_F[u][i] + mtii_t2[i]) % P; } void multi_index_interp(const int *X, const int *Y, int n, int *f) { mtii_dfs1(1, 1, n, X); memcpy(mtii_M1, mtii_M[1], sizeof(int) * mtii_LM[1]); assert(mtii_LM[1] == n + 1); FOR(i, 0, n - 1) mtii_M1[i] = mtii_M1[i + 1] * 1ll * (i + 1) % P; mtii_L1 = n; multi_index_calc(mtii_M1, mtii_L1, X, n, mtii_V); mtii_dfs2(1, 1, n, Y); FOR(i, 0, mtii_LF[1] - 1) f[i] = mtii_F[1][i]; } // make sure to init_G before dft!