#include using namespace std; // https://github.com/hos-lyric/libra/blob/master/algebra/modint.h template struct ModInt { static constexpr unsigned M = M_; unsigned x; constexpr ModInt() : x(0U) {} constexpr ModInt(unsigned x_) : x(x_ % M) {} constexpr ModInt(unsigned long long x_) : x(x_ % M) {} constexpr ModInt(int x_) : x(((x_ %= static_cast(M)) < 0) ? (x_ + static_cast(M)) : x_) {} constexpr ModInt(long long x_) : x(((x_ %= static_cast(M)) < 0) ? (x_ + static_cast(M)) : x_) {} ModInt &operator+=(const ModInt &a) { x = ((x += a.x) >= M) ? (x - M) : x; return *this; } ModInt &operator-=(const ModInt &a) { x = ((x -= a.x) >= M) ? (x + M) : x; return *this; } ModInt &operator*=(const ModInt &a) { x = (static_cast(x) * a.x) % M; return *this; } ModInt &operator/=(const ModInt &a) { return (*this *= a.inv()); } ModInt pow(long long e) const { if (e < 0) return inv().pow(-e); ModInt a = *this, b = 1U; for (; e; e >>= 1) { if (e & 1) b *= a; a *= a; } return b; } ModInt inv() const { unsigned a = M, b = x; int y = 0, z = 1; for (; b; ) { const unsigned q = a / b; const unsigned c = a - q * b; a = b; b = c; const int w = y - static_cast(q) * z; y = z; z = w; } assert(a == 1U); return ModInt(y); } ModInt operator+() const { return *this; } ModInt operator-() const { ModInt a; a.x = x ? (M - x) : 0U; return a; } ModInt operator+(const ModInt &a) const { return (ModInt(*this) += a); } ModInt operator-(const ModInt &a) const { return (ModInt(*this) -= a); } ModInt operator*(const ModInt &a) const { return (ModInt(*this) *= a); } ModInt operator/(const ModInt &a) const { return (ModInt(*this) /= a); } template friend ModInt operator+(T a, const ModInt &b) { return (ModInt(a) += b); } template friend ModInt operator-(T a, const ModInt &b) { return (ModInt(a) -= b); } template friend ModInt operator*(T a, const ModInt &b) { return (ModInt(a) *= b); } template friend ModInt operator/(T a, const ModInt &b) { return (ModInt(a) /= b); } explicit operator bool() const { return x; } bool operator==(const ModInt &a) const { return (x == a.x); } bool operator!=(const ModInt &a) const { return (x != a.x); } friend std::ostream &operator<<(std::ostream &os, const ModInt &a) { return os << a.x; } }; // https://yukicoder.me/submissions/737488 // Editorial 解 #include #include #include #include using namespace std; // Library Checker Characteristic Polynomial https://judge.yosupo.jp/problem/characteristic_polynomial // Upper Hessenberg reduction of square matrices // Complexity: O(n^3) // Reference: // http://www.phys.uri.edu/nigh/NumRec/bookfpdf/f11-5.pdf template void hessenberg_reduction(std::vector> &M) { assert(M.size() == M[0].size()); const int N = M.size(); for (int r = 0; r < N - 2; r++) { int piv = -1; for (int h = r + 1; h < N; ++h) { if (M[h][r] != 0) { piv = h; break; } } if (piv < 0) continue; for (int i = 0; i < N; i++) std::swap(M[r + 1][i], M[piv][i]); for (int i = 0; i < N; i++) std::swap(M[i][r + 1], M[i][piv]); const auto rinv = Tp(1) / M[r + 1][r]; for (int i = r + 2; i < N; i++) { const auto n = M[i][r] * rinv; for (int j = 0; j < N; j++) M[i][j] -= M[r + 1][j] * n; for (int j = 0; j < N; j++) M[j][r + 1] += M[j][i] * n; } } } // Characteristic polynomial of matrix M (|xI - M|) // Complexity: O(n^3) // R. Rehman, I. C. Ipsen, "La Budde's Method for Computing Characteristic Polynomials," 2011. template std::vector characteristic_poly(std::vector> M) { hessenberg_reduction(M); const int N = M.size(); // p[i + 1] = (Characteristic polynomial of i-th leading principal minor) std::vector> p(N + 1); p[0] = {1}; for (int i = 0; i < N; i++) { p[i + 1].assign(i + 2, 0); for (int j = 0; j < i + 1; j++) p[i + 1][j + 1] += p[i][j]; for (int j = 0; j < i + 1; j++) p[i + 1][j] -= p[i][j] * M[i][i]; Tp betas = 1; for (int j = i - 1; j >= 0; j--) { betas *= M[j + 1][j]; Tp hb = -M[j][i] * betas; for (int k = 0; k < j + 1; k++) p[i + 1][k] += hb * p[j][k]; } } return p[N]; } // Library Checker ここまで template std::vector det_of_first_degree_mat(std::vector> M0, std::vector> M1) { const int N = M0.size(); int multiply_by_x = 0; // 「特定の列に x をかける」操作を行った回数 T detAdetBinv = 1; // 解説中の 1 / (det A det B) の値 for (int p = 0; p < N; ++p) { // M1[p][p] に nonzero を持ってきて、M1 の第 p 列を掃き出す int pivot = -1; for (int row = p; row < N; ++row) { if (M1[row][p] != T()) { pivot = row; break; } } if (pivot < 0) { ++multiply_by_x; if (multiply_by_x > N) return std::vector(N + 1); // M1 の第 p 列で pivot が見つからなかった場合、M0 + x M1 の第 p 列に x をかけたい // かける前に M1 の第 p 列を第 1 ~ (p - 1) 列を使って掃き出して、 // x をかけた後で x の二次の項が出てこないようにする for (int row = 0; row < p; ++row) { T v = M1[row][p]; M1[row][p] = 0; for (int i = 0; i < N; ++i) M0[i][p] -= v * M0[i][row]; } for (int i = 0; i < N; ++i) swap(M0[i][p], M1[i][p]); --p; // 第 p 列をもう一度やり直す この処理は高々 N 回しか走らないので全体の計算量は O(n^3) が保たれる continue; } if (pivot != p) { M1[pivot].swap(M1[p]); M0[pivot].swap(M0[p]); detAdetBinv *= -1; } // p 行目を定数倍して M1[p][p] == 1 にする T v = M1[p][p], vinv = v.inv(); detAdetBinv *= v; for (int col = 0; col < N; ++col) { M0[p][col] *= vinv; M1[p][col] *= vinv; } // p 行目を使用して M1 の p 列目を p 行目以外ゼロにする for (int row = 0; row < N; ++row) { if (row == p) continue; T v = M1[row][p]; for (int col = 0; col < N; ++col) { M0[row][col] -= M0[p][col] * v; M1[row][col] -= M1[p][col] * v; } } } // この時点で M1 = I なので det(xI + M0) を求める for (auto &vec : M0) { for (auto &x : vec) x = -x; } auto poly = characteristic_poly(M0); for (auto &x : poly) x *= detAdetBinv; poly.erase(poly.begin(), poly.begin() + multiply_by_x); poly.resize(N + 1); return poly; } using mint = ModInt<998244353>; int main() { int n; cin >> n; vector> M0(n, vector(n)), M1(n, vector(n)); for (auto& v : M0) { for (auto &x: v) { int x1; cin >> x1; x = x1; } } for (auto& v : M1) { for (auto& x : v) { int x1; cin >> x1; x = x1; } } for (auto v : det_of_first_degree_mat(M0, M1)) { cout << v << '\n'; } }