#include <algorithm>
#include <iostream>
#include <numeric>
#include <string>
#include <tuple>
#include <utility>
#include <vector>
#include <limits>

#define rep(i, a, b) for (int i = int(a); i < int(b); i++)
using namespace std;
using ll = long long int;  // NOLINT
using P = pair<ll, ll>;

// clang-format off
#ifdef _DEBUG_
#define dump(...) do{ cerr << __LINE__ << ":\t" << #__VA_ARGS__ << " = "; debug_print(__VA_ARGS__); } while(false)
template<typename T, typename... Ts> void debug_print(const T &t, const Ts &...ts) { cerr << t; ((cerr << ", " << ts), ...); cerr << endl; }
#else
#define dump(...) do{ } while(false)
#endif
template<typename T> vector<T> make_v(size_t a, T b) { return vector<T>(a, b); }
template<typename... Ts> auto make_v(size_t a, Ts... ts) { return vector<decltype(make_v(ts...))>(a, make_v(ts...)); }
template<typename T> bool chmin(T &a, const T& b) { if (a > b) {a = b; return true; } return false; }
template<typename T> bool chmax(T &a, const T& b) { if (a < b) {a = b; return true; } return false; }
template<typename T, typename... Ts> void print(const T& t, const Ts&... ts) { cout << t; ((cout << ' ' << ts), ...); cout << '\n'; }
constexpr static struct PositiveInfinity { template<typename T> constexpr operator T() const { return numeric_limits<T>::max() / 2; } constexpr auto operator-() const; } inf;  // NOLINT
constexpr static struct NegativeInfinity { template<typename T> constexpr operator T() const { return numeric_limits<T>::lowest() / 2; } constexpr auto operator-() const; } NegativeInfinityVal;
constexpr auto PositiveInfinity::operator-() const { return NegativeInfinityVal; }
constexpr auto NegativeInfinity::operator-() const { return inf; }
// clang-format on

template<ll MOD>
class ModInt {
    ll n;
    auto constexpr inverse() const {
        return this->pow(*this, this->mod - 2);
    }

public:
    constexpr static ll mod = MOD;
    using mint = ModInt<MOD>;

    constexpr ModInt() : n(0) {}
    constexpr ModInt(const ll &nn) : n(((nn % MOD) + MOD) % MOD) {}
    constexpr mint operator+=(const mint &m) {
        n += m.n;
        if (n >= mint::mod) n -= mint::mod;
        return *this;
    }
    constexpr mint operator-=(const mint &m) {
        n -= m.n;
        if (n < 0) n += mint::mod;
        return *this;
    }
    constexpr mint operator*=(const mint &m) {
        n *= m.n;
        if (n >= mint::mod) n %= mint::mod;
        return *this;
    }
    constexpr mint operator/=(const mint &m) {
        return (*this) *= m.inverse();
    }
    friend constexpr mint operator+(mint t, const mint &m) {
        return t += m;
    }
    friend constexpr mint operator-(mint t, const mint &m) {
        return t -= m;
    }
    friend constexpr mint operator*(mint t, const mint &m) {
        return t *= m;
    }
    friend constexpr mint operator/(mint t, const mint &m) {
        return t /= m;
    }
    constexpr mint operator=(const ll &l) {
        n = l % mint::mod;
        if (n < 0) n += mint::mod;
        return *this;
    }
    friend ostream &operator<<(ostream &out, const mint &m) {
        out << m.n;
        return out;
    }
    friend istream &operator>>(istream &in, mint &m) {
        ll l;
        in >> l;
        m = l;
        return in;
    }
    static constexpr auto pow(const mint &x, ll p) {
        mint ans = 1;
        for (auto m = x; p > 0; p /= 2, m *= m) {
            if (p % 2) ans *= m;
        }
        return ans;
    }
    constexpr ll get_raw() const {
        return n;
    }
};
using mint = ModInt<998244353>::mint;
constexpr mint operator"" _m(unsigned long long m) {
    return mint(m);
}

template<typename T>
class LagrangeInterpolation {
    std::vector<T> x, y;

public:
    LagrangeInterpolation(const std::vector<T> &xp, const std::vector<T> &yp) : x(xp), y(yp) {}

    T interpolate(const T &t) const {
        const int N = static_cast<int>(x.size());
        std::vector<T> al(N, 1), ar(N, 1);
        for (int i = 0; i < N - 1; i++) {
            al[i + 1] = al[i] * (t - x[i]);
            ar[N - i - 2] = ar[N - i - 1] * (t - x[N - i - 1]);
        }
        T b = T{1} / std::accumulate(next(x.begin()), x.end(), T{1}, [&](T acc, const T &xi) { return acc * (x[0] - xi); });
        T ans = 0;
        for (uint i = 0; i < x.size(); i++) {
            ans += y[i] * al[i] * ar[i] * b;
            if (i + 1 < x.size()) {
                b /= x[i + 1] - x.front();
                b *= x[i] - x.back();
            }
        }
        return ans;
    }
};

mint factrial[] = {
    1,        295201906, 160030060, 957629942, 545208507, 213689172, 760025067, 939830261, 506268060, 39806322,  808258749, 440133909, 686156489, 741797144, 390377694, 12629586,  544711799, 104121967, 495867250, 421290700, 117153405, 57084755,  202713771, 675932866, 79781699,  956276337, 652678397, 35212756, 655645460, 468129309, 761699708, 533047427, 287671032, 206068022, 50865043,  144980423, 111276893, 259415897, 444094191, 593907889, 573994984, 892454686, 566073550, 128761001, 888483202, 251718753, 548033568, 428105027, 742756734, 546182474,
    62402409, 102052166, 826426395, 159186619, 926316039, 176055335, 51568171,  414163604, 604947226, 681666415, 511621808, 924112080, 265769800, 955559118, 763148293, 472709375, 19536133,  860830935, 290471030, 851685235, 242726978, 169855231, 612759169, 599797734, 961628039, 953297493, 62806842,  37844313, 909741023, 689361523, 887890124, 380694152, 669317759, 367270918, 806951470, 843736533, 377403437, 945260111, 786127243, 80918046,  875880304, 364983542, 623250998, 598764068, 804930040, 24257676,  214821357, 791011898, 954947696, 183092975,
};

int main() {
    cin.tie(nullptr);
    ios::sync_with_stdio(false);
    ll n, k;
    cin >> n >> k;

    auto powsum = [&](ll kk) -> mint {
        vector<mint> x(kk + 2), y(kk + 2, 1);
        iota(x.begin(), x.end(), 1);
        rep(i, 0, kk + 1) {
            y[i + 1] = mint::pow(x[i + 1], kk) + y[i];
        }

        LagrangeInterpolation<mint> li(x, y);
        return li.interpolate(n);
    };

    constexpr int B = 10'000'000;

    mint ans = n * powsum(k) - powsum(k + 1);
    n--;
    ans *= 2;
    ans *= factrial[n / B];
    rep(i, 0, n % B) {
        ans *= n / B * B + i + 1;
    }

    print(ans);

    // mint t = 1;
    // for (int i = 1; i <= n + 1; i++) {
    //     t *= i;
    //     if (i % B == 0) {
    //         print(t);
    //     }
    // }

    return 0;
}