Bài toán được giải bằng tư tưởng Bao hàm loại trừ dựa trên số cặp quân mã tấn công nhau như sau:

Gọi ~H(x)~ là số cấu hình mà không quan tâm đến thứ tự vị trí của ~3~ người có số cặp quân mã tấn công nhau lớn hơn hoặc bằng ~x~. Đáp án bài toán là:

~ans = 6!\times \sum_{x=0}^{\infty}(-1)^x\times H(x)~

Nhận thấy ~H(x)=0,\,x\geq 3~ nên công thức đáp án có thể rút gọn thành:

~ans = 6!\times \left[H(0) - H(1) + H(2)\right]~

Cách tính ~H(0)~

Dễ dàng nhận được bằng cách tính số cách đặt ~3~ quân mã ở vị trí bất kỳ vào ~n\times m~ ô bàn cờ theo công thức:

~H(0)=\binom{n\times m}{3}~

~\binom{n}{k}~ là tổ hợp chập ~k~ của ~n~. Xem thêm định nghĩa tại đây.

Cách tính ~H(1)~

Để tính được ~H(1)~, ta cần cố định ~1~ cặp quân mã tấn công nhau và quân mã còn lại sẽ được chọn vị trí tự do, Gọi ~G(1)~ là số cấu hình đặt ~1~ cặp quân mã tấn công nhau trên bảng, ta có công thức:

~H(1)=\binom{n\times m - 2}{1}\times G(1)=(n\times m-2) \times G(1)~

Để tính ~G(1)~ ta sẽ xét các trường hợp:

~(x;y)\leftrightarrow (x+1; y+2)~

~(x;y)\leftrightarrow (x+2; y+1)~

~(x;y)\leftrightarrow (x-1; y-2)~

~(x;y)\leftrightarrow (x-2; y-1)~

Với ~(x;y)~ là vị trí của quân mã trên bàn cờ.

Trường hợp thứ ~i~ mình sẽ xem nó là một hình chữ nhật với độ lớn ~a_i\times b_i~ sao cho hình chữ nhật là nhỏ nhất và bao gọn toàn bộ các quân mã đang xét trong từng trường hợp. Cụ thể:

~a_1\times b_1 = 2\times 3~

~a_2\times b_2 = 3\times 2~

~a_3\times b_3 = 2\times 3~

~a_4\times b_4 = 3\times 2~

Với mỗi trường hợp thì số cấu hình trong trường hợp đó là ~(n-a_i+1)\times (m-b_i+1)~ (Nhận được từ việc tịnh tiến hình chữ nhật đến mọi vị trí có thể của bàn cờ). Vậy ta có công thức của ~G(1)~:

~G(1)=\sum_{i}\left[(n-a_i+1)\times (m-b_i+1)\right]~

Suy ra công thức ~H(1)~:

~H(1)=(n\times m-2) \times \sum_{i}\left[(n-a_i+1)\times (m-b_i+1)\right]~

Cách tính ~H(2)~

Tương tự tư tưởng tính ~H(1)~ nhưng lần này có nhiều trường hợp hơn không thể xét bằng tay được nên mình sẽ sử dụng thuật trâu với độ phức tạp phù hợp để tìm ra các trường hợp. Tổng cộng có ~28~ trường hợp:

~\#(1\times 4) = 2~

~\#(2\times 2) = 8~

~\#(2\times 3) = 4~

~\#(2\times 4) = 2~

~\#(3\times 2) = 4~

~\#(3\times 3) = 4~

~\#(4\times 1) = 2~

~\#(4\times 2) = 2~

Với ~\#(a_i\times b_i)~ là số trường hợp có hình chữ nhật kích thước ~a_i\times b_i~.

Công thức của ~H(2)~:

~H(2) = \sum_i \left[(n-a_i+1)\times (m-b_i+1)\right]~

Code mẫu

/**
 *    Author :  Nguyen Tuan Vu
 *    Created : 12.11.2022
 **/

