すべての行に少なくとも 1 つの 1 を配置する必要があるため、行ごとに処理し、その行にいくつかの 1 を追加することを考えます。現在の状態を表すために、dp[i][j] を「最初の i 行を埋めた時点で、1 がすでに存在する列が j 列ある場合の埋め方の総数」と定義します。遷移では、現在の行で新たに 1 を追加する位置の数を k とし、残りの位置には任意の値を入れます。

しかし、この遷移は k = 0 の場合に破綻します。なぜなら、この式は行内に 1 が少なくとも 1 つあることを保証できないからです。そのため、k = 0 の場合は特別に扱います。この場合、j 個の位置の中に少なくとも 1 つ 1 がある場合の数を包除原理で計算し、dp[i][j] = (m ^ j \times (m - 1) ^ {n - j} - (m - 1) ^ n) dp[i - 1][j] とします。

#include <bits/stdc++.h>
using namespace std;
#define rep(i, l, r) for(int i = l; i <= r; ++i)
const int N = 250 + 5;
const int Mod = 1000000000 + 7;
int n, m, powM[N], powM1[N], C[N][N], dp[N][N];

int read() {
    char c; int x = 0, f = 1;
    c = getchar();
    while (c > '9' || c < '0') { if (c == '-') f = -1; c = getchar(); }
    while (c >= '0' && c <= '9') x = x * 10 + c - '0', c = getchar();
    return x * f;
}
int add(int a, int b) { return (a += b) >= Mod ? a - Mod : a; }
int sub(int a, int b) { return (a -= b) < 0 ? a + Mod : a; }
int mul(int a, int b) { return (long long)a * b % Mod; }
int powMod(int a, int b) {
    int res = 1;
    while (b) {
        if (b & 1) res = mul(res, a);
        a = mul(a, a);
        b >>= 1;
    }
    return res;
}
int main() {
    n = read(); m = read();
    powM[0] = powM1[0] = 1;
    rep(i, 1, n) {
        powM[i] = mul(powM[i - 1], m);
        powM1[i] = mul(powM1[i - 1], m - 1);
    }
    rep(i, 0, n) C[i][0] = 1;
    rep(i, 1, n) rep(j, 1, i) C[i][j] = add(C[i - 1][j], C[i - 1][j - 1]);
    dp[0][0] = 1;
    rep(i, 1, n) rep(j, 0, n) {
        rep(k, 0, j - 1) {
            dp[i][j] = add(dp[i][j], mul(dp[i - 1][k], mul(C[n - k][j - k], mul(powM[k], powM1[n - j]))));
        }
        dp[i][j] = add(dp[i][j], mul(dp[i - 1][j], mul(sub(powM[j], powM1[j]), powM1[n - j])));
    }
    printf("%d\n", dp[n][n]);
    return 0;
}

上記の DP は O(n^3) で動作し、この問題を解くには十分ですが、より効率的な方法はないでしょうか？制約の難しい点は、すべての行と列に少なくとも 1 つの 1 が存在しなければならないことです。これを別の視点から考えます。特定の行や列に 1 を配置しないことを「強制」するのはどうでしょうか？この方が計算が容易です。f[i][j] を「正確に i 行と j 列に 1 を配置しないことが強制され、残りのセルは任意に埋める場合の総数」と定義します。すると次の式が成り立ちます。

次に g[i][j] を「正確に i 行と j 列に 1 が存在しない場合の総数」とします。すると f と g の間には次の関係があります:

高次元の二項係数包除原理を用いると:

したがって、求めたい g[0][0] は次のようになります。

ここで式を整理します。n^2 - n(i + j) + i \times j = (n - i)(n - j) および n(i + j) - i \times j = (n - j)i + nj という因数分解を利用します。

共通の因数でまとめます:

n - j と j の指数の和が n であり、\binom{n}{j} が現れることから、これは二項定理の形になっています。したがって:

最後に powM[i] = m^i および powM1[i] = (m-1)^i を事前計算しておけば、O(n log n) で解を求めることができます。

#include <bits/stdc++.h>
using namespace std;
#define rep(i, l, r) for(int i = l; i <= r; ++i)
const int N = 250 + 5;
const int Mod = 1000000000 + 7;
int n, m, ans, powM[N], powM1[N], fac[N], inv[N];

int read() {
    char c; int x = 0, f = 1;
    c = getchar();
    while (c > '9' || c < '0') { if (c == '-') f = -1; c = getchar(); }
    while (c >= '0' && c <= '9') x = x * 10 + c - '0', c = getchar();
    return x * f;
}
int add(int a, int b) { return (a += b) >= Mod ? a - Mod : a; }
int sub(int a, int b) { return (a -= b) < 0 ? a + Mod : a; }
int mul(int a, int b) { return (long long)a * b % Mod; }
int powMod(int a, int b) {
    int res = 1;
    while (b) {
        if (b & 1) res = mul(res, a);
        a = mul(a, a);
        b >>= 1;
    }
    return res;
}
int comb(int n, int m) {
    if (m > n) return 0;
    return mul(fac[n], mul(inv[m], inv[n - m]));
}
int main() {
    n = read(); m = read();
    fac[0] = inv[0] = powM[0] = powM1[0] = 1;
    rep(i, 1, n) {
        fac[i] = mul(fac[i - 1], i);
        inv[i] = powMod(fac[i], Mod - 2);
        powM[i] = mul(powM[i - 1], m);
        powM1[i] = mul(powM1[i - 1], m - 1);
    }
    rep(i, 0, n) {
        int term = mul(comb(n, i), powMod(sub(mul(powM[n - i], powM1[i]), powM1[n]), n));
        if (i & 1) ans = sub(ans, term);
        else ans = add(ans, term);
    }
    printf("%d\n", ans);
    return 0;
}

この問題は [JSOI2015] 染色問題 と同様の技法で解くことができます。