#line 1 "playspace/main.cpp" #include #line 8 "library/gandalfr/other/io_supporter.hpp" template std::ostream &operator<<(std::ostream &os, const std::vector &v) { for (int i = 0; i < (int)v.size(); i++) os << v[i] << (i + 1 != (int)v.size() ? " " : ""); return os; } template std::ostream &operator<<(std::ostream &os, const std::set &st) { for (const T &x : st) { std::cout << x << " "; } return os; } template std::ostream &operator<<(std::ostream &os, const std::multiset &st) { for (const T &x : st) { std::cout << x << " "; } return os; } template std::ostream &operator<<(std::ostream &os, const std::deque &dq) { for (const T &x : dq) { std::cout << x << " "; } return os; } template std::ostream &operator<<(std::ostream &os, const std::pair &p) { os << p.first << ' ' << p.second; return os; } template std::ostream &operator<<(std::ostream &os, std::queue &q) { int sz = q.size(); while (--sz) { os << q.front() << ' '; q.push(q.front()); q.pop(); } os << q.front(); q.push(q.front()); q.pop(); return os; } template std::istream &operator>>(std::istream &is, std::vector &v) { for (T &in : v) is >> in; return is; } template std::istream &operator>>(std::istream &is, std::pair &p) { is >> p.first >> p.second; return is; } #line 3 "library/gandalfr/math/matrix.hpp" #line 8 "library/gandalfr/math/matrix.hpp" template class matrix { private: int H, W; std::valarray> table; enum rowtrans_operation_name { SCALE, SWAP, ADD }; struct rowtrans_operation { int op, tar, res; T scl; }; using operations_history = std::vector; public: matrix() = default; matrix(int _H, int _W, T val = 0) : H(_H), W(_W), table(std::valarray(val, _W), _H) {} matrix(const std::vector> &vv) : H(vv.size()), W(vv[0].size()), table(std::valarray(W), H) { for (int i = 0; i < H; i++) for (int j = 0; j < W; j++) table[i][j] = vv[i][j]; } matrix(const std::valarray> &vv) : H(vv.size()), W(vv[0].size()), table(vv) {} /** * @brief 行列をリサイズする。 * @param val 拡張部分の値 */ void resize(int _H, int _W, T val = 0) { H = _H, W = _W; table.resize(_H, std::valarray(val, _H)); } int size_H() const { return H; } int size_W() const { return W; } void transpose() { matrix ret(W, H); for (int i = 0; i < H; i++) for (int j = 0; j < W; j++) ret.table[j][i] = table[i][j]; *this = std::move(ret); } void row_assign(int i, const std::valarray &row) { assert(W == (int)row.size()); table[i] = std::move(row); } void row_swap(int i, int j) { assert(0 <= i && i < H); assert(0 <= j && j < H); table[i].swap(table[j]); } /** * @attention O(n^3) * @attention 整数型では正しく計算できない。double や fraction を使うこと。 * @attention 枢軸選びをしていないので double では誤差が出るかも。 */ operations_history sweep_method() { operations_history hist; T ret = 1; for (int h = 0, w = 0; h < H && w < W; w++) { if (table[h][w] == 0) { for (int piv = h + 1; piv < H; piv++) { if (table[piv][w] != 0) { hist.push_back({SWAP, h, piv, 0}); row_swap(h, piv); break; } } if (table[h][w] == 0) { continue; } } T inv = 1 / table[h][w]; hist.push_back({SCALE, -1, w, inv}); table[h] *= inv; for (int j = h + 1; j < H; j++) { hist.push_back({ADD, h, j, -table[j][w]}); table[j] -= table[h] * table[j][w]; } h++; } return hist; } int rank() { auto U(*this); U.sweep_method(); int r = 0; for (int i = 0; i < H; ++i) { for (int j = i; j < W; ++j) { if (U.table[i][j] != 0) { r++; break; } } } return r; } T determinant() const { assert(H == W); matrix U(*this); T det = 1; auto hist = U.sweep_method(); if (U.table[H - 1][H - 1] == 0) return 0; for (auto &[op, tar, res, scl] : hist) { switch (op) { case SCALE: det /= scl; break; case SWAP: det *= -1; break; } } return det; } std::vector solve_system_of_equations(const std::vector &y) { assert(H == W); std::vector x(y); matrix U(*this); auto hist = U.sweep_method(); if (U.table[H - 1][H - 1] == 0) return {}; for (auto &[op, tar, res, scl] : hist) { switch (op) { case SCALE: x[res] *= scl; break; case SWAP: std::swap(x[tar], x[res]); break; case ADD: x[res] += x[tar] * scl; break; } } for (int i = H - 1; i >= 0; --i) { for (int j = 0; j < i; ++j) { x[j] -= U.table[j][i] * x[i]; } } return x; } matrix inverse() { assert(H == W); matrix INV(matrix::E(H)), U(*this); auto hist = U.sweep_method(); if (U.table[H - 1][H - 1] == 0) return matrix(0, 0); for (auto &[op, tar, res, scl] : hist) { switch (op) { case SCALE: INV.table[res] *= scl; break; case SWAP: std::swap(INV.table[tar], INV.table[res]); break; case ADD: INV.table[res] += INV.table[tar] * scl; break; } } for (int i = H - 1; i >= 0; --i) { for (int j = 0; j < i; ++j) { INV.table[j] -= INV.table[i] * U.table[j][i]; } } return INV; } void print() const { for (int i = 0; i < H; i++) { for (int j = 0; j < W; j++) { std::cout << table[i][j] << (j == W - 1 ? "" : " "); } std::cout << std::endl; } } matrix &operator+=(const matrix &a) { this->table += a.table; return *this; } matrix &operator-=(const matrix &a) { this->table -= a.table; return *this; } matrix &operator*=(const T &a) { this->table *= a; return *this; } matrix &operator*=(const matrix &a) { assert(W == a.H); matrix a_t(a), ret(H, a.W); a_t.transpose(); for (int i = 0; i < H; i++) { for (int j = 0; j < a_t.H; j++) { ret.table[i][j] = (table[i] * a_t.table[j]).sum(); } } *this = std::move(ret); return *this; } matrix &operator/=(const T &a) { this->table /= a; return *this; } /** * @brief 行列の冪乗。 * @param n 指数 * @attention n が 0 なら単位行列。 * @attention 演算子の優先度に注意。 */ matrix operator^=(long long n) { assert(H == W); if (n == 0) return *this = E(H); n--; matrix x(*this); while (n) { if (n & 1) *this *= x; x *= x; n >>= 1; } return *this; } matrix operator+() { return *this; } matrix operator-() { return matrix(*this) *= -1; } matrix operator+(const matrix &a) { return matrix(*this) += a; } matrix operator-(const matrix &a) { return matrix(*this) -= a; } template matrix operator*(const S &a) { return matrix(*this) *= a; } matrix operator/(const T &a) { return matrix(*this) /= a; } matrix operator^(long long n) { return matrix(*this) ^= n; } friend std::istream &operator>>(std::istream &is, matrix &mt) { for (auto &arr : mt.table) for (auto &x : arr) is >> x; return is; } T &operator()(int h, int w) { assert(0 <= h && h < H && 0 <= w && w <= W); return table[h][w]; } /** * @brief サイズ n の単位行列。 */ static matrix E(int N) { matrix ret(N, N); for (int i = 0; i < N; i++) ret.table[i][i] = 1; return ret; } }; #line 4 "playspace/main.cpp" using namespace std; using ll = long long; const int INF = 1001001001; const int MAXINT = std::numeric_limits::max(); const int MININT = std::numeric_limits::min(); const ll INFLL = 1001001001001001001; const ll MAXLL = std::numeric_limits::max(); const ll MINLL = std::numeric_limits::min(); const ll MOD = 1000000007; const ll _MOD = 998244353; #define rep(i, j, n) for(ll i = (ll)(j); i < (ll)(n); i++) #define rrep(i, j, n) for(ll i = (ll)(n-1); i >= (ll)(j); i--) #define all(a) (a).begin(),(a).end() #define debug(a) std::cerr << #a << ": " << a << std::endl #define LF cout << endl #ifdef ENABLE_MULTI_FOR #define mrep(it, mr) for(multi_iter it(mr); !it.fin(); ++it) #endif template inline bool chmax(T1 &a, T2 b) { return a < b && (a = b, true); } template inline bool chmin(T1 &a, T2 b) { return a > b && (a = b, true); } void Yes(bool ok){ std::cout << (ok ? "Yes" : "No") << std::endl; } int main(void) { matrix mt(2, 2); cin >> mt; (mt ^ 3).print(); }