#define SINGLE_TESTCASE // #define MULTI_TESTCASE // #define AOJ_TESTCASE #define FAST_IO // #define FAST_CIO // #define INTERACTIVE #define INF 4'000'000'000'000'000'037LL #define EPS 1e-11 // https://github.com/miscalculation53/library/tree/wip/template/template_all.hpp // https://github.com/miscalculation53/library/tree/wip/template/template_all_but_modint.hpp // https://github.com/miscalculation53/library/tree/wip/template/template_types.hpp #include namespace { using namespace std; using ll = long long; using uint = unsigned int; using ull = unsigned long long; using pll = pair; #define vc vector template using vvc = vc>; using vstr = vc; using i128 = __int128_t; using u128 = __uint128_t; i128 stoi128(const string &s) { const bool neg = s.front() == '-'; u128 res = 0; for (int i = neg; i < (int)s.size(); i++) res = 10 * res + s[i] - '0'; if (neg) return -i128(res - 1) - 1; return i128(res); } string i128tos(i128 x) { if (x == 0) return "0"; string sign = "", res = ""; u128 ux; if (x < 0) ux = u128(-(x + 1)) + 1, sign = "-"; else ux = x; while (ux > 0) { res += '0' + ux % 10; ux /= 10; } reverse(res.begin(), res.end()); return sign + res; } istream &operator>>(istream &is, i128 &a) { string s; is >> s; a = stoi128(s); return is; } ostream &operator<<(ostream &os, const i128 &a) { os << i128tos(a); return os; } #define cauto const auto // https://github.com/miscalculation53/library/tree/wip/template/template_rep.hpp #define overload4(_1,_2,_3,_4,name,...) name #define rep1(i,n) for (ll i = 0, nnnnn = ll(n); i < nnnnn; i++) #define rep2(i,l,r) for (ll i = ll(l), rrrrr = ll(r); i < rrrrr; i++) #define rep3(i,l,r,d) for (ll i = ll(l), rrrrr = ll(r), ddddd = ll(d); ddddd > 0 ? i < rrrrr : i > rrrrr; i += d) #define rep(...) overload4(__VA_ARGS__, rep3, rep2, rep1)(__VA_ARGS__) #define repi1(i,n) for (int i = 0, nnnnn = int(n); i < nnnnn; i++) #define repi2(i,l,r) for (int i = int(l), rrrrr = int(r); i < rrrrr; i++) #define repi3(i,l,r,d) for (int i = int(l), rrrrr = int(r), ddddd = int(d); ddddd > 0 ? i < rrrrr : i > rrrrr; i += d) #define repi(...) overload4(__VA_ARGS__, repi3, repi2, repi1)(__VA_ARGS__) #define fe(...) for (auto __VA_ARGS__) #define fec(...) for (cauto &__VA_ARGS__) #define fem(...) for (auto &__VA_ARGS__) // https://github.com/miscalculation53/library/tree/wip/template/template_math.hpp // https://github.com/miscalculation53/library/tree/wip/utils/is_integral_ext.hpp template constexpr bool is_integral_ext = is_integral_v || is_same_v || is_same_v; template constexpr bool is_signed_ext = is_signed_v || is_same_v; template constexpr bool is_unsigned_ext = is_unsigned_v || is_same_v; // https://github.com/miscalculation53/library/tree/wip/utils/default_infty.hpp namespace default_infty_detail { template inline constexpr bool unsupported = false; } template && is_integral_ext>> inline constexpr T divfloor(U a, V b) { return T(a) / T(b) - (T(a) % T(b) && (T(a) ^ T(b)) < 0); } template && is_integral_ext>> inline constexpr T safemod(U a, V b) { return T(a) - T(b) * divfloor(a, b); } // https://github.com/miscalculation53/library/tree/wip/template/template_vector.hpp #define ALL(a) (a).begin(), (a).end() template inline T SZ(const V &x) { return x.size(); } #define eb emplace_back #define LMD(x,fx) ([&](const auto &x) { return fx; }) #define GEN_VEC(n,i,fi) (gen_vec(n, LMD(i, fi))) // https://github.com/miscalculation53/library/tree/wip/template/template_algo.hpp // https://github.com/miscalculation53/library/tree/wip/utils/resolved_infty.hpp // https://github.com/miscalculation53/library/tree/wip/utils/resolved_value.hpp template vvc top(const vvc &a) { if (a.empty()) return {}; const int n = a.size(), m = a[0].size(); vvc b(m, vc(n)); repi(i, n) { assert(SZ(a[i]) == m); repi(j, m) b[j][i] = a[i][j]; } return b; } vstr top(const vstr &a) { vvc a_(a.size()); repi(i, SZ(a)) a_[i] = {ALL(a[i])}; vvc b_ = top(a_); vstr b(b_.size()); repi(i, SZ(b)) b[i] = {ALL(b_[i])}; return b; } template struct has_e0 : false_type {}; template struct has_e0> : true_type {}; template inline constexpr bool has_e0_v = has_e0::value; template struct MonoidAdd { using S = T; }; namespace internal { } constexpr array DRULgrid = {{{1, 0}, {0, 1}, {-1, 0}, {0, -1}}}; constexpr array DRULplane = {{{0, -1}, {1, 0}, {0, 1}, {-1, 0}}}; // https://github.com/miscalculation53/library/tree/wip/template/template_binsearch.hpp template struct is_random_access_iterator { static constexpr bool value = is_same_v< typename iterator_traits::iterator_category, random_access_iterator_tag >; }; template constexpr bool is_random_access_iterator_v = is_random_access_iterator::value; #define DEFAULT_COMP ranges::less namespace internal { }; // https://github.com/miscalculation53/library/tree/wip/template/template_bit.hpp template inline constexpr ull MASK(T k) { return (1ULL << k) - 1ULL; } inline constexpr ll bit_width(ll x) { return std::bit_width((ull)x); } inline constexpr ll bit_floor(ll x) { return std::bit_floor((ull)x); } inline constexpr ll bit_ceil(ll x) { return std::bit_ceil((ull)x); } inline constexpr ll countr_zero(ll x) { assert(x != 0); return std::countr_zero((ull)x); } inline constexpr ll popcount(ll x) { return std::popcount((ull)x); } inline constexpr bool has_single_bit(ll x) { return std::has_single_bit((ull)x); } inline constexpr ull lsb_pos(ull x) { assert(x != 0); return countr_zero(x); } inline constexpr ull msb_pos(ull x) { assert(x != 0); return bit_width(x) - 1; } inline constexpr ull lsb_mask(ull x) { assert(x != 0); return x & -x; } inline constexpr ull msb_mask(ull x) { assert(x != 0); return bit_floor(x); } inline constexpr bool btest(ull x, uint k) { return (x >> k) & 1; } inline constexpr bool bsubset(ull x, ull y) { return (x & y) == x; } inline constexpr bool bsupset(ull x, ull y) { return (x & y) == y; } inline constexpr ull bsetminus(ull x, ull y) { return x & ~y; } // https://github.com/miscalculation53/library/tree/wip/template/template_inout.hpp // https://github.com/miscalculation53/library/tree/wip/template/template_dump.hpp // https://github.com/miscalculation53/library/tree/wip/template/template_dump_map.hpp #define CPP_DUMP_DEFINE_DATA(...) #define dump(...) #define local(...) #define oj(...) __VA_ARGS__ #define local_oj(a,b) (b) namespace fastio { template struct unsigned_integer { using type = make_unsigned_t; }; template <> struct unsigned_integer { using type = u128; }; template <> struct unsigned_integer { using type = u128; }; template using unsigned_integer_t = typename unsigned_integer::type; static constexpr uint32_t SIZ = 1 << 17; char ibuf[SIZ]; char obuf[SIZ]; char out[100]; uint32_t pil = 0, pir = 0, por = 0; struct Pre { char num[10000][4]; constexpr Pre() : num() { for (int i = 0; i < 10000; i++) { int n = i; for (int j = 3; j >= 0; j--) { num[i][j] = n % 10 | '0'; n /= 10; } } } } constexpr pre; inline void load() { memcpy(ibuf, ibuf + pil, pir - pil); pir = pir - pil + fread(ibuf + pir - pil, 1, SIZ - pir + pil, stdin); pil = 0; if (pir < SIZ) ibuf[pir++] = '\n'; } inline void flush() { fwrite(obuf, 1, por, stdout); por = 0; } void rd1(char &c) { do { if (pil + 1 > pir) load(); c = ibuf[pil++]; } while (c <= ' '); } void rd1(string &x) { x.clear(); while (true) { if (pil == pir) load(); while (pil < pir && ibuf[pil] <= ' ') ++pil; if (pil < pir) break; } while (true) { uint32_t p = pil; while (pil < pir && ibuf[pil] > ' ') ++pil; x.append(ibuf + p, pil - p); if (pil < pir) { ++pil; return; } load(); } } template void rd1_real(T &x) { string s; rd1(s); if constexpr (!is_same_v) { auto [p, ec] = from_chars(s.data(), s.data() + s.size(), x); if (ec == errc{} && p == s.data() + s.size()) return; } if constexpr (is_same_v) x = stold(s); else x = stod(s); } template void rd1_integer(T &x) { using U = unsigned_integer_t; bool minus = false; U val = 0; if constexpr (check_buffer) if (pil + 100 > pir) load(); uint32_t p = pil; while (ibuf[p] < '-') ++p; if constexpr (is_signed::value || is_same_v) { if (ibuf[p] == '-') minus = true, ++p; } while ('0' <= ibuf[p]) val = val * 10 + (ibuf[p++] & 15); pil = p; if constexpr (is_signed::value || is_same_v) { if (minus) { const U min_abs = U(numeric_limits::max()) + 1; x = val == min_abs ? numeric_limits::lowest() : -T(val); } else x = T(val); } else x = T(val); } void rd1(int &x) { rd1_integer(x); } void rd1(ll &x) { rd1_integer(x); } void rd1(i128 &x) { rd1_integer(x); } void rd1(uint &x) { rd1_integer(x); } void rd1(ull &x) { rd1_integer(x); } void rd1(u128 &x) { rd1_integer(x); } void rd1(double &x) { rd1_real(x); } void rd1(long double &x) { rd1_real(x); } template void rd1(vc &x) { for (auto &d: x) rd1(d); } template void read(T &...x) { if constexpr (sizeof...(T) <= SIZ / 100 && ((!is_same_v && (is_integral_v || is_same_v || is_same_v)) && ...)) { if (pil + 100 * sizeof...(T) > pir) load(); (rd1_integer(x), ...); } else (rd1(x), ...); } void wt1(const char c) { if (por == SIZ) flush(); obuf[por++] = c; } void wt1(string_view s) { while (!s.empty()) { if (por == SIZ) flush(); size_t n = min(s.size(), SIZ - por); memcpy(obuf + por, s.data(), n); por += n; s.remove_prefix(n); } } template void wt1_integer(T x) { if (por > SIZ - 100) flush(); using U = unsigned_integer_t; U ux; if constexpr (is_signed::value || is_same_v) { if (x < 0) obuf[por++] = '-', ux = U(0) - U(x); else ux = U(x); } else ux = x; int outi; for (outi = 96; ux >= 10000; outi -= 4) { memcpy(out + outi, pre.num[ux % 10000], 4); ux /= 10000; } if (ux >= 1000) { memcpy(obuf + por, pre.num[ux], 4); por += 4; } else if (ux >= 100) { memcpy(obuf + por, pre.num[ux] + 1, 3); por += 3; } else if (ux >= 10) { int q = (ux * 103) >> 10; obuf[por] = q | '0'; obuf[por + 1] = (ux - q * 10) | '0'; por += 2; } else obuf[por++] = ux | '0'; memcpy(obuf + por, out + outi + 4, 96 - outi); por += 96 - outi; } template void wt1_real(T x) { if constexpr (!is_same_v) { auto [p, ec] = to_chars(out, out + sizeof(out), x, chars_format::fixed, 15); if (ec == errc{}) { wt1(string_view(out, p)); return; } } ostringstream oss; oss << fixed << setprecision(15) << x; wt1(oss.str()); } void wt1(int x) { wt1_integer(x); } template , int> = 0> void wt1(T x) { wt1_integer(x); } void wt1(i128 x) { wt1_integer(x); } void wt1(u128 x) { wt1_integer(x); } void wt1(double x) { wt1_real(x); } void wt1(long double x) { wt1_real(x); } template void wt1(const vector &val) { auto n = val.size(); for (size_t i = 0; i < n; i++) { if (i) wt1(' '); wt1(val[i]); } } template void write(T &&...x) { (wt1(std::forward(x)), ...); } template void print(T &&...x) { if constexpr (sizeof...(T)) { int i = 0; ((i++ ? wt1(' ') : void(), wt1(std::forward(x))), ...); } wt1('\n'); } } struct Dummy { Dummy() { atexit(fastio::flush); } } dummy; namespace internal { template void READnodump(Ts &...a) { fastio::read(a...); } template void READVEC2nodump(int n, int m, vvc &...v) { (v.assign(n, vc(m)), ...); READnodump(v...); } }; #define READ(...) internal::READnodump(__VA_ARGS__); dump(__VA_ARGS__) #define IN(T,...) T __VA_ARGS__; READ(__VA_ARGS__) #define CHAR(...) IN(char, __VA_ARGS__) #define INT(...) IN(int, __VA_ARGS__) #define LL(...) IN(ll, __VA_ARGS__) #define STR(...) IN(string, __VA_ARGS__) #define ARR(T,n,...) array __VA_ARGS__; READ(__VA_ARGS__) #define READVEC(...) internal::READVECnodump(__VA_ARGS__); dump(__VA_ARGS__) #define READVEC2(...) internal::READVEC2nodump(__VA_ARGS__); dump(__VA_ARGS__) #define VEC(T,n,...) vc __VA_ARGS__; READVEC(n, __VA_ARGS__) #define VEC2(T,n,m,...) vvc __VA_ARGS__; READVEC2(n, m, __VA_ARGS__) #define READJAG(...) internal::READJAGnodump(__VA_ARGS__); dump(__VA_ARGS__) #define JAG(T,n,...) vvc __VA_ARGS__; READJAG(n, __VA_ARGS__) #define ENDL '\n' #define WRITE fastio::write #define PRINT fastio::print #define PRINTEXIT(...) do { PRINT(__VA_ARGS__); exit(0); } while (false) #define PRINTRETURN(...) do { PRINT(__VA_ARGS__); return; } while (false) template void PRINTV(const vc &v) { for (auto &vi : v) PRINT(vi); } #define PRINTVEXIT(...) do { PRINTV(__VA_ARGS__); exit(0); } while (false) #define PRINTVRETURN(...) do { PRINTV(__VA_ARGS__); return; } while (false) namespace internal { }; namespace internal { }; #define UNZIP(vt,...) auto [__VA_ARGS__] = unzip(vt) #define ZIP(vt,...) auto vt = zip(tuple{__VA_ARGS__}) // https://github.com/miscalculation53/library/tree/wip/template/template_random.hpp mt19937_64 mt; bool randbool(double p) { assert(0 <= p && p <= 1); return bernoulli_distribution(p)(mt); } namespace internal { }; // https://github.com/miscalculation53/library/tree/wip/math/modint/template_modint.hpp // https://github.com/miscalculation53/library/tree/wip/math/modint/modint.hpp // https://github.com/miscalculation53/library/tree/wip/math/modint/modint_internal_static.hpp // https://github.com/miscalculation53/library/tree/wip/utils/larger_int.hpp namespace larger_int_detail { } template struct larger_int { private: static constexpr bool check(); static_assert(check()); public: using type = T; }; #define LARGER_INT(T,U) template <> struct larger_int { using type = U; }; LARGER_INT(signed char, short) LARGER_INT(short, int) LARGER_INT(int, long long) LARGER_INT(long, __int128_t) LARGER_INT(long long, __int128_t) LARGER_INT(unsigned char, unsigned short) LARGER_INT(unsigned short, unsigned int) LARGER_INT(unsigned int, unsigned long long) LARGER_INT(unsigned long, __uint128_t) LARGER_INT(unsigned long long, __uint128_t) #undef LARGER_INT template struct Rational; template struct larger_int> { using type = Rational::type>; }; template using larger_int_t = typename larger_int::type; // https://github.com/miscalculation53/library/tree/wip/math/modint/modint_internal_isprime.hpp namespace internal { template constexpr ll powmod_constexpr(ll x, ll n, T m) { if (m == 1) return 0; using U = make_unsigned_t; using L = larger_int_t; U r = 1, y = safemod(x, m); while (n) { if (n & 1) r = L(r) * y % m; y = L(y) * y % m; n >>= 1; } return r; } template constexpr bool isprime = isprime_constexpr(n); }; namespace internal { template struct policy_static { using mod_type = decltype(M); using value_type = make_unsigned_t; using calc_type = larger_int_t; static constexpr bool is_prime = isprime_constexpr(M); static constexpr mod_type mod() { return M; } static constexpr value_type umod() { return M; } static constexpr value_type init(value_type v) { return v; } static constexpr mod_type val(value_type v) { return v; } static constexpr value_type mul(value_type a, value_type b) { return (value_type)((calc_type(a) * b) % M); } }; }; // https://github.com/miscalculation53/library/tree/wip/math/modint/modint_internal_barrett32.hpp namespace internal { struct barrett32 { uint m; ull im; explicit barrett32(uint m) : m(m), im((ull)(-1) / m + 1) {} uint umod() const { return m; } uint mul(uint a, uint b) const { ull z = a; z *= b; ull x = ull((u128(z) * im) >> 64); ull y = x * m; return uint(z - y + (z < y ? m : 0)); } }; template struct policy_barrett32 { using value_type = uint; using calc_type = ull; using mod_type = int; static constexpr bool is_prime = false; static inline barrett32 reducer{998244353}; static void set_mod(mod_type m) { reducer = barrett32(m); } static mod_type mod() { return reducer.umod(); } static value_type umod() { return reducer.umod(); } static value_type init(value_type v) { return v; } static mod_type val(value_type v) { return v; } static value_type mul(value_type a, value_type b) { return reducer.mul(a, b); } }; }; // https://github.com/miscalculation53/library/tree/wip/math/modint/modint_internal_montgomery64.hpp namespace internal { inline constexpr ull inv64(ull a) { ull x = a; while (a * x != 1) x *= 2 - a * x; return x; } struct montgomery64odd { ull m, im, sq; explicit montgomery64odd(ull m) : m(m), im(inv64(m)), sq(-u128(m) % m) {} ull umod() const { return m; } ull reduce(u128 x) const { auto t = (x + u128(m) * (-im * ull(x))) >> 64; if (t >= m) t -= m; return (ull)t; } ull inv_reduce(i128 v) const { return reduce(u128(v % m + m) * sq); } }; struct montgomery64 { ull m, mx, imx, d, q; uint b; explicit montgomery64(ull m) : m(m) { b = countr_zero(m), mx = m >> b; imx = inv64(mx); d = powmod_constexpr((mx + 1) / 2, b, mx); u128 sq = -u128(mx) % mx; q = (1 + (((sq - 1) * d) << b)) % m; } ull umod() const { return m; } ull reduce(u128 x) const { if (b == 0) { auto t = (x + u128(mx) * (-imx * ull(x))) >> 64; if (t >= m) t -= m; return (ull)t; } ull p = x & MASK(b); x = (x >> b) + p * d; ull y = p << (64 - b); auto t = (x + u128(mx) * (imx * (y - ull(x)))) >> (64 - b); if (t >= m) { t -= m; if (t >= m) t -= m; } return (ull)t; } ull inv_reduce(i128 v) const { return reduce(u128(v % m + m) * q); } }; template struct policy_montgomery64_odd { using value_type = ull; using calc_type = u128; using mod_type = ll; static constexpr bool is_prime = false; static inline montgomery64odd reducer{(1LL << 61) - 1}; static mod_type mod() { return reducer.umod(); } static value_type umod() { return reducer.umod(); } static value_type init(value_type v) { return reducer.inv_reduce(v); } static mod_type val(value_type v) { return reducer.reduce(v); } static value_type mul(value_type a, value_type b) { return reducer.reduce((calc_type)a * b); } }; template struct policy_montgomery64 { using value_type = ull; using calc_type = u128; using mod_type = ll; static constexpr bool is_prime = false; static inline montgomery64 reducer{(1LL << 61) - 1}; static mod_type mod() { return reducer.umod(); } static value_type umod() { return reducer.umod(); } static value_type init(value_type v) { return reducer.inv_reduce(v); } static mod_type val(value_type v) { return reducer.reduce(v); } static value_type mul(value_type a, value_type b) { return reducer.reduce((calc_type)a * b); } }; }; // https://github.com/miscalculation53/library/tree/wip/math/extgcd.hpp template constexpr tuple extgcd(T a, T b) { if (a == 0 && b == 0) return {0, 0, 0}; T x1 = 1, y1 = 0, z1 = a; T x2 = 0, y2 = 1, z2 = b; while (z2 != 0) { T q = z1 / z2; tie(x1, x2) = make_pair(x2, x1 - q * x2); tie(y1, y2) = make_pair(y2, y1 - q * y2); tie(z1, z2) = make_pair(z2, z1 - q * z2); } if (z1 < 0) x1 = -x1, y1 = -y1, z1 = -z1; return {z1, x1, y1}; } namespace internal { template struct modint_impl { using V = typename Policy::value_type; using M = typename Policy::mod_type; using mint = modint_impl; private: V _v; public: static constexpr M mod() { return Policy::mod(); } template static auto set_mod(M m) -> decltype(T::set_mod(m)) { return T::set_mod(m); } static mint raw(V v) { mint x; x._v = v; return x; } modint_impl() : _v(0) {} template >> modint_impl(T v) { V rem; if constexpr (is_signed_ext) { using S = make_signed_t; S x = v % S(Policy::umod()); if (x < 0) x += Policy::umod(); rem = x; } else rem = V(v % Policy::umod()); _v = Policy::init(rem); }; M val() const { return Policy::val(_v); } mint &operator+=(const mint &rhs) { _v += rhs._v; if (_v >= Policy::umod()) _v -= Policy::umod(); return *this; } mint &operator-=(const mint &rhs) { _v -= rhs._v; if (_v >= Policy::umod()) _v += Policy::umod(); return *this; } mint &operator*=(const mint &rhs) { _v = Policy::mul(_v, rhs._v); return *this; } mint &operator/=(const mint &rhs) { return *this *= rhs.inv(); } mint operator-() const { return mint() - *this; } mint inv() const { if constexpr (Policy::is_prime) { return pow(mod() - 2); } else { auto [g, x, y] = extgcd(val(), mod()); assert(g == 1); return mint(x); } } friend mint operator+(const mint &lhs, const mint &rhs) { return mint(lhs) += rhs; } friend mint operator-(const mint &lhs, const mint &rhs) { return mint(lhs) -= rhs; } friend mint operator*(const mint &lhs, const mint &rhs) { return mint(lhs) *= rhs; } friend mint operator/(const mint &lhs, const mint &rhs) { return mint(lhs) /= rhs; } friend bool operator==(const mint &lhs, const mint &rhs) { return lhs._v == rhs._v; } friend bool operator!=(const mint &lhs, const mint &rhs) { return lhs._v != rhs._v; } friend M safe_hash_key(const mint &x) { return x.val(); } friend void rd1(mint &x) { long long a; fastio::rd1(a); x = a; } friend void wt1(const mint &x) { fastio::wt1(x.val()); } }; }; template using static_modint32 = internal::modint_impl>; template using dynamic_modint32 = internal::modint_impl>; template using static_modint64 = internal::modint_impl>; template using dynamic_modint64_odd = internal::modint_impl>; template using dynamic_modint64 = internal::modint_impl>; using modint = dynamic_modint32<-1>; template struct is_modint : std::false_type { }; template struct is_modint> : std::true_type { }; template inline constexpr bool is_modint_v = is_modint::value; template struct is_static_modint : false_type {}; template struct is_static_modint> : true_type {}; template struct is_static_modint> : true_type {}; template inline constexpr bool is_static_modint_v = is_static_modint::value; template struct is_dynamic_modint : false_type {}; template struct is_dynamic_modint> : true_type {}; template struct is_dynamic_modint> : true_type {}; template struct is_dynamic_modint> : true_type {}; template inline constexpr bool is_dynamic_modint_v = is_dynamic_modint::value; // https://github.com/miscalculation53/library/tree/wip/math/modint/power_table.hpp // https://github.com/miscalculation53/library/tree/wip/math/modint/binomial.hpp template struct Binomial { private: inline static decltype(T::mod()) mod; public: inline static vc fac_, finv_, inv_; }; // https://github.com/miscalculation53/library/tree/wip/math/modint/stom.hpp // https://github.com/miscalculation53/library/tree/wip/math/modint/to_rational.hpp // https://github.com/miscalculation53/library/tree/wip/math/svp2d.hpp namespace cpp_dump { template struct rat_value { T p, q; friend ostream &operator<<(ostream &os, const rat_value &x) { os << x.p; if (x.q != 1) os << '/' << x.q; return os; } }; template inline constexpr bool rat_is_map = false; template inline constexpr bool rat_is_map> = true; template inline constexpr bool rat_is_map> = true; template inline constexpr bool rat_is_map> = true; template inline constexpr bool rat_is_map> = true; struct rat_closure : std::ranges::range_adaptor_closure { }; constexpr rat_closure rat() { return rat_closure{}; } } // using mint = modint998244353; // using mint = modint1000000007; // using mint = static_modint<1000000000>; using mint = modint; void init() { oj(mt.seed(random_device()())); } // https://github.com/miscalculation53/library/tree/wip/math/linalg/matrix.hpp // https://github.com/miscalculation53/library/tree/wip/algebra/algebra_basic_ops.hpp // https://github.com/miscalculation53/library/tree/wip/algebra/algebra_base.hpp template struct SemiRingFromMonoidMonoid { static_assert(is_same_v, "Madd::S and Mmul::S must be identical"); using S = typename Madd::S; static constexpr auto add = Madd::op; static constexpr auto e0 = Madd::e; static constexpr auto mul = Mmul::op; static constexpr auto e1 = Mmul::e; }; template struct RingFromGroupMonoid { static_assert(is_same_v, "Gadd::S and Mmul::S must be identical"); using S = typename Gadd::S; static constexpr auto add = Gadd::op; static constexpr auto e0 = Gadd::e; static constexpr auto minus = Gadd::inv; static constexpr auto mul = Mmul::op; static constexpr auto e1 = Mmul::e; }; template struct FieldFromGroupGroup { static_assert(is_same_v, "Gadd::S and Gmul::S must be identical"); using S = typename Gadd::S; static constexpr auto add = Gadd::op; static constexpr auto e0 = Gadd::e; static constexpr auto minus = Gadd::inv; static constexpr auto mul = Gmul::op; static constexpr auto e1 = Gmul::e; static constexpr auto inv = Gmul::inv; }; template struct has_e1 : false_type {}; template struct has_e1> : true_type {}; template inline constexpr bool has_e1_v = has_e1::value; template struct MonoidMul { using S = T; }; template struct GroupAddSub { using S = T; static constexpr S op(S a, S b) { return a + b; } static constexpr S e() { if constexpr (has_e0_v) return S::e0(); else return S{}; } static constexpr S inv(S a) { return -a; } }; template struct GroupMulDiv { using S = T; static constexpr S op(S a, S b) { return a * b; } static constexpr S e() { if constexpr (has_e1_v) return S::e1(); else return S(1); } static constexpr S inv(S a) { return e() / a; } }; template using RingAddSubMul = RingFromGroupMonoid, MonoidMul>; template using FieldAddSubMulDiv = FieldFromGroupGroup, GroupMulDiv>; // https://github.com/miscalculation53/library/tree/wip/math/linalg/bbla.hpp // https://github.com/miscalculation53/library/tree/wip/math/fps/berlekamp_massey.hpp // https://github.com/miscalculation53/library/tree/wip/math/dot_product.hpp // https://github.com/miscalculation53/library/tree/wip/math/modint/internal_mod32.hpp namespace internal { template struct ordinary_mod32 : false_type {}; template struct ordinary_mod32> : true_type {}; template struct ordinary_mod32> : true_type {}; template struct ordinary_mod32_field : false_type {}; template struct ordinary_mod32_field> : ordinary_mod32 {}; } namespace internal { template using dot_product_mod32_value = ordinary_mod32; template struct dot_product_mod32 : false_type {}; template struct dot_product_mod32> : dot_product_mod32_value {}; template struct dot_product_mod32> : dot_product_mod32_value {}; template struct dot_product_mod32, MonoidMul>> : dot_product_mod32_value {}; template S dot_product_mod32_impl(int n, It1 a, It2 b) { const ull mod = S::mod(); ull sum = 0; for (; n >= block; n -= block) { repi(j, block) sum += ull((a++)->val()) * (b++)->val(); sum %= mod; } repi(j, n) sum += ull((a++)->val()) * (b++)->val(); return S::raw(sum % mod); } } template typename SR::S dot_product(int n, It1 a, It2 b) { using S = typename SR::S; assert(n >= 0); if constexpr (internal::dot_product_mod32::value) { if (S::mod() <= (1 << 30)) return internal::dot_product_mod32_impl<16, S>(n, a, b); return internal::dot_product_mod32_impl<4, S>(n, a, b); } else { S sum = SR::e0(); repi(i, n) sum = SR::add(sum, SR::mul(*a++, *b++)); return sum; } } namespace bbla { } // https://github.com/miscalculation53/library/tree/wip/ds/csr.hpp // https://github.com/miscalculation53/library/tree/wip/math/linalg/row_reduction_mod32.hpp namespace internal { template pair row_reduction_mod32(Rows &a, int n, int m, bool rref) { if (n == 0 || m == 0) return {0, S(1)}; if (none_of(a.begin(), a.end(), [](const auto &row) { return any_of(row.begin(), row.end(), [](S x) { return x.val() != 0; }); })) return {0, S(0)}; const ull mod = S::mod(); const int block = mod <= (1u << 30) ? 16 : 4; vvc w(n, vc(m)); vc used(n); vc first(n); repi(i, n) repi(j, m) w[i][j] = a[i][j].val(); int rk = 0; S de = 1; for (int i = 0, j = 0; i < n && j < m; ++j) { int k = i, next = m; for (; k < n; ++k) { while (first[k] < m && (w[k][first[k]] %= mod) == 0) ++first[k]; if (first[k] == j) break; next = min(next, first[k]); } if (k == n) { de = 0; if (next == m) break; j = next - 1; continue; } if (k != i) { swap(w[k], w[i]); swap(used[k], used[i]); swap(first[k], first[i]); de = -de; } S pivot = S::raw(w[i][j]); de *= pivot; S inv = pivot.inv(); repi(l, j, m) w[i][l] = (S::raw(w[i][l] % mod) * inv).val(); used[i] = 0; repi(k, rref ? 0 : i + 1, n) if (k != i) { if (first[k] > j) continue; ull c = w[k][j] % mod; w[k][j] = 0; if (k > i) first[k] = j + 1; if (c == 0) continue; if (used[k] == block) { repi(l, j + 1, m) w[k][l] %= mod; used[k] = 0; } c = mod - c; repi(l, j + 1, m) w[k][l] += c * w[i][l]; ++used[k]; } ++i; ++rk; } repi(i, n) repi(j, m) a[i][j] = S::raw(first[i] == m ? 0 : w[i][j] % mod); return {rk, de}; } } template struct Matrix : vvc { using S = typename F::S; using M = Matrix; using V = vc; using vvc::vector; using vvc::operator=; Matrix(int n, int m, const S &diag = F::e0(), const S &non_diag = F::e0()) { *this = vvc(n, vc(m, non_diag)); repi(i, min(n, m))(*this)[i][i] = diag; } Matrix(const vvc &a) { *this = a; } template pair shape() const { const int n = (*this).size(); if (n == 0) return {0, 0}; const int m = (*this)[0].size(); return {n, m}; } V operator*(const V &v) const { auto [n, m] = shape(); assert(SZ(v) == m); V res(n, F::e0()); repi(i, n) res[i] = dot_product(m, (*this)[i].begin(), v.begin()); return res; } M operator*(const M &b) const { auto [n, m] = shape(); auto [m_, p] = b.shape(); assert(m == m_); M res(n, p); if constexpr (internal::dot_product_mod32::value) { if (n >= 16 && m >= 16 && p >= 16) { ll nonzero = 0; for (const auto &row : *this) nonzero += m - count(row.begin(), row.end(), F::e0()); if (nonzero == 0) return res; if (nonzero * 2 >= ll(n) * m) { vvc bt(p, vc(m)); repi(k, m) repi(j, p) bt[j][k] = b[k][j]; repi(i, n) repi(j, p) res[i][j] = dot_product(m, (*this)[i].begin(), bt[j].begin()); return res; } } } repi(ii, 0, n, BS) repi(kk, 0, m, BS) repi(jj, 0, p, BS) { repi(i, ii, min(ii + BS, n)) repi(k, kk, min(kk + BS, m)) { S aik = (*this)[i][k]; if (aik == F::e0()) continue; repi(j, jj, min(jj + BS, p)) res[i][j] = F::add(res[i][j], F::mul(aik, b[k][j])); } } return res; } M operator*(const S &x) const { return M(*this) *= x; } template tuple row_reduction(bool rref = false) const { auto [n, m] = shape(); M a(*this); if constexpr (internal::ordinary_mod32_field::value) { if (n >= 32 && m >= 32) { auto [rk, de] = internal::row_reduction_mod32(a, n, m, rref); return {std::move(a), I(rk), de}; } } I rk = 0; S de = F::e1(); for (int i = 0, j = 0; i < n && j < m; j++) { repi(k, i, n) { if (a[k][j] != F::e0()) { swap(a[k], a[i]); if (k != i) de = F::minus(de); break; } } if (a[i][j] == F::e0()) { de = 0; continue; } de = F::mul(de, a[i][j]); S aij_inv = F::inv(a[i][j]); repi(l, j + 1, m) a[i][l] = F::mul(a[i][l], aij_inv); a[i][j] = F::e1(); if (rref) { repi(k, i) { S akj = a[k][j]; if (akj == F::e0()) continue; a[k][j] = F::e0(); repi(l, j + 1, m) a[k][l] = F::add(a[k][l], F::minus(F::mul(a[i][l], akj))); } } repi(k, i + 1, n) { S akj = a[k][j]; if (akj == F::e0()) continue; a[k][j] = F::e0(); repi(l, j + 1, m) a[k][l] = F::add(a[k][l], F::minus(F::mul(a[i][l], akj))); } i++; rk++; } return {a, rk, de}; } }; void main2() { LL(P, N, M); mint::set_mod(P); VEC2(mint, N, M, A_); Matrix> A(A_); auto B = get<0>(A.row_reduction(true)); dump(A, B); vvc C_(N); { ll prv = -1; rep(j, M) { ll tmp = -1; rep(i, N - 1, -1, -1) { if (B[i][j] != 0) { tmp = i; break; } } if (prv != tmp) { // is pivot rep(i, N) C_[i].eb(A[i][j]); } prv = tmp; } } if (C_[0].empty()) rep(i, N) C_[i].eb(0); vvc F_; rep(i, N) if (B[i] != vc(M, 0)) F_.eb(B[i]); if (F_.empty()) F_.eb(vc(M, 0)); decltype(A) C(C_), F(F_); dump(C, F); assert(C * F == A); ll score1 = N * M; auto [CN, CM] = C.shape(); auto [FN, FM] = F.shape(); ll score2 = CN * CM + FN * FM; if (score1 < score2) PRINT(1), PRINT(N, M), PRINTV(A); else PRINT(2), PRINT(CN, CM), PRINTV(C), PRINT(FN, FM), PRINTV(F); } void test() { } // https://github.com/miscalculation53/library/tree/wip/template/template_main.hpp template struct Main { Main() { cauto CERR = [](string val, string color) { string s = "\033[" + color + "m" + val + "\033[m"; }; CERR("\n[FAST_IO]\n\n", "32"); cout << fixed << setprecision(20); init(); CERR("\n[SINGLE_TESTCASE]\n\n", "36"); main2(); } }; Main main_dummy; } int main() {}