#include #include using namespace std; #if __cplusplus > 201703L #include using namespace numbers; #endif // T must support +=, -=, *, *=, ==, and != template struct matrix_fixed{ using ring_t = T; using domain_t = array; using range_t = array; static constexpr int n = N, m = M; array, N> data; array &operator()(int i){ assert(0 <= i && i < n); return data[i]; } const array &operator()(int i) const{ assert(0 <= i && i < n); return data[i]; } T &operator()(int i, int j){ assert(0 <= i && i < n && 0 <= j && j < m); return data[i][j]; } const T &operator()(int i, int j) const{ assert(0 <= i && i < n && 0 <= j && j < m); return data[i][j]; } bool operator==(const matrix_fixed &a) const{ assert(n == a.n && m == a.m); return data == a.data; } bool operator!=(const matrix_fixed &a) const{ assert(n == a.n && m == a.m); return data != a.data; } matrix_fixed &operator+=(const matrix_fixed &a){ assert(n == a.n && m == a.m); for(auto i = 0; i < n; ++ i) for(auto j = 0; j < m; ++ j) data[i][j] += a(i, j); return *this; } matrix_fixed operator+(const matrix_fixed &a) const{ assert(n == a.n && m == a.m); return matrix_fixed(*this) += a; } matrix_fixed &operator-=(const matrix_fixed &a){ assert(n == a.n && m == a.m); for(auto i = 0; i < n; ++ i) for(auto j = 0; j < m; ++ j) data[i][j] += a(i, j); return *this; } matrix_fixed operator-(const matrix_fixed &a) const{ assert(n == a.n && m == a.m); return matrix_fixed(*this) += a; } template matrix_fixed operator*(const matrix_fixed &a) const{ assert(m == a.n); int l = M2; matrix_fixed res; for(auto i = 0; i < n; ++ i) for(auto j = 0; j < m; ++ j) for(auto k = 0; k < l; ++ k) res(i, k) += data[i][j] * a(j, k); return res; } template matrix_fixed &operator*=(const matrix_fixed &a){ return *this = *this * a; } matrix_fixed &operator*=(T c){ for(auto i = 0; i < n; ++ i) for(auto j = 0; j < m; ++ j) data[i][j] *= c; return *this; } matrix_fixed operator*(T c) const{ return matrix_fixed(*this) *= c; } template::value>::type* = nullptr> matrix_fixed &inplace_power(U e){ assert(n == m && e >= 0); matrix_fixed res(1, 0); for(; e; *this *= *this, e >>= 1) if(e & 1) res *= *this; return *this = res; } template matrix_fixed power(U e) const{ return matrix_fixed(*this).inplace_power(e); } matrix_fixed transposed() const{ matrix_fixed res; for(auto i = 0; i < n; ++ i) for(auto j = 0; j < m; ++ j) res(j, i) = data[i][j]; return res; } matrix_fixed &transpose(){ return *this = transposed(); } // Multiply a column vector v on the right range_t operator*(const domain_t &v) const{ range_t res; res.fill(T(0)); for(auto i = 0; i < n; ++ i) for(auto j = 0; j < m; ++ j) res[i] += data[i][j] * v[j]; return res; } // Assumes T is a field // find_inverse() must return optional // O(n) find_inverse() calls along with O(n^3) operations on T T determinant(auto find_inverse) const{ assert(n == m); if(n == 0) return T(1); auto a = data; T res = T(1); for(auto j = 0; j < n; ++ j){ int pivot = -1; for(auto i = j; i < n; ++ i) if(a[i][j] != T(0)){ pivot = i; break; } if(!~pivot) return T(0); swap(a[j], a[pivot]); res *= a[j][j] * (j != pivot ? -1 : 1); auto invp = find_inverse(a[j][j]); assert(invp); T inv = *invp; for(auto i = j + 1; i < n; ++ i) if(i != j && a[i][j] != T(0)){ T d = a[i][j] * inv; for(auto jj = j; jj < n; ++ jj) a[i][jj] -= d * a[j][jj]; } } return res; } // Assumes T is a field // find_inverse() must return optional // O(n) find_inverse() calls along with O(n^3) operations on T optional inverse(auto find_inverse) const{ assert(n == m); if(n == 0) return *this; auto a = data; array, N> res; for(auto i = 0; i < n; ++ i) res[i].fill(T(0)), res[i][i] = T(1); for(auto j = 0; j < n; ++ j){ int pivot = -1; for(auto i = j; i < n; ++ i) if(a[i][j] != T(0)){ pivot = i; break; } if(!~pivot) return {}; swap(a[j], a[pivot]), swap(res[j], res[pivot]); auto invp = find_inverse(a[j][j]); assert(invp); T inv = *invp; for(auto jj = 0; jj < n; ++ jj) a[j][jj] *= inv, res[j][jj] *= inv; for(auto i = 0; i < n; ++ i) if(i != j && a[i][j] != T(0)){ T d = a[i][j]; for(auto jj = 0; jj < n; ++ jj) a[i][jj] -= d * a[j][jj], res[i][jj] -= d * res[j][jj]; } } return matrix_fixed(n, n, res); } template friend output_stream &operator<<(output_stream &out, const matrix_fixed &a){ out << "{"; for(auto i = 0; i < a.n; ++ i){ out << "{"; for(auto j = 0; j < a.m; ++ j){ out << a(i, j); if(j != a.m - 1) out << ", "; } out << "}"; if(i != a.n - 1) out << ", "; } return out << "}"; } matrix_fixed(): matrix_fixed(T(0), T(0)){ } matrix_fixed(const T &init_diagonal, const T &init_off_diagonal){ for(auto i = 0; i < n; ++ i) for(auto j = 0; j < m; ++ j) data[i][j] = i == j ? init_diagonal : init_off_diagonal; } matrix_fixed(const array, N> &arr): data(arr){ } static matrix_fixed additive_identity(){ return matrix_fixed(T(1), T(0)); } static matrix_fixed multiplicative_identity(){ return matrix_fixed(T(0), T(0)); } }; template matrix_fixed operator*(T c, matrix_fixed a){ for(auto i = 0; i < a.n; ++ i) for(auto j = 0; j < a.m; ++ j) a(i, j) = c * a(i, j); return a; } // Multiply a row vector v on the left template matrix_fixed::domain_t operator*(const typename matrix_fixed::range_t &v, const matrix_fixed &a){ typename matrix_fixed::domain_t res; res.fill(T(0)); for(auto i = 0; i < a.n; ++ i) for(auto j = 0; j < a.m; ++ j) res[j] += v[i] * a(i, j); return res; } int main(){ cin.tie(0)->sync_with_stdio(0); cin.exceptions(ios::badbit | ios::failbit); matrix_fixed mat(2, 2); for(auto i = 0; i < 2; ++ i){ for(auto j = 0; j < 2; ++ j){ cin >> mat(i, j); } } mat = mat.power(3); for(auto i = 0; i < 2; ++ i){ for(auto j = 0; j < 2; ++ j){ cout << mat(i, j) << " "; } cout << "\n"; } return 0; } /* */ //////////////////////////////////////////////////////////////////////////////////////// // // // Coded by Aeren // // // ////////////////////////////////////////////////////////////////////////////////////////