#define SINGLE_TESTCASE #define FAST_CIO #define INF 4'000'000'000'000'000'037LL #define EPS 1e-11 #include using namespace std; namespace { using ll = long long; using uint = unsigned int; using ull = unsigned long long; using pll = pair; #define vc vector template using vvc = vc>; using vpll = vc; using vstr = vc; #ifdef __SIZEOF_INT128__ using i128 = __int128_t; using u128 = __uint128_t; i128 stoi128(const string &s) { i128 res = 0; if (s.front() == '-') { for (int i = 1; i < (int)s.size(); i++) res = 10 * res + s[i] - '0'; res = -res; } else { for (auto &&c : s) res = 10 * res + c - '0'; } return res; } string i128tos(i128 x) { if (x == 0) return "0"; string sign = "", res = ""; if (x < 0) x = -x, sign = "-"; while (x > 0) { res += '0' + x % 10; x /= 10; } reverse(res.begin(), res.end()); return sign + res; } istream &operator>>(istream &is, i128 &a) { string s; is >> s; a = stoi128(s); return is; } ostream &operator<<(ostream &os, const i128 &a) { os << i128tos(a); return os; } #endif #define cauto const auto #define overload4(_1, _2, _3, _4, name, ...) name #define rep1(i, n) for (ll i = 0, nnnnn = ll(n); i < nnnnn; i++) #define rep2(i, l, r) for (ll i = ll(l), rrrrr = ll(r); i < rrrrr; i++) #define rep3(i, l, r, d) for (ll i = ll(l), rrrrr = ll(r), ddddd = ll(d); ddddd > 0 ? i < rrrrr : i > rrrrr; i += d) #define rep(...) overload4(__VA_ARGS__, rep3, rep2, rep1)(__VA_ARGS__) #define repi1(i, n) for (int i = 0, nnnnn = int(n); i < nnnnn; i++) #define repi2(i, l, r) for (int i = int(l), rrrrr = int(r); i < rrrrr; i++) #define repi3(i, l, r, d) for (int i = int(l), rrrrr = int(r), ddddd = int(d); ddddd > 0 ? i < rrrrr : i > rrrrr; i += d) #define repi(...) overload4(__VA_ARGS__, repi3, repi2, repi1)(__VA_ARGS__) #define fe(...) for (auto __VA_ARGS__) #define fec(...) for (cauto &__VA_ARGS__) #define fem(...) for (auto &__VA_ARGS__) template inline bool chmin(T &a, U b) { return a > b ? a = b, true : false; } template inline constexpr T divfloor(U a, V b) { return T(a) / T(b) - (T(a) % T(b) && (T(a) ^ T(b)) < 0); } template inline constexpr T divceil(U a, V b) { return T(a) / T(b) + (T(a) % T(b) && (T(a) ^ T(b)) >= 0); } template inline constexpr T safemod(U a, V b) { return T(a) - T(b) * divfloor(a, b); } 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, INF); } 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, INF); } #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) ([&](auto x) { return fx; }) #define GEN_VEC(n, i, fi) (gen_vec(n, LMD(i, fi))) template auto dvec(const V (&sz)[d], const T &init) { if constexpr (i < d) return vc(sz[i], dvec(sz, init)); else return init; } 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; } template vc permuted(const vc &p, const vc &q, const vc &...rs) { return permuted(permuted(p, q), rs...); } #if __cplusplus < 202002L #else #endif template void unique(V &v) { v.erase(std::unique(ALL(v)), v.end()); } template void rotate(V &v, U k) { const U n = v.size(); k = (k % n + n) % n; std::rotate(v.begin(), v.begin() + k, v.end()); } 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; } const vpll DRULgrid = {{1, 0}, {0, 1}, {-1, 0}, {0, -1}}; const vpll DRULplane = {{0, -1}, {1, 0}, {0, 1}, {-1, 0}}; 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; #if __cplusplus < 202002L namespace internal { }; #define DEFAULT_COMP less<> #else #define DEFAULT_COMP ranges::less #endif #if __cplusplus < 202002L #else inline constexpr ll bit_width(ll x) { return std::bit_width((ull)x); } inline constexpr ll bit_floor(ll x) { return std::bit_floor((ull)x); } inline constexpr ll bit_ceil(ll x) { return std::bit_ceil((ull)x); } inline constexpr ll countr_zero(ll x) { assert(x != 0); return std::countr_zero((ull)x); } inline constexpr ll popcount(ll x) { return std::popcount((ull)x); } inline constexpr bool has_single_bit(ll x) { return std::has_single_bit((ull)x); } #endif #define dump(...) #define local(...) #define oj(...) __VA_ARGS__ #define local_oj(a, b) (b) 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; } namespace internal { template void CIN(Ts &...a) { (cin >> ... >> a); } template void READnodump(Ts &...a) { CIN(a...); } template void READVECnodump(int n, vc &v) { v.resize(n); READnodump(v); } template void READVECnodump(int n, vc &v, vc &...vs) { READVECnodump(n, v), READVECnodump(n, vs...); } template void READVEC2nodump(int n, int m, vvc &v) { v.assign(n, vc(m)); READnodump(v); } template void READVEC2nodump(int n, int m, vvc &v, vvc &...vs) { READVEC2nodump(n, m, v), READVEC2nodump(n, m, vs...); } template void READJAGnodump(int n, vvc &v) { v.resize(n); repi(i, n) { int k; READnodump(k); READVECnodump(k, v[i]); } } template void READJAGnodump(int n, vvc &v, vvc &...vs) { READJAGnodump(n, v), READJAGnodump(n, 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__) #define ENDL '\n' template ostream &operator<<(ostream &os, const pair &p) { os << p.first << ' ' << p.second; return os; } namespace internal { template void cout_tuple(ostream &os, const T &t) { if constexpr (N < std::tuple_size::value) { if constexpr (N > 0) { os << ' '; } const auto x = std::get(t); os << x; cout_tuple(os, t); } } }; // namespace internal template ostream &operator<<(ostream &os, const tuple &t) { internal::cout_tuple(os, 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 COUTP() { cout << ENDL; } template void COUTP(const T &a) { cout << a << ENDL; } template void COUTP(const T &a, const Ts &...b) { cout << a; (cout << ... << (cout << ' ', b)); cout << ENDL; } }; // namespace internal #define WRITE internal::COUTW #define PRINT internal::COUTP #define PRINTEXIT(...) do { PRINT(__VA_ARGS__); exit(0); } while (false) #define PRINTRETURN(...) do { PRINT(__VA_ARGS__); return; } while (false) #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 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; } namespace internal { template auto tuple_add_impl(A &a, const B &b, const index_sequence) { ((get(a) += get(b)), ...); 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 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> top(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> top(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> top(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> top(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 top(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 top(const tuple...> &tv) { size_t n = get<0>(tv).size(); apply([&](auto &...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; } mt19937_64 mt; namespace internal { constexpr ll powmod32_constexpr(ll x, ll n, int m) { if (m == 1) return 0; uint _m = (uint)m; ull r = 1; ull y = safemod(x, m); while (n) { if (n & 1) r = (r * y) % _m; y = (y * y) % _m; n >>= 1; } return r; } constexpr bool isprime32_constexpr(int 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; constexpr ll bases[3] = {2, 7, 61}; for (ll a : bases) { ll t = d; ll y = powmod32_constexpr(a, t, n); while (t != n - 1 && y != 1 && y != n - 1) { y = y * y % n; t <<= 1; } if (y != n - 1 && t % 2 == 0) return false; } return true; } template constexpr bool isprime32 = isprime32_constexpr(n); struct barrett32 { uint m; ull im; explicit barrett32(uint m) : m(m), im((ull)(-1) / m + 1) {} uint umod() const { return m; } uint mul(uint a, uint b) const { ull z = a; z *= b; ull x = (ull)((u128(z)*im) >> 64); ull y = x * m; return (uint)(z - y + (z < y ? m : 0)); } }; } namespace internal { #define REF static_cast(*this) #define CREF static_cast(*this) #define VAL *static_cast(this) template struct modint_base { mint &operator+=(const mint &rhs) { mint &self = REF; self._v += rhs._v; if (self._v >= self.umod()) self._v -= self.umod(); return self; } mint &operator-=(const mint &rhs) { mint &self = REF; self._v -= rhs._v; if (self._v >= self.umod()) self._v += self.umod(); return self; } mint &operator/=(const mint &rhs) { mint &self = REF; return self = self * rhs.inv(); } mint &operator++() { mint &self = REF; self._v++; if (self._v == self.umod()) self._v = 0; return self; } mint &operator--() { mint &self = REF; if (self._v == 0) self._v = self.umod(); self._v--; return self; } mint operator++(int) { mint res = VAL; ++REF; return res; } mint operator--(int) { mint res = VAL; --REF; return res; } mint operator+() const { return VAL; } mint operator-() const { return mint() - VAL; } mint pow(ll n) const { assert(n >= 0); mint x = VAL, r = 1; while (n) { if (n & 1) r *= x; x *= x; n >>= 1; } return r; } 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 mint(lhs).eq(rhs); } friend bool operator!=(const mint &lhs, const mint &rhs) { return mint(lhs).neq(rhs); } private: bool eq(const mint &rhs) { return REF._v == rhs._v; } bool neq(const mint &rhs) { return REF._v != rhs._v; } }; } template , T>, int> = 0> istream &operator>>(istream &is, T &x) { ll a; is >> a; x = a; return is; } template , T>, int> = 0> ostream &operator<<(ostream &os, const T &x) { os << x.val(); return os; } template constexpr tuple extgcd(T a, T b) { if (a == 0 && b == 0) return {0, 0, 0}; T x1 = 1, y1 = 0, z1 = a; T x2 = 0, y2 = 1, z2 = b; while (z2 != 0) { T q = z1 / z2; tie(x1, x2) = make_pair(x2, x1 - q * x2); tie(y1, y2) = make_pair(y2, y1 - q * y2); tie(z1, z2) = make_pair(z2, z1 - q * z2); } if (z1 < 0) x1 = -x1, y1 = -y1, z1 = -z1; return {z1, x1, y1}; } template struct static_modint : internal::modint_base> { using mint = static_modint; private: friend struct internal::modint_base>; uint _v; static constexpr uint umod() { return m; } static constexpr bool prime = internal::isprime32; public: static constexpr int mod() { return m; } static mint raw(int v) { mint x; x._v = v; return x; } static_modint() : _v(0) {} template ::value>> static_modint(T v) { if constexpr (is_signed_v) { ll x = (ll)(v % (ll)(umod())); if (x < 0) x += umod(); _v = (uint)x; } else { _v = (uint)(v % umod()); } } int val() const { return (int)_v; } mint& operator*=(const mint &rhs) { ull z = _v; z *= rhs._v; _v = (uint)(z % umod()); return *this; } mint inv() const { if (prime) { assert(_v != 0); return CREF.pow(umod() - 2); } else { auto [g, x, y] = extgcd(_v, m); assert(g == 1); return x; } } }; template struct dynamic_modint : internal::modint_base> { using mint = dynamic_modint; private: friend struct internal::modint_base>; uint _v; static internal::barrett32 bt; static uint umod() { return bt.umod(); } public: static int mod() { return (int)(bt.umod()); } static mint raw(int v) { mint x; x._v = v; return x; } dynamic_modint() : _v(0) {} template ::value>> dynamic_modint(T v) { if constexpr (is_signed_v) { ll x = (ll)(v % (ll)(umod())); if (x < 0) x += umod(); _v = (uint)x; } else { _v = (uint)(v % umod()); } } int val() const { return (int)_v; } mint& operator*=(const mint &rhs) { _v = bt.mul(_v, rhs._v); return *this; } mint inv() const { auto [g, x, y] = extgcd(_v, mod()); assert(g == 1); return x; } }; template internal::barrett32 dynamic_modint::bt(998244353); using modint998244353 = static_modint<998244353>; template struct is_static_modint : false_type {}; template struct is_static_modint> : true_type {}; template inline constexpr bool is_static_modint_v = is_static_modint::value; template struct is_dynamic_modint : false_type {}; template struct is_dynamic_modint> : true_type {}; template inline constexpr bool is_dynamic_modint_v = is_dynamic_modint::value; template inline constexpr bool is_modint_v = is_static_modint_v || is_dynamic_modint_v; template struct has_mod : false_type {}; template struct has_mod().mod)>> : true_type {}; 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) { reserve(n); return pw[n]; } }; template struct Binomial { private: static decltype(T::mod()) mod; static vc fac_, finv_, inv_; public: static void reserve(int n) { if (mod != T::mod()) { mod = T::mod(); fac_ = {1, 1}, finv_ = {1, 1}, inv_ = {0, 1}; } int i = fac_.size(); chmin(n, T::mod() - 1); if (n < i) return; fac_.resize(n + 1), finv_.resize(n + 1), inv_.resize(n + 1); for (; i <= n; i++) { 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 inv(T n) { assert(n != 0); reserve(n.val()); return inv_[n.val()]; } }; template decltype(T::mod()) Binomial::mod{}; template vc Binomial::fac_{}; template vc Binomial::finv_{}; template vc Binomial::inv_{}; using mint = modint998244353; using bi = Binomial; void init() { oj(mt.seed(random_device()())); } template struct FenwickTree { private: int N; vector bit; public: FenwickTree(int n = 0) { N = n; bit.assign(N + 1, 0); } T sum(int i) const { assert(0 <= i && i <= N); T s = 0; while (i > 0) { s += bit[i]; i -= i & -i; } return s; } T sum(int l, int r) const { assert(l <= r); return sum(r) - sum(l); } T get(int i) const { return sum(i, i + 1); } void add(int i, T x) { assert(0 <= i && i < N); i++; while (i <= N) { bit[i] += x; i += i & -i; } } }; template struct FactorialSystem : vector { using vector::vector; using vector::operator=; using F = FactorialSystem; void carry() { if ((*this).empty()) return; for (int i = 0; i < (int)(*this).size() - 1; i++) { (*this)[i + 1] += divfloor((*this)[i], i + 1); (*this)[i] = safemod((*this)[i], i + 1); } for (int i = (int)(*this).size() - 1; (*this)[i] > i; i++) { (*this).emplace_back(divfloor((*this)[i], i + 1)); (*this)[i] = safemod((*this)[i], i + 1); } } void zero_suppress() { while (!(*this).empty() && (*this).back() == 0) (*this).pop_back(); } bool operator==(const F &b) const { for (int i = max((*this).size(), b.size()) - 1; i >= 0; i--) { T ai = i < (int)(*this).size() ? (*this)[i] : 0; T bi = i < (int)b.size() ? b[i] : 0; if (ai != bi) return false; } return true; } bool operator<(const F &b) const { for (int i = max((*this).size(), b.size()) - 1; i >= 0; i--) { T ai = i < (int)(*this).size() ? (*this)[i] : 0; T bi = i < (int)b.size() ? b[i] : 0; if (ai < bi) return true; if (ai > bi) return false; } return false; } bool operator>(const F &b) const { for (int i = max((*this).size(), b.size()) - 1; i >= 0; i--) { T ai = i < (int)(*this).size() ? (*this)[i] : 0; T bi = i < (int)b.size() ? b[i] : 0; if (ai > bi) return true; if (ai < bi) return false; } return false; } bool operator!=(const F &b) const { return !((*this) == b); } bool operator<=(const F &b) const { return !((*this) > b); } bool operator>=(const F &b) const { return !((*this) < b); } F &operator+=(ll k) { if ((*this).empty()) (*this).resize(1); (*this).front() += k; carry(); return *this; } F operator+(ll k) const { return F(*this) += k; } F &operator-=(ll k) { return (*this) += (-k); } // assert は面倒なのでつけてません F operator-(ll k) const { return F(*this) -= k; } F &operator+=(const F &b) { if ((*this).size() < b.size()) (*this).resize(b.size()); for (int i = 0; i < (int)b.size(); i++) (*this)[i] += b[i]; carry(); return *this; } F operator+(const F &b) const { return F(*this) += b; } F &operator-=(const F &b) { assert((*this) >= b); if ((*this).size() < b.size()) (*this).resize(b.size()); for (int i = 0; i < (int)b.size(); i++) (*this)[i] -= b[i]; carry(); return *this; } F operator-(const F &b) const { return F(*this) -= b; } F &operator*=(ll k) { for (T &a : (*this)) a *= k; carry(); return *this; } F operator*(ll k) const { return F(*this) *= k; } pair divmod(ll k) const { F a(*this), res((*this).size(), 0); for (int i = (int)a.size() - 1; i > 0; i--) { res[i] = a[i] / k; a[i - 1] += a[i] % k * i; } return make_pair(res, a.front()); } F operator/(ll k) const { return (*this).divmod(k).first; } F &operator/=(ll k) { return (*this) = (*this) / k; } ll operator%(ll k) { return (*this).divmod(k).second; } }; template FactorialSystem perm_to_fac(const vector &P) { int N = P.size(); FenwickTree fw(N); for (int i = 0; i < N; i++) fw.add(i, 1); FactorialSystem fac(N); for (int i = 0; i < N; i++) { fac[N - 1 - i] = fw.sum(P[i]); fw.add(P[i], -1); } return fac; } #include template struct bigint { private: static constexpr int pw(const int &a, const int &b) { int res = 1; for (int i = 0; i < b; i++) res *= a; return res; } static constexpr int BASE = pw(base, digit); vector vec; bool is_nega = false; static constexpr int MOD1 = 1711276033, MOD2 = 1811939329; static constexpr int INV12 = atcoder::internal::inv_gcd(MOD1, MOD2).second; static constexpr ll MOD12 = (ll)MOD1 * (ll)MOD2; static constexpr ll mycrt(const ll &r1, const ll &r2) { ll t = safemod((r2 - r1) * INV12, MOD2); return r1 + t * MOD1; } static constexpr vector mycrt(const vector &vec1, const vector &vec2) { vector res(vec1.size()); for (int i = 0; i < (int)vec1.size(); i++) res[i] = mycrt(vec1[i], vec2[i]); return res; } void zero_suppress() { while (!vec.empty() && vec.back() == 0) vec.pop_back(); } void carry(int d = -1) { const int n = vec.size(); if (n == 0) return; if (d < 0) d = n - 1; for (int i = 0; i < d; i++) { vec[i + 1] += divfloor(vec[i], BASE); vec[i] = safemod(vec[i], BASE); } for (int i = d; !(0 <= vec[i] && vec[i] < BASE); i++) { vec.emplace_back(divfloor(vec[i], BASE)); vec[i] = safemod(vec[i], BASE); } } bigint &operator<<=(size_t k) { vec.insert(vec.begin(), k, 0); return *this; } bigint &operator>>=(size_t k) { vec.erase(vec.begin(), vec.begin() + min(vec.size(), k)); return *this; } bigint operator<<(size_t k) const { return bigint(*this) <<= k; } bigint operator>>(size_t k) const { return bigint(*this) >>= k; } public: bigint() {} bigint(string s) { if (!s.empty() && s[0] == '-') { is_nega = true; s.erase(s.begin()); } else is_nega = false; vec.resize(divceil((int)s.size(), digit)); for (int i = 0; i < (int)vec.size(); i++) { int r = (int)s.size() - i * digit; int l = max(0, r - digit); string t = s.substr(l, r - l); if (!t.empty()) vec[i] = stoi(t, 0, base); } zero_suppress(); } template bigint(const T &x, typename enable_if::value>::type* = 0) { *this = bigint(std::to_string(x)); } string to_string() const { const int n = vec.size(); if (n == 0) return "0"; string s = is_nega ? "-" : ""; auto my_to_string = [&](int x) -> string { if (base == 10) return std::to_string(x); string res = ""; while (x > 0) { int d = x % base; res += d >= 10 ? d - 10 + 'A' : d + '0'; x /= base; } reverse(res.begin(), res.end()); return res; }; s += my_to_string(vec.back()); for (int i = (int)vec.size() - 2; i >= 0; i--) { string t = my_to_string(vec[i]); s += string(digit - (int)t.size(), '0') + t; } return s; } auto operator<=>(const bigint &b) const { if (is_nega ^ b.is_nega) return is_nega ? -1 : 1; if (vec.size() != b.vec.size()) return ((vec.size() < b.vec.size()) ^ is_nega) ? -1 : 1; for (int i = (int)vec.size() - 1; i >= 0; i--) { if (vec[i] != b.vec[i]) return ((vec[i] < b.vec[i]) ^ is_nega) ? -1 : 1; } return 0; } bool operator==(const bigint &b) const { return (*this).vec == b.vec && (*this).is_nega == b.is_nega; } bigint operator-() const { bigint res(*this); res.is_nega ^= 1; return res; } bigint abs() const { bigint res(*this); res.is_nega = false; return res; } bigint &operator+=(const bigint &b) { if (is_nega ^ b.is_nega) return (*this) -= (-b); if (vec.size() < b.vec.size()) vec.resize(b.vec.size()); for (int i = 0; i < (int)b.vec.size(); i++) vec[i] += b.vec[i]; carry(); return *this; } bigint &operator-=(const bigint &b) { if (is_nega ^ b.is_nega) return (*this) += (-b); if (vec.size() < b.vec.size()) vec.resize(b.vec.size()); bool rev = abs() < b.abs(); for (int i = 0; i < (int)b.vec.size(); i++) vec[i] = rev ? b.vec[i] - vec[i] : vec[i] - b.vec[i]; carry(); zero_suppress(); is_nega ^= rev; return *this; } bigint &operator*=(const bigint &b) { vector vec1 = atcoder::convolution(vec, b.vec); vector vec2 = atcoder::convolution(vec, b.vec); vec = mycrt(vec1, vec2); carry(); zero_suppress(); is_nega ^= b.is_nega; return *this; } bigint operator+(const bigint &b) const { return bigint(*this) += b; } bigint operator-(const bigint &b) const { return bigint(*this) -= b; } bigint operator*(const bigint &b) const { return bigint(*this) *= b; } bigint inv(int d) const { assert(!vec.empty()); bigint a(abs()), b; { int l = 0, r = BASE + 1; while (r - l > 1) { int m = l + (r - l) / 2; bigint c = a * m; if (c.vec.size() > a.vec.size()) r = m; else l = m; } b = l; } bigint c, c2; const bigint TWO(2); for (int k = 1;; k = min(2 * k, d)) { c = a * b; if ((int)b.vec.size() >= d + 1 && (int)c.vec.size() == (int)a.vec.size() + (int)b.vec.size() - 1) { c2 = c + a; if (c2.vec.size() == a.vec.size() + b.vec.size()) { if (c2 == bigint(1) << ((int)c2.vec.size() - 1)) b += 1; break; } } b *= (TWO << ((int)a.vec.size() + (int)b.vec.size() - 1)) - c; if ((int)b.vec.size() >= k + 1) b >>= ((int)b.vec.size() - k - 1); } b >>= 1; if (is_nega) b = -b; return b; } bigint operator/(const bigint &b) const { assert(!b.vec.empty()); bigint a = abs(); if (a < b) return bigint(0); bigint ib = b.abs().inv((int)a.vec.size()); bigint c = (a * ib) >> ((int)a.vec.size() + (int)b.vec.size() - 1); if (a >= b.abs() * (c + 1)) c += 1; c.is_nega = is_nega ^ b.is_nega; return c; } pair divmod(const bigint &b) const { bigint q = *this / b; bigint r = *this - q * b; return make_pair(q, r); } bigint operator%(const bigint &b) const { return divmod(b).second; } bigint &operator/=(const bigint &b) { return *this = *this / b; } bigint &operator%=(const bigint &b) { return *this = *this % b; } }; template istream &operator>>(istream &is, bigint &a) { string s; is >> s; a = bigint(s); return is; } template ostream &operator<<(ostream &os, const bigint &a) { os << a.to_string(); return os; } void main2() { LL(N); VEC(int, N, P); offset(P, -1); auto fac = perm_to_fac(P); dump(fac); bigint ans = 0; rep(i, SZ(fac) - 1, 0, -1) ans = (ans + fac[i]) * i; PRINT(ans + 1); } void test() { } template struct Main { Main() { cauto CERR = [](string val, string color) { string s = "\033[" + color + "m" + val + "\033[m"; /* コードテストで確認する際にコメントアウトを外す cerr << val; //*/ }; CERR("\n[FAST_CIO]\n\n", "32"); cin.tie(0); ios::sync_with_stdio(false); cout << fixed << setprecision(20); init(); CERR("\n[SINGLE_TESTCASE]\n\n", "36"); main2(); } }; Main main_dummy; } int main() {}