#define SINGLE_TESTCASE // #define MULTI_TESTCASE // #define AOJ_TESTCASE #ifndef LOCAL #define FAST_IO // #define FAST_CIO // #define INTERACTIVE #endif #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 /** * @brief テンプレート(型) * @docs docs/template/template_types.md */ #include using namespace std; // https://github.com/miscalculation53/library/tree/wip/template/template_float128.hpp using namespace std; /** * @brief 四倍精度浮動小数点数 * @docs docs/template/template_types.md */ #if __cplusplus > 202002L && defined(__STDCPP_FLOAT128_T__) #include using f128 = std::float128_t; #elif defined(__SIZEOF_FLOAT128__) using f128 = __float128; #else #error "f128 requires a compiler with IEEE binary128 support, such as GCC on x86_64." #endif #define LIBRARY_HAS_FLOAT128 1 #if __cplusplus > 202002L && defined(__STDCPP_FLOAT128_T__) && \ (defined(_GLIBCXX_HAVE_FLOAT128_MATH) || defined(_GLIBCXX_LDOUBLE_IS_IEEE_BINARY128)) #define LIBRARY_USE_QUADMATH 0 #else #include #define LIBRARY_USE_QUADMATH 1 //! https://gcc.gnu.org/onlinedocs/libquadmath/Math-Library-Routines.html #define _p_F128_UNARY(name) \ inline f128 name(f128 x) { return f128(name##q(__float128(x))); } _p_F128_UNARY(acos) _p_F128_UNARY(acosh) _p_F128_UNARY(asin) _p_F128_UNARY(asinh) _p_F128_UNARY(atan) _p_F128_UNARY(atanh) _p_F128_UNARY(cbrt) _p_F128_UNARY(ceil) _p_F128_UNARY(cos) _p_F128_UNARY(cosh) _p_F128_UNARY(erf) _p_F128_UNARY(erfc) _p_F128_UNARY(exp) _p_F128_UNARY(exp2) _p_F128_UNARY(expm1) _p_F128_UNARY(fabs) _p_F128_UNARY(floor) _p_F128_UNARY(lgamma) _p_F128_UNARY(log) _p_F128_UNARY(log10) _p_F128_UNARY(log1p) _p_F128_UNARY(log2) _p_F128_UNARY(logb) _p_F128_UNARY(nearbyint) _p_F128_UNARY(rint) _p_F128_UNARY(round) _p_F128_UNARY(sin) _p_F128_UNARY(sinh) _p_F128_UNARY(sqrt) _p_F128_UNARY(tan) _p_F128_UNARY(tanh) _p_F128_UNARY(tgamma) _p_F128_UNARY(trunc) #undef _p_F128_UNARY template , int> = 0> inline f128 abs(T x) { return f128(fabsq(__float128(x))); } #define _p_F128_BINARY(name) \ inline f128 name(f128 x, f128 y) { return f128(name##q(__float128(x), __float128(y))); } \ template , int> = 0> \ inline f128 name(f128 x, T y) { return name(x, f128(y)); } \ template , int> = 0> \ inline f128 name(T x, f128 y) { return name(f128(x), y); } _p_F128_BINARY(atan2) _p_F128_BINARY(copysign) _p_F128_BINARY(fdim) _p_F128_BINARY(fmax) _p_F128_BINARY(fmin) _p_F128_BINARY(fmod) _p_F128_BINARY(hypot) _p_F128_BINARY(nextafter) _p_F128_BINARY(pow) _p_F128_BINARY(remainder) #undef _p_F128_BINARY inline f128 fma(f128 x, f128 y, f128 z) { return f128(fmaq(__float128(x), __float128(y), __float128(z))); } inline f128 hypot(f128 x, f128 y, f128 z) { return hypot(hypot(x, y), z); } namespace template_float128_detail { template constexpr bool arithmetic = is_arithmetic_v || is_same_v; } #define _p_F128_TERNARY_MIXED(name) \ template && template_float128_detail::arithmetic, int> = 0> \ inline f128 name(f128 x, T y, U z) { return name(x, f128(y), f128(z)); } \ template && !is_same_v && template_float128_detail::arithmetic, int> = 0> \ inline f128 name(T x, f128 y, U z) { return name(f128(x), y, f128(z)); } \ template && !is_same_v && template_float128_detail::arithmetic && !is_same_v, int> = 0> \ inline f128 name(T x, U y, f128 z) { return name(f128(x), f128(y), z); } _p_F128_TERNARY_MIXED(fma) _p_F128_TERNARY_MIXED(hypot) #undef _p_F128_TERNARY_MIXED inline f128 frexp(f128 x, int *exponent) { return f128(frexpq(__float128(x), exponent)); } inline f128 ldexp(f128 x, int exponent) { return f128(ldexpq(__float128(x), exponent)); } inline f128 scalbn(f128 x, int exponent) { return f128(scalbnq(__float128(x), exponent)); } inline f128 scalbln(f128 x, long exponent) { return f128(scalblnq(__float128(x), exponent)); } inline f128 modf(f128 x, f128 *integer) { __float128 value; const auto fraction = modfq(__float128(x), &value); *integer = f128(value); return f128(fraction); } inline f128 remquo(f128 x, f128 y, int *quotient) { return f128(remquoq(__float128(x), __float128(y), quotient)); } inline long lrint(f128 x) { return lrintq(__float128(x)); } inline long lround(f128 x) { return lroundq(__float128(x)); } inline long long llrint(f128 x) { return llrintq(__float128(x)); } inline long long llround(f128 x) { return llroundq(__float128(x)); } inline int ilogb(f128 x) { return ilogbq(__float128(x)); } template , int> = 0> inline bool isfinite(T x) { return finiteq(__float128(x)); } template , int> = 0> inline bool isinf(T x) { return isinfq(__float128(x)); } template , int> = 0> inline bool isnan(T x) { return isnanq(__float128(x)); } template , int> = 0> inline bool signbit(T x) { return signbitq(__float128(x)); } template , int> = 0> inline bool isnormal(T x) { return isfinite(x) && abs(x) >= numeric_limits::min(); } template , int> = 0> inline int fpclassify(T x) { if (isnan(x)) return FP_NAN; if (isinf(x)) return FP_INFINITE; if (x == 0) return FP_ZERO; return isnormal(x) ? FP_NORMAL : FP_SUBNORMAL; } #endif inline f128 nexttoward(f128 x, long double y) { return nextafter(x, f128(y)); } namespace template_float128_detail { inline bool parse(string_view s, f128 &x) { if (!s.empty() && s.front() == '+') { s.remove_prefix(1); if (s.empty() || s.front() == '-') return false; } #if LIBRARY_USE_QUADMATH string token(s); char *end; errno = 0; const __float128 value = strtoflt128(token.c_str(), &end); if (end == token.data() || end != token.data() + token.size()) return false; if (errno == ERANGE && (value == 0 || isinfq(value))) return false; x = f128(value); return true; #else auto [p, ec] = from_chars(s.data(), s.data() + s.size(), x); return ec == errc{} && p == s.data() + s.size(); #endif } inline string format(f128 x, chars_format style, int precision) { string s(128, '\0'); while (true) { #if LIBRARY_USE_QUADMATH const char *spec = style == chars_format::fixed ? "%.*Qf" : style == chars_format::scientific ? "%.*Qe" : "%.*Qg"; const int n = quadmath_snprintf(s.data(), s.size(), spec, precision, __float128(x)); if (n < 0) abort(); if (size_t(n) < s.size()) { s.resize(size_t(n)); return s; } s.resize(size_t(n) + 1); #else auto [p, ec] = to_chars(s.data(), s.data() + s.size(), x, style, precision); if (ec == errc{}) { s.resize(size_t(p - s.data())); return s; } s.resize(s.size() * 2); #endif } } } // namespace template_float128_detail #ifndef EPS #define EPS 1e-11 #endif using ld = decltype(EPS); using ll = long long; using uint = unsigned int; using ull = unsigned long long; using pll = pair; using tlll = tuple; using tllll = tuple; #define vc vector template using vvc = vc>; template using vvvc = vc>>; using vb = vc; using vl = vc; using vpll = vc; using vtlll = vc; using vtllll = vc; using vstr = vc; using vvb = vvc; using vvl = vvc; template using pql = priority_queue, greater>; template using pqg = priority_queue; #ifdef __SIZEOF_INT128__ using i128 = __int128_t; using u128 = __uint128_t; i128 stoi128(const string &s) { const bool neg = s.front() == '-'; u128 res = 0; for (size_t i = size_t(neg); i < s.size(); i++) res = 10 * res + unsigned(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 = u128(x); while (ux > 0) { res += char('0' + ux % 10); ux /= 10; } reverse(res.begin(), res.end()); return sign + res; } istream &operator>>(istream &is, i128 &a) { string s; if (!(is >> s)) return is; a = stoi128(s); return is; } ostream &operator<<(ostream &os, const i128 &a) { os << i128tos(a); return os; } #endif #define cauto const auto // https://github.com/miscalculation53/library/tree/wip/template/template_rep.hpp /** * @brief テンプレート(rep) * @docs docs/template/template_rep.md */ #define overload4(_1, _2, _3, _4, name, ...) name #define LIBRARY_REPEAT_LOOP(T, i, l, r, d) \ for (auto [i, library_repeat_end, library_repeat_step] = array{T(l), T(r), T(d)}; \ library_repeat_step > 0 ? i < library_repeat_end : i > library_repeat_end; \ i += library_repeat_step) #define rep1(i, n) LIBRARY_REPEAT_LOOP(ll, i, 0, n, 1) #define rep2(i, l, r) LIBRARY_REPEAT_LOOP(ll, i, l, r, 1) #define rep3(i, l, r, d) LIBRARY_REPEAT_LOOP(ll, i, l, r, d) #define rep(...) overload4(__VA_ARGS__, rep3, rep2, rep1)(__VA_ARGS__) #define repi1(i, n) LIBRARY_REPEAT_LOOP(int, i, 0, n, 1) #define repi2(i, l, r) LIBRARY_REPEAT_LOOP(int, i, l, r, 1) #define repi3(i, l, r, d) LIBRARY_REPEAT_LOOP(int, i, l, r, 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 #ifndef INF #define INF 4'000'000'000'000'000'037LL #endif #ifndef EPS #define EPS 1e-11 #endif // https://github.com/miscalculation53/library/tree/wip/utils/is_integral_ext.hpp /** * @brief $128$ ビット整数を含めた整数判定 * @docs docs/utils/is_integral_ext.md */ 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/integer_arithmetic.hpp // https://github.com/miscalculation53/library/tree/wip/utils/make_unsigned_ext.hpp /** * @brief $128$ ビット整数に対応した make_singed と make_unsigned * @docs docs/utils/make_unsigned_ext.md */ template struct make_unsigned_ext { using type = make_unsigned_t; }; template <> struct make_unsigned_ext { using type = u128; }; template <> struct make_unsigned_ext { using type = u128; }; template using make_unsigned_ext_t = typename make_unsigned_ext::type; template struct make_signed_ext { using type = make_signed_t; }; template <> struct make_signed_ext { using type = i128; }; template using make_signed_ext_t = typename make_signed_ext::type; /** * @brief 整数の符号・絶対値・剰余を保つ型変換 * @docs docs/utils/integer_arithmetic.md */ template constexpr bool is_negative(T x) { static_assert(is_integral_ext); if constexpr (is_signed_ext) return x < 0; return false; } template constexpr make_unsigned_ext_t unsigned_abs(T x) { static_assert(is_integral_ext && !is_same_v); using U = make_unsigned_ext_t; return is_negative(x) ? U(U(0) - U(x)) : U(x); } // 整数は整数昇格後も unsigned となる型、それ以外は元の型。 template > struct unsigned_arithmetic { using type = T; }; template struct unsigned_arithmetic { static_assert(!is_same_v); using type = common_type_t>; }; template using unsigned_arithmetic_t = typename unsigned_arithmetic::type; // mod 2^w (w は R の幅) で合同な値を返す。 template constexpr R wrapping_cast(T x) { static_assert(is_integral_ext && !is_same_v); static_assert(is_integral_ext && !is_same_v); if constexpr (is_signed_ext) { using U = make_unsigned_ext_t; const U bits = U(x); if (bits <= U(numeric_limits::max())) return R(bits); return R(-1) - R(U(~bits)); } else return R(x); } // https://github.com/miscalculation53/library/tree/wip/utils/default_infty.hpp #ifndef INF #define INF 4'000'000'000'000'000'037LL #endif /** * @brief 型ごとの既定の無限大 * @docs docs/utils/default_infty.md */ namespace default_infty_detail { template struct has_infty : false_type {}; template struct has_infty> : true_type {}; template const T &custom_value() { static const T value = T::infty(); return value; } template inline constexpr bool unsupported = false; } template constexpr decltype(auto) default_infty() { if constexpr (default_infty_detail::has_infty::value) return default_infty_detail::custom_value(); else if constexpr (is_same_v || is_same_v) return T(INF) * T(INF); else if constexpr (is_integral_ext && !is_same_v) { if constexpr (sizeof(T) >= sizeof(ll)) return T(INF); else if constexpr (numeric_limits::digits >= 31) return (T(1) << 30) - 1; else return T(numeric_limits::max() / 2); } else if constexpr (numeric_limits::has_infinity) return numeric_limits::infinity(); else static_assert(default_infty_detail::unsupported, "No default infinity for this type; specify infty explicitly or define T::infty()."); } /** * @brief テンプレート(演算) * @docs docs/template/template_math.md */ template inline bool chmin(T &a, U b) { return a > b ? a = b, true : false; } template inline bool chmax(T &a, U b) { return a < b ? a = b, true : false; } template || is_same_v, int> = 0> inline constexpr T abs(T x) { return is_negative(x) ? -x : x; } namespace template_math_detail { template constexpr common_type_t gcd(T a, U b) { using R = common_type_t; using W = make_unsigned_ext_t; W x = W(unsigned_abs(a)), y = W(unsigned_abs(b)); while (y != 0) { W r = x % y; x = y, y = r; } return R(x); } template constexpr common_type_t lcm(T a, U b) { using R = common_type_t; using W = make_unsigned_ext_t; W x = W(unsigned_abs(a)), y = W(unsigned_abs(b)); return x == 0 || y == 0 ? R(0) : R(x / template_math_detail::gcd(x, y) * y); } } // namespace template_math_detail template && !is_same_v, int> = 0> constexpr common_type_t gcd(i128 a, T b) { return template_math_detail::gcd(a, b); } template && !is_same_v && !is_same_v && !is_same_v, int> = 0> constexpr common_type_t gcd(T a, i128 b) { return template_math_detail::gcd(a, b); } template && !is_same_v, int> = 0> constexpr common_type_t gcd(u128 a, T b) { return template_math_detail::gcd(a, b); } template && !is_same_v && !is_same_v && !is_same_v, int> = 0> constexpr common_type_t gcd(T a, u128 b) { return template_math_detail::gcd(a, b); } template && !is_same_v, int> = 0> constexpr common_type_t lcm(i128 a, T b) { return template_math_detail::lcm(a, b); } template && !is_same_v && !is_same_v && !is_same_v, int> = 0> constexpr common_type_t lcm(T a, i128 b) { return template_math_detail::lcm(a, b); } template && !is_same_v, int> = 0> constexpr common_type_t lcm(u128 a, T b) { return template_math_detail::lcm(a, b); } template && !is_same_v && !is_same_v && !is_same_v, int> = 0> constexpr common_type_t lcm(T a, u128 b) { return template_math_detail::lcm(a, b); } 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 divceil(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) { const T x = T(a), y = T(b); if constexpr (is_signed_ext) if (y == T(-1)) return T(0); const T r = x % y; return r != T(0) && is_negative(x) != is_negative(y) ? T(r + y) : r; } template constexpr T ipow(U a, V b) { assert(b >= 0); if (b == 0) return 1; if (a == 0 || a == 1) return a; if (a < 0 && a == -1) return b & 1 ? -1 : 1; T res = 1, tmp = a; while (true) { if (b & 1) res *= tmp; b >>= 1; if (b == 0) break; tmp *= tmp; } return res; } template T mul_limited(A a, B b, M m) { assert(a >= 0 && b >= 0 && m >= 0); if (b == 0) return 0; return T(a) > T(m) / T(b) ? T(m) : T(a) * T(b); } template T mul_limited(A a, B b) { return mul_limited(a, b, default_infty()); } template T pow_limited(A a, B b, M m) { assert(a >= 0 && b >= 0 && m >= 0); if (a <= 1 || b == 0) return min(ipow(a, b), T(m)); T res = 1, tmp = a; while (true) { if (b & 1) { if (res > T(m) / tmp) return m; res *= tmp; } b >>= 1; if (b == 0) break; if (tmp > T(m) / tmp) return m; tmp *= tmp; } return res; } template T pow_limited(A a, B b) { return pow_limited(a, b, default_infty()); } template constexpr T iroot(A a, K k) { assert(a >= 0 && k >= 1); if (a <= 1 || k == 1) return a; if (k == 2) { const T aa = T(a); T x = T(sqrtl((long double)a)); while (x > aa / x) x--; while (x < numeric_limits::max()) { const T y = x + 1; if (y > aa / y) break; x = y; } return x; } auto isok = [&](T x) -> bool { if (x == 0) return true; T res = 1, k2 = k; while (true) { if (k2 & 1) { if (res > T(a) / x) return false; res *= x; } k2 >>= 1; if (k2 == 0) break; if (x > T(a) / x) return false; x *= x; } return res <= T(a); }; T x = pow(a, 1.0 / k); bool up = true; while (!isok(x)) up = false, x--; if (up) { while (x < numeric_limits::max() && isok(x + 1)) x++; } return x; } template constexpr T iroot_ceil(A a, K k) { T x = iroot(a, k); return ipow(x, k) == a ? x : x + 1; } //! https://misawa.github.io/others/avoid_errors/techniques_to_avoid_errors.html template int SGN(A a, D eps = EPS) { return int(a > eps) - int(a < -eps); } // 位取り記数法と同じ順番(下位桁が後ろ) // 0 に対しては {0} が返る template vc base_repr(U val, V base) { assert(val >= 0); assert(base >= 2); if (val == 0) return {0}; vc a; while (val > 0) { a.emplace_back(val % base); val /= base; } reverse(a.begin(), a.end()); return a; } // 位取り記数法と同じ順番(下位桁が後ろ) template vc base_repr(U val, V base, int n) { assert(val >= 0); assert(base >= 2); assert(n >= 0); vc a(n); repi(i, n) { a[i] = val % base; val /= base; } reverse(a.begin(), a.end()); return a; } template string base_repr_str(U val, int base) { assert(val >= 0); assert(2 <= base && base <= 36); auto a = base_repr(val, base); string s = ""; for (cauto &ai : a) s += (ai < 10 ? '0' + ai : (use_upper ? 'A' : 'a') + (ai - 10)); return s; } template string base_repr_str(U val, int base, int n) { assert(val >= 0); assert(2 <= base && base <= 36); assert(n >= 0); auto a = base_repr(val, base, n); string s = ""; for (cauto &ai : a) s += (ai < 10 ? '0' + ai : (use_upper ? 'A' : 'a') + (ai - 10)); return s; } // https://github.com/miscalculation53/library/tree/wip/template/template_vector.hpp /** * @brief テンプレート(vector) * @docs docs/template/template_vector.md */ #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; }) template auto gen_vec(int n, const F &f) { assert(n >= 0); vc res; res.reserve(n); repi(i, n) res.eb(f(i)); return res; } #define GEN_VEC(n, i, fi) (gen_vec(n, LMD(i, fi))) namespace template_vector_detail { template struct nested_vector { using type = vc::type>; }; template struct nested_vector { using type = T; }; template void populate_dimensions(typename nested_vector::type &result, const V *sizes, const T &init) { if constexpr (depth == 1) { typename nested_vector::type replacement(size_t(*sizes), init); result.swap(replacement); } else { result.resize(size_t(*sizes)); for (auto &child : result) populate_dimensions(child, sizes + 1, init); } } template typename nested_vector::type build_dimensions(const V *sizes, const T &init) { if constexpr (depth == 0) return init; else { typename nested_vector::type result; populate_dimensions(result, sizes, init); return result; } } } template auto dvec(const V (&sz)[d], const T &init) { if constexpr (i >= d) return init; else return template_vector_detail::build_dimensions(sz + i, init); } template auto dvec(const ll (&sz)[d], const T &init) { return dvec(sz, init); } template void fill(V &v, const T &value) { if constexpr (is_assignable_v) std::fill(std::begin(v), std::end(v), value); else for (auto &row : v) fill(row, value); } template optional ctol(const char &c, const string &s) { repi(i, SZ(s)) if (s[i] == c) return T(i); return nullopt; } template vc stov(const string &s, char first) { return gen_vec(SZ(s), [&](int i) -> T { return s[i] - first; }); } template optional> stov(const string &s, const string &t) { vc res; res.reserve(s.size()); for (char c : s) { auto index = ctol(c, t); if (!index) return nullopt; res.push_back(*index); } return res; } template string vtos(const vc &v, char first) { string res = ""; fe(vi : v) res += vi + first; return res; } template string vtos(const vc &v, const string &t) { string res = ""; fe(vi : v) res += t[vi]; return res; } template vc concat(const vvc &vs) { vc res; for (cauto &v : vs) res.insert(res.end(), ALL(v)); return res; } template vc concat(const vc &v) { return v; } template vc concat(vc v, const vc &...vs) { (v.insert(v.end(), ALL(vs)), ...); return v; } template vc merged(const vc &a, const vc &b) { vc res; merge(ALL(a), ALL(b), back_inserter(res)); return res; } template T vecget(const vc &v, I i, const T &dflt_negative = -default_infty(), const T &dflt_positive = default_infty()) { if (i < 0) return dflt_negative; if (i >= SZ(v)) return dflt_positive; return v[i]; } // https://github.com/miscalculation53/library/tree/wip/template/template_algo.hpp #ifndef INF #define INF 4'000'000'000'000'000'037LL #endif // https://github.com/miscalculation53/library/tree/wip/utils/resolved_infty.hpp // https://github.com/miscalculation53/library/tree/wip/utils/resolved_value.hpp /** * @brief 値または関数から値を取得 * @docs docs/utils/resolved_value.md */ template , int> = 0> constexpr T resolved_value() { return T(x); } template , long> = 0> const T &resolved_value() { static const T value = T(x()); return value; } /** * @brief 既定値・値・関数から無限大を取得 * @docs docs/utils/resolved_infty.md */ template constexpr decltype(auto) resolved_infty() { if constexpr (is_same_v) return default_infty(); else return resolved_value(); } /** * @brief テンプレート(アルゴリズム) * @docs docs/template/template_algo.md */ template auto SUM(const V &v) { typename V::value_type s{}; fec(vi : v) s += vi; return s; } template T SUM(const V &v) { T s{}; fec(vi : v) s += vi; return s; } template auto MAX(const V &v) { return *max_element(ALL(v)); } template auto MIN(const V &v) { return *min_element(ALL(v)); } template optional ARGMAX(const V &v) { if (v.empty()) return nullopt; return I(max_element(ALL(v)) - v.begin()); } template optional ARGMIN(const V &v) { if (v.empty()) return nullopt; return I(min_element(ALL(v)) - v.begin()); } template T mex(const V &a) { int n = a.size(); vector exists(n, false); repi(i, n) if (0 <= a[i] && a[i] < n) exists[a[i]] = true; repi(x, n) if (!exists[x]) return x; return n; } // (0, 1. ..., n-1) の順列か判定 template bool is_permutation(const vc &p) { const int n = p.size(); vc b(n, false); repi(i, n) { if (!(0 <= p[i] && p[i] < n)) return false; b[p[i]] = true; } return all_of(ALL(b), [](bool bi) { return bi; }); } template vc permid(const int &n, const int &base_index = 0) { vc p(n); repi(i, n) p[i] = i + base_index; return p; } template vc perminv(const vc &p) { if (p.empty()) return {}; const int n = p.size(); vc q(MAX(p) + 1); repi(i, n) if (p[i] >= 0) q[p[i]] = i; return q; } // a[p[i]] for all i template vc permuted(const vc &a, const vc &p) { const int n = p.size(); vc res(n); repi(i, n) { assert(0 <= p[i] && p[i] < U(a.size())); res[i] = a[p[i]]; } return res; } // p[q[r[i]]] for all i など template vc permuted(const vc &p, const vc &q, const vc &...rs) { return permuted(permuted(p, q), rs...); } template V reversed(const V &v) { return V(v.rbegin(), v.rend()); } #if __cplusplus < 202002L template V sorted(V v, Args&&... args) { sort(ALL(v), forward(args)...); return v; } #else template V sorted(V v, Args&&... args) { ranges::sort(v, forward(args)...); return v; } #endif template > void unique(V &v, Equal equal = {}) { v.erase(std::unique(ALL(v), equal), v.end()); } template > V uniqued(V v, Equal equal = {}) { unique(v, equal); return v; } template , class Equal = equal_to<>> void sortunique(V &v, Compare comp = {}, Equal equal = {}) { sort(ALL(v), comp); unique(v, equal); } template , class Equal = equal_to<>> V sortuniqued(V v, Compare comp = {}, Equal equal = {}) { sortunique(v, comp, equal); return v; } // 01234 -> 12340 template void rotate(V &v, U k) { const U n = v.size(); if (n == 0) return; k = (k % n + n) % n; std::rotate(v.begin(), v.begin() + k, v.end()); } // 01234 -> 12340 template V rotated(V v, U k) { rotate(v, k); return v; } 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; 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 {}; } template () * declval())> static constexpr S pow(const S &a, I k) { return a * k; } }; template struct MonoidMin { using S = T; static constexpr S op(S a, S b) { return min(a, b); } static constexpr decltype(auto) e() { return resolved_infty(); } template static constexpr S pow(const S &a, I k) { return k == 0 ? e() : a; } }; template struct MonoidMax { using S = T; static constexpr S op(S a, S b) { return max(a, b); } static constexpr S e() { return -resolved_infty(); } template static constexpr S pow(const S &a, I k) { return k == 0 ? e() : a; } }; namespace internal { template struct HasMonoidPow : false_type { }; template struct HasMonoidPow(), declval()))>> : true_type { }; } template typename M::S pow_monoid(typename M::S a, I k) { if constexpr (is_signed_ext) assert(k >= 0); if constexpr (internal::HasMonoidPow::value) return M::pow(a, k); else { typename M::S c = M::e(); for (; k; k >>= 1) { if (k & 1) c = M::op(c, a); a = M::op(a, a); } return c; } } template typename G::S pow_group(typename G::S a, I k) { if constexpr (is_signed_ext) { if (k < 0) { a = G::inv(a); return G::op(pow_monoid(a, -(k + 1)), a); } } return pow_monoid(a, k); } // left_index が 0 なら、長さ n+1 で a.front() が e() // left_index が 1 なら、長さ n で e() がない template vc cuml(const vc &v, int left_index = 0) { const int n = v.size(); vc res(n + 1); res[0] = M::e(); repi(i, n) res[i + 1] = M::op(res[i], v[i]); res.erase(res.begin(), res.begin() + left_index); return res; } // right_index が 0 なら、長さ n+1 で a.back() が e() // right_index が 1 なら、長さ n で e() がない template vc cumr(const vc &v, int right_index = 0) { return reversed(cuml(reversed(v), right_index)); } template vc cumlsum(const vc &v, int left_index = 0) { return cuml>(v, left_index); } template vc cumrsum(const vc &v, int right_index = 0) { return cumr>(v, right_index); } template vc cumlmin(const vc &v, int left_index = 0) { return cuml>(v, left_index); } template vc cumrmin(const vc &v, int right_index = 0) { return cumr>(v, right_index); } template vc cumlmax(const vc &v, int left_index = 0) { return cuml>(v, left_index); } template vc cumrmax(const vc &v, int right_index = 0) { return cumr>(v, right_index); } // デフォルトでは長さ n+1 // left_index, right_index をそれぞれ 1 にすると、左右が削除される template vc adjd(const vc &v, int left_index = 0, int right_index = 0) { int n = v.size(); assert(0 <= left_index && 0 <= right_index && left_index + right_index <= n + 1); vc res(n + 1); if (n == 0) { res[0] = T{}; res.erase(res.end() - right_index, res.end()); res.erase(res.begin(), res.begin() + left_index); return res; } res[0] = v[0]; repi(i, 1, n) res[i] = v[i] - v[i - 1]; res[n] = -v[n - 1]; res.erase(res.end() - right_index, res.end()); res.erase(res.begin(), res.begin() + left_index); return res; } 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 /** * @brief テンプレート(二分探索) * @docs docs/template/template_binsearch.md */ 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; // --- LB, UB --- #if __cplusplus < 202002L struct identity { template constexpr T &&operator()(T &&t) const noexcept { return forward(t); } }; namespace internal { template inline T bound_helper(const V &v, Judge judge) { int l = -1, r = v.size(); while (r - l > 1) { int m = (l + r) / 2; if (judge(m)) l = m; else r = m; } return r; } }; // val <= v[i] となる最小の i (val 未満の値の個数) template , class Proj = identity> inline T LB(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) { return internal::bound_helper(v, [&](int i) -> bool { return comp(proj(*(v.begin() + i)), val); }); } // val < v[i] となる最小の i (val 以下の値の個数) template , class Proj = identity> inline T UB(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) { return internal::bound_helper(v, [&](int i) -> bool { return !comp(val, proj(*(v.begin() + i))); }); } #define DEFAULT_COMP less<> #else // val <= v[i] となる最小の i (val 未満の値の個数) template inline T LB(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) { return ranges::lower_bound(v, val, comp, proj) - v.begin(); } // val < v[i] となる最小の i (val 以下の値の個数) template inline T UB(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) { return ranges::upper_bound(v, val, comp, proj) - v.begin(); } #define DEFAULT_COMP ranges::less #endif // --- vector 等の lt, leq, gt, geq --- // 条件を満たす添字を返す。左側は -1、右側は v.size() を番兵とする。 template inline auto lt_max(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) -> enable_if_t, T> { static_assert(numeric_limits::is_signed, "index type must represent -1"); return LB(v, val, comp, proj) - 1; } template inline auto leq_max(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) -> enable_if_t, T> { static_assert(numeric_limits::is_signed, "index type must represent -1"); return UB(v, val, comp, proj) - 1; } template inline auto gt_min(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) -> enable_if_t, T> { return UB(v, val, comp, proj); } template inline auto geq_min(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) -> enable_if_t, T> { return LB(v, val, comp, proj); } // v[i] < val となる i の個数 template inline auto lt_cnt(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) -> enable_if_t, T> { return LB(v, val, comp, proj); } // v[i] <= val となる i の個数 template inline auto leq_cnt(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) -> enable_if_t, T> { return UB(v, val, comp, proj); } // val < v[i] となる i の個数 template inline auto gt_cnt(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) -> enable_if_t, T> { return SZ(v) - UB(v, val, comp, proj); } // val <= v[i] となる i の個数 template inline auto geq_cnt(const V &v, const Value &val, Comp comp = {}, Proj proj = {}) -> enable_if_t, T> { return SZ(v) - LB(v, val, comp, proj); } // l <= v[i] < r となる i の個数 template inline auto in_cnt(const V &v, L l, R r, Comp comp = {}, Proj proj = {}) -> enable_if_t, T> { if (l > r) return 0; return lt_cnt(v, r, comp, proj) - lt_cnt(v, l, comp, proj); } // --- set 等の lt, leq, gt, geq --- // 条件を満たす要素へのイテレータを返す。存在しなければ nullopt。 template inline auto lt_max(const V &v, const Value &val) -> enable_if_t, optional> { auto it = v.lower_bound(val); if (it == v.begin()) return nullopt; return prev(it); } template inline auto leq_max(const V &v, const Value &val) -> enable_if_t, optional> { auto it = v.upper_bound(val); if (it == v.begin()) return nullopt; return prev(it); } template inline auto gt_min(const V &v, const Value &val) -> enable_if_t, optional> { auto it = v.upper_bound(val); if (it == v.end()) return nullopt; return it; } template inline auto geq_min(const V &v, const Value &val) -> enable_if_t, optional> { auto it = v.lower_bound(val); if (it == v.end()) return nullopt; return it; } // --- 自作二分探索 --- namespace internal { template bool binsearch_adjacent(T a, T b) { if (a < b) return a + 1 == b; if (b < a) return b + 1 == a; return false; } }; // (ok, ng) template pair binsearch(const Judge &judge, InitOk init_ok, InitNg init_ng, bool check_ok = true, bool check_ng = true) { T ok(init_ok), ng(init_ng); if (check_ok) assert(judge(ok)); if (check_ng) assert(!judge(ng)); while (!internal::binsearch_adjacent(ok, ng)) { T mid = (ok & ng) + ((ok ^ ng) >> 1); (judge(mid) ? ok : ng) = mid; } return {ok, ng}; } template T binsearch_real(const Judge &judge, InitOk init_ok, InitNg init_ng, int iteration_count = 100, bool check_ok = true, bool check_ng = true) { T ok(init_ok), ng(init_ng); if (check_ok) assert(judge(ok)); if (check_ng) assert(!judge(ng)); repi(_, iteration_count) { T mid = (ok + ng) / 2; (judge(mid) ? ok : ng) = mid; } return ok; } // (ok, ng) template pair expsearch(const Judge &judge, InitVal init_val, bool positive = true) { T cur(init_val), step = 1; const bool cur_ok = judge(cur); auto advance = [&](T x, T d, bool pos) -> T { if (pos) return x > numeric_limits::max() - d ? numeric_limits::max() : x + d; else return x < numeric_limits::lowest() + d ? numeric_limits::lowest() : x - d; }; T prv = advance(cur, 1, !positive); if (prv != cur && judge(prv) != cur_ok) { if (cur_ok) return {cur, prv}; else return {prv, cur}; } while (true) { T nxt = advance(cur, step, positive); assert(nxt != cur && "the boundary must exist in the searched direction"); if (nxt == cur || judge(nxt) != cur_ok) { T ok = cur_ok ? cur : nxt; T ng = cur_ok ? nxt : cur; return binsearch(judge, ok, ng, false, false); } cur = nxt; if (step > numeric_limits::max() / 2) step = numeric_limits::max(); else step *= 2; } } // https://github.com/miscalculation53/library/tree/wip/template/template_bit.hpp /** * @brief テンプレート(ビット演算) * @docs docs/template/template_bit.md */ template inline constexpr ull pow2(T k) { return 1ULL << k; } template inline constexpr ull MASK(T k) { return (1ULL << k) - 1ULL; } #if __cplusplus < 202002L // x == 0 ならば 0、そうでなければ 1 + floor(log2(x)) // 0, 1, 2, 2, 3, 3, 3, 3, 4, 4, ... inline constexpr ull bit_width(ull x) { return x == 0 ? 0 : 64 - __builtin_clzll(x); } // 0, 1, 2, 2, 4, 4, 4, 4, 8, 8, ... inline constexpr ull bit_floor(ull x) { return x == 0 ? 0ULL : 1ULL << (bit_width(x) - 1); } // 1, 1, 2, 4, 4, 8, 8, 8, 8, 16, ... inline constexpr ull bit_ceil(ull x) { return x == 0 ? 1ULL : 1ULL << bit_width(x - 1); } inline constexpr ull countr_zero(ull x) { assert(x != 0); return __builtin_ctzll(x); } inline constexpr ull popcount(ull x) { return __builtin_popcountll(x); } inline constexpr bool has_single_bit(ull x) { return popcount(x) == 1; } #else // 0, 1, 2, 2, 3, 3, 3, 3, 4, 4, ... inline constexpr ll bit_width(ll x) { return std::bit_width((ull)x); } // 0, 1, 2, 2, 4, 4, 4, 4, 8, 8, ... inline constexpr ll bit_floor(ll x) { return std::bit_floor((ull)x); } // 1, 1, 2, 4, 4, 8, 8, 8, 8, 16, ... 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); } #endif 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; } template inline void bset(T &x, uint k, bool b = 1) { b ? x |= (1ULL << k) : x &= ~(1ULL << k); } template inline void bflip(T &x, uint k) { x ^= (1ULL << k); } 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 #if defined(LOCAL) && LIBRARY_HAS_FLOAT128 namespace cpp_dump::_detail { struct export_command; inline string export_var(const f128 &, const string &, size_t, size_t, bool, const export_command &); } #endif // https://github.com/miscalculation53/library/tree/wip/template/template_dump_map.hpp /*! * Portions adapted from https://github.com/sassan-dev/cpp-dump/tree/aabd98dd8a4b5ae3f62cf12387a68a9f6591197e * * MIT License * * Copyright (c) 2023 Ryota Sasaki * * Permission is hereby granted, free of charge, to any person obtaining a copy * of this software and associated documentation files (the "Software"), to deal * in the Software without restriction, including without limitation the rights * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell * copies of the Software, and to permit persons to whom the Software is * furnished to do so, subject to the following conditions: * * The above copyright notice and this permission notice shall be included in all * copies or substantial portions of the Software. * * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE * SOFTWARE. */ #ifdef LOCAL // cpp-dump 内の再帰的な出力からも、map 用のオーバーロードを参照できるようにする。 namespace cpp_dump::_detail { struct export_command; #define _p_LIBRARY_DECLARE_MAP_EXPORT(Map) \ template \ string export_var(const std::Map &, const string &, size_t, \ size_t, bool, const export_command &); _p_LIBRARY_DECLARE_MAP_EXPORT(map) _p_LIBRARY_DECLARE_MAP_EXPORT(multimap) _p_LIBRARY_DECLARE_MAP_EXPORT(unordered_map) _p_LIBRARY_DECLARE_MAP_EXPORT(unordered_multimap) #undef _p_LIBRARY_DECLARE_MAP_EXPORT } #include //! https://github.com/philip82148/cpp-dump namespace cpp_dump::options { inline string map_key_color = "\x1b[36m"; } namespace cpp_dump::_detail { template struct map_values_for_dump { It first, last; const Format &format; map_values_for_dump(It first, It last, const Format &format) : first(first), last(last), format(format) {} struct iterator { It it; const Format &format; decltype(auto) operator*() const { return format(it->second); } iterator &operator++() { ++it; return *this; } bool operator!=(const iterator &other) const { return it != other.it; } }; iterator begin() const { return {first, format}; } iterator end() const { return {last, format}; } }; inline string color_map_key(const string &text) { if (!use_es() || options::map_key_color.empty()) return text; string plain; for (size_t i = 0; i < text.size();) { if (text[i] == '\x1b' && i + 1 < text.size() && text[i + 1] == '[') { size_t end = i + 2; while (end < text.size() && (isdigit(static_cast(text[end])) || text[end] == ';')) ++end; if (end < text.size() && text[end] == 'm') { i = end + 1; continue; } } plain += text[i++]; } return es::apply(options::map_key_color, plain); } template string export_library_map( const Map &map, const string &indent, size_t last_line_length, size_t current_depth, bool fail_on_newline, const export_command &command, const Format &format) { if (map.empty()) return es::bracket("{ }", current_depth); if (current_depth >= options::max_depth) return es::bracket("{ ", current_depth) + es::op("...") + es::bracket(" }", current_depth); const auto &key_command = command.next_for_map_key(); const auto &value_command = command.next_for_map_value(); const size_t next_depth = current_depth + 1; auto map_wrapper = [&]() { if constexpr (is_multimap) return _multimap_wrapper(map); else return _map_dummy_wrapper(map); }(); bool multiline = options::cont_indent_style == types::cont_indent_style_t::always; if (options::cont_indent_style == types::cont_indent_style_t::when_nested) multiline = is_multimap || is_iterable_like || is_iterable_like; if (options::cont_indent_style == types::cont_indent_style_t::when_non_tuples_nested) multiline = is_multimap || (is_iterable_like && !is_tuple) || (is_iterable_like && !is_tuple); for (;; multiline = true) { if (multiline && fail_on_newline) return "\n"; const string element_indent = multiline ? indent + " " : indent; string output = es::bracket(multiline ? "{" : "{ ", current_depth); bool first = true, retry = false; auto skipped_map = command.create_skip_container(map_wrapper); for (const auto &[ellipsis, it, index] : skipped_map) { if (!first) output += es::op(multiline ? "," : ", "); first = false; if (multiline) output += "\n" + element_indent; if (ellipsis) output += es::op("..."); else { const auto &[key, value] = *it; auto render = [&](const auto &x, const export_command &child_command) { return export_var(format(x), element_indent, get_last_line_length(output, last_line_length), next_depth, !multiline, child_command); }; output += color_map_key(render(key, key_command)); if constexpr (is_multimap) output += es::member(" (" + to_string(map.count(key)) + ")"); output += es::op(": "); if constexpr (is_multimap) { auto [begin, end] = map.equal_range(key); // この view は所有しない。描画中だけ元の map を参照する。 auto values = map_values_for_dump(begin, end, format); output += export_var(values, element_indent, get_last_line_length(output, last_line_length), next_depth, !multiline, value_command); } else output += render(value, value_command); } if (!multiline && (has_newline(output) || last_line_length + get_length(output) + 2 > options::max_line_width)) { retry = true; break; } } if (retry) continue; return output + (multiline ? "\n" + indent : " ") + es::bracket("}", current_depth); } } #define _p_LIBRARY_DEFINE_MAP_EXPORT(Map) \ template \ string export_var(const std::Map &value, const string &indent, \ size_t last_line_length, size_t current_depth, \ bool fail_on_newline, const export_command &command) \ { \ return export_library_map(value, indent, last_line_length, current_depth, \ fail_on_newline, command, \ [](const auto &x) -> const auto & { return x; }); \ } _p_LIBRARY_DEFINE_MAP_EXPORT(map) _p_LIBRARY_DEFINE_MAP_EXPORT(multimap) _p_LIBRARY_DEFINE_MAP_EXPORT(unordered_map) _p_LIBRARY_DEFINE_MAP_EXPORT(unordered_multimap) #undef _p_LIBRARY_DEFINE_MAP_EXPORT } #endif /** * @brief テンプレート(dump) * @docs docs/template/template_dump.md */ #ifdef LOCAL namespace cpp_dump::_detail { #if LIBRARY_HAS_FLOAT128 inline string export_var( const f128 &x, const string &, size_t, size_t, bool, const export_command & ) { return es::signed_number(template_float128_detail::format( x, chars_format::general, numeric_limits::max_digits10)); } #endif inline string export_var( const i128 &x, const string &, size_t, size_t, bool, const export_command & ) { return es::signed_number(i128tos(x)); } template inline auto export_object_generic( const T &value, const string &indent, size_t last_line_length, size_t current_depth, bool fail_on_newline, const export_command &command ) -> decltype(value.dump_data(), string()) { string class_name = es::class_name(get_typename()); _p_CPP_DUMP_DEFINE_EXPORT_OBJECT_COMMON1; apply( [&](const auto &...member) { (append_output(member.first, member.second.get()), ...); }, value.dump_data() ); _p_CPP_DUMP_DEFINE_EXPORT_OBJECT_COMMON2; } } // namespace cpp_dump::_detail #define _p_LIBRARY_CPP_DUMP_MEMBER(member) pair{string_view(#member), cref(member)} #define CPP_DUMP_DEFINE_DATA(...) \ auto dump_data() const \ { \ return make_tuple(_p_CPP_DUMP_EXPAND_VA(_p_LIBRARY_CPP_DUMP_MEMBER, __VA_ARGS__)); \ } #define dump(...) cpp_dump(__VA_ARGS__) namespace cp = cpp_dump; CPP_DUMP_SET_OPTION_GLOBAL(log_label_func, cp::log_label::line()); CPP_DUMP_SET_OPTION_GLOBAL(max_iteration_count, 100); #define local(...) __VA_ARGS__ #define oj(...) #define local_oj(a, b) (a) CPP_DUMP_DEFINE_EXPORT_OBJECT_GENERIC(content()); #else #define CPP_DUMP_DEFINE_DATA(...) #define dump(...) ((void)0) #define local(...) #define oj(...) __VA_ARGS__ #define local_oj(a, b) (b) #endif template vc content(queue que) { vc res; while (!que.empty()) { res.eb(que.front()); que.pop(); } return res; } template vc content(priority_queue pque) { vc res; while (!pque.empty()) { res.eb(pque.top()); pque.pop(); } return res; } template auto content(const T &obj) { return obj.content(); } /** * @brief テンプレート(入出力) * @docs docs/template/template_inout.md */ #if LIBRARY_HAS_FLOAT128 inline istream &operator>>(istream &is, f128 &x) { string s; if (is >> s) { f128 value; if (template_float128_detail::parse(s, value)) x = value; else is.setstate(ios::failbit); } return is; } inline ostream &operator<<(ostream &os, f128 x) { const auto field = os.flags() & ios::floatfield; const auto format = field == ios::fixed ? chars_format::fixed : field == ios::scientific ? chars_format::scientific : chars_format::general; string s = template_float128_detail::format(x, format, int(os.precision())); if ((os.flags() & ios::showpos) && !signbit(x)) s.insert(s.begin(), '+'); if (os.flags() & ios::uppercase) for (char &c : s) c = char(toupper(static_cast(c))); return os << s; } #endif // 十進変換とバッファ処理は公開 API の仕様から独立に設計。 #if defined FAST_IO and not defined LOCAL namespace fastio { constexpr uint32_t SIZ = 1 << 17; inline char ibuf[SIZ + 1], obuf[SIZ]; inline uint32_t pil = 0, pir = 0, por = 0; __attribute__((noinline)) inline void load() { const uint32_t remaining = pir - pil; if (remaining && pil) memmove(ibuf, ibuf + pil, remaining); pil = 0; pir = remaining + uint32_t(fread(ibuf + remaining, 1, SIZ - remaining, stdin)); ibuf[pir] = '\x7f'; } inline bool ascii_space(unsigned char c) { return c == ' ' || unsigned(c - '\t') < 5; } inline void find_token() { for (;;) { const char *p = ibuf + pil; const char *end = ibuf + pir; while (ascii_space(*p)) ++p; pil = uint32_t(p - ibuf); if (p != end) return; load(); if (!pir) abort(); } } __attribute__((noinline)) inline void flush() { size_t written = 0; while (written < por) { const size_t n = fwrite(obuf + written, 1, por - written, stdout); if (!n) abort(); written += n; } por = 0; } inline void append_bytes(const char *data, size_t size) { while (size) { const size_t n = min(size, size_t(SIZ - por)); memcpy(obuf + por, data, n); por += uint32_t(n); data += n; size -= n; if (por == SIZ) flush(); } } inline void append_char(char c) { if (por == SIZ) flush(); obuf[por++] = c; } template using integer_magnitude = conditional_t<(sizeof(T) > sizeof(ull)), u128, conditional_t<(sizeof(T) > sizeof(uint)), ull, uint>>; // 各桁対、4 桁、8 桁の順に十進の重みをまとめる。 inline uint decode_eight(uint64_t columns) { columns &= 0x0f0f0f0f0f0f0f0fULL; const uint64_t pairs = (columns & 0x00ff00ff00ff00ffULL) * 10 + ((columns >> 8) & 0x00ff00ff00ff00ffULL); const uint64_t quads = (pairs & 0x0000ffff0000ffffULL) * 100 + ((pairs >> 16) & 0x0000ffff0000ffffULL); return uint(quads) * 10000 + uint(quads >> 32); } inline constexpr uint powers_of_ten[8] = {1, 10, 100, 1000, 10000, 100000, 1000000, 10000000}; template __attribute__((always_inline)) inline void rd1_integer(T &x) { using U = integer_magnitude; const char *p = ibuf + pil; for (;;) { while (ascii_space(*p)) ++p; if (p != ibuf + pir) break; pil = pir; load(); if (!pir) abort(); p = ibuf; } bool negative = false; if (*p == '+') ++p; else if constexpr (numeric_limits::is_signed) { if (*p == '-') negative = true, ++p; } U value = 0; for (;;) { while (size_t(ibuf + pir - p) >= 8) { uint64_t columns; memcpy(&columns, p, 8); #if __BYTE_ORDER__ == __ORDER_BIG_ENDIAN__ columns = __builtin_bswap64(columns); #endif const uint64_t stop = ((columns + 0x4646464646464646ULL) | (columns - 0x3030303030303030ULL)) & 0x8080808080808080ULL; if (stop) { if constexpr (sizeof(U) > sizeof(ull)) { const unsigned digits = unsigned(__builtin_ctzll(stop)) / 8; if (digits) { uint tail; if (digits <= 3) { tail = unsigned(p[0] - '0'); if (digits >= 2) tail = tail * 10 + unsigned(p[1] - '0'); if (digits == 3) tail = tail * 10 + unsigned(p[2] - '0'); } else tail = decode_eight(columns << (8 * (8 - digits))); value = value * powers_of_ten[digits] + tail; p += digits; } goto integer_complete; } else { if (stop & 0x80808080ULL) { if (!(stop & 0x808080ULL)) { value = value * 1000 + uint(columns & 15) * 100 + uint(columns >> 8 & 15) * 10 + uint(columns >> 16 & 15); p += 3; } else if (!(stop & 0x8080ULL)) { value = value * 100 + uint(columns & 15) * 10 + uint(columns >> 8 & 15); p += 2; } else if (!(stop & 0x80ULL)) { value = value * 10 + uint(columns & 15); ++p; } goto integer_complete; } break; } } columns &= 0x0f0f0f0f0f0f0f0fULL; const uint64_t pairs = (columns & 0x00ff00ff00ff00ffULL) * 10 + ((columns >> 8) & 0x00ff00ff00ff00ffULL); const uint64_t quads = (pairs & 0x0000ffff0000ffffULL) * 100 + ((pairs >> 16) & 0x0000ffff0000ffffULL); value = value * 100000000 + uint(quads) * 10000 + uint(quads >> 32); p += 8; } while (unsigned(*p - '0') < 10) { value = value * 10 + unsigned(*p - '0'); ++p; } pil = uint32_t(p - ibuf); if (pil != pir) break; load(); p = ibuf; if (!pir) break; } integer_complete: pil = uint32_t(p - ibuf); if constexpr (numeric_limits::is_signed) { if (negative && value) x = T(-T(value - 1) - 1); else x = T(value); } else x = T(value); } template ::is_integer, int> = 0> __attribute__((always_inline)) inline void rd1(T &x) { rd1_integer(x); } inline void rd1(char &x) { find_token(); x = ibuf[pil++]; } inline void rd1(string &x) { find_token(); x.clear(); for (;;) { const char *p = ibuf + pil; const char *end = ibuf + pir; const char *first = p; while (p != end && !ascii_space(*p)) ++p; x.append(first, size_t(p - first)); pil = uint32_t(p - ibuf); if (p != end) return; load(); if (!pir) return; } } template void rd1_real(T &x) { string token; rd1(token); const char *first = token.data(); const char *last = first + token.size(); if (first != last && *first == '+') ++first; auto result = from_chars(first, last, x, chars_format::general); if (result.ec != errc{} || result.ptr != last) { char *end; if constexpr (is_same_v) x = strtof(token.c_str(), &end); else if constexpr (is_same_v) x = strtod(token.c_str(), &end); else x = strtold(token.c_str(), &end); if (end != token.data() + token.size()) abort(); } } inline void rd1(float &x) { rd1_real(x); } inline void rd1(double &x) { rd1_real(x); } inline void rd1(long double &x) { rd1_real(x); } #if LIBRARY_HAS_FLOAT128 inline void rd1(f128 &x) { string token; rd1(token); if (!template_float128_detail::parse(token, x)) abort(); } #endif template void rd1(pair &x); template void rd1(tuple &x); template void rd1(array &x); template void rd1(vc &x); template void rd1(pair &x) { rd1(x.first), rd1(x.second); } template void rd1(tuple &x) { apply([](auto &...element) { (rd1(element), ...); }, x); } template void rd1(array &x) { for (auto &element : x) rd1(element); } template void rd1(vc &x) { if constexpr (is_same_v) for (size_t i = 0; i < x.size(); ++i) { bool value; rd1(value); x[i] = value; } else for (auto &element : x) rd1(element); } inline constexpr auto decimal_triplets = [] { array, 1000> result{}; for (unsigned i = 0; i < 1000; ++i) { result[i][0] = char('0' + i / 100); result[i][1] = char('0' + i / 10 % 10); result[i][2] = char('0' + i % 10); } return result; }(); inline void store_nine(char *out, uint x) { const uint q = x / 1000, r = q / 1000; memcpy(out, decimal_triplets[r].data(), 3); memcpy(out + 3, decimal_triplets[q - r * 1000].data(), 3); memcpy(out + 6, decimal_triplets[x - q * 1000].data(), 3); } inline char *prepend_small(char *end, uint x) { while (x >= 1000) { const uint q = x / 1000; end -= 3; memcpy(end, decimal_triplets[x - q * 1000].data(), 3); x = q; } const unsigned digits = 1 + unsigned(x >= 10) + unsigned(x >= 100); end -= digits; memcpy(end, decimal_triplets[x].data() + 3 - digits, digits); return end; } inline char *prepend_wide(char *end, ull x) { while (x >= 1000000000ULL) { const ull q = x / 1000000000ULL; end -= 9; store_nine(end, uint(x - q * 1000000000ULL)); x = q; } return prepend_small(end, uint(x)); } inline char *store_leading_triplet(char *out, uint x) { if (x >= 100) { memcpy(out, decimal_triplets[x].data(), 3); return out + 3; } if (x >= 10) { memcpy(out, decimal_triplets[x].data() + 1, 2); return out + 2; } *out = char('0' + x); return out + 1; } inline char *store_short(char *out, uint x) { if (x >= 1000000) { const uint q = x / 1000000; const uint r = x - q * 1000000; out = store_leading_triplet(out, q); memcpy(out, decimal_triplets[r / 1000].data(), 3); memcpy(out + 3, decimal_triplets[r % 1000].data(), 3); return out + 6; } if (x >= 1000) { const uint q = x / 1000; out = store_leading_triplet(out, q); memcpy(out, decimal_triplets[x - q * 1000].data(), 3); return out + 3; } return store_leading_triplet(out, x); } template inline char *store_ordinary(char *out, U x) { if constexpr (sizeof(U) > sizeof(uint)) { if (x >= 1000000000000000000ULL) { const uint hi = uint(x / 1000000000000000000ULL); const ull tail = x - ull(hi) * 1000000000000000000ULL; const uint middle = uint(tail / 1000000000ULL); out = store_leading_triplet(out, hi); store_nine(out, middle); store_nine(out + 9, uint(tail - ull(middle) * 1000000000ULL)); return out + 18; } } if (x >= 1000000000) { const uint hi = uint(x / 1000000000); out = store_short(out, hi); store_nine(out, uint(x - U(hi) * 1000000000)); return out + 9; } return store_short(out, uint(x)); } template __attribute__((always_inline)) inline void wt1_integer(T x) { using U = integer_magnitude; U magnitude = U(x); bool negative = false; if constexpr (numeric_limits::is_signed) { negative = x < 0; const U sign_mask = U(0) - U(negative); magnitude = (magnitude ^ sign_mask) + U(negative); } if constexpr (sizeof(U) <= sizeof(ull)) { if (por > SIZ - 40) flush(); char *out = obuf + por; // 負号の位置だけを進め、非負なら数字で上書きする。 *out = '-'; out += negative; out = store_ordinary(out, magnitude); if constexpr (newline) *out++ = '\n'; por = uint32_t(out - obuf); return; } char text[48]; char *const end = text + sizeof(text); char *first = end; if constexpr (sizeof(U) > sizeof(ull)) { constexpr ull base = 1000000000000000000ULL; while (magnitude >= base) { const U q = magnitude / base; const ull tail = ull(magnitude - q * base); const uint hi = uint(tail / 1000000000ULL); first -= 18; store_nine(first, hi); store_nine(first + 9, uint(tail - ull(hi) * 1000000000ULL)); magnitude = q; } first = prepend_wide(first, ull(magnitude)); } else if constexpr (sizeof(U) > sizeof(uint)) first = prepend_wide(first, magnitude); else first = prepend_small(first, magnitude); if (negative) *--first = '-'; append_bytes(first, size_t(end - first)); if constexpr (newline) append_char('\n'); } template ::is_integer, int> = 0> __attribute__((always_inline)) inline void wt1(T x) { wt1_integer(x); } inline void wt1(char x) { append_char(x); } inline void wt1(const char *x) { append_bytes(x, strlen(x)); } inline void wt1(string_view x) { append_bytes(x.data(), x.size()); } inline void wt1(const string &x) { append_bytes(x.data(), x.size()); } template void wt1_real(T x) { char text[numeric_limits::max_exponent10 + 64]; auto result = to_chars(text, text + sizeof(text), x, chars_format::fixed, 15); if (result.ec != errc{}) abort(); append_bytes(text, size_t(result.ptr - text)); } inline void wt1(float x) { wt1_real(x); } inline void wt1(double x) { wt1_real(x); } inline void wt1(long double x) { wt1_real(x); } #if LIBRARY_HAS_FLOAT128 inline void wt1(f128 x) { wt1(template_float128_detail::format(x, chars_format::fixed, 35)); } #endif template void wt1(const pair &x); template void wt1(const tuple &x); template void wt1(const array &x); template void wt1(const vc &x); template void wt1(const pair &x) { wt1(x.first), append_char(' '), wt1(x.second); } template void wt1(const tuple &x) { apply([](const auto &...element) { size_t index = 0; ((index++ ? void(append_char(' ')) : void(), wt1(element)), ...); }, x); } template void wt1(const array &x) { for (size_t i = 0; i < N; ++i) { if (i) append_char(' '); wt1(x[i]); } } template void wt1(const vc &x) { for (size_t i = 0; i < x.size(); ++i) { if (i) append_char(' '); wt1(x[i]); } } template void read(T &...x) { (rd1(x), ...); } template void write(T &&...x) { (wt1(forward(x)), ...); } template void print(T &&...x) { if constexpr (sizeof...(T) == 1 && ((is_integral_v> && !is_same_v, char>) && ...)) (wt1_integer, true>(forward(x)), ...); else { size_t index = 0; ((index++ ? void(append_char(' ')) : void(), wt1(forward(x))), ...); append_char('\n'); } } } #endif #if defined FAST_IO and not defined LOCAL struct Dummy { Dummy() { atexit(fastio::flush); } } dummy; #endif //! https://trap.jp/post/1224/ // ---- 入力 ---- #if defined LOCAL or not defined FAST_IO template istream &operator>>(istream &is, pair &p) { is >> p.first >> p.second; return is; } template istream &operator>>(istream &is, tuple &t) { apply([&](auto &...a) { (is >> ... >> a); }, t); return is; } template istream &operator>>(istream &is, array &a) { for (size_t i = 0; i < n; i++) is >> a[i]; return is; } template istream &operator>>(istream &is, vc &a) { const size_t n = a.size(); for (size_t i = 0; i < n; i++) is >> a[i]; return is; } #endif namespace internal { #if defined LOCAL or not defined FAST_IO template void CIN(Ts &...a) { (cin >> ... >> a); } #endif #if defined FAST_IO and not defined LOCAL template void READnodump(Ts &...a) { fastio::read(a...); } #else template void READnodump(Ts &...a) { CIN(a...); } #endif template void READVECnodump(int n, vc &...v) { (v.resize(size_t(n)), ...); READnodump(v...); } template void READVEC2nodump(int n, int m, vvc &...v) { (v.assign(size_t(n), vc(size_t(m))), ...); READnodump(v...); } template void READJAGnodump(int n, vvc &...vs) { auto read_one = [&](auto &v) { v.resize(size_t(n)); for (auto &row : v) { int k; READnodump(k); row.resize(size_t(k)); READnodump(row); } }; (read_one(vs), ...); } }; // namespace internal #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__) // ---------- // ----- 出力 ----- #ifdef INTERACTIVE #define ENDL endl #else #define ENDL '\n' #endif #if defined LOCAL or not defined FAST_IO template ostream &operator<<(ostream &os, const pair &p) { os << p.first << ' ' << p.second; return os; } template ostream &operator<<(ostream &os, const tuple &t) { if constexpr (sizeof...(Ts)) { apply([&](const auto &...x) { int i = 0; ((os << (i++ ? " " : "") << x), ...); }, t); } return os; } template ostream &operator<<(ostream &os, const array &a) { for (size_t i = 0; i < n; i++) { if (i) os << ' '; os << a[i]; } return os; } template ostream &operator<<(ostream &os, const vc &v) { const size_t n = v.size(); for (size_t i = 0; i < n; i++) { if (i) os << ' '; os << v[i]; } return os; } namespace internal { template void COUTW(const Ts &...a) { if constexpr (sizeof...(Ts)) (cout << ... << a); } template void COUTP(const Ts &...a) { if constexpr (sizeof...(Ts)) { int i = 0; ((cout << (i++ ? " " : "") << a), ...); } cout << ENDL; } }; // namespace internal #endif #if defined FAST_IO and not defined LOCAL #define WRITE fastio::write #define PRINT fastio::print #else #define WRITE internal::COUTW #define PRINT internal::COUTP #endif #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) // ---------- // ----- 基準ずらし ----- template pair &operator+=(pair &a, const P &b) { a.first += b.first; a.second += b.second; return a; } template pair operator+(pair a, const P &b) { return a += b; } template pair &operator-=(pair &a, const P &b) { a.first -= b.first; a.second -= b.second; return a; } template pair operator-(pair a, const P &b) { return a -= b; } template pair operator-(pair a) { a.first = -a.first; a.second = -a.second; return a; } template array &operator+=(array &a, const A &b) { for (size_t i = 0; i < n; i++) a[i] += b[i]; return a; } template array operator+(array a, const A &b) { return a += b; } template array &operator-=(array &a, const A &b) { for (size_t i = 0; i < n; i++) a[i] -= b[i]; return a; } template array operator-(array a, const A &b) { return a -= b; } template array operator-(array a) { for (auto &ai : a) ai = -ai; return a; } namespace internal { template auto &tuple_add_impl(A &a, const B &b, const index_sequence) { ((get(a) += get(b)), ...); return a; } template auto &tuple_sub_impl(A &a, const B &b, const index_sequence) { ((get(a) -= get(b)), ...); return a; } template auto &tuple_neg_impl(A &a, const index_sequence) { ((get(a) = -get(a)), ...); return a; } }; // namespace internal template tuple &operator+=(tuple &a, const Tp &b) { return internal::tuple_add_impl(a, b, make_index_sequence>>{}); } template tuple operator+(tuple a, const Tp &b) { return a += b; } template tuple &operator-=(tuple &a, const Tp &b) { return internal::tuple_sub_impl(a, b, make_index_sequence>>{}); } template tuple operator-(tuple a, const Tp &b) { return a -= b; } template tuple operator-(tuple a) { internal::tuple_neg_impl(a, make_index_sequence{}); return a; } template void offset(vc &v, const Add &add) { for (auto &vi : v) vi += add; } template void offset(vvc &v, const Add &add) { for (auto &vi : v) for (auto &vij : vi) vij += add; } // ---------- // ----- 転置 ----- template array, m> unzip(const vc> &vt) { const size_t n = vt.size(); array, m> tv; tv.fill(vc(n)); for (size_t i = 0; i < n; i++) for (size_t j = 0; j < m; j++) tv[j][i] = vt[i][j]; return tv; } template vc> zip(const array, m> &tv) { if (tv.empty()) return {}; const size_t n = tv[0].size(); vc> vt(n); for (size_t j = 0; j < m; j++) { assert(tv[j].size() == n); for (size_t i = 0; i < n; i++) vt[i][j] = tv[j][i]; } return vt; } template pair, vc> unzip(const vc> &vt) { const size_t n = vt.size(); pair, vc> tv; tv.first.resize(n), tv.second.resize(n); for (size_t i = 0; i < n; i++) tie(tv.first[i], tv.second[i]) = vt[i]; return tv; } template vc> zip(const pair, vc> &tv) { const size_t n = tv.first.size(); assert(n == tv.second.size()); vc> vt(n); for (size_t i = 0; i < n; i++) vt[i] = make_pair(tv.first[i], tv.second[i]); return vt; } namespace internal { template auto vt_to_tv_impl(V &tv, const Tp &t, index_sequence, size_t index) { ((get(tv)[index] = get(t)), ...); } template auto tv_to_vt_impl(const Tp &tv, index_sequence, size_t index) { return make_tuple(get(tv)[index]...); } }; template auto unzip(const vc> &vt) { const size_t n = vt.size(); tuple...> tv; apply([&](auto &...v) { ((v.resize(n)), ...); }, tv); for (size_t i = 0; i < n; i++) internal::vt_to_tv_impl(tv, vt[i], make_index_sequence>{}, i); return tv; } template auto zip(const tuple...> &tv) { size_t n = get<0>(tv).size(); apply([&](auto &...v) { ((void(v), assert(v.size() == n)), ...); }, tv); vc> vt(n); for (size_t i = 0; i < n; i++) vt[i] = internal::tv_to_vt_impl(tv, index_sequence_for{}, i); return vt; } namespace internal { template auto zip_vectors(const vc &...v) { if constexpr (sizeof...(Ts) == 2) return zip(pair{v...}); else return zip(tuple{v...}); } }; #define UNZIP(vt, ...) auto [__VA_ARGS__] = unzip(vt) #define ZIP(vt, ...) auto vt = internal::zip_vectors(__VA_ARGS__) // ---------- // https://github.com/miscalculation53/library/tree/wip/template/template_random.hpp /** * @brief テンプレート(ランダム生成) * @docs docs/template/template_random.md */ mt19937_64 mt; // [l, r] から等確率 template T randint(U1 l, U2 r) { assert(T(l) <= T(r)); return uniform_int_distribution(T(l), T(r))(mt); } // [l, r) から等確率 template T randrange(U1 l, U2 r) { assert(T(l) < T(r)); return uniform_int_distribution(T(l), T(r) - 1)(mt); } // [l, r) から一様ランダムな実数を返す template T randreal(U1 l, U2 r) { assert(T(l) < T(r)); return uniform_real_distribution(T(l), T(r))(mt); } // 確率 p で true を返す bool randbool(double p) { assert(0 <= p && p <= 1); return bernoulli_distribution(p)(mt); } namespace internal { template void random_sample_range(V &res, T l, T r) { int k = res.size(); T n = r - l; if (k <= 256) { repi(i, k) { T j = n - T(k) + T(i), x = randint(0, j); if (find(res.begin(), res.begin() + i, x) != res.begin() + i) x = j; res[i] = x; } } else { unordered_set used; used.reserve(2 * size_t(k)); repi(i, k) { T j = n - T(k) + T(i), x = randint(0, j); if (!used.insert(x).second) x = j, used.insert(x); res[i] = x; } } for (T &x : res) x += l; if constexpr (does_sort) sort(res.begin(), res.end()); else shuffle(res.begin(), res.end(), mt); } }; // namespace internal // [l, r) から相異なる k 個を選ぶ // does_sort: ソートするかどうか template array random_sample_range_array(U1 l, U2 r) { assert(T(r) - T(l) >= T(k)); array res; internal::random_sample_range(res, T(l), T(r)); return res; } // [l, r) から相異なる k 個を選ぶ // does_sort: ソートするかどうか template vc random_sample_range_vector(U1 l, U2 r, int k) { assert(k >= 0); assert(T(r) - T(l) >= T(k)); vc res(k); internal::random_sample_range(res, T(l), T(r)); return res; } // https://github.com/miscalculation53/library/tree/wip/math/modint/template_modint.hpp // https://github.com/miscalculation53/library/tree/wip/math/modint/modint.hpp /*! * Portions corresponding to AtCoder Library (CC0-1.0): modular arithmetic operators. * Upstream: atcoder/modint.hpp * https://github.com/atcoder/ac-library/tree/864245a00b00dd008d1abfdc239618fdb7d139da * License: https://github.com/atcoder/ac-library/blob/864245a00b00dd008d1abfdc239618fdb7d139da/LICENSE */ // https://github.com/miscalculation53/library/tree/wip/math/modint/modint_internal_static.hpp // https://github.com/miscalculation53/library/tree/wip/utils/larger_int.hpp /** * @brief 積の計算に使う拡張型 * @docs docs/utils/larger_int.md */ namespace larger_int_detail { template [[deprecated("larger_int: 128-bit integer is not widened; intermediate arithmetic may overflow")]] constexpr bool warn_no_wider_integer() { return true; } } template struct larger_int { private: static constexpr bool check() { // 型エイリアス経由でも、128 ビット整数を入力にしたときだけ警告する。 if constexpr (is_same_v || is_same_v) return larger_int_detail::warn_no_wider_integer(); return true; } 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 /*! * Portions corresponding to AtCoder Library (CC0-1.0): primality testing. * Upstream: atcoder/internal_math.hpp * https://github.com/atcoder/ac-library/tree/864245a00b00dd008d1abfdc239618fdb7d139da * License: https://github.com/atcoder/ac-library/blob/864245a00b00dd008d1abfdc239618fdb7d139da/LICENSE */ 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_constexpr(T n) { if constexpr (sizeof(T) > 4) { if (n <= INT_MAX) return isprime_constexpr(n); } if (n <= 1) return false; if (n == 2 || n == 7 || n == 61) return true; if (n % 2 == 0) return false; ll d = n - 1; while (d % 2 == 0) d /= 2; using U = make_unsigned_t; using L = larger_int_t; auto miller_rabin = [&](const auto &bases) constexpr { for (ll a : bases) { ll t = d, y = powmod_constexpr(a, t, n); while (t != n - 1 && y != 1 && y != n - 1) { y = L(y) * y % n; t <<= 1; } if (y != n - 1 && t % 2 == 0) return false; } return true; }; if constexpr(sizeof(T) <= 4) { constexpr ll bases[3] = {2, 7, 61}; return miller_rabin(bases); } else { constexpr ll bases[7] = {2, 325, 9375, 28178, 450775, 9780504, 1795265022}; return miller_rabin(bases); } } 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 = conditional_t<(M <= UINT_MAX), ull, 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 /*! * Portions corresponding to AtCoder Library (CC0-1.0): Barrett reduction. * Upstream: atcoder/internal_math.hpp * https://github.com/atcoder/ac-library/tree/864245a00b00dd008d1abfdc239618fdb7d139da * License: https://github.com/atcoder/ac-library/blob/864245a00b00dd008d1abfdc239618fdb7d139da/LICENSE */ namespace internal { struct barrett32 { uint m = 0; ull im = 0; barrett32() = default; 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{}; static void set_mod(mod_type m) { assert(m >= 1 && "modint: set_mod(m) requires m >= 1"); 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 = 0, im = 0, sq = 0; montgomery64odd() = default; 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); } }; //! https://www.mathenachia.blog/even-mod-montgomery-impl/ struct montgomery64 { ull m = 0, mx = 0, imx = 0, d = 0, q = 0; uint b = 0; montgomery64() = default; 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{}; static void set_mod(mod_type m) { assert(m >= 1 && m % 2 == 1 && "modint: set_mod(m) requires a positive odd modulus"); reducer = montgomery64odd(m); } 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{}; static void set_mod(mod_type m) { assert(m >= 1 && "modint: set_mod(m) requires m >= 1"); reducer = montgomery64(m); } 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 /** * @brief 拡張ユークリッド互除法 (extgcd) * @docs docs/math/extgcd.md */ // g == gcd(x, y) >= 0, ax + by == g を満たす (g, x, y) // max(|x|, |y|) <= max(|a|, |b|) template constexpr tuple extgcd(T a, T b) { if (a == 0 && b == 0) return {0, 0, 0}; // a*x1 + b*y1 == z1 ...(1) // a*x2 + b*y2 == z2 ...(2) T x1 = 1, y1 = 0, z1 = a; T x2 = 0, y2 = 1, z2 = b; while (z2 != 0) { // (1)' = (2) // (2)' = (1) - q*(2) 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() { M m = Policy::mod(); assert(m > 0 && "modint: modulus is not set; call mint::set_mod(m) before use"); return m; } template static auto set_mod(M m) -> decltype(T::set_mod(m)) { return T::set_mod(m); } static void precompute() { if constexpr (!is_void_v) (void)PowInv::get(); } static mint raw(V v) { mint x; x._v = v; return x; } modint_impl() : _v(0) { (void)mod(); } template >> modint_impl(T v) { const V m = mod(); V rem; if constexpr (is_signed_ext) { using S = make_signed_t; S x = v % S(m); if (x < 0) x += m; rem = x; } else rem = V(v % m); _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++() { _v++; if (_v == Policy::umod()) _v = 0; return *this; } mint &operator--() { if (_v == 0) _v = Policy::umod(); _v--; return *this; } mint operator++(int) { mint res = *this; ++(*this); return res; } mint operator--(int) { mint res = *this; --(*this); return res; } mint operator+() const { return *this; } mint operator-() const { return mint() - *this; } template mint pow(T n) const { assert(n >= 0); if constexpr (!is_void_v) return raw(PowInv::get().pow(val(), n)); mint x = *this, r = 1; while (n) { if (n & 1) r *= x; x *= x; n >>= 1; } return r; } mint inv() const { if constexpr (!is_void_v) { return raw(PowInv::get().inv(val())); } else 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(); } #if defined LOCAL or not defined FAST_IO friend std::istream &operator>>(std::istream &is, mint &x) { long long a; is >> a; x = a; return is; } friend std::ostream &operator<<(std::ostream &os, const mint &x) { os << x.val(); return os; } #else friend void rd1(mint &x) { long long a; fastio::rd1(a); x = a; } friend void wt1(const mint &x) { fastio::wt1(x.val()); } #endif }; }; 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 modint998244353 = static_modint32<998244353>; using modint1000000007 = static_modint32<1000000007>; using modint = dynamic_modint32<-1>; using modint61 = static_modint64<(1LL << 61) - 1>; using modint64 = dynamic_modint64<-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, PowInv>> : true_type {}; template struct is_static_modint, PowInv>> : 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; template struct has_mod : std::false_type { }; template struct has_mod> : std::true_type { }; template struct modint_less { bool operator()(const mint &a, const mint &b) const { if constexpr (is_modint_v) return a.val() < b.val(); else return a < b; } }; template struct modint_hash { auto operator()(const mint &x) const { if constexpr (is_modint_v) return std::hash{}(x.val()); else return std::hash{}(x); } }; // https://github.com/miscalculation53/library/tree/wip/math/modint/power_table.hpp /** * @brief 累乗テーブル * @docs docs/math/modint/power_table.md */ template struct PowerTable { private: decltype(mint::mod()) mod; mint base; vc pw; public: PowerTable() {} PowerTable(const mint &base) : mod(mint::mod()), base(base), pw(1, 1) {} void reserve(int n) { if (mod != mint::mod()) { mod = mint::mod(); pw = {1}; } int i = pw.size(); if (n < i) return; pw.resize(n + 1); for (; i <= n; i++) pw[i] = pw[i - 1] * base; } mint pow(int n) { assert(n >= 0); reserve(n); return pw[n]; } }; // https://github.com/miscalculation53/library/tree/wip/math/modint/binomial.hpp /** * @brief 二項係数 * @docs docs/math/modint/binomial.md */ template struct Binomial { private: inline static decltype(T::mod()) mod; public: inline static vc fac_, finv_, inv_; static void reserve(int n) { if constexpr (is_dynamic_modint_v) { if (mod != T::mod()) { mod = T::mod(); fac_ = {1, 1}, finv_ = {1, 1}, inv_ = {0, 1}; } } else { if (fac_.empty()) fac_ = {1, 1}, finv_ = {1, 1}, inv_ = {0, 1}; } if (n < SZ(fac_)) return; chmin(n, T::mod() - 1); int si = fac_.size(); fac_.resize(n + 1), finv_.resize(n + 1), inv_.resize(n + 1); repi(i, si, n + 1) { fac_[i] = fac_[i - 1] * T::raw(i); inv_[i] = -inv_[T::mod() % i] * T::raw(T::mod() / i); finv_[i] = finv_[i - 1] * inv_[i]; } } static T fac(int n) { assert(n >= 0); if (n >= T::mod()) return 0; reserve(n); return fac_[n]; } static T finv(int n) { assert(n < T::mod()); if (n < 0) return 0; reserve(n); return finv_[n]; } static T inv(int n) { n %= T::mod(); if (n < 0) n += T::mod(); assert(n != 0); reserve(n); return inv_[n]; } static T P(int n, int k) { if (n < k) return 0; if (n < 0 || k < 0) return 0; if (n >= T::mod()) return 0; reserve(n); return fac_[n] * finv_[n - k]; } static T C(int n, int k) { if (n < k) return 0; if (n < 0 || k < 0) return 0; if (n >= T::mod()) return 0; reserve(n); return fac_[n] * finv_[k] * finv_[n - k]; } static T H(int n, int k) { if (n == 0 && k == 0) return 1; return C(n + k - 1, k); } }; // https://github.com/miscalculation53/library/tree/wip/math/modint/stom.hpp /** * @brief string を mint に * @docs docs/math/modint/stom.md */ template mint stom(string s) { mint res = 0; fec(c : s) { res *= 10; res += c - '0'; } return res; } // https://github.com/miscalculation53/library/tree/wip/math/modint/to_rational.hpp // https://github.com/miscalculation53/library/tree/wip/math/svp2d.hpp /** * @brief 最短ベクトル問題($2$ 次元) * @docs docs/math/svp2d.md */ // 入力: a = (a1, a2), b = (b1, b2) // 出力: 整数 x, y を用いて v = (v1, v2) = xa + yb と表されるベクトルのうち、 // ノルムの 2 乗 v1^2 + v2^2 が最小のもののひとつ (符号は不定) // T は整数、U は T の 2 倍のビット数 // 初期値の √2 倍が T から溢れないなら正しく動く template > pair svp2d(const pair &a, const pair &b) { assert((a != pair{0, 0} && b != pair{0, 0})); auto [a1, a2] = a; auto [b1, b2] = b; if constexpr (is_same_v && is_integral_v && is_signed_v && sizeof(T) <= sizeof(ll)) { // 各入力が [-2^30, 2^30) なら、中間値も ll に収まる。 constexpr ull offset = 1ULL << 30; if (((ull(a1) + offset) | (ull(a2) + offset) | (ull(b1) + offset) | (ull(b2) + offset)) < 2 * offset) return svp2d(a, b); } if ((U)a1 * a1 + (U)a2 * a2 < (U)b1 * b1 + (U)b2 * b2) swap(a1, b1), swap(a2, b2); while ((U)a1 * a1 + (U)a2 * a2 > (U)b1 * b1 + (U)b2 * b2) { swap(a1, b1), swap(a2, b2); const U dot = U(a1) * b1 + U(a2) * b2, norm = U(a1) * a1 + U(a2) * a2; U k = dot / norm, r = dot % norm; // 最も近い整数を選び、同距離なら大きい方を選ぶ。 if (r > 0 && r >= norm - r) ++k; if (r < 0 && -r > norm + r) --k; b1 = T(U(b1) - k * a1), b2 = T(U(b2) - k * a2); if (b1 == 0 && b2 == 0) return {a1, a2}; } return {a1, a2}; } /** * @brief modint を有理数に * @docs docs/math/modint/to_rational.md */ template pair mint_to_rat(const mint &x) { auto [p, q] = svp2d({x.val(), 1}, {mint::mod(), 0}); if (q < 0) p = -p, q = -q; return {p, q}; } #if __cplusplus >= 202302L 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; } }; struct mint_to_rat_fn { template requires requires(const T &x) { { T::mod() } -> std::convertible_to; { x.val() } -> std::convertible_to; } constexpr auto operator()(const T &x) const -> rat_value { auto [p, q] = mint_to_rat(x); return {p, q}; } }; template struct rat_object { T value; }; 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 { template static constexpr bool can_convert() { #ifdef LOCAL if constexpr (_detail::is_value_with_command) return can_convert().value)>(); else #endif if constexpr (std::ranges::range) { using Ref = std::ranges::range_reference_t; if constexpr (std::same_as, std::remove_cvref_t>) return false; else return can_convert(); } else if constexpr (requires { std::tuple_size>::value; }) { return [](std::index_sequence) { if constexpr ((requires { std::get(std::declval()); } && ...)) return (can_convert(std::declval()))>() || ...); else return false; }(std::make_index_sequence>>{}); } else return std::invocable; } // 範囲アダプタ同士の合成は標準ライブラリの operator| に任せる。 template requires(!requires { std::views::all | std::declval(); }) constexpr decltype(auto) operator()(T &&t) const { if constexpr (!std::ranges::range && std::invocable) { return mint_to_rat_fn{}(std::forward(t)); } else if constexpr (can_convert()) { #ifdef LOCAL // 表示設定の一時オブジェクトを保持し、元の値は参照したままにする。 if constexpr (_detail::is_value_with_command) return rat_object>{std::forward(t)}; else #endif if constexpr (!std::ranges::range || rat_is_map>) return rat_object{std::forward(t)}; else return std::forward(t) | std::views::transform(*this); } else { return std::forward(t); } } }; #ifdef LOCAL template string export_object_generic( const rat_object &object, const string &indent, size_t last_line_length, size_t current_depth, bool fail_on_newline, const _detail::export_command &command) { using U = std::remove_cvref_t; const auto &value = object.value; auto render = [&](const auto &x) { return _detail::export_var(x, indent, last_line_length, current_depth, fail_on_newline, command); }; if constexpr (_detail::is_value_with_command) return _detail::export_var(rat_closure{}(value.value), indent, last_line_length, current_depth, fail_on_newline, value.command); else if constexpr (!rat_closure::can_convert()) return render(value); else if constexpr (_detail::is_map) return _detail::export_library_map(value, indent, last_line_length, current_depth, fail_on_newline, command, rat_closure{}); else if constexpr (!std::ranges::range && requires { std::tuple_size::value; }) return std::apply([&](const auto &...member) { return render(std::tuple{rat_object &>{member}...}); }, value); else return render(rat_closure{}(value)); } // rat() の前後にある表示設定を、cpp-dump と同じ順序で合成する。 template requires _detail::is_value_with_command auto operator|(rat_object object, const _detail::export_command &command) { return rat_object{object.value | command}; } template requires _detail::is_value_with_command auto operator|(rat_object object, _detail::export_command &&command) { return rat_object{object.value | std::move(command)}; } #endif template requires(!requires { std::views::all | std::declval(); }) constexpr decltype(auto) operator|(T &&t, rat_closure c) { return c(std::forward(t)); } constexpr rat_closure rat() { return rat_closure{}; } } #endif using mint = modint998244353; // using mint = modint1000000007; // using mint = static_modint<1000000000>; // using mint = modint; using bi = Binomial; void init() { oj(mt.seed(random_device()())); } // https://github.com/miscalculation53/library/tree/wip/ds/fenwick_tree/fenwick_tree.hpp // https://github.com/miscalculation53/library/tree/wip/algebra/algebra_basic_ops.hpp // https://github.com/miscalculation53/library/tree/wip/algebra/algebra_base.hpp /** * @brief 代数的構造の struct(基本) * @docs docs/algebra/algebra_base.md */ template struct Monoid { using S = S_; static constexpr auto op = op_; static constexpr auto e = e_; }; template struct Group { using S = S_; static constexpr auto op = op_; static constexpr auto e = e_; static constexpr auto inv = inv_; }; template struct SemiRing { using S = S_; static constexpr auto add = add_; static constexpr auto e0 = e0_; static constexpr auto mul = mul_; static constexpr auto e1 = e1_; }; template struct Ring { using S = S_; static constexpr auto add = add_; static constexpr auto e0 = e0_; static constexpr auto minus = minus_; static constexpr auto mul = mul_; static constexpr auto e1 = e1_; }; template struct Field { using S = S_; static constexpr auto add = add_; static constexpr auto e0 = e0_; static constexpr auto minus = minus_; static constexpr auto mul = mul_; static constexpr auto e1 = e1_; static constexpr auto inv = inv_; }; template struct HasInverse : false_type {}; template struct HasInverse()))>> : true_type {}; template struct OppositeMonoid { using S = typename M::S; static constexpr S op(const S &a, const S &b) { return M::op(b, a); } static constexpr auto e = M::e; }; template struct OppositeGroup { using S = typename G::S; static constexpr S op(const S &a, const S &b) { return G::op(b, a); } static constexpr auto e = G::e; static constexpr auto inv = G::inv; }; template struct NormalAndOppositeMonoid { struct S { typename M::S normal; typename M::S opposite; S() {} template , std::nullptr_t> = nullptr> S(Args &&...args) : normal(std::forward(args)...), opposite(normal) {} S rev() { return {opposite, normal}; } S(const typename M::S &normal, const typename M::S &opposite) : normal(normal), opposite(opposite) {} }; static constexpr S op(const S &a, const S &b) { return {M::op(a.normal, b.normal), M::op(b.opposite, a.opposite)}; } static constexpr S e() { return {M::e(), M::e()}; } }; template struct NormalAndOppositeGroup { struct S { typename G::S normal; typename G::S opposite; S() {} template , std::nullptr_t> = nullptr> S(Args &&...args) : normal(std::forward(args)...), opposite(normal) {} S rev() { return {opposite, normal}; } S(const typename G::S &normal, const typename G::S &opposite) : normal(normal), opposite(opposite) {} }; static constexpr S op(const S &a, const S &b) { return {G::op(a.normal, b.normal), G::op(b.opposite, a.opposite)}; } static constexpr S e() { return {G::e(), G::e()}; } static constexpr S inv(const S &a) { return {G::inv(a.normal), G::inv(a.opposite)}; } }; template using MonoidOfSemiRingAdd = Monoid; template using MonoidOfSemiRingMul = Monoid; template using GroupOfRingAdd = Group; template using GroupOfFieldMul = Group; // Madd は可換 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; }; // Gadd は可換 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; }; // Gadd, Gmul は可換 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; }; #ifdef LOCAL CPP_DUMP_DEFINE_EXPORT_OBJECT_GENERIC(normal, opposite); #endif /** * @brief 代数的構造(四則演算と min, max) * @docs docs/algebra/algebra_basic_ops.md */ 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; 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 1; } }; 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 () * declval())> static constexpr S pow(const S &a, I k) { return a * k; } }; 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 SemiRingMinPlus = SemiRingFromMonoidMonoid, MonoidAdd>; template using SemiRingMaxPlus = SemiRingFromMonoidMonoid, MonoidAdd>; template using RingAddSubMul = RingFromGroupMonoid, MonoidMul>; template using FieldAddSubMulDiv = FieldFromGroupGroup, GroupMulDiv>; /** * @brief Fenwick Tree * @docs docs/ds/fenwick_tree/fenwick_tree.md */ // G は可換群 (prefix だけなどであれば可換モノイドでも OK) template struct FenwickTree { using S = typename G::S; private: int n; vc dat; public: FenwickTree() {} FenwickTree(int n) : n(n), dat(n + 1, G::e()) {} FenwickTree(const vc &v) : FenwickTree(v.size()) { repi(i, n) dat[i + 1] = v[i]; repi(i, 1, n + 1) { int p = i + (i & -i); if (p <= n) dat[p] = G::op(dat[p], dat[i]); } } template I size() const { return n; } // [0, r) // 逆元は必要ない S sum(int r) const { assert(0 <= r && r <= n); S s = G::e(); while (r > 0) { s = G::op(s, dat[r]); r -= r & -r; } return s; } // [l, r) // 逆元が必要 S sum(int l, int r) const { assert(0 <= l && l <= r && r <= n); // 両端の prefix に共通する節点までで打ち切る。 S right = G::e(), left = G::e(); while (l < r) { right = G::op(right, dat[r]); r -= r & -r; } while (r < l) { left = G::op(left, dat[l]); l -= l & -l; } return G::op(G::inv(left), right); } // 逆元が必要 S get(int i) const { assert(0 <= i && i < n); return sum(i, i + 1); } // 逆元は必要ない void add(int i, S x) { assert(0 <= i && i < n); i++; while (i <= n) { dat[i] = G::op(dat[i], x); i += i & -i; } } // 逆元が必要 void set(int i, S x) { add(i, G::op(G::inv(get(i)), x)); } private: pair prefix_lt_raw(S w) const { if (w <= G::e()) return {-1, G::e()}; int k = bit_floor(n); int x = 0; S v = G::e(); while (k > 0) { if (x + k <= n) { S nv = G::op(v, dat[x + k]); if (nv < w) v = nv, x += k; } k >>= 1; } return {x, v}; } public: // 非負の通常の整数和を想定する。該当する接頭辞がなければ nullopt。 template optional> lt_max_id_sum(S w) const { auto result = prefix_lt_raw(w); return result.first < 0 ? nullopt : optional>(pair{I(result.first), std::move(result.second)}); } template I lt_max(S w) const { static_assert(numeric_limits::is_signed, "index type must represent -1"); return prefix_lt_raw(w).first; } template I geq_min(S w) const { return prefix_lt_raw(w).first + 1; } template I leq_max(S w) const { return lt_max(w + 1); } template I gt_min(S w) const { return geq_min(w + 1); } // 要素が [0, size()) の多重集合を管理するのに使ったとき、多重集合のサイズを返す inline S size_of_multiset() const { return sum(n); } // 要素が [0, size()) の多重集合を管理するのに使ったとき、値 x は何番目から何番目か [l, r) inline pair order_in_multiset(int x) const { return {sum(x), sum(x + 1)}; } // 要素が [0, size()) の多重集合を管理するのに使ったとき、k 番目の値 // 範囲外の k は nullopt template inline optional kth_in_multiset(S k) const { if (k < 0) return nullopt; int key = leq_max(k); return key < n ? optional(I(key)) : nullopt; } // 要素が [0, size()) の多重集合を管理するのに使ったとき、x 未満で最大の要素の順位 (左番兵 -1) template inline I lt_max_in_multiset(T x) const { static_assert(numeric_limits::is_signed, "index type must represent -1"); S count = sum(clamp(x, T(0), T(n))); return I(count) - 1; } // x 以下で最大の要素の順位 (左番兵 -1) template inline I leq_max_in_multiset(T x) const { return lt_max_in_multiset(x + 1); } // x 以上で最小の要素の順位 (右番兵 size_of_multiset()) template inline I geq_min_in_multiset(T x) const { return sum(clamp(x, T(0), T(n))); } // x 超過で最小の要素の順位 (右番兵 size_of_multiset()) template inline I gt_min_in_multiset(T x) const { return geq_min_in_multiset(x + 1); } vc content() const { vc res(dat.begin() + 1, dat.end()); // 線形構築で親へ加えた値を、逆順に取り除く。 for (int i = n; i > 0; i--) { const int p = i + (i & -i); if (p <= n) res[p - 1] = G::op(res[p - 1], G::inv(res[i - 1])); } return res; } }; void main2() { LL(N); VEC(ll, N, A); FenwickTree> fw(A); set st{-1, N}; rep(i, N) if (A[i] == 0) st.insert(i); LL(Q); rep(_, Q) { LL(x, b); x--; fw.set(x, b); st.erase(x); ll l = *lt_max(st, x).value(), r = *gt_min(st, x).value(); dump(fw.content(), st, l, r); PRINT(fw.sum(l + 1, r)); fw.set(x, A[x]); if (A[x] == 0) st.insert(x); } } 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"; #ifdef LOCAL cerr << s; #endif /* コードテストで確認する際にコメントアウトを外す cerr << val; //*/ }; #if defined FAST_IO and not defined LOCAL CERR("\n[FAST_IO]\n\n", "32"); #endif #if defined FAST_CIO and not defined LOCAL CERR("\n[FAST_CIO]\n\n", "32"); #endif #if (defined FAST_CIO or defined INTERACTIVE) and not defined LOCAL cin.tie(nullptr); ios::sync_with_stdio(false); #endif cout << fixed << setprecision(20); init(); #ifdef LOCAL test(); #endif #if defined AOJ_TESTCASE or (defined LOCAL and defined SINGLE_TESTCASE) CERR("\n[AOJ_TESTCASE]\n\n", "35"); while (true) { dump("new testcase"); main2(); } #elif defined SINGLE_TESTCASE CERR("\n[SINGLE_TESTCASE]\n\n", "36"); main2(); #elif defined MULTI_TESTCASE CERR("\n[MULTI_TESTCASE]\n\n", "33"); local(while (true)) { dump("T"); IN(uint, T); while (T--) { dump("new testcase"); main2(); } } #endif } }; Main main_dummy; int main() {} // submit-code: hash=9ddc5a74935f source=9ddc5a74935fabdd5eea58ffda60d9925c26c5e118b679b98481a2f3b0c803ab input=fd06333d752761ecc7960434c8e4661027f827110aad04745e2e54772b1aa162 body=fd06333d752761ecc7960434c8e4661027f827110aad04745e2e54772b1aa162