#pragma GCC optimize("O2")
#pragma GCC target("avx,avx2,fma")
#include <bits/stdc++.h>
#define MASK(x) ((1) << (x))
#define BIT(x, i) (((x) >> (i)) & (1))
#define ALL(v) (v).begin(), (v).end()
#define REP(i, n) for (int i = 0, _n = (n); i < _n; ++i)
#define FOR(i, a, b) for (int i = (a), _b = (b); i <= _b; ++i)
#define FORD(i, b, a) for (int i = (b), _a = (a); i >= _a; --i)
#define db(val) "[" #val " = " << (val) << "] "
#define TIME (1.0 * clock() / CLOCKS_PER_SEC)

template <class X, class Y>
bool minimize(X &a, Y b) {
  if (a > b) return a = b, true;
  return false;
}
template <class X, class Y>
bool maximize(X &a, Y b) {
  if (a < b) return a = b, true;
  return false;
}

using namespace std;

mt19937 rng(chrono::steady_clock::now().time_since_epoch().count());
int Rand(int l, int r) { return l + rng() % (r - l + 1); }

const int Dx[] = {-2, -2, -1, -1, 1, 1, 2, 2};
const int Dy[] = {-1, 1, -2, 2, -2, 2, -1, 1};
const int Mod = 1e9 + 7;
int n, m;

int main() {
  ios_base::sync_with_stdio(0);
  cin.tie(0);

  int ntest;
  cin >> ntest;
  multiset<pair<int, int>> mp[3];
  set<vector<pair<int, int>>> s2, s3;
  REP(i, 8) {
    int x = Dx[i];
    int y = Dy[i];
    vector<pair<int, int>> tmp;
    tmp.push_back({-min(0, x), -min(0, y)});
    tmp.push_back({x - min(0, x), y - min(0, y)});
    sort(ALL(tmp));
    s2.insert(tmp);

    REP(j, 8) {
      if (i != j) {
        tmp.clear();
        int x2 = Dx[j];
        int y2 = Dy[j];
        tmp.push_back({-min({0, x, x2}), -min({0, y, y2})});
        tmp.push_back({x - min({0, x, x2}), y - min({0, y, y2})});
        tmp.push_back({x2 - min({0, x, x2}), y2 - min({0, y, y2})});
        sort(tmp.begin(), tmp.end());
        s3.insert(tmp);
      }

      if (-Dx[i] != Dx[j] || -Dy[i] != Dy[j]) {
        tmp.clear();
        int x2 = x + Dx[j];
        int y2 = y + Dy[j];
        tmp.push_back({-min({0, x, x2}), -min({0, y, y2})});
        tmp.push_back({x - min({0, x, x2}), y - min({0, y, y2})});
        tmp.push_back({x2 - min({0, x, x2}), y2 - min({0, y, y2})});
        sort(tmp.begin(), tmp.end());
        s3.insert(tmp);
      }
    }
  }

  for (auto x : s2) {
    mp[1].insert({abs(x[1].first - x[0].first) + 1, abs(x[1].second - x[0].second) + 1});
  }

  for (auto x : s3) {
    mp[2].insert({max({x[0].first, x[1].first, x[2].first}) - min({x[0].first, x[1].first, x[2].first}) + 1,
                  max({x[0].second, x[1].second, x[2].second}) - min({x[0].second, x[1].second, x[2].second}) + 1});
  }

  while (ntest--) {
    cin >> n >> m;
    int num = 1ll * n * m % Mod;

    // H[0]
    int ans = 1ll * num * (num - 1 + Mod) % Mod * (num - 2 + Mod) % Mod;

    // H[1]
    for (auto x : mp[1]) {
      ans = (ans - 1ll * max(0, n - x.first + 1) * max(0, m - x.second + 1) % Mod * (1ll * n * m % Mod - 2 + Mod) % Mod * 6 % Mod + Mod) % Mod;
    }

    // H[2]
    for (auto x : mp[2]) {
      ans = (ans + 1ll * max(0, n - x.first + 1) * max(0, m - x.second + 1) % Mod * 6 % Mod) % Mod;
    }
    cout << ans << '\n';
  }
  return 0;
}