#include using namespace std; using ll = long long; static const int MOD = 998244353; ll modpow(ll a, ll e) { ll r = 1; while (e > 0) { if (e & 1) r = r * a % MOD; a = a * a % MOD; e >>= 1; } return r; } struct Comb { vector fact, ifact; void init(int n) { fact.assign(n + 1, 1); ifact.assign(n + 1, 1); for (int i = 1; i <= n; i++) fact[i] = (ll)fact[i - 1] * i % MOD; ifact[n] = (int)modpow(fact[n], MOD - 2); for (int i = n; i >= 1; i--) ifact[i - 1] = (ll)ifact[i] * i % MOD; } int C(int n, int r) const { if (r < 0 || r > n) return 0; return (ll)fact[n] * ifact[r] % MOD * ifact[n - r] % MOD; } }; vector mul_qs_minus_1_full(const vector &p, int s, int E) { vector res(E + 1, 0); // p(q) * q^s if (s <= E) { for (int i = 0; i + s <= E; i++) { res[i + s] = p[i]; } } // - p(q) for (int i = 0; i <= E; i++) { res[i] -= p[i]; if (res[i] < 0) res[i] += MOD; } return res; } void add_shift_scaled_full(vector &dst, const vector &src, int shift, int scale, int E) { if (scale == 0 || shift > E) return; for (int i = 0; i + shift <= E; i++) { if (src[i] == 0) continue; dst[i + shift] = (dst[i + shift] + (ll)src[i] * scale) % MOD; } } int main() { ios::sync_with_stdio(false); cin.tie(nullptr); int N, M, K; cin >> N >> M >> K; int E = N * (N - 1) / 2; int R = E - M; int S = N - 2; Comb comb; comb.init(E); // dp[r][s] is a polynomial in q of fixed length E+1. // active[r][s] indicates whether the state has appeared. vector>> dp(S + 1, vector>(S + 2, vector(E + 1, 0))); vector> active(S + 1, vector(S + 2, 0)); dp[S][1][0] = 1; active[S][1] = 1; // Build layers 1,2,...,K-1 without reaching vertex N. for (int step = 0; step < K - 1; step++) { vector>> ndp(S + 1, vector>(S + 2, vector(E + 1, 0))); vector> nactive(S + 1, vector(S + 2, 0)); for (int r = 0; r <= S; r++) { for (int s = 0; s <= S + 1; s++) { if (!active[r][s]) continue; vector cur = dp[r][s]; // cur * (q^s - 1)^t for (int t = 0; t <= r; t++) { int shift = t * (t - 1) / 2; int ways = comb.C(r, t); add_shift_scaled_full(ndp[r - t][t], cur, shift, ways, E); nactive[r - t][t] = 1; if (t != r) { cur = mul_qs_minus_1_full(cur, s, E); } } } } dp.swap(ndp); active.swap(nactive); } // H(q): generating polynomial for dist(1,N)=K in q=1+x. vector H(E + 1, 0); // N must have at least one edge to the current frontier: (q^s - 1). // All edges among the remaining r vertices, between them and N, // and between them and the current frontier are free. for (int r = 0; r <= S; r++) { for (int s = 0; s <= S + 1; s++) { if (!active[r][s]) continue; const vector &p = dp[r][s]; int base = s * r + r + r * (r - 1) / 2; for (int i = 0; i <= E; i++) { int val = p[i]; if (val == 0) continue; if (i + base + s <= E) { H[i + base + s] += val; if (H[i + base + s] >= MOD) H[i + base + s] -= MOD; } if (i + base <= E) { H[i + base] -= val; if (H[i + base] < 0) H[i + base] += MOD; } } } } // We need [x^R] H(1+x). If H(q)=sum_c h_c q^c, then // [x^R] H(1+x) = sum_c h_c * C(c,R). ll ans = 0; for (int c = R; c <= E; c++) { if (H[c] == 0) continue; ans += (ll)H[c] * comb.C(c, R) % MOD; ans %= MOD; } cout << ans << '\n'; return 0; }