#include #include using namespace std; #if __cplusplus > 201703L #include using namespace numbers; #endif template struct modular_unfixed_base{ static unsigned int _mod; static unsigned long long _inverse_mod; static unsigned int &mod(){ return _mod; } static void precalc_barrett(){ _inverse_mod = (unsigned long long)-1 / _mod + 1; } static void setup(unsigned int mod = 0){ if(!mod) cin >> mod; _mod = mod; assert(_mod >= 1); precalc_barrett(); } template static vector precalc_power(T base, int SZ){ vector res(SZ + 1, 1); for(auto i = 1; i <= SZ; ++ i) res[i] = res[i - 1] * base; return res; } static vector _INV; static void precalc_inverse(int SZ){ if(_INV.empty()) _INV.assign(2, 1); for(auto x = _INV.size(); x <= SZ; ++ x) _INV.push_back(_mod / x * -_INV[_mod % x]); } // _mod must be a prime static modular_unfixed_base _primitive_root; static modular_unfixed_base primitive_root(){ if(_primitive_root) return _primitive_root; if(_mod == 2) return _primitive_root = 1; if(_mod == 998244353) return _primitive_root = 3; unsigned int divs[20] = {}; divs[0] = 2; int cnt = 1; unsigned int x = (_mod - 1) / 2; while(x % 2 == 0) x /= 2; for(auto i = 3; 1LL * i * i <= x; i += 2){ if(x % i == 0){ divs[cnt ++] = i; while(x % i == 0) x /= i; } } if(x > 1) divs[cnt ++] = x; for(auto g = 2; ; ++ g){ bool ok = true; for(auto i = 0; i < cnt; ++ i){ if((modular_unfixed_base(g).power((_mod - 1) / divs[i])) == 1){ ok = false; break; } } if(ok) return _primitive_root = g; } } constexpr modular_unfixed_base(): data(){ } modular_unfixed_base(const double &x){ data = normalize(llround(x)); } modular_unfixed_base(const long double &x){ data = normalize(llround(x)); } template::value>::type* = nullptr> modular_unfixed_base(const T &x){ data = normalize(x); } template::value>::type* = nullptr> static unsigned int normalize(const T &x){ if(_mod == 1) return 0; assert(_inverse_mod); int sign = x >= 0 ? 1 : -1; unsigned int v = _mod <= sign * x ? sign * x - ((__uint128_t)(sign * x) * _inverse_mod >> 64) * _mod : sign * x; if(v >= _mod) v += _mod; if(sign == -1 && v) v = _mod - v; return v; } const unsigned int &operator()() const{ return data; } template operator T() const{ return data; } modular_unfixed_base &operator+=(const modular_unfixed_base &otr){ if((data += otr.data) >= _mod) data -= _mod; return *this; } modular_unfixed_base &operator-=(const modular_unfixed_base &otr){ if((data += _mod - otr.data) >= _mod) data -= _mod; return *this; } template::value>::type* = nullptr> modular_unfixed_base &operator+=(const T &otr){ return *this += modular_unfixed_base(otr); } template::value>::type* = nullptr> modular_unfixed_base &operator-=(const T &otr){ return *this -= modular_unfixed_base(otr); } modular_unfixed_base &operator++(){ return *this += 1; } modular_unfixed_base &operator--(){ return *this += _mod - 1; } modular_unfixed_base operator++(int){ modular_unfixed_base result(*this); *this += 1; return result; } modular_unfixed_base operator--(int){ modular_unfixed_base result(*this); *this += _mod - 1; return result; } modular_unfixed_base operator-() const{ return modular_unfixed_base(_mod - data); } modular_unfixed_base &operator*=(const modular_unfixed_base &rhs){ data = normalize((unsigned long long)data * rhs.data); return *this; } template::value>::type* = nullptr> modular_unfixed_base &inplace_power(T e){ if(e < 0) *this = 1 / *this, e = -e; modular_unfixed_base res = 1; for(; e; *this *= *this, e >>= 1) if(e & 1) res *= *this; return *this = res; } template::value>::type* = nullptr> modular_unfixed_base power(T e) const{ return modular_unfixed_base(*this).inplace_power(e); } modular_unfixed_base &operator/=(const modular_unfixed_base &otr){ int a = otr.data, m = _mod, u = 0, v = 1; if(a < _INV.size()) return *this *= _INV[a]; while(a){ int t = m / a; m -= t * a; swap(a, m); u -= t * v; swap(u, v); } assert(m == 1); return *this *= u; } unsigned int data; }; template unsigned int modular_unfixed_base::_mod; template unsigned long long modular_unfixed_base::_inverse_mod; template vector> modular_unfixed_base::_INV; template modular_unfixed_base modular_unfixed_base::_primitive_root; template bool operator==(const modular_unfixed_base &lhs, const modular_unfixed_base &rhs){ return lhs.data == rhs.data; } template::value>::type* = nullptr> bool operator==(const modular_unfixed_base &lhs, T rhs){ return lhs == modular_unfixed_base(rhs); } template::value>::type* = nullptr> bool operator==(T lhs, const modular_unfixed_base &rhs){ return modular_unfixed_base(lhs) == rhs; } template bool operator!=(const modular_unfixed_base &lhs, const modular_unfixed_base &rhs){ return !(lhs == rhs); } template::value>::type* = nullptr> bool operator!=(const modular_unfixed_base &lhs, T rhs){ return !(lhs == rhs); } template::value>::type* = nullptr> bool operator!=(T lhs, const modular_unfixed_base &rhs){ return !(lhs == rhs); } template bool operator<(const modular_unfixed_base &lhs, const modular_unfixed_base &rhs){ return lhs.data < rhs.data; } template bool operator>(const modular_unfixed_base &lhs, const modular_unfixed_base &rhs){ return lhs.data > rhs.data; } template bool operator<=(const modular_unfixed_base &lhs, const modular_unfixed_base &rhs){ return lhs.data <= rhs.data; } template bool operator>=(const modular_unfixed_base &lhs, const modular_unfixed_base &rhs){ return lhs.data >= rhs.data; } template modular_unfixed_base operator+(const modular_unfixed_base &lhs, const modular_unfixed_base &rhs){ return modular_unfixed_base(lhs) += rhs; } template::value>::type* = nullptr> modular_unfixed_base operator+(const modular_unfixed_base &lhs, T rhs){ return modular_unfixed_base(lhs) += rhs; } template::value>::type* = nullptr> modular_unfixed_base operator+(T lhs, const modular_unfixed_base &rhs){ return modular_unfixed_base(lhs) += rhs; } template modular_unfixed_base operator-(const modular_unfixed_base &lhs, const modular_unfixed_base &rhs){ return modular_unfixed_base(lhs) -= rhs; } template::value>::type* = nullptr> modular_unfixed_base operator-(const modular_unfixed_base &lhs, T rhs){ return modular_unfixed_base(lhs) -= rhs; } template::value>::type* = nullptr> modular_unfixed_base operator-(T lhs, const modular_unfixed_base &rhs){ return modular_unfixed_base(lhs) -= rhs; } template modular_unfixed_base operator*(const modular_unfixed_base &lhs, const modular_unfixed_base &rhs){ return modular_unfixed_base(lhs) *= rhs; } template::value>::type* = nullptr> modular_unfixed_base operator*(const modular_unfixed_base &lhs, T rhs){ return modular_unfixed_base(lhs) *= rhs; } template::value>::type* = nullptr> modular_unfixed_base operator*(T lhs, const modular_unfixed_base &rhs){ return modular_unfixed_base(lhs) *= rhs; } template modular_unfixed_base operator/(const modular_unfixed_base &lhs, const modular_unfixed_base &rhs) { return modular_unfixed_base(lhs) /= rhs; } template::value>::type* = nullptr> modular_unfixed_base operator/(const modular_unfixed_base &lhs, T rhs) { return modular_unfixed_base(lhs) /= rhs; } template::value>::type* = nullptr> modular_unfixed_base operator/(T lhs, const modular_unfixed_base &rhs) { return modular_unfixed_base(lhs) /= rhs; } template istream &operator>>(istream &in, modular_unfixed_base &number){ long long x; in >> x; number.data = modular_unfixed_base::normalize(x); return in; } // #define _PRINT_AS_FRACTION template ostream &operator<<(ostream &out, const modular_unfixed_base &number){ #ifdef LOCAL #ifdef _PRINT_AS_FRACTION out << number(); cerr << "("; for(auto d = 1; ; ++ d){ if((number * d).data <= 1000000){ cerr << (number * d).data << "/" << d; break; } else if((-number * d).data <= 1000000){ cerr << "-" << (-number * d).data << "/" << d; break; } } cerr << ")"; return out; #else return out << number(); #endif #else return out << number(); #endif } #undef _PRINT_AS_FRACTION using modular = modular_unfixed_base<0>; // T must support +=, -=, *, *=, ==, and != template struct matrix{ using ring_t = T; using domain_t = vector; using range_t = vector; int n, m; vector> data; vector &operator()(int i){ assert(0 <= i && i < n); return data[i]; } const vector &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 &a) const{ assert(n == a.n && m == a.m); return data == a.data; } bool operator!=(const matrix &a) const{ assert(n == a.n && m == a.m); return data != a.data; } matrix &operator+=(const matrix &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 operator+(const matrix &a) const{ return matrix(*this) += a; } matrix &operator-=(const matrix &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 operator-(const matrix &a) const{ return matrix(*this) -= a; } matrix operator*=(const matrix &a){ assert(m == a.n); int l = a.m; matrix res(n, l); 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 *this = res; } matrix operator*(const matrix &a) const{ return matrix(*this) *= a; } matrix &operator*=(T c){ for(auto i = 0; i < n; ++ i) for(auto j = 0; j < m; ++ j) data[i][j] *= c; return *this; } matrix operator*(T c) const{ return matrix(*this) *= c; } template::value>::type* = nullptr> matrix &inplace_power(U e){ assert(n == m && e >= 0); matrix res(n, n, T(1)); for(; e; *this *= *this, e >>= 1) if(e & 1) res *= *this; return *this = res; } template matrix power(U e) const{ return matrix(*this).inplace_power(e); } matrix transposed() const{ matrix res(m, n); for(auto i = 0; i < n; ++ i) for(auto j = 0; j < m; ++ j) res(j, i) = data[i][j]; return res; } range_t operator*(const domain_t &v) const{ assert(m == (int)v.size()); vector res(n, 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 either a floating, integral, or a modular type. // If T is a floating type, O(n) divisions with O(n^3) additions, subtractions, and multiplications. // Otherwise, O(n^2 * log(size)) divisions with O(n^3) additions, subtractions, and multiplications. // Returns {REF matrix, determinant, rank} tuple inplace_REF(){ if(n == 0) return {*this, T(1), 0}; T det = 1; int rank = 0; for(auto j = 0; j < m; ++ j){ if constexpr(is_floating_point_v){ static const T eps = 1e-9; int pivot = rank; for(auto i = rank + 1; i < n; ++ i) if(abs(data[pivot][j]) < abs(data[i][j])) pivot = i; if(rank != pivot){ swap(data[rank], data[pivot]); det *= -1; } if(abs(data[rank][j]) <= eps) continue; det *= data[rank][j]; T inv = 1 / data[rank][j]; for(auto i = rank + 1; i < n; ++ i) if(abs(data[i][j]) > eps){ T coef = data[i][j] * inv; for(auto k = j; k < m; ++ k) data[i][k] -= coef * data[rank][k]; } } else{ for(auto i = rank + 1; i < n; ++ i) while(data[i][j]){ T q; if constexpr(is_integral_v) q = data[rank][j] / data[i][j]; else q = data[rank][j].data / data[i][j].data; if(q) for(auto k = j; k < m; ++ k) data[rank][k] -= q * data[i][k]; swap(data[rank], data[i]); det *= -1; } if(rank == j) det *= data[rank][j]; else det = T(0); if(data[rank][j]) ++ rank; } } return {*this, det, rank}; } // Assumes T is either a floating, integral, or a modular type. // If T is a floating type, O(n) divisions with O(n^3) additions, subtractions, and multiplications. // Otherwise, O(n^2 * log(size)) divisions with O(n^3) additions, subtractions, and multiplications. // Returns {REF matrix, determinant, rank} tuple REF() const{ return matrix(*this).inplace_REF(); } // Assumes T is a field. // O(n) divisions with O(n^3) additions, subtractions, and multiplications. // Returns {REF matrix, determinant, rank} tuple inplace_REF_field(){ if(n == 0) return {*this, T(1), 0}; T det = T(1); int rank = 0; for(auto j = 0; j < m; ++ j){ int pivot = -1; for(auto i = rank; i < n; ++ i) if(data[i][j] != T(0)){ pivot = i; break; } if(!~pivot){ det = T(0); continue; } if(rank != pivot){ swap(data[rank], data[pivot]); det *= -1; } det *= data[rank][j]; T inv = 1 / data[rank][j]; for(auto i = rank + 1; i < n; ++ i) if(data[i][j] != T(0)){ T coef = data[i][j] * inv; for(auto k = j; k < m; ++ k) data[i][k] -= coef * data[j][k]; } ++ rank; } return {*this, det, rank}; } // Assumes T is a field. // O(n) divisions with O(n^3) additions, subtractions, and multiplications. // Returns {REF matrix, determinant, rank} tuple REF_field() const{ return matrix(*this).inplace_REF_field(); } // Assumes T is a field. // O(n) divisions with O(n^3) additions, subtractions, and multiplications. optional inverse(auto find_inverse) const{ assert(n == m); if(n == 0) return *this; auto a = data; vector> res(n, vector(n, T(0))); for(auto i = 0; i < n; ++ i) res[i][i] = T(1); for(auto j = 0; j < n; ++ j){ int rank = j, pivot = -1; for(auto i = rank; i < n; ++ i) if(a[i][j] != T(0)){ pivot = i; break; } if(!~pivot) return {}; swap(a[rank], a[pivot]), swap(res[rank], res[pivot]); T inv = 1 / a[rank][j]; for(auto k = j; k < n; ++ k) a[rank][k] *= inv, res[rank][k] *= inv; for(auto i = 0; i < n; ++ i) if(i != j && a[i][j] != T(0)){ T d = a[i][j]; for(auto k = j; k < n; ++ k) a[i][k] -= d * a[j][k], res[i][k] -= d * res[j][k]; } } return res; } // Assumes T is either a floating, integral, or a modular type. // If T is a floating type, O(n) divisions with O(n^3) additions, subtractions, and multiplications. // Otherwise, O(n^2 * log(size)) divisions with O(n^3) additions, subtractions, and multiplications. T determinant() const{ assert(n == m); return get<1>(REF()); } // Assumes T is a field. // O(n) divisions with O(n^3) additions, subtractions, and multiplications. T determinant_field() const{ assert(n == m); return get<1>(REF_field()); } // Assumes T is either a floating, integral, or a modular type. // If T is a floating type, O(n) divisions with O(n^3) additions, subtractions, and multiplications. // Otherwise, O(n^2 * log(size)) divisions with O(n^3) additions, subtractions, and multiplications. int rank() const{ return get<2>(REF()); } // Assumes T is a field. // O(n) divisions with O(n^3) additions, subtractions, and multiplications. int rank_field() const{ return get<2>(REF_field()); } // O(n * 2^n) T permanent() const{ assert(n <= 30 && n == m); T perm = n ? 0 : 1; vector sum(n); for(auto order = 1; order < 1 << n; ++ order){ int j = __lg(order ^ order >> 1 ^ order - 1 ^ order - 1 >> 1), sign = (order ^ order >> 1) & 1 << j ? 1 : -1; T prod = order & 1 ? -1 : 1; if((order ^ order >> 1) & 1 << j) for(auto i = 0; i < n; ++ i) prod *= sum[i] += data[i][j]; else for(auto i = 0; i < n; ++ i) prod *= sum[i] -= data[i][j]; perm += prod; } return perm * (n & 1 ? -1 : 1); } template friend output_stream &operator<<(output_stream &out, const matrix &a){ out << "\n"; for(auto i = 0; i < a.n; ++ i){ for(auto j = 0; j < a.m; ++ j){ out << a(i, j) << " "; } out << "\n"; } return out; } matrix(int n, int m, T init_diagonal = T(0), T init_off_diagonal = T(0)): n(n), m(m){ data.assign(n, vector(m, T(0))); 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(const vector> &arr, int _n = -1, int _m = -1): n(~_n ? _n : (int)arr.size()), m(~_m ? _m : (int)arr[0].size()), data(arr){ } }; template matrix operator*(T c, matrix 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 typename matrix::domain_t operator*(const typename matrix::range_t &v, const matrix &a){ assert(a.n == (int)size(v)); typename matrix::domain_t res(a.m, 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); int n, mod; cin >> n >> mod, mod <<= 1; modular::setup(mod); matrix a(n, n); for(auto i = 0; i < n; ++ i){ for(auto j = 0; j < n; ++ j){ cin >> a(i, j); } } cout << (a.permanent() - a.determinant()).data / 2 << "\n"; return 0; } /* */ //////////////////////////////////////////////////////////////////////////////////////// // // // Coded by Aeren // // // ////////////////////////////////////////////////////////////////////////////////////////