// BEGIN: main.cpp #line 1 "main.cpp" // BEGIN: my_template.hpp #line 1 "my_template.hpp" #if defined(USE_PCH) #include #else #if defined(__GNUC__) #include #pragma GCC optimize("Ofast,unroll-loops") // 環境によってはコンパイル成功かつ実行時エラー #pragma GCC target("avx2,popcnt") #endif #include #include using namespace std; using ll = long long; using u8 = uint8_t; using u16 = uint16_t; using u32 = uint32_t; using u64 = uint64_t; using i128 = __int128; using u128 = unsigned __int128; using f128 = __float128; template constexpr bool dependent_false = false; template constexpr T infty = [] { static_assert(dependent_false, "infty is not defined"); return T{}; }(); template <> constexpr int infty = 1'010'000'000; template <> constexpr ll infty = 2'020'000'000'000'000'000; template <> constexpr u32 infty = infty; template <> constexpr u64 infty = infty; template <> constexpr i128 infty = i128(infty) * 2'000'000'000'000'000'000; template <> constexpr double infty = infty; template <> constexpr long double infty = infty; using pi = pair; using vi = vector; template using vc = vector; template using vvc = vector>; template using vvvc = vector>; template using vvvvc = vector>; template using pq_max = priority_queue; template using pq_min = priority_queue, greater>; #define vv(type, name, h, ...) \ vector> name(h, vector(__VA_ARGS__)) #define vvv(type, name, h, w, ...) \ vector>> name( \ h, vector>(w, vector(__VA_ARGS__))) #define vvvv(type, name, a, b, c, ...) \ vector>>> name( \ a, vector>>( \ b, vector>(c, vector(__VA_ARGS__)))) // https://trap.jp/post/1224/ #define FOR1(a) for (ll _ = 0; _ < ll(a); ++_) #define FOR2(i, a) for (ll i = 0; i < ll(a); ++i) #define FOR3(i, a, b) for (ll i = a; i < ll(b); ++i) #define FOR4(i, a, b, c) for (ll i = a; i < ll(b); i += (c)) #define FOR1_R(a) for (ll i = ll(a) - 1; i >= ll(0); --i) #define FOR2_R(i, a) for (ll i = ll(a) - 1; i >= ll(0); --i) #define FOR3_R(i, a, b) for (ll i = ll(b) - 1; i >= ll(a); --i) #define overload4(a, b, c, d, e, ...) e #define overload3(a, b, c, d, ...) d #define FOR(...) overload4(__VA_ARGS__, FOR4, FOR3, FOR2, FOR1)(__VA_ARGS__) #define FOR_R(...) overload3(__VA_ARGS__, FOR3_R, FOR2_R, FOR1_R)(__VA_ARGS__) #define all(x) (x).begin(), (x).end() #define len(x) ll(x.size()) #define elif else if #define eb emplace_back #define mp make_pair #define mt make_tuple #define fi first #define se second #define stoi stoll // require y > 0 template T floor(T x, T y) { return x / y - (x % y < 0); } // require y > 0 template T ceil(T x, T y) { return (x / y) + (x % y > 0); } // require y > 0 template T bmod(T x, T y) { T r = x % y; return (r < 0 ? r + y : r); } // require y > 0 template pair divmod(T x, T y) { T q = x / y, r = x % y; if (r < 0) --q, r += y; return {q, r}; } constexpr auto TEN = [] { array A{}; A[0] = 1; for (int i = 1; i < 20; ++i) A[i] = 10 * A[i - 1]; return A; }(); template T SUM(const U &A) { return std::accumulate(A.begin(), A.end(), T{}); } #define MIN(v) *min_element(all(v)) #define MAX(v) *max_element(all(v)) template inline long long LB(const C &c, const T &x) { return lower_bound(c.begin(), c.end(), x) - c.begin(); } template inline long long UB(const C &c, const T &x) { return upper_bound(c.begin(), c.end(), x) - c.begin(); } #define UNIQUE(x) sort(all(x)), x.erase(unique(all(x)), x.end()) template T POP(deque &que) { T a = que.front(); que.pop_front(); return a; } template T POP(priority_queue &que) { T a = que.top(); que.pop(); return a; } template T POP(vc &que) { T a = que.back(); que.pop_back(); return a; } template i128 binary_search(F check, i128 ok, i128 ng, bool check_ok = true) { if (check_ok) assert(check(ok)); while (1) { i128 x = (ok + ng) / 2; if (x == ok || x == ng) break; (check(x) ? ok : ng) = x; } return ok; } template double binary_search_real(F check, double ok, double ng, int iter = 100) { FOR(iter) { double x = (ok + ng) / 2; (check(x) ? ok : ng) = x; } return (ok + ng) / 2; } template inline bool chmax(T &a, const S &b) { T c = max(a, b); bool changed = (c != a); a = c; return changed; } template inline bool chmin(T &a, const S &b) { T c = min(a, b); bool changed = (c != a); a = c; return changed; } // ? は -1 vc s_to_vi(const string &S, char first_char) { vc A(S.size()); FOR(i, S.size()) { A[i] = (S[i] != '?' ? S[i] - first_char : -1); } return A; } template vc cumsum(const vc &A, int off = 1) { int N = A.size(); vc B(N + 1); FOR(i, N) { B[i + 1] = B[i] + A[i]; } if (off == 0) B.erase(B.begin()); return B; } // stable sort template vc argsort(const vc &A) { vc ids(len(A)); iota(all(ids), 0); sort(all(ids), [&](int i, int j) { return (A[i] == A[j] ? i < j : A[i] < A[j]); }); return ids; } // A[I[0]], A[I[1]], ... template vc rearrange(const vc &A, const vc &I) { vc B(len(I)); FOR(i, len(I)) B[i] = A[I[i]]; return B; } template void concat(vc &first, const Vectors &...others) { first.reserve(first.size() + (others.size() + ... + 0)); (first.insert(first.end(), others.begin(), others.end()), ...); } // i128 template , int> = 0> constexpr i128 abs(T x) { return x < 0 ? -x : x; } constexpr i128 gcd(i128 a, i128 b) { while (b != 0) { i128 c = a % b; a = b, b = c; } return abs(a); } #endif // END: my_template.hpp #line 2 "main.cpp" // BEGIN: other/io.hpp #line 1 "other/io.hpp" #define FASTIO // https://judge.yosupo.jp/submission/21623 namespace fastio { static constexpr uint32_t SZ = 1 << 17; char ibuf[SZ]; char obuf[SZ]; char out[100]; // pointer of ibuf, obuf uint32_t pil = 0, pir = 0, por = 0; bool input_eof = false; template constexpr bool is_signed_integer_v = is_signed_v || is_same_v; 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; [[noreturn]] inline void input_error(const char *message) { fputs(message, stderr); fputc('\n', stderr); exit(EXIT_FAILURE); } 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() { uint32_t n = pir - pil; memmove(ibuf, ibuf + pil, n); pil = 0; pir = n; if (input_eof) return; pir += fread(ibuf + pir, 1, SZ - pir, stdin); if (ferror(stdin)) input_error("fastio: input error"); if (feof(stdin)) { input_eof = true; // Allows the last token to end exactly at EOF without a trailing // whitespace. if (pir < SZ) ibuf[pir++] = '\n'; } } inline char get_char() { if (pil == pir) { load(); if (pil == pir) input_error("fastio: unexpected EOF"); } return ibuf[pil++]; } inline void flush() { fwrite(obuf, 1, por, stdout); por = 0; } void rd(char &c) { do c = get_char(); while (isspace(static_cast(c))); } void rd(string &x) { x.clear(); char c; do c = get_char(); while (isspace(static_cast(c))); do { x += c; c = get_char(); } while (!isspace(static_cast(c))); } template void rd_real(T &x) { string s; rd(s); x = stod(s); } template void rd_integer_slow(T &x) { char c; do c = get_char(); while (c < '-'); bool minus = 0; if constexpr (is_signed_integer_v) { if (c == '-') { minus = 1, c = get_char(); } } x = 0; assert('0' <= c && c <= '9'); while ('0' <= c && c <= '9') { x = x * 10 + (c & 15), c = get_char(); } assert(isspace(static_cast(c))); if constexpr (is_signed_integer_v) { if (minus) x = -x; } } template void rd_integer(T &x) { if (pil + 100 > pir) { load(); if (pil + 100 > pir) { rd_integer_slow(x); return; } } char c; do c = ibuf[pil++]; while (c < '-'); bool minus = 0; if constexpr (is_signed_integer_v) { if (c == '-') { minus = 1, c = ibuf[pil++]; } } x = 0; assert('0' <= c && c <= '9'); while ('0' <= c && c <= '9') { x = x * 10 + (c & 15), c = ibuf[pil++]; } assert(isspace(static_cast(c))); if constexpr (is_signed_integer_v) { if (minus) x = -x; } } template enable_if_t || is_same_v || is_same_v> rd( T &x) { rd_integer(x); } template enable_if_t || is_same_v> rd(T &x) { rd_real(x); } template void rd(pair &p) { rd(p.first), rd(p.second); } template void rd_tuple(T &t) { if constexpr (N < tuple_size::value) { auto &x = get(t); rd(x); rd_tuple(t); } } template void rd(tuple &tpl) { rd_tuple(tpl); } template void rd(array &x) { for (auto &d : x) rd(d); } template void rd(vc &x) { for (auto &d : x) rd(d); } template void read(T &...x) { (rd(x), ...); } inline void wt_range(const char *s, size_t n) { size_t i = 0; while (i < n) { if (por == SZ) flush(); size_t chunk = min(n - i, (size_t)(SZ - por)); memcpy(obuf + por, s + i, chunk); por += chunk; i += chunk; } } void wt(const char c) { if (por == SZ) flush(); obuf[por++] = c; } void wt(const char *s) { wt_range(s, strlen(s)); } void wt(const string &s) { wt_range(s.data(), s.size()); } template void wt_integer(T x) { if (por > SZ - 100) flush(); using U = unsigned_integer_t; U y = static_cast(x); if constexpr (is_signed_integer_v) { if (x < 0) { obuf[por++] = '-'; y = U(0) - y; } } int outi; for (outi = 96; y >= 10000; outi -= 4) { memcpy(out + outi, pre.num[y % 10000], 4); y /= 10000; } if (y >= 1000) { memcpy(obuf + por, pre.num[y], 4); por += 4; } else if (y >= 100) { memcpy(obuf + por, pre.num[y] + 1, 3); por += 3; } else if (y >= 10) { int q = (y * 103) >> 10; obuf[por] = q | '0'; obuf[por + 1] = (y - q * 10) | '0'; por += 2; } else obuf[por++] = y | '0'; memcpy(obuf + por, out + outi + 4, 96 - outi); por += 96 - outi; } template inline void wt_real(T x) { static char buf[1000]; int n = std::snprintf(buf, sizeof(buf), "%.15f", (double)x); wt_range(buf, (size_t)n); } template enable_if_t || is_same_v || is_same_v> wt( T x) { wt_integer(x); } template enable_if_t || is_same_v> wt(T x) { wt_real(x); } inline void wt(bool b) { wt(static_cast('0' + (b ? 1 : 0))); } template void wt(const pair &val) { wt(val.first); wt(' '); wt(val.second); } template void wt_tuple(const T &t) { if constexpr (N < tuple_size::value) { if constexpr (N > 0) wt(' '); wt(get(t)); wt_tuple(t); } } template void wt(const tuple &tpl) { wt_tuple(tpl); } template void wt(const array &val) { auto n = val.size(); for (size_t i = 0; i < n; i++) { if (i) wt(' '); wt(val[i]); } } template void wt(const vector &val) { auto n = val.size(); for (size_t i = 0; i < n; i++) { if (i) wt(' '); wt(val[i]); } } void print() { wt('\n'); } template void print(Head &&head, Tail &&...tail) { wt(forward(head)); ((wt(' '), wt(forward(tail))), ...); wt('\n'); } // gcc expansion. called automaticall after main. void __attribute__((destructor)) _d() { flush(); } } // namespace fastio using fastio::flush; using fastio::print; using fastio::read; #if defined(LOCAL) #define HDR "[DEBUG:", __func__, __LINE__, "]" #define SHOW(...) \ SHOW_IMPL(__VA_ARGS__, SHOW8, SHOW7, SHOW6, SHOW5, SHOW4, SHOW3, SHOW2, \ SHOW1) \ (__VA_ARGS__) #define SHOW_IMPL(_1, _2, _3, _4, _5, _6, _7, _8, NAME, ...) NAME #define SHOW1(x) print(HDR, #x, "=", (x)), flush() #define SHOW2(x, y) print(HDR, #x, "=", (x), #y, "=", (y)), flush() #define SHOW3(x, y, z) \ print(HDR, #x, "=", (x), #y, "=", (y), #z, "=", (z)), flush() #define SHOW4(x, y, z, w) \ print(HDR, #x, "=", (x), #y, "=", (y), #z, "=", (z), #w, "=", (w)), flush() #define SHOW5(x, y, z, w, v) \ print(HDR, #x, "=", (x), #y, "=", (y), #z, "=", (z), #w, "=", (w), #v, "=", \ (v)), \ flush() #define SHOW6(x, y, z, w, v, u) \ print(HDR, #x, "=", (x), #y, "=", (y), #z, "=", (z), #w, "=", (w), #v, "=", \ (v), #u, "=", (u)), \ flush() #define SHOW7(x, y, z, w, v, u, t) \ print(HDR, #x, "=", (x), #y, "=", (y), #z, "=", (z), #w, "=", (w), #v, "=", \ (v), #u, "=", (u), #t, "=", (t)), \ flush() #define SHOW8(x, y, z, w, v, u, t, s) \ print(HDR, #x, "=", (x), #y, "=", (y), #z, "=", (z), #w, "=", (w), #v, "=", \ (v), #u, "=", (u), #t, "=", (t), #s, "=", (s)), \ flush() #else #define SHOW(...) #endif #define INT(...) \ int __VA_ARGS__; \ read(__VA_ARGS__) #define LL(...) \ ll __VA_ARGS__; \ read(__VA_ARGS__) #define U32(...) \ u32 __VA_ARGS__; \ read(__VA_ARGS__) #define U64(...) \ u64 __VA_ARGS__; \ read(__VA_ARGS__) #define STR(...) \ string __VA_ARGS__; \ read(__VA_ARGS__) #define CHAR(...) \ char __VA_ARGS__; \ read(__VA_ARGS__) #define DBL(...) \ double __VA_ARGS__; \ read(__VA_ARGS__) #define VEC(type, name, size) \ vector name(size); \ read(name) #define VV(type, name, h, w) \ vector> name(h, vector(w)); \ read(name) void YES(bool t = 1) { print(t ? "YES" : "NO"); } void NO(bool t = 1) { YES(!t); } void Yes(bool t = 1) { print(t ? "Yes" : "No"); } void No(bool t = 1) { Yes(!t); } void yes(bool t = 1) { print(t ? "yes" : "no"); } void no(bool t = 1) { yes(!t); } void YA(bool t = 1) { print(t ? "YA" : "TIDAK"); } void TIDAK(bool t = 1) { YA(!t); } void Alice(bool t = 1) { print(t ? "Alice" : "Bob"); } void Bob(bool t = 1) { Alice(!t); }// END: other/io.hpp #line 3 "main.cpp" // BEGIN: nt/prime_table.hpp #line 1 "nt/prime_table.hpp" template vc prime_table(int LIM) { ++LIM; const int S = 32768; static int done = 2; static vc primes = {2}, sieve(S + 1); if (done < LIM) { done = LIM; primes = {2}, sieve.assign(S + 1, 0); const int R = LIM / 2; primes.reserve(int(LIM / log(LIM) * 1.1)); vc> cp; for (int i = 3; i <= S; i += 2) { if (!sieve[i]) { cp.eb(i, i * i / 2); for (int j = i * i; j <= S; j += 2 * i) sieve[j] = 1; } } for (int L = 1; L <= R; L += S) { array block{}; for (auto& [p, idx] : cp) for (int i = idx; i < S + L; idx = (i += p)) block[i - L] = 1; FOR(i, min(S, R - L)) if (!block[i]) primes.eb((L + i) * 2 + 1); } } int k = LB(primes, LIM); return {primes.begin(), primes.begin() + k}; } // END: nt/prime_table.hpp #line 5 "main.cpp" // BEGIN: bigint/base.hpp #line 1 "bigint/base.hpp" // BEGIN: poly/convolution.hpp #line 1 "poly/convolution.hpp" // BEGIN: mod/modint.hpp #line 1 "mod/modint.hpp" // BEGIN: mod/modint_common.hpp #line 1 "mod/modint_common.hpp" // BEGIN: other/bit.hpp #line 1 "other/bit.hpp" int popcnt(int x) { return __builtin_popcount(x); } int popcnt(u32 x) { return __builtin_popcount(x); } int popcnt(ll x) { return __builtin_popcountll(x); } int popcnt(u64 x) { return __builtin_popcountll(x); } int popcnt_sgn(int x) { return (__builtin_parity(unsigned(x)) & 1 ? -1 : 1); } int popcnt_sgn(u32 x) { return (__builtin_parity(x) & 1 ? -1 : 1); } int popcnt_sgn(ll x) { return (__builtin_parityll(x) & 1 ? -1 : 1); } int popcnt_sgn(u64 x) { return (__builtin_parityll(x) & 1 ? -1 : 1); } // (0, 1, 2, 3, 4) -> (-1, 0, 1, 1, 2) int topbit(int x) { return (x == 0 ? -1 : 31 - __builtin_clz(x)); } int topbit(u32 x) { return (x == 0 ? -1 : 31 - __builtin_clz(x)); } int topbit(ll x) { return (x == 0 ? -1 : 63 - __builtin_clzll(x)); } int topbit(u64 x) { return (x == 0 ? -1 : 63 - __builtin_clzll(x)); } // (0, 1, 2, 3, 4) -> (-1, 0, 1, 0, 2) int lowbit(int x) { return (x == 0 ? -1 : __builtin_ctz(x)); } int lowbit(u32 x) { return (x == 0 ? -1 : __builtin_ctz(x)); } int lowbit(ll x) { return (x == 0 ? -1 : __builtin_ctzll(x)); } int lowbit(u64 x) { return (x == 0 ? -1 : __builtin_ctzll(x)); } template T kth_bit(int k) { assert(0 <= k && k < int(8 * sizeof(T))); return T(1) << k; } template bool has_kth_bit(T x, int k) { assert(0 <= k && k < int(8 * sizeof(T))); return x >> k & 1; } template struct all_bit { static_assert(is_unsigned::value); UINT s; all_bit(UINT s) : s(s) {} struct iter { UINT s; int operator*() const { return lowbit(s); } void operator++() { s &= s - 1; } bool operator!=(nullptr_t) const { return s; } }; iter begin() const { return {s}; } nullptr_t end() const { return nullptr; } }; template struct all_subset { static_assert(is_unsigned::value); UINT s; all_subset(UINT s) : s(s) {} struct iter { UINT s, t; bool done = false; UINT operator*() const { return t; } void operator++() { done = (t == 0); t = (t - 1) & s; } bool operator!=(nullptr_t) const { return !done; } }; iter begin() const { return {s, s}; } nullptr_t end() const { return nullptr; } }; constexpr u64 full_mask(int n) { assert(0 <= n && n <= 64); return n == 64 ? -1ULL : (1ULL << n) - 1; } u64 bit_reverse(u64 x) { x = ((x & 0x5555555555555555ULL) << 1) | ((x >> 1) & 0x5555555555555555ULL); x = ((x & 0x3333333333333333ULL) << 2) | ((x >> 2) & 0x3333333333333333ULL); x = ((x & 0x0f0f0f0f0f0f0f0fULL) << 4) | ((x >> 4) & 0x0f0f0f0f0f0f0f0fULL); x = ((x & 0x00ff00ff00ff00ffULL) << 8) | ((x >> 8) & 0x00ff00ff00ff00ffULL); x = ((x & 0x0000ffff0000ffffULL) << 16) | ((x >> 16) & 0x0000ffff0000ffffULL); x = (x << 32) | (x >> 32); return x; }// END: other/bit.hpp #line 3 "mod/modint_common.hpp" struct has_mod_impl { template static auto check(T &&x) -> decltype(x.get_mod(), std::true_type{}); template static auto check(...) -> std::false_type; }; template class has_mod : public decltype(has_mod_impl::check(std::declval())) {}; template mint fact(int n) { static const int mod = mint::get_mod(); assert(0 <= n && n < mod); static vector dat = {1, 1}; if (len(dat) <= n) { int now = len(dat); int m = min(mod, 1 << (topbit(n) + 1)); dat.resize(m); FOR(i, now, m) dat[i] = dat[i - 1] * mint::raw(i); } return dat[n]; } template mint fact_inv(int n) { static const int mod = mint::get_mod(); static vector dat = {1, 1}; if (n < 0) return mint(0); if (len(dat) <= n) { int now = len(dat); int m = min(mod, 1 << (topbit(n) + 1)); dat.resize(m); dat[m - 1] = fact(m - 1).inverse(); FOR_R(i, now, m - 1) dat[i] = dat[i + 1] * mint::raw(i + 1); } return dat[n]; } template mint fact_invs(Ts... xs) { return (mint(1) * ... * fact_inv(xs)); } template mint inv(int n) { static const int mod = mint::get_mod(); assert(1 <= n && n < mod); return fact(n - 1) * fact_inv(n); } template <> double inv(int n) { assert(n != 0); return 1.0 / n; } template mint multinomial(Head &&head, Tail &&...tail) { return fact(head) * fact_invs(std::forward(tail)...); } template mint C_dense(int n, int k) { assert(n >= 0); if (k < 0 || n < k) return 0; static vvc C; static int H = 0, W = 0; auto calc = [&](int i, int j) -> mint { if (i == 0) return (j == 0 ? mint(1) : mint(0)); return C[i - 1][j] + (j ? C[i - 1][j - 1] : 0); }; if (W <= k) { FOR(i, H) { C[i].resize(k + 1); FOR(j, W, k + 1) { C[i][j] = calc(i, j); } } W = k + 1; } if (H <= n) { C.resize(n + 1); FOR(i, H, n + 1) { C[i].resize(W); FOR(j, W) { C[i][j] = calc(i, j); } } H = n + 1; } return C[n][k]; } template mint C(ll n, ll k) { assert(n >= 0); if (k < 0 || n < k) return 0; if constexpr (dense) return C_dense(n, k); if constexpr (!large) return multinomial(n, k, n - k); k = min(k, n - k); mint x(1); FOR(i, k) x *= mint(n - i); return x * fact_inv(k); } template mint C_inv(ll n, ll k) { assert(n >= 0); assert(0 <= k && k <= n); if (!large) return fact_inv(n) * fact(k) * fact(n - k); return mint(1) / C(n, k); } // [x^d](1-x)^{-n} template mint C_negative(ll n, ll d) { assert(n >= 0); if (d < 0) return mint(0); if (n == 0) { return (d == 0 ? mint(1) : mint(0)); } return C(n + d - 1, d); }// END: mod/modint_common.hpp #line 2 "mod/modint.hpp" template struct modint { static constexpr u32 umod = u32(mod); static_assert(0 < umod && umod < u32(1) << 31); u32 val; static modint raw(u32 v) { modint x; x.val = v; return x; } constexpr modint() : val(0) {} constexpr modint(u32 x) : val(x % umod) {} constexpr modint(u64 x) : val(x % umod) {} constexpr modint(u128 x) : val(x % umod) {} constexpr modint(int x) : val((x %= mod) < 0 ? x + mod : x){}; constexpr modint(ll x) : val((x %= mod) < 0 ? x + mod : x){}; constexpr modint(i128 x) : val((x %= mod) < 0 ? x + mod : x){}; bool operator<(const modint &other) const { return val < other.val; } modint &operator+=(const modint &p) { if ((val += p.val) >= umod) val -= umod; return *this; } modint &operator-=(const modint &p) { if ((val += umod - p.val) >= umod) val -= umod; return *this; } modint &operator*=(const modint &p) { val = u64(val) * p.val % umod; return *this; } modint &operator/=(const modint &p) { *this *= p.inverse(); return *this; } modint operator-() const { return modint::raw(val ? mod - val : u32(0)); } modint operator+(const modint &p) const { return modint(*this) += p; } modint operator-(const modint &p) const { return modint(*this) -= p; } modint operator*(const modint &p) const { return modint(*this) *= p; } modint operator/(const modint &p) const { return modint(*this) /= p; } bool operator==(const modint &p) const { return val == p.val; } bool operator!=(const modint &p) const { return val != p.val; } modint inverse() const { int a = val, b = mod, u = 1, v = 0, t; while (b > 0) { t = a / b; swap(a -= t * b, b), swap(u -= t * v, v); } return modint(u); } modint pow(ll n) const { if (n < 0) return inverse().pow(-n); assert(n >= 0); modint ret(1), mul(val); while (n > 0) { if (n & 1) ret *= mul; mul *= mul; n >>= 1; } return ret; } static constexpr int get_mod() { return mod; } // (n, r), r は 1 の 2^n 乗根 static constexpr pair ntt_info() { if (mod == 120586241) return {20, 74066978}; if (mod == 167772161) return {25, 17}; if (mod == 469762049) return {26, 30}; if (mod == 754974721) return {24, 362}; if (mod == 880803841) return {23, 211}; if (mod == 943718401) return {22, 663003469}; if (mod == 998244353) return {23, 31}; if (mod == 1004535809) return {21, 582313106}; if (mod == 1012924417) return {21, 368093570}; if (mod == 1224736769) return {24, 1191450770}; if (mod == 2013265921) return {27, 244035102}; return {-1, -1}; } static constexpr bool can_ntt() { return ntt_info().fi != -1; } }; #ifdef FASTIO template void rd(modint &x) { fastio::rd(x.val); x.val %= mod; // assert(0 <= x.val && x.val < mod); } template void wt(modint x) { fastio::wt(x.val); } #endif using modint107 = modint<1000000007>; using modint998 = modint<998244353>; // END: mod/modint.hpp #line 2 "poly/convolution.hpp" // BEGIN: mod/mod_inv.hpp #line 1 "mod/mod_inv.hpp" // long でも大丈夫 // (val * x - 1) が mod の倍数になるようにする // 特に mod=0 なら x=0 が満たす ll mod_inv(ll val, ll mod) { if (mod == 0) return 0; mod = abs(mod); val %= mod; if (val < 0) val += mod; ll a = val, b = mod, u = 1, v = 0, t; while (b > 0) { t = a / b; swap(a -= t * b, b), swap(u -= t * v, v); } if (u < 0) u += mod; return u; } // END: mod/mod_inv.hpp #line 3 "poly/convolution.hpp" // BEGIN: mod/crt3.hpp #line 1 "mod/crt3.hpp" constexpr u32 mod_pow_constexpr(u64 a, u64 n, u32 mod) { a %= mod; u64 res = 1; FOR(32) { if (n & 1) res = res * a % mod; a = a * a % mod, n /= 2; } return res; } template T CRT2(u64 a0, u64 a1) { static_assert(p0 < p1); static constexpr u64 x0_1 = mod_pow_constexpr(p0, p1 - 2, p1); u64 c = (a1 - a0 + p1) * x0_1 % p1; return a0 + c * p0; } template T CRT3(u64 a0, u64 a1, u64 a2) { static_assert(p0 < p1 && p1 < p2); static constexpr u64 x1 = mod_pow_constexpr(p0, p1 - 2, p1); static constexpr u64 x2 = mod_pow_constexpr(u64(p0) * p1 % p2, p2 - 2, p2); static constexpr u64 p01 = u64(p0) * p1; u64 c = (a1 - a0 + p1) * x1 % p1; u64 ans_1 = a0 + c * p0; c = (a2 - ans_1 % p2 + p2) * x2 % p2; return T(ans_1) + T(c) * T(p01); } template T CRT4(u64 a0, u64 a1, u64 a2, u64 a3) { static_assert(p0 < p1 && p1 < p2 && p2 < p3); static constexpr u64 x1 = mod_pow_constexpr(p0, p1 - 2, p1); static constexpr u64 x2 = mod_pow_constexpr(u64(p0) * p1 % p2, p2 - 2, p2); static constexpr u64 x3 = mod_pow_constexpr(u64(p0) * p1 % p3 * p2 % p3, p3 - 2, p3); static constexpr u64 p01 = u64(p0) * p1; u64 c = (a1 - a0 + p1) * x1 % p1; u64 ans_1 = a0 + c * p0; c = (a2 - ans_1 % p2 + p2) * x2 % p2; u128 ans_2 = ans_1 + c * static_cast(p01); c = (a3 - ans_2 % p3 + p3) * x3 % p3; return T(ans_2) + T(c) * T(p01) * T(p2); } template T CRT5(u64 a0, u64 a1, u64 a2, u64 a3, u64 a4) { static_assert(p0 < p1 && p1 < p2 && p2 < p3 && p3 < p4); static constexpr u64 x1 = mod_pow_constexpr(p0, p1 - 2, p1); static constexpr u64 x2 = mod_pow_constexpr(u64(p0) * p1 % p2, p2 - 2, p2); static constexpr u64 x3 = mod_pow_constexpr(u64(p0) * p1 % p3 * p2 % p3, p3 - 2, p3); static constexpr u64 x4 = mod_pow_constexpr(u64(p0) * p1 % p4 * p2 % p4 * p3 % p4, p4 - 2, p4); static constexpr u64 p01 = u64(p0) * p1; static constexpr u64 p23 = u64(p2) * p3; u64 c = (a1 - a0 + p1) * x1 % p1; u64 ans_1 = a0 + c * p0; c = (a2 - ans_1 % p2 + p2) * x2 % p2; u128 ans_2 = ans_1 + c * static_cast(p01); c = static_cast(a3 - ans_2 % p3 + p3) * x3 % p3; u128 ans_3 = ans_2 + static_cast(c * p2) * p01; c = static_cast(a4 - ans_3 % p4 + p4) * x4 % p4; return T(ans_3) + T(c) * T(p01) * T(p23); } // END: mod/crt3.hpp #line 4 "poly/convolution.hpp" // BEGIN: poly/convolution_naive.hpp #line 1 "poly/convolution_naive.hpp" template ::value>::type* = nullptr> vc convolution_naive(const vc& a, const vc& b) { int n = int(a.size()), m = int(b.size()); if (n > m) return convolution_naive(b, a); if (n == 0) return {}; vector ans(n + m - 1); FOR(i, n) FOR(j, m) ans[i + j] += a[i] * b[j]; return ans; } template ::value>::type* = nullptr> vc convolution_naive(const vc& a, const vc& b) { int n = int(a.size()), m = int(b.size()); if (n > m) return convolution_naive(b, a); if (n == 0) return {}; vc ans(n + m - 1); if (n <= 16 && (T::get_mod() < (1 << 30))) { for (int k = 0; k < n + m - 1; ++k) { int s = max(0, k - m + 1); int t = min(n, k + 1); u64 sm = 0; for (int i = s; i < t; ++i) { sm += u64(a[i].val) * (b[k - i].val); } ans[k] = sm; } } else { for (int k = 0; k < n + m - 1; ++k) { int s = max(0, k - m + 1); int t = min(n, k + 1); u128 sm = 0; for (int i = s; i < t; ++i) { sm += u64(a[i].val) * (b[k - i].val); } ans[k] = T::raw(sm % T::get_mod()); } } return ans; } // END: poly/convolution_naive.hpp #line 5 "poly/convolution.hpp" // BEGIN: poly/convolution_karatsuba.hpp #line 1 "poly/convolution_karatsuba.hpp" #line 2 "poly/convolution_karatsuba.hpp" // 任意の環でできる template vc convolution_karatsuba(const vc& f, const vc& g) { const int thresh = 30; if (min(len(f), len(g)) <= thresh) return convolution_naive(f, g); int n = max(len(f), len(g)); int m = ceil(n, 2); vc f1, f2, g1, g2; if (len(f) < m) f1 = f; if (len(f) >= m) f1 = {f.begin(), f.begin() + m}; if (len(f) >= m) f2 = {f.begin() + m, f.end()}; if (len(g) < m) g1 = g; if (len(g) >= m) g1 = {g.begin(), g.begin() + m}; if (len(g) >= m) g2 = {g.begin() + m, g.end()}; vc a = convolution_karatsuba(f1, g1); vc b = convolution_karatsuba(f2, g2); FOR(i, len(f2)) f1[i] += f2[i]; FOR(i, len(g2)) g1[i] += g2[i]; vc c = convolution_karatsuba(f1, g1); vc F(len(f) + len(g) - 1); assert(2 * m + len(b) <= len(F)); FOR(i, len(a)) F[i] += a[i], c[i] -= a[i]; FOR(i, len(b)) F[2 * m + i] += b[i], c[i] -= b[i]; if (c.back() == T(0)) c.pop_back(); FOR(i, len(c)) if (c[i] != T(0)) F[m + i] += c[i]; return F; }// END: poly/convolution_karatsuba.hpp #line 6 "poly/convolution.hpp" // BEGIN: poly/ntt.hpp #line 1 "poly/ntt.hpp" #line 2 "poly/ntt.hpp" template void ntt(vector& a, bool inverse) { assert(mint::can_ntt()); const int rank2 = mint::ntt_info().fi; const u32 mod = mint::get_mod(); static array root, iroot; static array rate2, irate2; static array rate3, irate3; assert(rank2 != -1 && len(a) <= (1 << max(0, rank2))); static bool prepared = 0; if (!prepared) { prepared = 1; root[rank2] = mint::ntt_info().se; iroot[rank2] = mint(1) / root[rank2]; FOR_R(i, rank2) { root[i] = root[i + 1] * root[i + 1]; iroot[i] = iroot[i + 1] * iroot[i + 1]; } mint prod = 1, iprod = 1; for (int i = 0; i <= rank2 - 2; i++) { rate2[i] = root[i + 2] * prod; irate2[i] = iroot[i + 2] * iprod; prod *= iroot[i + 2]; iprod *= root[i + 2]; } prod = 1, iprod = 1; for (int i = 0; i <= rank2 - 3; i++) { rate3[i] = root[i + 3] * prod; irate3[i] = iroot[i + 3] * iprod; prod *= iroot[i + 3]; iprod *= root[i + 3]; } } int n = int(a.size()); int h = topbit(n); assert(n == 1 << h); if (!inverse) { int len = 0; while (len < h) { if (h - len == 1) { int p = 1 << (h - len - 1); mint rot = 1; FOR(s, 1 << len) { int offset = s << (h - len); FOR(i, p) { auto l = a[i + offset]; auto r = a[i + offset + p] * rot; a[i + offset] = l + r; a[i + offset + p] = l - r; } rot *= rate2[topbit(~s & -~s)]; } len++; } else { int p = 1 << (h - len - 2); mint rot = 1, imag = root[2]; for (int s = 0; s < (1 << len); s++) { mint rot2 = rot * rot; mint rot3 = rot2 * rot; int offset = s << (h - len); for (int i = 0; i < p; i++) { u64 mod2 = u64(mod) * mod; u64 a0 = a[i + offset].val; u64 a1 = u64(a[i + offset + p].val) * rot.val; u64 a2 = u64(a[i + offset + 2 * p].val) * rot2.val; u64 a3 = u64(a[i + offset + 3 * p].val) * rot3.val; u64 a1na3imag = (a1 + mod2 - a3) % mod * imag.val; u64 na2 = mod2 - a2; a[i + offset] = a0 + a2 + a1 + a3; a[i + offset + 1 * p] = a0 + a2 + (2 * mod2 - (a1 + a3)); a[i + offset + 2 * p] = a0 + na2 + a1na3imag; a[i + offset + 3 * p] = a0 + na2 + (mod2 - a1na3imag); } rot *= rate3[topbit(~s & -~s)]; } len += 2; } } } else { mint coef = mint(1) / mint(len(a)); FOR(i, len(a)) a[i] *= coef; int len = h; while (len) { if (len == 1) { int p = 1 << (h - len); mint irot = 1; FOR(s, 1 << (len - 1)) { int offset = s << (h - len + 1); FOR(i, p) { u64 l = a[i + offset].val; u64 r = a[i + offset + p].val; a[i + offset] = l + r; a[i + offset + p] = (mod + l - r) * irot.val; } irot *= irate2[topbit(~s & -~s)]; } len--; } else { int p = 1 << (h - len); mint irot = 1, iimag = iroot[2]; FOR(s, (1 << (len - 2))) { mint irot2 = irot * irot; mint irot3 = irot2 * irot; int offset = s << (h - len + 2); for (int i = 0; i < p; i++) { u64 a0 = a[i + offset + 0 * p].val; u64 a1 = a[i + offset + 1 * p].val; u64 a2 = a[i + offset + 2 * p].val; u64 a3 = a[i + offset + 3 * p].val; u64 x = (mod + a2 - a3) * iimag.val % mod; a[i + offset] = a0 + a1 + a2 + a3; a[i + offset + 1 * p] = (a0 + mod - a1 + x) * irot.val; a[i + offset + 2 * p] = (a0 + a1 + 2 * mod - a2 - a3) * irot2.val; a[i + offset + 3 * p] = (a0 + 2 * mod - a1 - x) * irot3.val; } irot *= irate3[topbit(~s & -~s)]; } len -= 2; } } } } // END: poly/ntt.hpp #line 7 "poly/convolution.hpp" template vector convolution_ntt(vector a, vector b) { assert(mint::can_ntt()); if (a.empty() || b.empty()) return {}; int n = int(a.size()), m = int(b.size()); int sz = 1; while (sz < n + m - 1) sz *= 2; // sz = 2^k のときの高速化。分割統治的なやつで損しまくるので。 if ((n + m - 3) <= sz / 2) { auto a_last = a.back(), b_last = b.back(); a.pop_back(), b.pop_back(); auto c = convolution(a, b); c.resize(n + m - 1); c[n + m - 2] = a_last * b_last; FOR(i, len(a)) c[i + len(b)] += a[i] * b_last; FOR(i, len(b)) c[i + len(a)] += b[i] * a_last; return c; } a.resize(sz), b.resize(sz); bool same = a == b; ntt(a, 0); if (same) { b = a; } else { ntt(b, 0); } FOR(i, sz) a[i] *= b[i]; ntt(a, 1); a.resize(n + m - 1); return a; } template vector convolution_garner(const vector& a, const vector& b) { int n = len(a), m = len(b); if (!n || !m) return {}; static constexpr int p0 = 167772161; static constexpr int p1 = 469762049; static constexpr int p2 = 754974721; using mint0 = modint; using mint1 = modint; using mint2 = modint; vc a0(n), b0(m); vc a1(n), b1(m); vc a2(n), b2(m); FOR(i, n) a0[i] = a[i].val, a1[i] = a[i].val, a2[i] = a[i].val; FOR(i, m) b0[i] = b[i].val, b1[i] = b[i].val, b2[i] = b[i].val; auto c0 = convolution_ntt(a0, b0); auto c1 = convolution_ntt(a1, b1); auto c2 = convolution_ntt(a2, b2); vc c(len(c0)); FOR(i, n + m - 1) { c[i] = CRT3(c0[i].val, c1[i].val, c2[i].val); } return c; } vector convolution(vector a, vector b) { int n = len(a), m = len(b); if (!n || !m) return {}; if (min(n, m) <= 2500) return convolution_naive(a, b); ll mi_a = MIN(a), mi_b = MIN(b); for (auto& x : a) x -= mi_a; for (auto& x : b) x -= mi_b; assert(MAX(a) * MAX(b) <= 1e18); auto Ac = cumsum(a), Bc = cumsum(b); vi res(n + m - 1); for (int k = 0; k < n + m - 1; ++k) { int s = max(0, k - m + 1); int t = min(n, k + 1); res[k] += (t - s) * mi_a * mi_b; res[k] += mi_a * (Bc[k - s + 1] - Bc[k - t + 1]); res[k] += mi_b * (Ac[t] - Ac[s]); } static constexpr u32 MOD1 = 1004535809; static constexpr u32 MOD2 = 1012924417; using mint1 = modint; using mint2 = modint; vc a1(n), b1(m); vc a2(n), b2(m); FOR(i, n) a1[i] = a[i], a2[i] = a[i]; FOR(i, m) b1[i] = b[i], b2[i] = b[i]; auto c1 = convolution_ntt(a1, b1); auto c2 = convolution_ntt(a2, b2); FOR(i, n + m - 1) { res[i] += CRT2(c1[i].val, c2[i].val); } return res; } template vc convolution(const vc& a, const vc& b) { // static_assert(!is_same_v>, "use Bit_Array version for mod // 2"); int n = len(a), m = len(b); if (!n || !m) return {}; if (mint::can_ntt()) { if (min(n, m) <= 50) return convolution_karatsuba(a, b); return convolution_ntt(a, b); } if (min(n, m) <= 200) return convolution_karatsuba(a, b); return convolution_garner(a, b); }// END: poly/convolution.hpp #line 2 "bigint/base.hpp" // BEGIN: nt/digit_sum.hpp #line 1 "nt/digit_sum.hpp" int digit_sum(u64 x) { const int K = 100'000; static vc dp(K); if (dp[1] == 0) { FOR(x, 1, K) dp[x] = dp[x / 10] + (x % 10); } int res = 0; while (x) { res += dp[x % K]; x /= K; } return res; }// END: nt/digit_sum.hpp #line 3 "bigint/base.hpp" // 15桁の4素数畳み込みにしたけど高速にならなかった. // https://judge.yosupo.jp/submission/311757 struct BigInteger { static constexpr int LOG = 9; static constexpr int MOD = TEN[LOG]; using bint = BigInteger; int sgn; vc dat; BigInteger() : sgn(0) {} BigInteger(i128 val) { if (val == 0) { sgn = 0; return; } sgn = 1; if (val != 0) { if (val < 0) sgn = -1, val = -val; while (val > 0) { dat.eb(val % MOD), val /= MOD; } } } BigInteger(string s) { assert(!s.empty()); sgn = 1; if (s[0] == '-') { sgn = -1; s.erase(s.begin()); assert(!s.empty()); } { reverse(all(s)); while (len(s) >= 2 && s.back() == '0') s.pop_back(); reverse(all(s)); } if (s[0] == '0') { sgn = 0; return; } reverse(all(s)); int n = len(s); int m = ceil(n, LOG); dat.assign(m, 0); FOR(i, n) { dat[i / LOG] += TEN[i % LOG] * (s[i] - '0'); } } bint &operator=(const bint &p) { sgn = p.sgn, dat = p.dat; return *this; } bool operator<(const bint &p) const { if (sgn != p.sgn) { return sgn < p.sgn; } if (sgn == 0) return false; if (len(dat) != len(p.dat)) { if (sgn == 1) return len(dat) < len(p.dat); if (sgn == -1) return len(dat) > len(p.dat); } FOR_R(i, len(dat)) { if (dat[i] == p.dat[i]) continue; if (sgn == 1) return dat[i] < p.dat[i]; if (sgn == -1) return dat[i] > p.dat[i]; } return false; } bool operator>(const bint &p) const { return p < *this; } bool operator<=(const bint &p) const { return !(*this > p); } bool operator>=(const bint &p) const { return !(*this < p); } bint &operator+=(const bint p) { if (sgn == 0) { return *this = p; } if (p.sgn == 0) return *this; if (sgn != p.sgn) { *this -= (-p); return *this; } int n = max(len(dat), len(p.dat)); dat.resize(n + 1); FOR(i, n) { if (i < len(p.dat)) dat[i] += p.dat[i]; if (dat[i] >= MOD) dat[i] -= MOD, dat[i + 1] += 1; } while (len(dat) && dat.back() == 0) dat.pop_back(); if (dat.empty()) sgn = 0; return *this; } bint &operator-=(const bint p) { if (p.sgn == 0) return *this; if (sgn == 0) return *this = (-p); if (sgn != p.sgn) { *this += (-p); return *this; } if ((sgn == 1 && *this < p) || (sgn == -1 && *this > p)) { *this = p - *this; sgn = -sgn; return *this; } FOR(i, len(p.dat)) { dat[i] -= p.dat[i]; } FOR(i, len(dat) - 1) { if (dat[i] < 0) dat[i] += MOD, dat[i + 1] -= 1; } while (len(dat) && dat.back() == 0) { dat.pop_back(); } if (dat.empty()) sgn = 0; return *this; } bint &operator*=(const bint &p) { sgn *= p.sgn; if (sgn == 0) { dat.clear(); } else { dat = convolve(dat, p.dat); } return *this; } // bint &operator/=(const bint &p) { return *this; } bint operator-() const { bint p = *this; p.sgn *= -1; return p; } bint operator+(const bint &p) const { return bint(*this) += p; } bint operator-(const bint &p) const { return bint(*this) -= p; } bint operator*(const bint &p) const { return bint(*this) *= p; } // bint operator/(const modint &p) const { return modint(*this) /= p; } bool operator==(const bint &p) const { return (sgn == p.sgn && dat == p.dat); } bool operator!=(const bint &p) const { return !((*this) == p); } vc convolve(const vc &a, const vc &b) { int n = len(a), m = len(b); if (!n || !m) return {}; if (min(n, m) <= 500) { vc c(n + m - 1); u128 x = 0; FOR(k, n + m - 1) { int s = max(0, k + 1 - m), t = min(k, n - 1); FOR(i, s, t + 1) { x += u64(a[i]) * b[k - i]; } c[k] = x % MOD, x = x / MOD; } while (x > 0) { c.eb(x % MOD), x = x / MOD; } return c; } static constexpr int p0 = 167772161; static constexpr int p1 = 469762049; static constexpr int p2 = 754974721; using mint0 = modint; using mint1 = modint; using mint2 = modint; vc a0(all(a)), b0(all(b)); vc a1(all(a)), b1(all(b)); vc a2(all(a)), b2(all(b)); auto c0 = convolution_ntt(a0, b0); auto c1 = convolution_ntt(a1, b1); auto c2 = convolution_ntt(a2, b2); vc c(len(c0)); u128 x = 0; FOR(i, n + m - 1) { x += CRT3(c0[i].val, c1[i].val, c2[i].val); c[i] = x % MOD, x = x / MOD; } while (x) { c.eb(x % MOD), x = x / MOD; } return c; } string to_string() { if (dat.empty()) return "0"; string s; for (int x : dat) { FOR(LOG) { s += '0' + (x % 10); x = x / 10; } } while (s.back() == '0') s.pop_back(); if (sgn == -1) s += '-'; reverse(all(s)); return s; } // https://codeforces.com/contest/504/problem/D string to_binary_string() { assert(sgn >= 0); vc A(all(dat)); string ANS; while (1) { while (len(A) && A.back() == u32(0)) POP(A); if (A.empty()) break; u64 rem = 0; FOR_R(i, len(A)) { rem = rem * MOD + A[i]; A[i] = rem >> 32; rem &= u32(-1); } FOR(i, 32) { ANS += '0' + (rem >> i & 1); } } while (len(ANS) && ANS.back() == '0') ANS.pop_back(); reverse(all(ANS)); if (ANS.empty()) ANS += '0'; return ANS; } // https://codeforces.com/contest/759/problem/E pair divmod(int p) { vc after; ll rm = 0; FOR_R(i, len(dat)) { rm = rm * MOD + dat[i]; after.eb(rm / p); rm = rm % p; } reverse(all(after)); while (len(after) && after.back() == 0) POP(after); bint q; q.sgn = sgn; q.dat = after; rm *= sgn; if (rm < 0) { rm += p; q -= 1; } return {q, rm}; } int bmod(int p) { ll rm = 0; FOR_R(i, len(dat)) { rm = (rm * MOD + dat[i]) % p; } rm *= sgn; if (rm < 0) { rm += p; } return rm; } // https://codeforces.com/problemset/problem/582/D vc base_p_representation(int p) { vc A(all(dat)); vc res; while (1) { while (len(A) && A.back() == u32(0)) POP(A); if (A.empty()) break; u64 rm = 0; FOR_R(i, len(A)) { rm = rm * MOD + A[i]; A[i] = rm / p; rm %= p; } res.eb(rm); } reverse(all(res)); return res; } // overflow 無視して計算 ll to_ll() { ll x = 0; FOR_R(i, len(dat)) x = MOD * x + dat[i]; return sgn * x; } // https://codeforces.com/contest/986/problem/D bint pow(ll n) { assert(n >= 0); auto dfs = [&](auto &dfs, ll n) -> bint { if (n == 1) return (*this); bint x = dfs(dfs, n / 2); x *= x; if (n & 1) x *= (*this); return x; }; if (n == 0) return bint(1); return dfs(dfs, n); } // https://codeforces.com/contest/986/problem/D double log10() { assert(!dat.empty() && sgn == 1); if (len(dat) <= 3) { double x = 0; FOR_R(i, len(dat)) x = MOD * x + dat[i]; return std::log10(x); } double x = 0; FOR(i, 4) x = MOD * x + dat[len(dat) - 1 - i]; x = std::log10(x); x += double(LOG) * (len(dat) - 4); return x; } int digit_sum() { int ans = 0; for (auto &x : dat) ans += ::digit_sum(x); // global にある digit_sum return ans; } }; // END: bigint/base.hpp #line 6 "main.cpp" // BEGIN: graph/base.hpp #line 1 "graph/base.hpp" // BEGIN: ds/hashmap.hpp #line 1 "ds/hashmap.hpp" // u64 -> Val template struct HashMap { // n は入れたいものの個数で ok HashMap(u32 n = 0) { build(n); } void build(u32 n) { u32 k = 8; while (k < n * 2) k *= 2; cap = k / 2, mask = k - 1; key.resize(k), val.resize(k), used.assign(k, 0); } // size を保ったまま. size=0 にするときは build すること. void clear() { used.assign(len(used), 0); cap = (mask + 1) / 2; } int size() { return len(used) / 2 - cap; } int index(const u64& k) { int i = 0; for (i = hash(k); used[i] && key[i] != k; i = (i + 1) & mask) { } return i; } Val& operator[](const u64& k) { int i = index(k); if (used[i]) return val[i]; if (cap == 0) extend(), i = index(k); used[i] = 1, key[i] = k, val[i] = Val{}, --cap; return val[i]; } Val get(const u64& k, Val default_value) { int i = index(k); return (used[i] ? val[i] : default_value); } bool count(const u64& k) { int i = index(k); return used[i] && key[i] == k; } // f(key, val) template void enumerate_all(F f) { FOR(i, len(used)) if (used[i]) f(key[i], val[i]); } private: u32 cap, mask; vc key; vc val; vc used; u64 hash(u64 x) { static const u64 FIXED_RANDOM = std::chrono::steady_clock::now().time_since_epoch().count(); x += FIXED_RANDOM; x = (x ^ (x >> 30)) * 0xbf58476d1ce4e5b9; x = (x ^ (x >> 27)) * 0x94d049bb133111eb; return (x ^ (x >> 31)) & mask; } void extend() { vc> dat; dat.reserve(len(used) / 2 - cap); FOR(i, len(used)) { if (used[i]) dat.eb(key[i], val[i]); } build(2 * len(dat)); for (auto& [a, b] : dat) (*this)[a] = b; } };// END: ds/hashmap.hpp #line 2 "graph/base.hpp" template struct Edge { int frm, to; T cost; int id; }; template struct Graph { static constexpr bool is_directed = directed; int N, M; using cost_type = T; using edge_type = Edge; vector edges; vector indptr; vector csr_edges; vc vc_deg, vc_indeg, vc_outdeg; HashMap MP_FOR_EID; bool prepared; class OutgoingEdges { public: OutgoingEdges(const Graph* G, int l, int r) : G(G), l(l), r(r) {} const edge_type* begin() const { if (l == r) { return 0; } return &G->csr_edges[l]; } const edge_type* end() const { if (l == r) { return 0; } return &G->csr_edges[r]; } private: const Graph* G; int l, r; }; bool is_prepared() { return prepared; } Graph() : N(0), M(0), prepared(0) {} Graph(int N) : N(N), M(0), prepared(0) {} void build(int n) { N = n, M = 0; prepared = 0; edges.clear(); indptr.clear(); csr_edges.clear(); vc_deg.clear(); vc_indeg.clear(); vc_outdeg.clear(); MP_FOR_EID.clear(); } void add(int frm, int to, T cost = 1, int i = -1) { assert(!prepared); assert(0 <= frm && frm < N && 0 <= to && to < N); if (i == -1) i = M; auto e = edge_type({frm, to, cost, i}); edges.eb(e); ++M; } #ifdef FASTIO // wt, off void read_tree(bool wt = false, int off = 1) { read_graph(N - 1, wt, off); } void read_graph(int M, bool wt = false, int off = 1) { for (int m = 0; m < M; ++m) { INT(a, b); a -= off, b -= off; if (!wt) { add(a, b); } else { T c; read(c); add(a, b, c); } } build(); } #endif void build() { assert(!prepared); prepared = true; indptr.assign(N + 1, 0); for (auto&& e : edges) { indptr[e.frm + 1]++; if (!directed) indptr[e.to + 1]++; } for (int v = 0; v < N; ++v) { indptr[v + 1] += indptr[v]; } auto counter = indptr; csr_edges.resize(indptr.back() + 1); for (auto&& e : edges) { csr_edges[counter[e.frm]++] = e; if (!directed) csr_edges[counter[e.to]++] = edge_type({e.to, e.frm, e.cost, e.id}); } } OutgoingEdges operator[](int v) const { assert(prepared); return {this, indptr[v], indptr[v + 1]}; } vc deg_array() { if (vc_deg.empty()) calc_deg(); return vc_deg; } pair, vc> deg_array_inout() { if (vc_indeg.empty()) calc_deg_inout(); return {vc_indeg, vc_outdeg}; } int deg(int v) { if (vc_deg.empty()) calc_deg(); return vc_deg[v]; } int in_deg(int v) { if (vc_indeg.empty()) calc_deg_inout(); return vc_indeg[v]; } int out_deg(int v) { if (vc_outdeg.empty()) calc_deg_inout(); return vc_outdeg[v]; } #ifdef FASTIO void debug() { #ifdef LOCAL print("Graph"); if (!prepared) { print("frm to cost id"); for (auto&& e : edges) print(e.frm, e.to, e.cost, e.id); } else { print("indptr", indptr); print("frm to cost id"); FOR(v, N) for (auto&& e : (*this)[v]) print(e.frm, e.to, e.cost, e.id); } flush(); #endif } #endif vc new_idx; vc used_e; // G における頂点 V[i] が、新しいグラフで i になるようにする // {G, es} // sum(deg(v)) の計算量になっていて、 // 新しいグラフの n+m より大きい可能性があるので注意 Graph rearrange(vc V, bool keep_eid = 0) { if (len(new_idx) != N) new_idx.assign(N, -1); int n = len(V); FOR(i, n) new_idx[V[i]] = i; Graph G(n); vc history; FOR(i, n) { for (auto&& e : (*this)[V[i]]) { if (len(used_e) <= e.id) used_e.resize(e.id + 1); if (used_e[e.id]) continue; int a = e.frm, b = e.to; if (new_idx[a] != -1 && new_idx[b] != -1) { history.eb(e.id); used_e[e.id] = 1; int eid = (keep_eid ? e.id : -1); G.add(new_idx[a], new_idx[b], e.cost, eid); } } } FOR(i, n) new_idx[V[i]] = -1; for (auto&& eid : history) used_e[eid] = 0; G.build(); return G; } Graph to_directed_tree(int root = -1) { if (root == -1) root = 0; assert(!is_directed && prepared && M == N - 1); Graph G1(N); vc par(N, -1); auto dfs = [&](auto& dfs, int v) -> void { for (auto& e : (*this)[v]) { if (e.to == par[v]) continue; par[e.to] = v, dfs(dfs, e.to); } }; dfs(dfs, root); for (auto& e : edges) { int a = e.frm, b = e.to; if (par[a] == b) swap(a, b); assert(par[b] == a); G1.add(a, b, e.cost); } G1.build(); return G1; } int get_eid(u64 a, u64 b) { if (len(MP_FOR_EID) == 0) { MP_FOR_EID.build(N - 1); for (auto& e : edges) { u64 a = e.frm, b = e.to; u64 k = to_eid_key(a, b); MP_FOR_EID[k] = e.id; } } return MP_FOR_EID.get(to_eid_key(a, b), -1); } u64 to_eid_key(u64 a, u64 b) { if (!directed && a > b) swap(a, b); return N * a + b; } private: void calc_deg() { assert(vc_deg.empty()); vc_deg.resize(N); for (auto&& e : edges) vc_deg[e.frm]++, vc_deg[e.to]++; } void calc_deg_inout() { assert(vc_indeg.empty()); vc_indeg.resize(N); vc_outdeg.resize(N); for (auto&& e : edges) { vc_indeg[e.to]++, vc_outdeg[e.frm]++; } } }; // END: graph/base.hpp #line 7 "main.cpp" using B = BigInteger; void solve() { B LCM = 1; auto primes = prime_table(300); for (int p : primes) { int q = 1; while (q * p <= 300) { q *= p; LCM = LCM * p; } } // print(len(LCM.to_string())); // unit frac vc frac(301); FOR(b, 1, 301) { frac[b] = LCM.divmod(b).fi; } INT(N, M); vc vis(N); vc dp(N); pq_min> que; auto upd = [&](int v, const B& x) -> void { if (!vis[v]) { vis[v] = 1; dp[v] = x; que.emplace(x, v); return; } if (dp[v] > x) { dp[v] = x; vis[v] = 1; que.emplace(x, v); } }; vc wt(M); Graph G(N); FOR(i, M) { INT(u, v, a, b); wt[i] = LCM.divmod(b).fi * a; --u, --v; G.add(u, v); } G.build(); upd(0, B(0)); while (len(que)) { auto [x, v] = POP(que); if (dp[v] != x) continue; for (auto& e : G[v]) { upd(e.to, x + wt[e.id]); } } FOR(v, 1, N) { B a = dp[v]; B den = 1; for (int p : primes) { int q = 1; while (q * p <= 300) { q *= p; auto [aa, bb] = a.divmod(p); if (bb == 0) { a = aa; } else { den *= p; } } } print(a.to_string(), den.to_string()); } } signed main() { solve(); } // END: main.cpp