#include using namespace std; template struct modint { modint() : x(0) {} modint(long long v) : x(v % MOD) { if (x < 0) x += MOD; } long long x; long long val() const { return x; } static constexpr long long mod() noexcept { return MOD; } friend modint operator+(modint a, modint b) { return a.x + b.x; } friend modint operator-(modint a, modint b) { return a.x - b.x; } friend modint operator*(modint a, modint b) { return a.x * b.x; } friend modint operator/(modint a, modint b) { return a.x * b.inv(); } friend modint operator+=(modint &a, modint b) { return a = a + b; } friend modint operator-=(modint &a, modint b) { return a = a - b; } friend modint operator*=(modint &a, modint b) { return a = a * b; } friend modint operator/=(modint &a, modint b) { return a = a / b; } friend bool operator==(modint a, modint b) {return a.x == b.x; } friend bool operator!=(modint a, modint b) {return a.x != b.x; } modint operator+() const { return *this; } modint operator-() const { return modint() - *this; } // ここまでがテンプレ // これ以降は必要に応じて modint pow(long long k) const { modint a = x, res = 1; while (k > 0) { if (k & 1) res *= a; a *= a; k >>= 1; } return res; } modint inv() const { long long a = x, b = MOD; long long u = 1, v = 0; while (b > 0) { long long t = a / b; a -= t * b; swap(a, b); u -= t * v; swap(u, v); } return u; } // ++a とか a-- とかを使うなら書く modint& operator++() { x++; if (x == MOD ) x = 0; return *this; } modint& operator--() { if (x == 0) x = MOD; x--; return *this; } modint operator++(int) { modint res = *this; ++*this; return res; } modint operator--(int) { modint res = *this; --*this; return res; } }; using mint = modint<998244353>; #include vector convolution(vector f, vector g) { int n = f.size(), m = g.size(); vector h(n + m - 1); for (int i = 0; i < n; i++) { for (int j = 0; j < m; j++) { h[i + j] += f[i] * g[j]; } } return h; } vector div(vector f, vector g) { if (f.size() < g.size()) return {}; int n = f.size() - g.size(), m = g.size() - 1; vector h(n + 1); for (int i = n; i >= 0; i--) { h[i] = f[i + m] / g[m]; for (int j = m; j >= 0; j--) { f[i + j] -= h[i] * g[j]; } } return h; } vector Berlekamp_Massey(const vector &a) { int n = a.size(); vector b = {1}, c = {1}; mint y = 1; for (int d = 1; d <= n; d++) { int k = b.size(), l = c.size(); mint x = 0; for (int i = 0; i < l; i++) { x += c[i] * a[d - l + i]; } b.push_back(0); k++; if (x == 0) continue; mint buf = x / y; if (l < k) { vector tmp = c; c.insert(c.begin(), k - l, 0); for (int i = 0; i < k; i++) { c[k - i - 1] -= buf * b[k - i -1]; } b = tmp; y = x; } else { for (int i = 0; i < k; i++) { c[l - i - 1] -= buf * b[k - i - 1]; } } } reverse(c.begin(), c.end()); for (mint &x : c) x = -x; return c; } mint Bostan_Mori(vector p, vector q, long long k) { mint res = 0; if (p.size() >= q.size()) { vector r = div(p, q); auto qr = convolution(q, r); for (int i = 0; i < qr.size(); i++) p[i] -= qr[i]; while (p.size() && p.back() == 0) p.pop_back(); if (k < r.size()) res += r[k]; } if (p.empty()) return res; p.resize( q.size() - 1 ); auto sub = [&](const vector &f, bool odd = 0) -> vector { int n = f.size(); if (!odd) n++; vector g(n / 2); for (int i = odd; i < f.size(); i += 2) g[i / 2] = f[i]; return g; }; while (k) { auto q2 = q; for (int i = 1; i < q2.size(); i += 2) q2[i] = -q2[i]; p = sub(convolution(p, q2), k & 1); q = sub(convolution(q, q2)); k /= 2; } return res + p[0]; } mint linear_recurrence(vector a, vector c, long long k) { vector c2(c.size() + 1); for (int i = 0; i < c.size(); i++) c2[i + 1] = -c[i]; c2[0] = 1; auto c3 = convolution(a, c2); c3.resize(a.size()); return Bostan_Mori(c3, c2, k); } mint BMBM(const vector x, long long k) { auto tmp = Berlekamp_Massey(x); int n = tmp.size() - 1; vector a(n), c(n); for (int i = 0; i < n; i++) { a[i] = x[i]; c[i] = tmp[i + 1]; } return linear_recurrence(a, c, k); } #include int main() { ios::sync_with_stdio(false); cin.tie(nullptr); long long N; int K; cin >> N >> K; int M = 4 * K; vector F(M), A(M); F[0] = F[1] = 1; A[0] = 1, A[1] = 2; for (int i = 2; i < M; i++) { F[i] = F[i - 1] + F[i - 2]; A[i] = A[i - 1] + F[i].pow(K); } cout << BMBM(A, N - 1).val() << endl; }