#line 1 "No_62_\u30ea\u30d9\u30ea\u30aa\u30f3_Extra.cpp" #define YRSD #line 2 "YRS/all.hpp" #line 2 "YRS/aa/head.hpp" #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #define TE template #define TES template #define Z auto #define ep emplace_back #define eb emplace #define fi first #define se second #define all(x) (x).begin(), (x).end() #define OV4(a, b, c, d, e, ...) e #define FOR1(a) for (int _ = 0; _ < (a); ++_) #define FOR2(i, a) for (int i = 0; i < (a); ++i) #define FOR3(i, a, b) for (int i = (a); i < (b); ++i) #define FOR4(i, a, b, c) for (int i = (a); i < (b); i += (c)) #define FOR(...) OV4(__VA_ARGS__, FOR4, FOR3, FOR2, FOR1)(__VA_ARGS__) #define FOR1_R(a) for (int _ = (a) - 1; _ >= 0; --_) #define FOR2_R(i, a) for (int i = (a) - 1; i >= 0; --i) #define FOR3_R(i, a, b) for (int i = (b) - 1; i >= (a); --i) #define FOR4_R(i, a, b, c) for (int i = (b) - 1; i >= (a); i -= (c)) #define FOR_R(...) OV4(__VA_ARGS__, FOR4_R, FOR3_R, FOR2_R, FOR1_R)(__VA_ARGS__) #define FOR_subset(t, s) for (int t = (s); t > -1; t = (t == 0 ? -1 : (t - 1) & s)) #define sort ranges::sort using namespace std; TE using vc = vector; TE using vvc = vc>; TE using T1 = tuple; TE using T2 = tuple; TE using T3 = tuple; TE using T4 = tuple; TE using max_heap = priority_queue; TE using min_heap = priority_queue, greater>; using u8 = unsigned char; using uint = unsigned int; using ll = long long; using ull = unsigned long long; using ld = long double; using i128 = __int128; using u128 = __uint128_t; using f128 = __float128; using u16 = uint16_t; using PII = pair; using PLL = pair; #ifdef YRSD constexpr bool dbg = 1; #else constexpr bool dbg = 0; #endif #line 2 "YRS/IO/IO.hpp" istream &operator>>(istream &I, i128 &x) { static string s; I >> s; int f = s[0] == '-'; x = 0; const int N = (int)s.size(); FOR(i, f, N) x = x * 10 + s[i] - '0'; if (f) x = -x; return I; } ostream &operator<<(ostream &O, i128 x) { static string s; s.clear(); bool f = x < 0; if (f) x = -x; while (x) s += '0' + x % 10, x /= 10; if (s.empty()) s += '0'; if (f) s += '-'; reverse(all(s)); return O << s; } istream &operator>>(istream &I, f128 &x) { static string s; I >> s, x = stold(s); return I; } ostream &operator<<(ostream &O, const f128 x) { return O << ld(x); } template istream &operator>>(istream &I, tuple &t) { return apply([&I](Z &...s) { ((I >> s), ...); }, t), I; } template istream &operator>>(istream &I, pair &x) { return I >> x.fi >> x.se; } template ostream &operator<<(ostream &O, const pair &x) { return O << x.fi << ' ' << x.se; } TE requires requires(T &c) { begin(c); end(c); } and (not is_same_v, string>) istream &operator>>(istream &I, T &c) { for (Z &e : c) I >> e; return I; } TE requires requires(const T &c) { begin(c); end(c); } and (not is_same_v, const char*>) and (not is_same_v, string>) and (not is_array_v> or not is_same_v>, char>) ostream &operator<<(ostream &O, const T &a) { if (a.empty()) return O; Z i = a.begin(); O << *i++; for (; i != a.end(); ++i) O << ' ' << *i; return O; } void IN() {} TE void IN(T &x, Z &...s) { cin >> x, IN(s...); } void print() { cout << '\n'; } TES void print(T &&x, S &&...y) { cout << x; if constexpr (sizeof...(S)) cout << ' '; print(forward(y)...); } void put() { cout << ' '; } TES void put(T &&x, S &&...y) { cout << x; if constexpr (sizeof...(S)) cout << ' '; put(forward(y)...); } #define INT(...) int __VA_ARGS__; IN(__VA_ARGS__) #define UINT(...) uint __VA_ARGS__; IN(__VA_ARGS__) #define LL(...) ll __VA_ARGS__; IN(__VA_ARGS__) #define ULL(...) ull __VA_ARGS__; IN(__VA_ARGS__) #define I128(...) i128 __VA_ARGS__; IN(__VA_ARGS__) #define STR(...) string __VA_ARGS__; IN(__VA_ARGS__) #define CH(...) char __VA_ARGS__; IN(__VA_ARGS__) #define REAL(...) re __VA_ARGS__; IN(__VA_ARGS__) #define VEC(T, a, n) vc a(n); IN(a) void YES(bool o = 1) { print(o ? "YES" : "NO"); } void Yes(bool o = 1) { print(o ? "Yes" : "No"); } void yes(bool o = 1) { print(o ? "yes" : "no"); } void NO(bool o = 1) { YES(not o); } void No(bool o = 1) { Yes(not o); } void no(bool o = 1) { yes(not o); } void ALICE(bool o = 1) { print(o ? "ALICE" : "BOB"); } void Alice(bool o = 1) { print(o ? "Alice" : "Bob"); } void alice(bool o = 1) { print(o ? "alice" : "bob"); } void BOB(bool o = 1) { ALICE(not o); } void Bob(bool o = 1) { Alice(not o); } void bob(bool o = 1) { alice(not o); } void POSSIBLE(bool o = 1) { print(o ? "POSSIBLE" : "IMPOSSIBLE"); } void Possible(bool o = 1) { print(o ? "Possible" : "Impossible"); } void possible(bool o = 1) { print(o ? "possible" : "impossible"); } void IMPOSSIBLE(bool o = 1) { POSSIBLE(not o); } void Impossible(bool o = 1) { Possible(not o); } void impossible(bool o = 1) { possible(not o); } void TAK(bool o = 1) { print(o ? "TAK" : "NIE"); } void NIE(bool o = 1) { TAK(not o); } #line 5 "YRS/all.hpp" #if (__cplusplus >= 202002L) #include constexpr ld pi = numbers::pi; #endif TE constexpr T inf = numeric_limits::max(); template <> constexpr i128 inf = i128(inf) * 2'000'000'000'000'000'000; template constexpr pair inf> = {inf, inf}; TE constexpr static inline int pc(T x) { return popcount(make_unsigned_t(x)); } constexpr static inline ll len(const Z &a) { return a.size(); } void reverse(Z &a) { reverse(all(a)); } void unique(Z &a) { sort(a); a.erase(unique(all(a)), a.end()); } TE vc inverse(const vc &a) { int N = len(a); vc b(N, -1); FOR(i, N) if (a[i] != -1) b[a[i]] = i; return b; } Z QMAX(const Z &a) { return *max_element(all(a)); } Z QMIN(const Z &a) { return *min_element(all(a)); } constexpr bool chmax(Z &a, const Z &b) { return (a < b ? a = b, 1 : 0); } constexpr bool chmin(Z &a, const Z &b) { return (a > b ? a = b, 1 : 0); } vc argsort(const Z &a) { vc I(len(a)); iota(all(I), 0); sort(I, [&](int i, int k) { return a[i] < a[k] or (a[i] == a[k] and i < k); }); return I; } TE vc rearrange(const vc &a, const vc &I) { int N = len(I); vc b(N); FOR(i, N) b[i] = a[I[i]]; return b; } template vc pre_sum(const vc &a) { int N = len(a); vc c(N + 1); FOR(i, N) c[i + 1] = c[i] + a[i]; if (of == 0) c.erase(c.begin()); return c; } TE constexpr static int topbit(T x) { if (x == 0) return - 1; if constexpr (sizeof(T) <= 4) return 31 - __builtin_clz(x); else return 63 - __builtin_clzll(x); } TE constexpr static int lowbit(T x) { if (x == 0) return -1; if constexpr (sizeof(T) <= 4) return __builtin_ctz(x); else return __builtin_ctzll(x); } TE constexpr T floor(T x, T y) { return x / y - (x % y and (x ^ y) < 0); } TE constexpr T ceil(T x, T y) { return floor(x + y - 1, y); } TE constexpr T bmod(T x, T y) { return x - floor(x, y) * y; } TE constexpr pair divmod(T x, T y) { T q = floor(x, y); return pair{q, x - q * y}; } template T SUM(const Z &v) { return accumulate(all(v), T(0)); } int lb(const Z &a, Z x) { return lower_bound(all(a), x) - a.begin(); } int lb(Z l, Z r, Z x) requires(is_same_v) { return lower_bound(l, r, x) - l; } int ub(const Z &a, Z x) { return upper_bound(all(a), x) - a.begin(); } int ub(Z l, Z r, Z x) requires(is_same_v) { return upper_bound(l, r, x) - l; } template ll bina(Z f, ll l, ll r) { if constexpr (ck) assert(f(l)); while (abs(l - r) > 1) { ll x = (r + l) >> 1; (f(x) ? l : r) = x; } return l; } TE T bina_real(Z f, T l, T r, int c = 100) { while (c--) { T x = (l + r) / 2; (f(x) ? l : r) = x; } return (l + r) / 2; } Z pop(Z &s) { if constexpr (requires { s.pop_back(); }) { Z x = s.back(); return s.pop_back(), x; } else if constexpr (requires { s.top(); }) { Z x = s.top(); return s.pop(), x; } else { Z x = s.front(); return s.pop(), x; } } void setp(int x) { cout << fixed << setprecision(x); } TE inline void sh(vc &a, int N) { a.resize(N, T(0)); } #line 1 "YRS/debug.hpp" #ifdef YRSD void DBG() { cerr << "]" << endl; } TES void DBG(T &&x, S &&...y) { cerr << x; if constexpr (sizeof...(S)) cerr << ", "; DBG(forward(y)...); } #define debug(...) cerr << "[" << __LINE__ << "]: [" #__VA_ARGS__ "] = [", DBG(__VA_ARGS__) void ERR() { cerr << endl; } TES void ERR(T &&x, S &&...y) { cerr << x; if constexpr (sizeof...(S)) cerr << ", "; ERR(forward(y)...); } #define err(...) cerr << "[" << __LINE__ << "]: ", ERR(__VA_ARGS__) #define asser assert #else #define debug(...) void(0721) #define err(...) void(0721) #define asser(...) void(0721) #endif #line 4 "No_62_\u30ea\u30d9\u30ea\u30aa\u30f3_Extra.cpp" // #include "YRS/IO/fast_io.hpp" // #include "YRS/random/rng.hpp" #line 2 "YRS/pr/crt.hpp" #line 2 "YRS/mod/barrett.hpp" struct Barrett { uint m; ull im; explicit Barrett(uint m = 1) : m(m), im(ull(-1) / m + 1) {} uint modulo(ull z) const { if (m == 1) return 0; ull x = ull((u128(z) * im) >> 64); ull y = x * m; return (z - y + (z < y ? m : 0)); } inline uint mul(uint a, uint b) const { return modulo(ull(a) * b); } uint pow(uint a, uint b) const { uint s = 1; for (; b; b >>= 1, a = mul(a, a)) { if (b & 1) s = mul(s, a); } return s - (s >= m ? m : 0); } ull floor(ull z) const { if (m == 1) return z; ull x = (ull)(((u128)z * im) >> 64); ull y = x * m; return (z < y ? x - 1 : x); } pair divmod(ull z) const { if (m == 1) return {z, 0}; ull x = ull((u128(z) * im) >> 64), y = x * m; if (z < y) return {x - 1, z - y + m}; return {x, z - y}; } uint umod() const { return m; } }; struct Barrett_64 { u128 m, mh, ml; explicit Barrett_64(ull mod = 1) : m(mod) { u128 m = u128(-1) / mod; if (m * mod + mod == u128(0)) ++m; mh = m >> 64; ml = m & ull(-1); } ull modulo(u128 x) const { u128 z = (x & ull(-1)) * ml; z = (x & ull(-1)) * mh + (x >> 64) * ml + (z >> 64); z = (x >> 64) * mh + (z >> 64); x -= z * m; return x < m ? x : x - m; } ull mul(ull a, ull b) const { return modulo(u128(a) * b); } ull pow(ull a, ull b) const { ull s = 1; for (; b; b >>= 1, a = mul(a, a)) { if (b & 1) s = mul(s, a); } return s - (s >= m ? m : 0); } ull umod() const { return m; } }; #line 2 "YRS/mod/mod_inv.hpp" TE constexpr T mod_inv(T x, T m) { if (m == 0) return 0; m = abs(m); x %= m; if (x < 0) x += m; T a = x, b = m, 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 += m; return u; } #line 2 "YRS/pr/copr_bas.hpp" // 分解出互质基底 return {bas, dec(a)} dec:= {{id, c}} TE pair, vc>> copr_bas(const vc &a) { vc s; for (T e : a) { vc f; for (T x : s) { if (e == 1) { f.ep(x); continue; } vc dat{e, x}; FOR(p, 1, len(dat)) { FOR(i, p) { while (1) { if (dat[p] > 1 and dat[i] % dat[p] == 0) dat[i] /= dat[p]; else if (dat[i] > 1 and dat[p] % dat[i] == 0) dat[p] /= dat[i]; else break; } T g = gcd(dat[i], dat[p]); if (g == 1 or g == dat[i] or g == dat[p]) continue; dat[i] /= g, dat[p] /= g, dat.ep(g); } } e = dat[0]; int sz = len(dat); FOR(i, 1, sz) if (dat[i] != 1) f.ep(dat[i]); } if (e > 1) f.ep(e); s.swap(f); } sort(s); int N = len(a), sz = len(s); vc> res(N); FOR(i, N) { T x = a[i]; FOR(k, sz) { int t = 0; while (x % s[k] == 0) x /= s[k], ++t; if (t) res[i].ep(k, t); } } return {s, res}; } #line 2 "YRS/pr/factors.hpp" #line 2 "YRS/pr/prims_test.hpp" struct MM { using uu = unsigned __int128; inline static ull m, r, nn; static void set_mod(ull m) { MM::m = m; nn = -uu(m) % m; r = m; FOR(5) r *= 2 - m * r; r = -r; } static ull reduce(uu x) { return (x + uu(ull(x) * r) * m) >> 64; } ull x; MM() : x(0) {} MM(ull x) : x(reduce(uu(x) * nn)) {} ull val() const { ull y = reduce(x); return y >= m ? y - m : y; } MM &operator+=(MM y) { x += y.x - (m << 1); x = (ll(x) < 0 ? x + (m << 1) : x); return *this; } MM &operator-=(MM y) { x -= y.x; x = (ll(x) < 0 ? x + (m << 1) : x); return *this; } MM &operator*=(MM y) { x = reduce(uu(x) * y.x); return *this; } MM operator+(MM y) const { return MM(*this) += y; } MM operator-(MM y) const { return MM(*this) -= y; } MM operator*(MM y) const { return MM(*this) *= y; } bool operator==(MM y) const { return (x >= m ? x - m : x) == (y.x >= m ? y.x - m : y.x); } bool operator!=(MM y) const { return not operator==(y); } MM pow(ull k) const { MM r = 1, a = *this; for (; k; k >>= 1, a *= a) if (k & 1) r *= a; return r; } }; bool primetest(const ull x) { if (x == 2 or x == 3 or x == 5 or x == 7) return 1; if (x % 2 == 0 or x % 3 == 0 or x % 5 == 0 or x % 7 == 0) return 0; if (x < 121) return x > 1; const ull d = (x - 1) >> __builtin_ctzll(x - 1); MM::set_mod(x); const MM o(1), mo(x - 1); Z f = [&](ull a) -> bool { MM y = MM(a).pow(d); ull t = d; while (y != o and y != mo and t != x - 1) y *= y, t <<= 1; if (y != mo and t % 2 == 0) return 1; return 0; }; if (x < (1ull << 32)) { for (ull a : {2, 7, 61}) if (f(a)) return 0; } else { for (ull a : {2, 325, 9'375, 281'78, 450'775, 978'050'4, 179'526'502'2}) { if (x <= a) return 1; if (f(a)) return 0; } } return 1; } ll rho(ll n, ll c) { MM::set_mod(n); const MM cc(c); Z f = [&](MM x) { return x * x + cc; }; MM x = 1, y = 2, z = 1, q = 1; ll g = 1; const ll m = 1ll << (__lg(n) / 5); for (ll r = 1; g == 1; r <<= 1) { x = y; FOR(r) y = f(y); for (ll k = 0; k < r and g == 1; k += m) { z = y; FOR(i, min(m, r - k)) y = f(y), q *= x - y; g = gcd(q.val(), n); } } if (g == n) do { z = f(z); g = gcd((x - z).val(), n); } while (g == 1); return g; } #line 2 "YRS/random/rng.hpp" #include #ifdef MeIoN std::mt19937 rg(0); std::mt19937_64 rd_64(0); #else std::mt19937 rg(std::chrono::steady_clock::now().time_since_epoch().count()); std::mt19937_64 rd_64(std::chrono::steady_clock::now().time_since_epoch().count()); #endif uint rng() { return rg(); } uint rng(uint lim) { return rg() % lim; } int rng(int l, int r) { return l + rg() % (r - l); } ull rng_64() { return rd_64(); } ull rng_64(ull lim) { return rd_64() % lim; } ll rng_64(ll l, ll r) { return l + rd_64() % (r - l); } template void shuffle(vector &v) { const int N = len(v); FOR(i, 1, N) { int k = rng(0, i + 1); if (i != k) swap(v[i], v[k]); } } #line 5 "YRS/pr/factors.hpp" // https://yukicoder.me/problems/no/36 => factor(); ll find_pr_e(ll x) { assert(x > 1); if (primetest(x)) return x; FOR(100) { ll e = rho(x, rng_64(x)); if (primetest(e)) return e; x = e; } err("failed"); assert(0); return -1; } vc> factor(ll x) { assert(x >= 1); vc> r; for (int e = 2; e < 100; ++e) { if (e * e > x) break; if (x % e == 0) { int c = 0; do { x /= e, c += 1; } while (x % e == 0); r.ep(e, c); } } while (x > 1) { ll e = find_pr_e(x); int c = 0; do { x /= e, c += 1; } while (x % e == 0); r.ep(e, c); } return sort(r), r; } vc> factor_by_lpf(ll n, vc &lpf) { vc> s; while (n > 1) { int p = lpf[n], e = 0; while (n % p == 0) n /= p, ++e; s.ep(p, e); } return s; } #line 7 "YRS/pr/crt.hpp" // x = a_i (mod b_i) return min x >= 0 TE bool reduce_by_fac(vc &a, vc &mods) { int N = len(a); unordered_map> mp; FOR(i, N) { for (Z [e, c] : factor(mods[i])) { T mod = 1; FOR(c) mod *= e; T x = a[i] % mod; if (not mp.contains(e)) { mp[e] = {mod, x}; continue; } Z &[m, xx] = mp[e]; if (mod > m) swap(m, mod), swap(x, xx); if (xx % mod != x) return 0; } } mods.clear(), a.clear(); for (Z [p, x] : mp) mods.ep(x.fi), a.ep(x.se); return 1; } TE bool reduce_by_copr(vc &a, vc &mods) { Z [bs, pfs] = copr_bas(mods); int k = len(bs), N = len(a); vc> dat(k, {1, 0}); FOR(i, N) { for (Z [id, c] : pfs[i]) { T mod = 1; FOR(c) mod *= bs[id]; T val = a[i] % mod; Z &[mm, xx] = dat[id]; if (mod > mm) swap(mod, mm), swap(val, xx); if (xx % mod != val) return 0; } } mods.clear(), a.clear(); for (Z x : dat) mods.ep(x.fi), a.ep(x.se); return 1; } i128 crt_reduce(vc a, vc mods, ll mod = -1) { int N = len(a); FOR(i, N) a[i] = ((a[i] %= mods[i]) >= 0 ? a[i] : a[i] + mods[i]); bool ok = (N <= 10 ? reduce_by_copr(a, mods) : reduce_by_fac(a, mods)); N = len(a); if (not ok) return -1; if (N == 0) return 0; vc s(N); FOR(i, N) { Barrett X(mods[i]); ll x = a[i], pr = 1; FOR(k, i) { x = X.modulo(x + s[k] * (mods[i] - pr)); pr = X.mul(pr, mods[k]); } s[i] = X.mul(mod_inv(pr, mods[i]), x); } i128 res = 0, pr = 1; FOR(i, N) { res += pr * s[i], pr *= mods[i]; if (mod != -1) res %= mod, pr %= mod; } return res; } i128 crt_reduce(vc a, vc mods, ll mod = -1) { int N = len(a); FOR(i, N) a[i] = ((a[i] %= mods[i]) >= 0 ? a[i] : a[i] + mods[i]); bool ok = (N <= 10 ? reduce_by_copr(a, mods) : reduce_by_fac(a, mods)); N = len(a); if (not ok) return -1; if (N == 0) return 0; vc s(N); FOR(i, N) { Barrett_64 X(mods[i]); ll x = a[i], pr = 1; FOR(k, i) { x = X.modulo(x + s[k] * (mods[i] - pr)); pr = X.mul(pr, mods[k]); } s[i] = X.mul(mod_inv(pr, mods[i]), x); } i128 res = 0, pr = 1; FOR(i, N) { res += pr * s[i], pr *= mods[i]; if (mod != -1) res %= mod, pr %= mod; } return res; } TE requires (not is_integral_v) T crt_reduce(vc a, vc mods, T mod = -1) { int N = len(a); assert(N <= 10); FOR(i, N) a[i] = ((a[i] %= mods[i]) >= 0 ? a[i] : a[i] + mods[i]); bool ok = reduce_by_copr(a, mods); N = len(a); if (not ok) return -1; if (N == 0) return 0; vc s(N); FOR(i, N) { T x = a[i], pr = 1; FOR(k, i) { x = (x + s[k] * (mods[i] - pr)) % mods[i]; if (x < 0) x += mods[i]; pr = pr * mods[k] % mods[i]; } s[i] = mod_inv(pr, mods[i]) * x % mods[i]; } T res = 0, pr = 1; FOR(i, N) { res += pr * s[i], pr *= mods[i]; if (mod != -1) res %= mod, pr %= mod; } return res; } i128 crt_copr(vc a, vc mods, int mod = -1) { int N = len(a); FOR(i, N) a[i] = ((a[i] %= mods[i]) >= 0 ? a[i] : a[i] + mods[i]); vc s(N); FOR(i, N) { Barrett X(mods[i]); ll x = a[i], pr = 1; FOR(k, i) { x = X.modulo(x + s[k] * (mods[i] - pr)); pr = X.mul(pr, mods[k]); } s[i] = X.mul(mod_inv(pr, mods[i]), x); } i128 res = 0, pr = 1; FOR(i, N) { res += pr * s[i], pr *= mods[i]; if (mod != -1) res %= mod, pr %= pr; } return res; } i128 crt_copr(vc a, vc mods, ll mod = -1) { int N = len(a); FOR(i, N) a[i] = ((a[i] %= mods[i]) >= 0 ? a[i] : a[i] + mods[i]); vc s(N); FOR(i, N) { Barrett_64 X(mods[i]); ll x = a[i], pr = 1; FOR(k, i) { x = X.modulo(x + s[k] * (mods[i] - pr)); pr = X.mul(pr, mods[k]); } s[i] = X.mul(mod_inv(pr, mods[i]), x); } i128 res = 0, pr = 1; FOR(i, N) { res += pr * s[i], pr *= mods[i]; if (mod != -1) res %= mod, pr %= pr; } return res; } #line 2 "YRS/pr/f/line_con_1.hpp" #line 4 "YRS/pr/f/line_con_1.hpp" // ax = b (mod) return t(=== x), mod (reduced) TE pair line_con_1(T a, T b, T M) { assert(M != 0); // a = (a %= M) < 0 ? a + M : a; b = (b %= M) < 0 ? b + M : b; T g = gcd(abs(a), M); if (b % g != 0) return {-1, -1}; a /= g, b /= g, M /= g; T ia = mod_inv(a, M); return {b * ia % M, M}; } #line 2 "YRS/nt/bigint/big.hpp" #line 2 "YRS/mod/mint.hpp" #line 2 "YRS/mod/modint_common.hpp" TE concept is_mint = requires(T x) { { T::get_mod() }; { T::gen(0ull) } -> same_as; x.val; }; TE concept has_const_mod = requires { integral_constant {}; }; TE static vc &invs() { static vc a{0, 1}; return a; } TE static vc &fac() { static vc a{1, 1}; return a; } TE static vc &ifac() { static vc a{1, 1}; return a; } TE static int Set_inv(int N) { static vc &inv = invs(); if (len(inv) >= N) return N; inv.resize(N + 1); inv[0] = 1, inv[1] = 1; FOR(i, 1, N) inv[i + 1] = inv[i] * i; T t = pop(inv).inv(); FOR_R(i, N) inv[i] *= t, t *= i; return N; } TE static int Set_comb(int N) { static vc &fa = fac(), &ifa = ifac(); if (len(fa) >= N) return N; fa.resize(N); ifa.resize(N); FOR(i, 1, N) fa[i] = fa[i - 1] * i; ifa[N - 1] = fa[N - 1].inv(); FOR_R(i, N - 1) ifa[i] = ifa[i + 1] * (i + 1); return N; } template mint inv(int n) { static const int mod = mint::get_mod(); static vc &a = invs(); assert(0 <= n); while (len(a) <= n) { int k = len(a); int q = (mod + k - 1) / k; int r = k * q - mod; a.ep(a[r] * mint(q)); } return a[n]; } template mint fact(int n) { static const int mod = mint::get_mod(); static vc &a = fac(); assert(0 <= n); if (n >= mod) return 0; while (len(a) <= n) { int k = len(a); a.ep(a[k - 1] * mint(k)); } return a[n]; } template mint fact_inv(int n) { static vc &a = ifac(); if (n < 0) return mint(0); while (len(a) <= n) a.ep(a[len(a) - 1] * inv(len(a))); return a[n]; } template mint fact_invs(Ts... xs) { return (mint(1) * ... * fact_inv(xs)); } template mint multinomial(X&& a, S&&... b) { return fact(a) * fact_invs(forward(b)...); } template mint C_dense(int n, int k) { assert(n >= 0); if (k < 0 or n < k) return 0; static vc> C; static int H = 0, W = 0; Z 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 (int i = 0; i < H; ++i) { C[i].resize(k + 1); for (int j = W; j < k + 1; ++j) { C[i][j] = calc(i, j); } } W = k + 1; } if (H <= n) { C.resize(n + 1); for (int i = H; i < n + 1; ++i) { C[i].resize(W); for (int j = 0; j < W; ++j) { C[i][j] = calc(i, j); } } H = n + 1; } return C[n][k]; } template mint C(int N, int K) { assert(N >= 0); if (K < 0 or N < K) return 0; return fact(N) * fact_inv(K) * fact_inv(N - K); } template mint lucas(ll N, ll K) { static constexpr int P = mint::get_mod(); if (K > N) return 0; if (K == 0) return 1; return C(N % P, K % P) * lucas(N / P, K / P); } template mint binom(ll n, ll k) { assert(n >= 0); if (k < 0 or n < k) return 0; if constexpr (dense) return C_dense(n, k); if constexpr (not 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 and k <= n); if (not large) return fact_inv(n) * fact(k) * fact(n - k); return mint(1) / binom(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 binom(n + d - 1, d); } #define CC C #define fac fact #define ifac fact_inv #define set_comb Set_comb #define set_inv Set_inv #line 4 "YRS/mod/mint.hpp" #define C constexpr template struct mint_t { using mint = mint_t; static C uint m = mod; uint x; C uint val() const { return x; } C mint_t() : x(0) {} C mint_t(uint x) : x(x % m) {} C mint_t(ull x) : x(x % m) {} C mint_t(u128 x) : x(x % m) {} C mint_t(int x) : x((x %= mod) < 0 ? x + mod : x) {} C mint_t(ll x) : x((x %= mod) < 0 ? x + mod : x) {} C mint_t(i128 x) : x((x %= mod) < 0 ? x + mod : x) {} C mint &operator+=(mint p) { if ((x += p.x) >= m) x -= m; return *this; } C mint &operator-=(mint p) { if ((x += m - p.x) >= m) x -= m; return *this; } C mint operator+(mint p) const { return mint(*this) += p; } C mint operator-(mint p) const { return mint(*this) -= p; } C mint &operator*=(mint p) { x = ull(x) * p.x % m; return *this; } C mint operator*(mint p) const { return mint(*this) *= p; } C mint &operator/=(mint p) { return *this *= p.inv(); } C mint operator/(mint p) const { return mint(*this) /= p; } C mint operator-() const { return mint::gen(x ? mod - x : 0); } C mint inv() const { int a = x, b = mod, x = 1, y = 0; while (b > 0) { int t = a / b; swap(a -= t * b, b); swap(x -= t * y, y); } return mint(x); } C mint pow(ll k) const { if (k < 0) return inv().pow(-k); mint s(1), a(x); for (; k; k >>= 1, a *= a) if (k & 1) s *= a; return s; } C bool operator<(mint p) const { return x < p.x; } C bool operator==(mint p) const { return x == p.x; } C bool operator!=(mint p) const { return x != p.x; } static C mint gen(uint x) { mint s; s.x = x; return s; } friend istream &operator>>(istream &cin, mint &p) { ll t; cin >> t; p = t; return cin; } friend ostream &operator<<(ostream &cout, mint p) { return cout << p.x; } static C int get_mod() { return mod; } static C PII ntt_info() { if (mod == 167772161) return {25, 17}; if (mod == 469762049) return {26, 30}; if (mod == 754974721) return {24, 362}; if (mod == 998244353) return {23, 31}; if (mod == 120586241) return {20, 74066978}; if (mod == 880803841) return {23, 211}; if (mod == 943718401) return {22, 663003469}; if (mod == 1004535809) return {21, 582313106}; if (mod == 1012924417) return {21, 368093570}; return {-1, -1}; } static C bool can_ntt() { return ntt_info().fi != -1; } }; #undef C using M99 = mint_t<998244353>; using M17 = mint_t<1000000007>; #ifdef FIO template void rd(mint_t &x) { LL(y); x = y; } template void wt(mint_t x) { wt(x.x); } #endif #line 2 "YRS/po/convolution.hpp" #line 2 "YRS/po/c/ntt.hpp" #line 4 "YRS/po/c/ntt.hpp" template void ntt(vc &a, bool in) { assert(mint::can_ntt()); const int p = mint::ntt_info().fi; const uint m = mint::get_mod(); static array r, ir, ra, ira, rat, irat; assert(p != -1 and len(a) <= (1 << max(0, p))); static bool ok = 0; if (not ok) { ok = 1; r[p] = mint::ntt_info().se; ir[p] = mint(1) / r[p]; FOR_R(i, p) { r[i] = r[i + 1] * r[i + 1]; ir[i] = ir[i + 1] * ir[i + 1]; } mint s = 1, in = 1; FOR(i, p - 1) { ra[i] = r[i + 2] * s; ira[i] = ir[i + 2] * in; s *= ir[i + 2]; in *= r[i + 2]; } s = 1, in = 1; FOR(i, p - 2) { rat[i] = r[i + 3] * s; irat[i] = ir[i + 3] * in; s *= ir[i + 3]; in *= r[i + 3]; } } int N = len(a), n = topbit(N); if (not in) { int sz = 0; while (sz < n) { if (n - sz == 1) { int p = 1 << (n - sz - 1); mint c = 1; FOR(s, 1 << sz) { int of = s << (n - sz); FOR(i, p) { mint l = a[i + of], r = a[i + of + p] * c; a[i + of] = l + r, a[i + of + p] = l - r; } c *= ra[topbit(~s & -~s)]; } ++sz; } else { int p = 1 << (n - sz - 2); mint c = 1, in = r[2]; FOR(s, 1 << sz) { mint r2 = c * c, r3 = r2 * c; int of = s << (n - sz); FOR(i, p) { const ull mm = ull(m) * m; ull a0 = a[i + of].val(), a1 = ull(a[i + of + p].val()) * c.val(); ull aa = ull(a[i + of + 2 * p].val()) * r2.val(); ull bb = ull(a[i + of + 3 * p].val()) * r3.val(); ull t = (a1 + mm - bb) % m * in.val(); ull na = mm - aa; a[i + of] = a0 + a1 + aa + bb; a[i + of + p] = a0 + aa + mm * 2 - a1 - bb; a[i + of + 2 * p] = a0 + na + t; a[i + of + 3 * p] = a0 + na + mm - t; } c *= rat[topbit(~s & -~s)]; } sz += 2; } } } else { mint c = mint(1) / mint(N); FOR(i, N) a[i] *= c; int sz = n; while (sz) { if (sz == 1) { int p = 1 << (n - sz); mint c = 1; FOR(s, 1 << (sz - 1)) { int of = s << (n - sz + 1); FOR(i, p) { ull l = a[i + of].val(), r = a[i + of + p].val(); a[i + of] = l + r; a[i + of + p] = (m + l - r) * c.val(); } c *= ira[topbit(~s & -~s)]; } --sz; } else { int p = 1 << (n - sz); mint c = 1, in = ir[2]; FOR(s, 1 << (sz - 2)) { mint r2 = c * c, r3 = r2 * c; int of = s << (n - sz + 2); FOR(i, p) { ull a0 = a[i + of].val(), a1 = a[i + of + p].val(); ull aa = a[i + of + 2 * p].val(); ull bb = a[i + of + 3 * p].val(); ull x = (m + aa - bb) * in.val() % m; a[i + of] = a0 + a1 + aa + bb; a[i + of + p] = (a0 + m - a1 + x) * c.val(); a[i + of + 2 * p] = (a0 + a1 + 2 * m - aa - bb) * r2.val(); a[i + of + 3 * p] = (a0 + 2 * m - a1 - x) * r3.val(); } c *= irat[topbit(~s & -~s)]; } sz -= 2; } } } } #line 2 "YRS/mod/crt3.hpp" constexpr uint pw_c(ull a, ull b, uint mod) { a %= mod; ull res = 1; FOR(32) { if (b & 1) res = res * a % mod; a = a * a % mod, b >>= 1; } return res; } template T crt(ull a0, ull a1) { static_assert(p0 < p1); static constexpr ull x0_1 = pw_c(p0, p1 - 2, p1); ull c = (a1 - a0 + p1) * x0_1 % p1; return a0 + c * p0; } template T crt(ull a0, ull a1, ull a2) { static_assert(p0 < p1 and p1 < p2); static constexpr ull x1 = pw_c(p0, p1 - 2, p1); static constexpr ull x2 = pw_c(ull(p0) * p1 % p2, p2 - 2, p2); static constexpr ull p01 = ull(p0) * p1; ull c = (a1 - a0 + p1) * x1 % p1; ull ans_1 = a0 + c * p0; c = (a2 - ans_1 % p2 + p2) * x2 % p2; return T(ans_1) + T(c) * T(p01); } template T crt(ull a0, ull a1, ull a2, ull a3) { static_assert(p0 < p1 and p1 < p2 and p2 < p3); static constexpr ull x1 = pw_c(p0, p1 - 2, p1); static constexpr ull x2 = pw_c(ull(p0) * p1 % p2, p2 - 2, p2); static constexpr ull x3 = pw_c(ull(p0) * p1 % p3 * p2 % p3, p3 - 2, p3); static constexpr ull p01 = ull(p0) * p1; ull c = (a1 - a0 + p1) * x1 % p1; ull ans_1 = a0 + c * p0; c = (a2 - ans_1 % p2 + p2) * x2 % p2; u128 ans_2 = ans_1 + c * u128(p01); c = (a3 - ans_2 % p3 + p3) * x3 % p3; return T(ans_2) + T(c) * T(p01) * T(p2); } template T crt(ull a0, ull a1, ull a2, ull a3, ull a4) { static_assert(p0 < p1 and p1 < p2 and p2 < p3 and p3 < p4); static constexpr ull x1 = pw_c(p0, p1 - 2, p1); static constexpr ull x2 = pw_c(ull(p0) * p1 % p2, p2 - 2, p2); static constexpr ull x3 = pw_c(ull(p0) * p1 % p3 * p2 % p3, p3 - 2, p3); static constexpr ull x4 = pw_c(ull(p0) * p1 % p4 * p2 % p4 * p3 % p4, p4 - 2, p4); static constexpr ull p01 = ull(p0) * p1; static constexpr ull p23 = ull(p2) * p3; ull c = (a1 - a0 + p1) * x1 % p1; ull ans_1 = a0 + c * p0; c = (a2 - ans_1 % p2 + p2) * x2 % p2; u128 ans_2 = ans_1 + c * u128(p01); c = ull(a3 - ans_2 % p3 + p3) * x3 % p3; u128 ans_3 = ans_2 + u128(c * p2) * p01; c = ull(a4 - ans_3 % p4 + p4) * x4 % p4; return T(ans_3) + T(c) * T(p01) * T(p23); } #line 5 "YRS/po/convolution.hpp" template vc conv_naive(const vc &a, const vc &b) { int N = len(a), M = len(b), sz = N + M - 1; if (not N or not M) return {}; if (N > M) return conv_naive(b, a); vc c(sz); FOR(i, N) FOR(k, M) c[i + k] += a[i] * b[k]; return c; } template vc conv_ntt(vc a, vc b) { assert(mint::can_ntt()); if (a.empty() or b.empty()) return {}; int N = len(a), M = len(b), sz = 1; while (sz < N + M - 1) sz <<= 1; sh(a, sz), sh(b, sz); bool ok = a == b; ntt(a, 0); if (ok) b = a; else ntt(b, 0); FOR(i, sz) a[i] *= b[i]; ntt(a, 1); sh(a, N + M - 1); return a; } template vc conv_mtt(const vc &a, const vc &b) { int N = len(a), M = len(b); if (not N or not M) return {}; static constexpr int p0 = 167772161; static constexpr int p1 = 469762049; static constexpr int p2 = 754974721; using M0 = mint_t; using M1 = mint_t; using M2 = mint_t; 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(); vc c0 = conv_ntt(a0, b0); vc c1 = conv_ntt(a1, b1); vc c2 = conv_ntt(a2, b2); vc c(len(c0)); FOR(i, N + M - 1) c[i] = crt(c0[i].val(), c1[i].val(), c2[i].val()); return c; } template vc convolution(const vc &a, const vc &b) { int N = len(a), M = len(b); if (not N or not M) return {}; if (min(N, M) <= 30) return conv_naive(a, b); if (mint::can_ntt()) return conv_ntt(a, b); return conv_mtt(a, b); } #line 5 "YRS/nt/bigint/big.hpp" // https://www.luogu.com.cn/problem/P2152 高精度gcd struct bigint { static constexpr int TEN[] {1, 10, 100, 1000, 10000, 100000, 1000000, 10000000, 100000000, 1000000000}; static constexpr int LOG = 9, mod = TEN[LOG]; using T = bigint; int op; vc a; bigint() : op(0), a() {} bigint(int op, const vc &a) : op(op), a(a) {} bigint(ll x) { if (x == 0) { op = 0; return; } op = 1; if (x < 0) op = -1, x = -x; while (x) a.ep(x % mod), x /= mod; } bigint(string s) { if (s[0] == '0') { op = 0; return; } op = 1; if (s[0] == '-') op = -1, s.erase(s.begin()); reverse(s); int N = len(s), M = ceil(N, LOG); a.assign(M, 0); FOR(i, N) a[i / LOG] += TEN[i % LOG] * (s[i] - '0'); } bool operator<(const T &p) const { if (op != p.op) return op < p.op; if (op == 0) return 0; if (op == 1) return less(a, p.a); else return less(p.a, a); } bool operator>(const T &p) const { return p < *this; } bool operator<=(const T &p) const { return not(*this > p); } bool operator>=(const T &p) const { return not(*this < p); } bool operator==(const T &p) const { return op == p.op and a == p.a; } bool operator!=(const T &p) const { return op != p.op or a != p.a; } T operator-() const { T p = *this; return p.op = -op, p; } T operator+() const { return *this; } T &operator+=(const T &p) { if (op == 0) return *this = p; if (p.op == 0) return *this; if (op == p.op) return a = add(a, p.a), *this; if (less(a, p.a)) return a = sub(p.a, a), op = -op, *this; a = sub(a, p.a); if (is_zero(a)) op = 0; return *this; } T &operator-=(const T &p) { if (p.op == 0) return *this; if (op == 0) return *this = -p; if (op != p.op) return a = add(a, p.a), *this; if (less(a, p.a)) return a = sub(p.a, a), op = -op, *this; a = sub(a, p.a); if (is_zero(a)) op = 0; return *this; } T &operator*=(const T &p) { op *= p.op; if (not op) a.clear(); else a = mul(a, p.a); return *this; } T &operator/=(const T &p) { return *this = divmod(p).fi; } T &operator%=(const T &p) { return *this = divmod(p).se; } T operator+(const T &p) const { return T(*this) += p; } T operator-(const T &p) const { return T(*this) -= p; } T operator*(const T &p) const { return T(*this) *= p; } T operator/(const T &p) const { return T(*this) /= p; } T operator%(const T &p) const { return T(*this) %= p; } pair divmod(const T &p) const { assert(p.op != 0); if (op == 0) return {T(), T()}; Z res = divmod_newton(a, p.a); int op1 = op * p.op, op2 = op; if (is_zero(res.fi)) op1 = 0; if (is_zero(res.se)) op2 = 0; return {{op1, res.fi}, {op2, res.se}}; } T pow(ll k) const { return pow(*this, k); } static T pow(T a, ll k) { if (k == 0) return 1; T ls = pow(a, k >> 1), res = ls * ls; return (k & 1) ? res * a : res; } string to_string() const { if (not op) return "0"; string s = to_string(a); if (op == -1) s += '-'; return reverse(s), s; } string to_binary_string() const { assert(op != -1); vc A(all(a)); string s; while (1) { while (not A.empty() and A.back() == uint(0)) pop(A); if (A.empty()) break; ull r = 0; int N = len(A); FOR_R(i, N) { r = r * mod + A[i]; A[i] = r >> 32; r &= uint(-1); } FOR(i, 32) s += '0' + (r >> i & 1); } while (not s.empty() and s.back() == '0') pop(s); if (s.empty()) s += '0'; return reverse(s), s; } ll to_ll() const { if (op == 0) return 0; ll x = to_ll(a); return op == -1 ? -x : x; } i128 to_i128() const { if (op == 0) return 0; i128 x = to_i128(a); return op == -1 ? -x : x; } friend ostream &operator<<(ostream &cout, const T &b) { return cout << b.to_string(); } friend istream &operator>>(istream &cin, T &b) { static string s; return cin >> s, b = s, cin; } bool is_zero() const { return op == 0; } bool is_one() const { return op == 1 and len(a) == 1 and a[0] == 1; } bool is_odd() const { return op != 0 and (a[0] & 1); } bool is_even() const { return not is_odd(); } T div2() const { T r = *this; int N = len(a); for (int i = N; i--; r.a[i] >>= 1) if ((r.a[i] & 1) and i) r.a[i - 1] += mod; sh(r.a); if (r.a.empty()) r.op = 0; return r; } T gcd(const T &x) const { T a = this->abs(), b = x.abs(); if (a < b) swap(a, b); if (b.is_zero()) return a; int t = 0; while (a.is_even() and b.is_even()) a = a.div2(), b = b.div2(), ++t; while (b > 0) { if (a.is_even()) a = a.div2(); else if (b.is_even()) b = b.div2(); else a -= b; if (a < b) swap(a, b); } while (t--) a += a; return a; } T lcm(const T &x) const { return *this / gcd(x) * x; } T abs() const { if (op == 0) return 0; if (op < 0) return -(*this); return *this; } private: using vec = vc; using PVV = pair; static vc mul(const vec &a, const vec &b) { int N = len(a), M = len(b); if (not N or not M) return {}; if (min(N, M) <= 500) { vec 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 += ull(a[i]) * b[k - i]; c[k] = x % mod; x /= mod; } while (x > 0) c.ep(x % mod), x /= mod; return c; } static constexpr int p0 = 167772161, p1 = 469762049, p2 = 754974721; vc> a0(all(a)), b0(all(b)); vc> a1(all(a)), b1(all(b)); vc> a2(all(a)), b2(all(b)); Z c0 = conv_ntt(a0, b0); Z c1 = conv_ntt(a1, b1); Z c2 = conv_ntt(a2, b2); vec c(len(c0)); u128 x = 0; FOR(i, N + M - 1) { x += crt(c0[i].val(), c1[i].val(), c2[i].val()); c[i] = x % mod, x = x / mod; } while (x) c.ep(x % mod), x /= mod; return c; } static bool is_zero(const vec &a) { return a.empty(); } static bool is_one(const vec &a) { return len(a) == 1 and a[0] == 1; } static bool eq(const vec &a, const vec &b) { return a == b; } static bool less(const vec &a, const vec &b) { if (len(a) != len(b)) return len(a) < len(b); int N = len(a); FOR_R(i, N) if (a[i] != b[i]) return a[i] < b[i]; return 0; } static bool greater(const vec &a, const vec &b) { return less(b, a); } static bool less_eq(const vec &a, const vec &b) { return not greater(a, b); } static bool greater_eq(const vec &a, const vec &b) { return not less(a, b); } static void sh(vec &a) { while (not a.empty() and a.back() == 0) pop(a); } static vec to_vec(ll x) { vec s; while (x) s.ep(x % mod), x /= mod; return s; } static ll to_ll(const vec &a) { ll s = 0; int N = len(a); FOR_R(i, N) s = s * mod + a[i]; return s; } static i128 to_i128(const vec &a) { i128 s = 0; int N = len(a); FOR_R(i, N) s = s * mod + a[i]; return s; } static string to_string(const vec &a) { string s; for (int x : a) FOR(LOG) s += '0' + x % 10, x /= 10; while (s.back() == '0') pop(s); return s; } static vec add(const vec &a, const vec &b) { vc c(all(a)); c.resize(max(len(a), len(b)) + 1); int N = len(b); FOR(i, N) c[i] += b[i]; N = len(c) - 1; FOR(i, N) if (c[i] >= mod) c[i] -= mod, ++c[i + 1]; return sh(c), c; } static vec sub(const vec &a, const vec &b) { vc c(all(a)); int N = len(b); FOR(i, N) c[i] -= b[i]; N = len(a) - 1; FOR(i, N) if (c[i] < 0) c[i] += mod, --c[i + 1]; return sh(c), c; } // 0 <= A < 1e18, 1 <= B > 1e9 static PVV divmod_ll_int(const vec &a, const vec &b) { assert(0 <= len(a) and len(a) <= 2); assert(len(b) == 1); ll x = to_ll(a); int y = b[0]; return {to_vec(x / y), to_vec(x % y)}; } // 0 <= A < 1e18, 1 <= B < 1e18 static PVV divmod_ll_ll(const vec &a, const vec &b) { assert(0 <= len(a) and len(a) <= 2); assert(1 <= len(b) and len(b) <= 2); ll x = to_ll(a), y = to_ll(b); return {to_vec(x / y), to_vec(x % y)}; } // 1 <= B < 1e9 static PVV divmod_1e9(const vec &a, const vec &b) { assert(len(b) == 1); if (len(a) <= 2) return divmod_ll_int(a, b); int N = len(a); vc s(N); ll d = 0; int bb = b[0]; FOR_R(i, N) { d = d * mod + a[i]; assert(d <= 1ll * mod * bb); int q = d / bb, r = d % bb; s[i] = q, d = r; } return sh(s), pair{s, d ? vc{int(d)} : vc{}}; } // 0 <= A, 1 <= B static PVV divmod_naive(const vec &a, const vec &b) { assert(not is_zero(b)); if (len(b) == 1) return divmod_1e9(a, b); if (max(len(a), len(b)) <= 2) return divmod_ll_ll(a, b); if (less(a, b)) return {{}, a}; // B >= 1e9, A >= B int norm = mod / (b.back() + 1); vec x = mul(a, {norm}), y = mul(b, {norm}); int yb = y.back(); vec s(len(x) - len(y) + 1); vec r(x.end() - len(y), x.end()); int N = len(s); FOR_R(i, N) { if (len(r) < len(y)); else if (len(r) == len(y)) { if (less_eq(y, r)) s[i] = 1, r = sub(r, y); } else { assert(len(y) + 1 == len(r)); ll rb = 1ll * r.back() * mod + r.end()[-2]; int q = rb / yb; vec yq = mul(y, {q}); while (less(r, yq)) --q, yq = sub(yq, y); r = sub(r, yq); while (less_eq(y, r)) ++q, r = sub(r, y); s[i] = q; } if (i) r.insert(r.begin(), x[i - 1]); } sh(s), sh(r); Z [ss, rr] = divmod_1e9(r, {norm}); assert(is_zero(rr)); return {s, ss}; }; // 1 / a を 絶対誤差 B^{-deg} で求める static vec keis_inv(const vec &a, int deg) { assert(not a.empty() and mod / 2 <= a.back() and a.back() < mod); int k = deg, N = len(a); while (k > 64) k = (k + 1) >> 1; vec b(N + k + 1); b.back() = 1; b = divmod_naive(b, a).fi; while (k < deg) { vec s = mul(b, b); s.insert(s.begin(), 0); int d = min(N, k * 2 + 1); vec t{a.end() - d, a.end()}, v = mul(s, t); v.erase(v.begin(), v.begin() + d); vec w(k + 1), ww = add(b, b); copy(all(ww), back_inserter(w)); b = sub(w, v); b.erase(b.begin()); k <<= 1; } return b.erase(b.begin(), b.begin() + k - deg), b; } static PVV divmod_newton(const vec &a, const vec &b) { assert(not is_zero(b)); if (len(b) <= 64) return divmod_naive(a, b); if (len(a) - len(b) <= 64) return divmod_naive(a, b); int norm = mod / (b.back() + 1); vec x = mul(a, {norm}), y = mul(b, {norm}); int N = len(x), M = len(y); int deg = N - M + 2; vec z = keis_inv(y, deg), q = mul(x, z); q.erase(q.begin(), q.begin() + M + deg); vec yq = mul(y, q); while (less(x, yq)) q = sub(q, {1}), yq = sub(yq, y); vec r = sub(x, yq); while (less_eq(y, r)) q = add(q, {1}), r = sub(r, y); sh(q), sh(r); Z [qq, rr] = divmod_1e9(r, {norm}); assert(is_zero(rr)); return {q, qq}; } }; bigint abs(const bigint &x) { return x.abs(); } bigint gcd(const bigint &a, const bigint &b) { return a.gcd(b); } bigint lcm(const bigint &a, const bigint &b) { return a.lcm(b); } #line 9 "No_62_\u30ea\u30d9\u30ea\u30aa\u30f3_Extra.cpp" #define tests 1 #define fl 0 #define DB 10 void Yorisou() { using ll = bigint; LL(N, M, D, tx, ty, sx, sy, vx, vy); ll d = gcd(vx, vy); vx /= d, vy /= d, D *= d; ll f[]{tx, tx, N * 2 - tx, N * 2 - tx}, g[]{ty, M * 2 - ty, ty, M * 2 - ty}; FOR(i, 4) { ll tx = f[i], ty = g[i]; vc a, mds; Z [x, mod] = line_con_1(vx, tx - sx, N * 2); if (x == -1) continue; a.ep(x); mds.ep(mod); tie(x, mod) = line_con_1(vy, ty - sy, M * 2); if (x == -1) continue; a.ep(x); mds.ep(mod); ll st = crt_reduce(a, mds); if (st != -1 and st <= D) return print("Hit"); } print("Miss"); } #line 1 "YRS/aa/main.hpp" int main() { cin.tie(nullptr)->sync_with_stdio(0); int T = 1; if (fl) cerr.tie(0); if (tests and not fl) IN(T); for (int i = 0; i < T or fl; ++i) { Yorisou(); if (fl and i % DB == 0) cerr << "Case: " << i << '\n'; } return 0; } #line 36 "No_62_\u30ea\u30d9\u30ea\u30aa\u30f3_Extra.cpp"