#include <bits/stdc++.h>
using namespace std;
#define F first
#define S second
#define rep(i, a, b) for(int i = a; i < (b); ++i)
#define per(i, a, b) for(int i = b-1; i>=a ; i--)
#define trav(a, x) for(auto& a : x)
#define allin(a , x) for(auto a : x)
#define all(x) begin(x), end(x)
#define sz(x) (int)(x).size()
typedef long long ll;
#define int ll
#define ff first
#define ss second
#define pb push_back
#define eb emplace_back
#define mkp make_pair
typedef pair<int, int> pii;
typedef vector<ll> vl;
typedef vector<pii> vpi;
typedef pair<ll,ll> pll;
typedef vector<string> vs;
typedef vector<pll> vpl;
typedef vector<int> vi;
std::mt19937 rng((int) std::chrono::steady_clock::now().time_since_epoch().count());
ll cdiv(ll a, ll b) { return a/b+((a^b)>0&&a%b); } // divide a by b rounded up
ll fdiv(ll a, ll b) { return a/b-((a^b)<0&&a%b); } // divide a by b rounded down
#define Unique(v) sort(all(v));v.erase(unique(all(v)),v.end());
typedef complex<double> C;
typedef vector<double> vd;
const ll inf=LLONG_MAX;
const int maxn=1e6+10;
const int md=998244353;
void fft(vector<C>& a) {
    int n = sz(a), L = 31 - __builtin_clz(n);
    static vector<complex<long double>> R(2, 1);
    static vector<C> rt(2, 1);  // (^ 10% faster if double)
    for (static int k = 2; k < n; k *= 2) {
        R.resize(n); rt.resize(n);
        auto x = polar(1.0L, acos(-1.0L) / k);
        rep(i,k,2*k) rt[i] = R[i] = i&1 ? R[i/2] * x : R[i/2];
    }
    vi rev(n);
    rep(i,0,n) rev[i] = (rev[i / 2] | (i & 1) << L) / 2;
    rep(i,0,n) if (i < rev[i]) swap(a[i], a[rev[i]]);
    for (int k = 1; k < n; k *= 2)
        for (int i = 0; i < n; i += 2 * k) rep(j,0,k) {
            // C z = rt[j+k] * a[i+j+k]; // (25% faster if hand-rolled)  /// include-line
            auto x = (double *)&rt[j+k], y = (double *)&a[i+j+k];        /// exclude-line
            C z(x[0]*y[0] - x[1]*y[1], x[0]*y[1] + x[1]*y[0]);           /// exclude-line
            a[i + j + k] = a[i + j] - z;
            a[i + j] += z;
        }
}
vd conv(const vd& a, const vd& b) {
    if (a.empty() || b.empty()) return {};
    vd res(sz(a) + sz(b) - 1);
    int L = 32 - __builtin_clz(sz(res)), n = 1 << L;
    vector<C> in(n), out(n);
    copy(all(a), begin(in));
    rep(i,0,sz(b)) in[i].imag(b[i]);
    fft(in);
    for (C& x : in) x *= x;
    rep(i,0,n) out[i] = in[-i & (n - 1)] - conj(in[i]);
    fft(out);
    rep(i,0,sz(res)) res[i] = imag(out[i]) / (4 * n);
    return res;
}
typedef vector<ll> vl;
vl convMod(const vl &a, const vl &b) {
    if (a.empty() || b.empty()) return {};
    vl res(sz(a) + sz(b) - 1);
    int B=32-__builtin_clz(sz(res)), n=1<<B, cut=int(sqrt(md));
    vector<C> L(n), R(n), outs(n), outl(n);
    rep(i,0,sz(a)) L[i] = C((int)a[i] / cut, (int)a[i] % cut);
    rep(i,0,sz(b)) R[i] = C((int)b[i] / cut, (int)b[i] % cut);
    fft(L), fft(R);
    rep(i,0,n) {
        int j = -i & (n - 1);
        outl[j] = (L[i] + conj(L[j])) * R[i] / (2.0 * n);
        outs[j] = (L[i] - conj(L[j])) * R[i] / (2.0 * n) / 1i;
    }
    fft(outl), fft(outs);
    rep(i,0,sz(res)) {
        ll av = ll(real(outl[i])+.5), cv = ll(imag(outs[i])+.5);
        ll bv = ll(imag(outl[i])+.5) + ll(real(outs[i])+.5);
        res[i] = ((av % md * cut + bv) % md * cut + cv) % md;
    }
    return res;
}

int fac[maxn], inv_fac[maxn];
int modInv(int a, int m)
{
    int m0=m;
    int y=0, x=1;
    if(m==1)
        return 0;
    while(a>1)
    {
        int q=a/m;
        int t=m;
        m=a%m, a=t;
        t=y;
        y=x-q*y;
        x=t;
    }
    if(x<0)
        x+=m0;
    return x;
}
int bc(int n, int k)
{
    if (n<k || k<0) return 0;
    return fac[n]*inv_fac[k]%md*inv_fac[n-k]%md;
}
ll pw(ll b, ll e)
{
    ll ans=1;
    for(; e; b=b*b%md, e/=2)
        if(e&1)
            ans=ans*b%md;
    return (ans+md)%md;
}

int n, m;

vector<int> polys[1010];

int32_t main()
{
    #ifndef ONLINE_JUDGE
    freopen("input.txt","r",stdin);
    freopen("output.txt","w",stdout); 
    #endif
    ios::sync_with_stdio(0);cin.tie(0);

    fac[0]=1, inv_fac[0]=1;
    for (int i=1; i<maxn; i++) {
        fac[i]=fac[i-1]*i%md;
    }
    inv_fac[maxn-1]=modInv(fac[maxn-1], md);
    for (int i=maxn-2; i>=1; i--) {
        inv_fac[i]=inv_fac[i+1]*(i+1)%md;
    }

    cin>> n>> m;
    vector<int> a(m+1);
    for(int i=0; i<=m; i++)
    {
        a[i]=bc(m, i)*bc(m, i)%md*fac[i]%md;
        // cout<< i<<" " << a[i]<< "\n";
    }
    // for(int x: a)
    //     cout<< x<< " ";
    // cout<< "\n";
    //quero b=a^n
    for(int i = 1; i <= n; ++i)
        polys[i] = a;

    for(int sz = 2; sz <= 2*n; sz *= 2)
        for(int i = 1; i+sz/2 <= n; i += sz)
            polys[i] = convMod(polys[i], polys[i + sz / 2]);

    vector<int> b = polys[1];
    // for(int x: b)
    //     cout<< x<<" ";
    // cout<< "\n";
    int ans=fac[m*n];
    for(int i=1; i<=m*n; i++)
    {   
        ans+=(i%2?md-1:1)*b[i]%md*fac[m*n-i]%md;
        // cout<< i<< " "<< b[i]*fac[m*n-i]%md<< "\n";
        if(ans>=md)
            ans-=md;
    }
    cout<< ans<< "\n";
 
    return 0;
}