#line 2 "library/template/template.hpp" #include using namespace std; #line 2 "library/template/macro.hpp" #define rep(i, a, b) for (int i = (a); i < (int)(b); i++) #define rrep(i, a, b) for (int i = (int)(b) - 1; i >= (a); i--) #define ALL(v) (v).begin(), (v).end() #define UNIQUE(v) sort(ALL(v)), (v).erase(unique(ALL(v)), (v).end()) #define SZ(v) (int)v.size() #define MIN(v) *min_element(ALL(v)) #define MAX(v) *max_element(ALL(v)) #define LB(v, x) int(lower_bound(ALL(v), (x)) - (v).begin()) #define UB(v, x) int(upper_bound(ALL(v), (x)) - (v).begin()) #define YN(b) cout << ((b) ? "YES" : "NO") << "\n"; #define Yn(b) cout << ((b) ? "Yes" : "No") << "\n"; #define yn(b) cout << ((b) ? "yes" : "no") << "\n"; #line 6 "library/template/template.hpp" #line 2 "library/template/util.hpp" using uint = unsigned int; using ll = long long int; using ull = unsigned long long; using i128 = __int128_t; using u128 = __uint128_t; template using priority_queue_asc = priority_queue, greater>; template S SUM(const vector& a) { return accumulate(ALL(a), S(0)); } template inline bool chmin(T1& a, T2 b) { if (a > b) { a = b; return true; } return false; } template inline bool chmax(T1& a, T2 b) { if (a < b) { a = b; return true; } return false; } template inline bool chmin_opt(optional& a, T2 b) { if (!a || a > b) { a = b; return true; } return false; } template inline bool chmax_opt(optional& a, T2 b) { if (!a || a < b) { a = b; return true; } return false; } template int popcnt(T x) { return __builtin_popcountll(x); } template int topbit(T x) { return (x == 0 ? -1 : 63 - __builtin_clzll(x)); } template int lowbit(T x) { return (x == 0 ? -1 : __builtin_ctzll(x)); } #line 8 "library/template/template.hpp" #line 2 "library/template/inout.hpp" struct Fast { Fast() { cin.tie(nullptr); ios_base::sync_with_stdio(false); cout << fixed << setprecision(15); } } fast; ostream& operator<<(ostream& os, __uint128_t x) { char buf[40]; size_t k = 0; while (x > 0) buf[k++] = (char)(x % 10 + '0'), x /= 10; if (k == 0) buf[k++] = '0'; while (k) os << buf[--k]; return os; } ostream& operator<<(ostream& os, __int128_t x) { return x < 0 ? (os << '-' << (__uint128_t)(-x)) : (os << (__uint128_t)x); } template ostream& operator<<(ostream& os, const array& a); template istream& operator>>(istream& is, pair& p) { return is >> p.first >> p.second; } template ostream& operator<<(ostream& os, const pair& p) { return os << p.first << " " << p.second; } template istream& operator>>(istream& is, vector& a) { for (auto& v : a) is >> v; return is; } template ostream& operator<<(ostream& os, const vector& a) { for (auto it = a.begin(); it != a.end();) { os << *it; if (++it != a.end()) os << " "; } return os; } template ostream& operator<<(ostream& os, const array& a) { for (auto it = a.begin(); it != a.end();) { os << *it; if (++it != a.end()) os << " "; } return os; } template ostream& operator<<(ostream& os, const set& st) { os << "{"; for (auto it = st.begin(); it != st.end();) { os << *it; if (++it != st.end()) os << ","; } os << "}"; return os; } template ostream& operator<<(ostream& os, const map& mp) { os << "{"; for (auto it = mp.begin(); it != mp.end();) { os << it->first << ":" << it->second; if (++it != mp.end()) os << ","; } os << "}"; return os; } void in() {} template void in(T& t, U&... u) { cin >> t; in(u...); } template void in_zip(int n, T&... t) { assert(n >= 0 && ((size(t) >= static_cast(n)) && ...)); for (int i = 0; i < n; i++) in(t[i]...); } void out() { cout << "\n"; } template void out(const T& t, const U&... u) { cout << t; if (sizeof...(u)) cout << sep; out(u...); } template void out_opt(const optional& opt, const U& fallback, ostream& os = cout) { if (opt.has_value()) os << opt.value(); else os << fallback; os << "\n"; } template void out_opt(const vector>& vec, const U& fallback, ostream& os = cout) { for (auto it = vec.begin(); it != vec.end();) { if ((*it).has_value()) os << (*it).value(); else os << fallback; if (++it != vec.end()) os << " "; } os << "\n"; } namespace IO { template T read(U&&... u) { T t = T(forward(u)...); in(t); return t; } namespace Graph { vector> unweighted(int n, int m, bool directed = false, int offset = 1) { vector> g(n); for (int i = 0; i < m; i++) { int u, v; cin >> u >> v; u -= offset, v -= offset; g[u].push_back(v); if (!directed) g[v].push_back(u); } return g; } template vector>> weighted(int n, int m, bool directed = false, int offset = 1) { vector>> g(n); for (int i = 0; i < m; i++) { int u, v; T w; cin >> u >> v >> w; u -= offset, v -= offset; g[u].push_back({v, w}); if (!directed) g[v].push_back({u, w}); } return g; } } // namespace Graph namespace Tree { vector> unweighted(int n, bool directed = false, int offset = 1) { return Graph::unweighted(n, n - 1, directed, offset); } template vector>> weighted(int n, bool directed = false, int offset = 1) { return Graph::weighted(n, n - 1, directed, offset); } vector> rooted(int n, bool to_root = true, bool to_leaf = true, int offset = 1) { vector> g(n); for (int i = 1; i < n; i++) { int p; cin >> p; p -= offset; if (to_root) g[i].push_back(p); if (to_leaf) g[p].push_back(i); } return g; } } // namespace Tree } // namespace IO #line 10 "library/template/template.hpp" #line 2 "library/template/debug.hpp" #ifdef LOCAL #define debug 1 #define show(...) _show(0, #__VA_ARGS__, __VA_ARGS__) #else #define debug 0 #define show(...) true #endif template void _show(int, T) { cerr << '\n'; } template void _show(int i, const T1& a, const T2& b, const T3&... c) { for (; a[i] != ',' && a[i] != '\0'; i++) cerr << a[i]; cerr << ":" << b << " "; _show(i + 1, a, c...); } #line 2 "main.cpp" #line 2 "library/math/util.hpp" namespace Math { template T safe_mod(T a, T b) { assert(b != 0); if (b < 0) a = -a, b = -b; a %= b; return a >= 0 ? a : a + b; } template T floor(T a, T b) { assert(b != 0); if (b < 0) a = -a, b = -b; return a >= 0 ? a / b : (a + 1) / b - 1; } template T ceil(T a, T b) { assert(b != 0); if (b < 0) a = -a, b = -b; return a > 0 ? (a - 1) / b + 1 : a / b; } long long isqrt(long long n) { if (n <= 0) return 0; long long x = sqrt(n); while ((__int128)(x + 1) * (x + 1) <= n) x++; while ((__int128)x * x > n) x--; return x; } long long floor_root(long long n, int k) { assert(n >= 0); if (n == 0) return 0; assert(k >= 1); if (k == 1) return n; if (k > 64) return 1; long long x = round(pow((long double)n, 1.0L / k)); auto check = [&](long long a) { if (a <= 0) return true; __int128_t p = 1; for (int i = 0; i < k; ++i) if ((p *= a) > n) return false; return true; }; while (check(x + 1)) x++; while (!check(x)) x--; return x; } unsigned long long floor_root_unsigned(unsigned long long n, int k) { assert(k >= 1); if (n <= 1 || k == 1) return n; if (k >= 64) return 1; int bits = (64 + k - 1) / k; unsigned long long ok = 1, ng = min(n, 1ULL << bits); auto check = [&](unsigned long long a) { __uint128_t p = 1; for (int i = 0; i < k; i++) { p *= a; if (p > n) return false; } return true; }; while (ok + 1 < ng) { unsigned long long mid = ok + (ng - ok) / 2; (check(mid) ? ok : ng) = mid; } return ok; } // return g=gcd(a,b) // a*x+b*y=g // - b!=0 -> 0<=x<|b|/g // - b=0 -> ax=g template T ext_gcd(T a, T b, T& x, T& y) { T a0 = a, b0 = b; bool sgn_a = a < 0, sgn_b = b < 0; if (sgn_a) a = -a; if (sgn_b) b = -b; if (b == 0) { x = sgn_a ? -1 : 1; y = 0; return a; } T x00 = 1, x01 = 0, x10 = 0, x11 = 1; while (b != 0) { T q = a / b, r = a - b * q; x00 -= q * x01; x10 -= q * x11; swap(x00, x01); swap(x10, x11); a = b, b = r; } x = x00, y = x10; if (sgn_a) x = -x; if (sgn_b) y = -y; if (b0 != 0) { a0 /= a, b0 /= a; if (b0 < 0) a0 = -a0, b0 = -b0; T q = x >= 0 ? x / b0 : (x + 1) / b0 - 1; x -= b0 * q; y += a0 * q; } return a; } constexpr long long inv_mod(long long x, long long m) { x %= m; if (x < 0) x += m; long long a = m, b = x; long long y0 = 0, y1 = 1; while (b > 0) { long long q = a / b; swap(a -= q * b, b); swap(y0 -= q * y1, y1); } if (y0 < 0) y0 += m / a; return y0; } long long pow_mod(long long x, long long n, long long m) { if (m == 1) return 0; x = (x % m + m) % m; long long y = 1; while (n) { if (n & 1) y = y * x % m; x = x * x % m; n >>= 1; } return y; } constexpr long long pow_mod_constexpr(long long x, long long n, int m) { if (m == 1) return 0; unsigned int _m = (unsigned int)(m); unsigned long long r = 1; unsigned long long y = x % m; if (y >= m) y += m; while (n) { if (n & 1) r = (r * y) % _m; y = (y * y) % _m; n >>= 1; } return r; } constexpr bool is_prime_constexpr(int n) { if (n <= 1) return false; if (n == 2 || n == 7 || n == 61) return true; if (n % 2 == 0) return false; long long d = n - 1; while (d % 2 == 0) d /= 2; constexpr long long bases[3] = {2, 7, 61}; for (long long a : bases) { long long t = d; long long y = pow_mod_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 is_prime = is_prime_constexpr(n); }; // namespace Math #line 3 "library/modint/modint.hpp" template struct ModInt { using mint = ModInt; static constexpr unsigned int get_mod() { return m; } static mint raw(int v) { mint x; x._v = v; return x; } ModInt() : _v(0) {} ModInt(int64_t v) { long long x = (long long)(v % (long long)(umod())); if (x < 0) x += umod(); _v = (unsigned int)(x); } unsigned int val() const { return _v; } mint& operator++() { _v++; if (_v == umod()) _v = 0; return *this; } mint& operator--() { if (_v == 0) _v = umod(); _v--; return *this; } mint operator++(int) { mint result = *this; ++*this; return result; } mint operator--(int) { mint result = *this; --*this; return result; } mint& operator+=(const mint& rhs) { _v += rhs._v; if (_v >= umod()) _v -= umod(); return *this; } mint& operator-=(const mint& rhs) { _v -= rhs._v; if (_v >= umod()) _v += umod(); return *this; } mint& operator*=(const mint& rhs) { unsigned long long z = _v; z *= rhs._v; _v = (unsigned int)(z % umod()); return *this; } mint& operator/=(const mint& rhs) { return *this *= rhs.inv(); } mint operator+() const { return *this; } mint operator-() const { return mint() - *this; } mint pow(long long n) const { assert(0 <= n); mint x = *this, r = 1; while (n) { if (n & 1) r *= x; x *= x; n >>= 1; } return r; } mint inv() const { if (is_prime) { assert(_v); return pow(umod() - 2); } else { auto inv = Math::inv_mod(_v, umod()); return raw(inv); } } 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 istream& operator>>(istream& is, mint& x) { int64_t v; is >> v; x = mint(v); return is; } friend ostream& operator<<(ostream& os, const mint& x) { return os << x.val(); } private: unsigned int _v; static constexpr unsigned int umod() { return m; } static constexpr bool is_prime = Math::is_prime; }; using ModInt998244353 = ModInt<998244353>; using ModInt1000000007 = ModInt<1000000007>; #line 4 "main.cpp" using mint = ModInt998244353; #line 2 "library/fps/fps-ntt-friendly.hpp" #line 2 "library/fft/ntt.hpp" template struct NTT { static constexpr unsigned int mod = mint::get_mod(); static constexpr unsigned long long pow_constexpr(unsigned long long x, unsigned long long n, unsigned long long m) { unsigned long long y = 1; while (n) { if (n & 1) y = y * x % m; x = x * x % m; n >>= 1; } return y; } static constexpr unsigned int get_g() { unsigned long long x = 2; while (pow_constexpr(x, (mod - 1) >> 1, mod) == 1) x += 1; return x; } static constexpr unsigned int g = get_g(); static constexpr int rank2 = __builtin_ctzll(mod - 1); array root; array iroot; array rate2; array irate2; array rate3; array irate3; NTT() { root[rank2] = mint(g).pow((mod - 1) >> rank2); iroot[rank2] = root[rank2].inv(); for (int i = rank2 - 1; i >= 0; i--) { root[i] = root[i + 1] * root[i + 1]; iroot[i] = iroot[i + 1] * iroot[i + 1]; } { mint prod = 1, iprod = 1; for (int i = 0; i <= rank2 - 2; i++) { rate2[i] = root[i + 2] * prod; irate2[i] = iroot[i + 2] * iprod; prod *= iroot[i + 2]; iprod *= root[i + 2]; } } { mint prod = 1, iprod = 1; for (int i = 0; i <= rank2 - 3; i++) { rate3[i] = root[i + 3] * prod; irate3[i] = iroot[i + 3] * iprod; prod *= iroot[i + 3]; iprod *= root[i + 3]; } } } void ntt(vector& a) { int n = int(a.size()); int h = __builtin_ctzll((unsigned int)n); assert(h <= rank2); a.resize(1 << h); int len = 0; // a[i, i+(n>>len), i+2*(n>>len), ..] is transformed while (len < h) { if (h - len == 1) { int p = 1 << (h - len - 1); mint rot = 1; for (int s = 0; s < (1 << len); s++) { int offset = s << (h - len); for (int i = 0; i < p; i++) { auto l = a[i + offset]; auto r = a[i + offset + p] * rot; a[i + offset] = l + r; a[i + offset + p] = l - r; } if (s + 1 != (1 << len)) rot *= rate2[__builtin_ctzll(~(unsigned int)(s))]; } len++; } else { // 4-base int p = 1 << (h - len - 2); mint rot = 1, imag = root[2]; for (int s = 0; s < (1 << len); s++) { mint rot2 = rot * rot; mint rot3 = rot2 * rot; int offset = s << (h - len); for (int i = 0; i < p; i++) { auto mod2 = 1ULL * mint::get_mod() * mint::get_mod(); auto a0 = 1ULL * a[i + offset].val(); auto a1 = 1ULL * a[i + offset + p].val() * rot.val(); auto a2 = 1ULL * a[i + offset + 2 * p].val() * rot2.val(); auto a3 = 1ULL * a[i + offset + 3 * p].val() * rot3.val(); auto a1na3imag = 1ULL * mint(a1 + mod2 - a3).val() * imag.val(); auto na2 = mod2 - a2; a[i + offset] = a0 + a2 + a1 + a3; a[i + offset + 1 * p] = a0 + a2 + (2 * mod2 - (a1 + a3)); a[i + offset + 2 * p] = a0 + na2 + a1na3imag; a[i + offset + 3 * p] = a0 + na2 + (mod2 - a1na3imag); } if (s + 1 != (1 << len)) rot *= rate3[__builtin_ctzll(~(unsigned int)(s))]; } len += 2; } } } void intt(vector& a) { int n = int(a.size()); int h = __builtin_ctzll((unsigned int)n); assert(h <= rank2); a.resize(1 << h); int len = h; // a[i, i+(n>>len), i+2*(n>>len), ..] is transformed while (len) { if (len == 1) { int p = 1 << (h - len); mint irot = 1; for (int s = 0; s < (1 << (len - 1)); s++) { int offset = s << (h - len + 1); for (int i = 0; i < p; i++) { auto l = a[i + offset]; auto r = a[i + offset + p]; a[i + offset] = l + r; a[i + offset + p] = (unsigned long long)(mint::get_mod() + l.val() - r.val()) * irot.val(); } if (s + 1 != (1 << (len - 1))) irot *= irate2[__builtin_ctzll(~(unsigned int)(s))]; } len--; } else { // 4-base int p = 1 << (h - len); mint irot = 1, iimag = iroot[2]; for (int s = 0; s < (1 << (len - 2)); s++) { mint irot2 = irot * irot; mint irot3 = irot2 * irot; int offset = s << (h - len + 2); for (int i = 0; i < p; i++) { auto a0 = 1ULL * a[i + offset + 0 * p].val(); auto a1 = 1ULL * a[i + offset + 1 * p].val(); auto a2 = 1ULL * a[i + offset + 2 * p].val(); auto a3 = 1ULL * a[i + offset + 3 * p].val(); auto a2na3iimag = 1ULL * mint((mint::get_mod() + a2 - a3) * iimag.val()).val(); a[i + offset] = a0 + a1 + a2 + a3; a[i + offset + 1 * p] = (a0 + (mint::get_mod() - a1) + a2na3iimag) * irot.val(); a[i + offset + 2 * p] = (a0 + a1 + (mint::get_mod() - a2) + (mint::get_mod() - a3)) * irot2.val(); a[i + offset + 3 * p] = (a0 + (mint::get_mod() - a1) + (mint::get_mod() - a2na3iimag)) * irot3.val(); } if (s + 1 != (1 << (len - 2))) irot *= irate3[__builtin_ctzll(~(unsigned int)(s))]; } len -= 2; } } mint e = mint(n).inv(); for (auto& x : a) x *= e; } vector multiply(const vector& a, const vector& b) { if (a.empty() || b.empty()) return vector(); int n = a.size(), m = b.size(); int sz = n + m - 1; if (n <= 30 || m <= 30) { if (n > 30) return multiply(b, a); vector res(sz); for (int i = 0; i < n; i++) for (int j = 0; j < m; j++) res[i + j] += a[i] * b[j]; return res; } int sz1 = 1; while (sz1 < sz) sz1 <<= 1; vector res(sz1); for (int i = 0; i < n; i++) res[i] = a[i]; ntt(res); if (a == b) for (int i = 0; i < sz1; i++) res[i] *= res[i]; else { vector c(sz1); for (int i = 0; i < m; i++) c[i] = b[i]; ntt(c); for (int i = 0; i < sz1; i++) res[i] *= c[i]; } intt(res); res.resize(sz); return res; } // c[i]=sum[j]a[j]b[i+j] vector middle_product(const vector& a, const vector& b) { if (b.empty() || a.size() > b.size()) return {}; int n = a.size(), m = b.size(); int sz = m - n + 1; if (n <= 30 || sz <= 30) { vector res(sz); for (int i = 0; i < sz; i++) for (int j = 0; j < n; j++) res[i] += a[j] * b[i + j]; return res; } int sz1 = 1; while (sz1 < m) sz1 <<= 1; vector res(sz1), b2(sz1); reverse_copy(a.begin(), a.end(), res.begin()); copy(b.begin(), b.end(), b2.begin()); ntt(res); ntt(b2); for (int i = 0; i < res.size(); i++) res[i] *= b2[i]; intt(res); res.resize(m); res.erase(res.begin(), res.begin() + n - 1); return res; } void ntt_doubling(vector& a) { int n = (int)a.size(); auto b = a; intt(b); mint r = 1, zeta = mint(g).pow((mint::get_mod() - 1) / (n << 1)); for (int i = 0; i < n; i++) b[i] *= r, r *= zeta; ntt(b); copy(b.begin(), b.end(), back_inserter(a)); } }; /** * @brief NTT (数論変換) * @docs docs/fft/ntt.md */ #line 2 "library/fps/formal-power-series.hpp" template struct FormalPowerSeries : vector { using vector::vector; using FPS = FormalPowerSeries; FormalPowerSeries(const vector& r) : vector(r) {} FormalPowerSeries(vector&& r) : vector(std::move(r)) {} FPS& operator=(const vector& r) { vector::operator=(r); return *this; } FPS& operator+=(const FPS& r) { if (r.size() > this->size()) this->resize(r.size()); for (int i = 0; i < (int)r.size(); i++) (*this)[i] += r[i]; return *this; } FPS& operator+=(const mint& r) { if (this->empty()) this->resize(1); (*this)[0] += r; return *this; } FPS& operator-=(const FPS& r) { if (r.size() > this->size()) this->resize(r.size()); for (int i = 0; i < (int)r.size(); i++) (*this)[i] -= r[i]; return *this; } FPS& operator-=(const mint& r) { if (this->empty()) this->resize(1); (*this)[0] -= r; return *this; } FPS& operator*=(const mint& v) { for (int k = 0; k < (int)this->size(); k++) (*this)[k] *= v; return *this; } FPS& operator/=(const FPS& r) { if (this->size() < r.size()) { this->clear(); return *this; } int n = this->size() - r.size() + 1; if ((int)r.size() <= 64) { FPS f(*this), g(r); g.shrink(); mint coeff = g.at(g.size() - 1).inv(); for (auto& x : g) x *= coeff; int deg = (int)f.size() - (int)g.size() + 1; int gs = g.size(); FPS quo(deg); for (int i = deg - 1; i >= 0; i--) { quo[i] = f[i + gs - 1]; for (int j = 0; j < gs; j++) f[i + j] -= quo[i] * g[j]; } *this = quo * coeff; this->resize(n, mint(0)); return *this; } return *this = ((*this).rev().pre(n) * r.rev().inv(n)).pre(n).rev(); } FPS& operator%=(const FPS& r) { *this -= *this / r * r; shrink(); return *this; } FPS operator+(const FPS& r) const { return FPS(*this) += r; } FPS operator+(const mint& v) const { return FPS(*this) += v; } FPS operator-(const FPS& r) const { return FPS(*this) -= r; } FPS operator-(const mint& v) const { return FPS(*this) -= v; } FPS operator*(const FPS& r) const { return FPS(*this) *= r; } FPS operator*(const mint& v) const { return FPS(*this) *= v; } FPS operator/(const FPS& r) const { return FPS(*this) /= r; } FPS operator%(const FPS& r) const { return FPS(*this) %= r; } FPS operator-() const { FPS ret(this->size()); for (int i = 0; i < (int)this->size(); i++) ret[i] = -(*this)[i]; return ret; } void shrink() { while (this->size() && this->back() == mint(0)) this->pop_back(); } FPS rev() const { FPS ret(*this); reverse(begin(ret), end(ret)); return ret; } FPS dot(FPS r) const { FPS ret(min(this->size(), r.size())); for (int i = 0; i < (int)ret.size(); i++) ret[i] = (*this)[i] * r[i]; return ret; } FPS pre(int sz) const { return FPS(begin(*this), begin(*this) + min((int)this->size(), sz)); } FPS operator>>=(int sz) { assert(sz >= 0); if ((int)this->size() <= sz) this->clear(); else this->erase(this->begin(), this->begin() + sz); return *this; } FPS operator>>(int sz) const { if ((int)this->size() <= sz) return {}; FPS ret(*this); ret.erase(ret.begin(), ret.begin() + sz); return ret; } FPS operator<<=(int sz) { assert(sz >= 0); this->insert(this->begin(), sz, mint(0)); return *this; } FPS operator<<(int sz) const { FPS ret(*this); ret.insert(ret.begin(), sz, mint(0)); return ret; } FPS diff() const { const int n = (int)this->size(); FPS ret(max(0, n - 1)); mint one(1), coeff(1); for (int i = 1; i < n; i++) { ret[i - 1] = (*this)[i] * coeff; coeff += one; } return ret; } FPS integral() const { const int n = (int)this->size(); FPS ret(n + 1); ret[0] = mint(0); if (n > 0) ret[1] = mint(1); auto mod = mint::get_mod(); for (int i = 2; i <= n; i++) ret[i] = (-ret[mod % i]) * (mod / i); for (int i = 0; i < n; i++) ret[i + 1] *= (*this)[i]; return ret; } mint eval(mint x) const { mint r = 0, w = 1; for (auto& v : *this) r += w * v, w *= x; return r; } FPS log(int deg = -1) const { assert((*this)[0] == mint(1)); if (deg == -1) deg = (int)this->size(); return (this->diff() * this->inv(deg)).pre(deg - 1).integral(); } FPS pow(int64_t k, int deg = -1) const { const int n = (int)this->size(); if (deg == -1) deg = n; if (k == 0) { FPS ret(deg); if (deg) ret[0] = 1; return ret; } for (int i = 0; i < n; i++) { if ((*this)[i] != mint(0)) { mint rev = mint(1) / (*this)[i]; FPS ret = (((*this * rev) >> i).log(deg) * k).exp(deg); ret *= (*this)[i].pow(k); ret = (ret << (i * k)).pre(deg); if ((int)ret.size() < deg) ret.resize(deg, mint(0)); return ret; } if (__int128_t(i + 1) * k >= deg) return FPS(deg, mint(0)); } return FPS(deg, mint(0)); } static void* ntt_ptr; static void set_ntt(); FPS& operator*=(const FPS& r); FPS middle_product(const FPS& r) const; void ntt(); void intt(); void ntt_doubling(); static int ntt_root(); FPS inv(int deg = -1) const; FPS exp(int deg = -1) const; }; template void* FormalPowerSeries::ntt_ptr = nullptr; #line 5 "library/fps/fps-ntt-friendly.hpp" template void FormalPowerSeries::set_ntt() { if (!ntt_ptr) ntt_ptr = new NTT; } template FormalPowerSeries& FormalPowerSeries::operator*=(const FormalPowerSeries& r) { if (this->empty() || r.empty()) { this->clear(); return *this; } set_ntt(); auto ret = static_cast*>(ntt_ptr)->multiply(*this, r); return *this = FormalPowerSeries(ret.begin(), ret.end()); } template FormalPowerSeries FormalPowerSeries::middle_product(const FormalPowerSeries& r) const { set_ntt(); auto ret = static_cast*>(ntt_ptr)->middle_product(*this, r); return FormalPowerSeries(ret.begin(), ret.end()); } template void FormalPowerSeries::ntt() { set_ntt(); static_cast*>(ntt_ptr)->ntt(*this); } template void FormalPowerSeries::intt() { set_ntt(); static_cast*>(ntt_ptr)->intt(*this); } template void FormalPowerSeries::ntt_doubling() { set_ntt(); static_cast*>(ntt_ptr)->ntt_doubling(*this); } template int FormalPowerSeries::ntt_root() { set_ntt(); return static_cast*>(ntt_ptr)->g; } template FormalPowerSeries FormalPowerSeries::inv(int deg) const { assert((*this)[0] != mint(0)); if (deg == -1) deg = (*this).size(); FPS ret{mint(1) / (*this)[0]}; for (int i = 1; i < deg; i <<= 1) ret = (ret + ret - ret * ret * (*this).pre(i << 1)).pre(i << 1); return ret.pre(deg); } template FormalPowerSeries FormalPowerSeries::exp(int deg) const { assert((*this)[0] == mint(0)); if (deg == -1) deg = (*this).size(); FPS ret{mint(1)}; for (int i = 1; i < deg; i <<= 1) ret = (ret * ((*this).pre(i << 1) - ret.log(i << 1) + 1)).pre(i << 1); return ret.pre(deg); } #line 6 "main.cpp" using fps = FormalPowerSeries; #line 2 "library/modint/factorial.hpp" template struct Factorial { static void reserve(int n) { inv(n); fact(n); fact_inv(n); } static mint inv(int n) { static long long mod = mint::get_mod(); static vector buf({0, 1}); assert(n != 0); if (mod != mint::get_mod()) { mod = mint::get_mod(); buf = vector({0, 1}); } while ((int)buf.capacity() <= n) buf.reserve(buf.capacity() * 2); while ((int)buf.size() <= n) { long long k = buf.size(), q = (mod + k - 1) / k; buf.push_back(q * buf[k * q - mod]); } return buf[n]; } static mint fact(int n) { static long long mod = mint::get_mod(); static vector buf({1, 1}); assert(n >= 0); if (mod != mint::get_mod()) { mod = mint::get_mod(); buf = vector({1, 1}); } while ((int)buf.capacity() <= n) buf.reserve(buf.capacity() * 2); while ((int)buf.size() <= n) { long long k = buf.size(); buf.push_back(buf.back() * k); } return buf[n]; } static mint fact_inv(int n) { static long long mod = mint::get_mod(); static vector buf({1, 1}); assert(n >= 0); if (mod != mint::get_mod()) { mod = mint::get_mod(); buf = vector({1, 1}); } if ((int)buf.size() <= n) inv(n); while ((int)buf.capacity() <= n) buf.reserve(buf.capacity() * 2); while ((int)buf.size() <= n) { long long k = buf.size(); buf.push_back(buf.back() * inv(k)); } return buf[n]; } static mint binom(int n, int r) { if (r < 0 || r > n) return 0; return fact(n) * fact_inv(r) * fact_inv(n - r); } static mint binom_naive(int n, int r) { if (r < 0 || r > n) return 0; mint res = fact_inv(r); for (int i = 0; i < r; i++) res *= n - i; return res; } static mint multinom(const vector& r) { int n = 0; for (auto& x : r) { if (x < 0) return 0; n += x; } mint res = fact(n); for (auto& x : r) res *= fact_inv(x); return res; } static mint P(int n, int r) { if (r < 0 || r > n) return 0; return fact(n) * fact_inv(n - r); } // partition n items to r groups (allow empty group) static mint H(int n, int r) { if (n < 0 || r < 0) return 0; return r == 0 ? 1 : binom(n + r - 1, r); } }; /** * @brief 階乗, 二項係数 */ #line 4 "library/fps/taylor-shift.hpp" // f(x+a) template FormalPowerSeries TaylorShift(FormalPowerSeries f, mint a) { using fps = FormalPowerSeries; if (f.size() <= 1) return f; int n = f.size(); using fact = Factorial; fact::reserve(n); for (int i = 0; i < n; i++) f[i] *= fact::fact(i); reverse(f.begin(), f.end()); fps g(n, mint(1)); for (int i = 1; i < n; i++) g[i] = g[i - 1] * a * fact::inv(i); f = (f * g).pre(n); reverse(f.begin(), f.end()); for (int i = 0; i < n; i++) f[i] *= fact::fact_inv(i); return f; } /** * @brief Taylor Shift * @docs docs/fps/taylor-shift.md */ #line 4 "library/fps/composition.hpp" template FormalPowerSeries Composition(FormalPowerSeries f, FormalPowerSeries g, int n = -1) { if (n == -1) n = f.size(); assert(n >= 0); if (n == 0 || f.size() == 0) return FormalPowerSeries(n); f = TaylorShift(f, g[0]); int k = 1; while (k < f.size()) k <<= 1; int l = 1, m = 2 * k * l; int log = __builtin_ctz((unsigned int)m); mint wi = mint(FormalPowerSeries::ntt_root()).inv().pow((mint::get_mod() - 1) >> (log + 1)); vector rev(m); for (int i = 0; i < rev.size(); i++) rev[i] = (rev[i / 2] / 2) | ((i & 1) << (log - 1)); vector pow(m, 1); for (int i = 1; i < m; i++) pow[rev[i]] = pow[rev[i - 1]] * wi; mint inv2 = mint(2).inv(); FormalPowerSeries P(m), Q(m); for (int i = 1; i < g.size() && i < k; i++) Q[i] = -g[i]; auto rec = [&](auto rec, int k, int l, FormalPowerSeries Q) { if (k <= 1) { for (int i = 0; i < l && i < f.size(); i++) P[(l - 1 - i) * 2] = f[i]; return; } Q.resize(2 * m), Q.ntt(); FormalPowerSeries Q1(m); for (int i = 0; i < m; i++) { mint b = (i >> (log - 1)) == 0 ? 1 : -1; int j = i << 1; Q1[i] = Q[j] * Q[j ^ 1] + (Q[j] + Q[j ^ 1]) * b; } Q1.intt(); k >>= 1, l <<= 1; for (int i = k; i < k * 2; i++) for (int j = 0; j < l; j++) Q1[i + j * k * 2] = 0; rec(rec, k, l, Q1); for (int i = k; i < k * 2; i++) for (int j = 0; j < l; j++) P[i + j * k * 2] = 0; P.resize(m); reverse(P.begin() + 1, P.end()); P.ntt(), P.resize(2 * m); for (int i = m - 1; i >= 0; i--) { mint b = (i >> (log - 1)) == 0 ? 1 : -1; P[2 * i + 1] = -(Q[2 * i] + b) * pow[i] * P[i]; P[2 * i] = (Q[2 * i + 1] + b) * pow[i] * P[i]; } P.intt(); reverse(P.begin() + 1, P.end()); }; rec(rec, k, l, Q); P.resize(k); P = P.rev(); P.resize(n); return P; } /** * @brief FPS 合成 * @docs docs/fps/composition.md */ #line 8 "main.cpp" void solve() { int n, m; in(n, m); fps f(n), g(n), h(n); in(f, g, h); auto ans = Composition(h, g * f[1].pow(m)); out(ans); } int main() { int T = 1; // in(T); while (T--) solve(); }