結果
問題 | No.2688 Cell Proliferation (Hard) |
ユーザー | hitonanode |
提出日時 | 2024-04-04 22:29:18 |
言語 | C++23 (gcc 12.3.0 + boost 1.83.0) |
結果 |
AC
|
実行時間 | 673 ms / 4,000 ms |
コード長 | 24,582 bytes |
コンパイル時間 | 3,485 ms |
コンパイル使用メモリ | 231,208 KB |
実行使用メモリ | 67,956 KB |
最終ジャッジ日時 | 2024-10-01 00:38:16 |
合計ジャッジ時間 | 16,006 ms |
ジャッジサーバーID (参考情報) |
judge5 / judge1 |
(要ログイン)
テストケース
テストケース表示入力 | 結果 | 実行時間 実行使用メモリ |
---|---|---|
testcase_00 | AC | 2 ms
5,248 KB |
testcase_01 | AC | 2 ms
5,248 KB |
testcase_02 | AC | 1 ms
5,248 KB |
testcase_03 | AC | 78 ms
15,708 KB |
testcase_04 | AC | 654 ms
66,916 KB |
testcase_05 | AC | 353 ms
52,284 KB |
testcase_06 | AC | 168 ms
27,972 KB |
testcase_07 | AC | 168 ms
27,796 KB |
testcase_08 | AC | 645 ms
66,644 KB |
testcase_09 | AC | 638 ms
66,076 KB |
testcase_10 | AC | 654 ms
67,368 KB |
testcase_11 | AC | 644 ms
66,868 KB |
testcase_12 | AC | 173 ms
28,168 KB |
testcase_13 | AC | 673 ms
67,940 KB |
testcase_14 | AC | 624 ms
65,252 KB |
testcase_15 | AC | 660 ms
67,956 KB |
testcase_16 | AC | 623 ms
65,216 KB |
testcase_17 | AC | 353 ms
52,736 KB |
testcase_18 | AC | 655 ms
67,516 KB |
testcase_19 | AC | 358 ms
52,652 KB |
testcase_20 | AC | 335 ms
51,140 KB |
testcase_21 | AC | 340 ms
51,004 KB |
testcase_22 | AC | 354 ms
52,472 KB |
testcase_23 | AC | 337 ms
51,448 KB |
testcase_24 | AC | 653 ms
67,348 KB |
testcase_25 | AC | 172 ms
28,344 KB |
testcase_26 | AC | 663 ms
67,552 KB |
testcase_27 | AC | 329 ms
52,832 KB |
testcase_28 | AC | 566 ms
65,364 KB |
ソースコード
#include <algorithm> #include <array> #include <bitset> #include <cassert> #include <chrono> #include <cmath> #include <complex> #include <deque> #include <forward_list> #include <fstream> #include <functional> #include <iomanip> #include <ios> #include <iostream> #include <limits> #include <list> #include <map> #include <memory> #include <numeric> #include <optional> #include <queue> #include <random> #include <set> #include <sstream> #include <stack> #include <string> #include <tuple> #include <type_traits> #include <unordered_map> #include <unordered_set> #include <utility> #include <vector> using namespace std; using lint = long long; using pint = pair<int, int>; using plint = pair<lint, lint>; struct fast_ios { fast_ios(){ cin.tie(nullptr), ios::sync_with_stdio(false), cout << fixed << setprecision(20); }; } fast_ios_; #define ALL(x) (x).begin(), (x).end() #define FOR(i, begin, end) for(int i=(begin),i##_end_=(end);i<i##_end_;i++) #define IFOR(i, begin, end) for(int i=(end)-1,i##_begin_=(begin);i>=i##_begin_;i--) #define REP(i, n) FOR(i,0,n) #define IREP(i, n) IFOR(i,0,n) template <typename T> bool chmax(T &m, const T q) { return m < q ? (m = q, true) : false; } template <typename T> bool chmin(T &m, const T q) { return m > q ? (m = q, true) : false; } const std::vector<std::pair<int, int>> grid_dxs{{1, 0}, {-1, 0}, {0, 1}, {0, -1}}; int floor_lg(long long x) { return x <= 0 ? -1 : 63 - __builtin_clzll(x); } template <class T1, class T2> T1 floor_div(T1 num, T2 den) { return (num > 0 ? num / den : -((-num + den - 1) / den)); } template <class T1, class T2> std::pair<T1, T2> operator+(const std::pair<T1, T2> &l, const std::pair<T1, T2> &r) { return std::make_pair(l.first + r.first, l.second + r.second); } template <class T1, class T2> std::pair<T1, T2> operator-(const std::pair<T1, T2> &l, const std::pair<T1, T2> &r) { return std::make_pair(l.first - r.first, l.second - r.second); } template <class T> std::vector<T> sort_unique(std::vector<T> vec) { sort(vec.begin(), vec.end()), vec.erase(unique(vec.begin(), vec.end()), vec.end()); return vec; } template <class T> int arglb(const std::vector<T> &v, const T &x) { return std::distance(v.begin(), std::lower_bound(v.begin(), v.end(), x)); } template <class T> int argub(const std::vector<T> &v, const T &x) { return std::distance(v.begin(), std::upper_bound(v.begin(), v.end(), x)); } template <class IStream, class T> IStream &operator>>(IStream &is, std::vector<T> &vec) { for (auto &v : vec) is >> v; return is; } template <class OStream, class T> OStream &operator<<(OStream &os, const std::vector<T> &vec); template <class OStream, class T, size_t sz> OStream &operator<<(OStream &os, const std::array<T, sz> &arr); template <class OStream, class T, class TH> OStream &operator<<(OStream &os, const std::unordered_set<T, TH> &vec); template <class OStream, class T, class U> OStream &operator<<(OStream &os, const pair<T, U> &pa); template <class OStream, class T> OStream &operator<<(OStream &os, const std::deque<T> &vec); template <class OStream, class T> OStream &operator<<(OStream &os, const std::set<T> &vec); template <class OStream, class T> OStream &operator<<(OStream &os, const std::multiset<T> &vec); template <class OStream, class T> OStream &operator<<(OStream &os, const std::unordered_multiset<T> &vec); template <class OStream, class T, class U> OStream &operator<<(OStream &os, const std::pair<T, U> &pa); template <class OStream, class TK, class TV> OStream &operator<<(OStream &os, const std::map<TK, TV> &mp); template <class OStream, class TK, class TV, class TH> OStream &operator<<(OStream &os, const std::unordered_map<TK, TV, TH> &mp); template <class OStream, class... T> OStream &operator<<(OStream &os, const std::tuple<T...> &tpl); template <class OStream, class T> OStream &operator<<(OStream &os, const std::vector<T> &vec) { os << '['; for (auto v : vec) os << v << ','; os << ']'; return os; } template <class OStream, class T, size_t sz> OStream &operator<<(OStream &os, const std::array<T, sz> &arr) { os << '['; for (auto v : arr) os << v << ','; os << ']'; return os; } template <class... T> std::istream &operator>>(std::istream &is, std::tuple<T...> &tpl) { std::apply([&is](auto &&... args) { ((is >> args), ...);}, tpl); return is; } template <class OStream, class... T> OStream &operator<<(OStream &os, const std::tuple<T...> &tpl) { os << '('; std::apply([&os](auto &&... args) { ((os << args << ','), ...);}, tpl); return os << ')'; } template <class OStream, class T, class TH> OStream &operator<<(OStream &os, const std::unordered_set<T, TH> &vec) { os << '{'; for (auto v : vec) os << v << ','; os << '}'; return os; } template <class OStream, class T> OStream &operator<<(OStream &os, const std::deque<T> &vec) { os << "deq["; for (auto v : vec) os << v << ','; os << ']'; return os; } template <class OStream, class T> OStream &operator<<(OStream &os, const std::set<T> &vec) { os << '{'; for (auto v : vec) os << v << ','; os << '}'; return os; } template <class OStream, class T> OStream &operator<<(OStream &os, const std::multiset<T> &vec) { os << '{'; for (auto v : vec) os << v << ','; os << '}'; return os; } template <class OStream, class T> OStream &operator<<(OStream &os, const std::unordered_multiset<T> &vec) { os << '{'; for (auto v : vec) os << v << ','; os << '}'; return os; } template <class OStream, class T, class U> OStream &operator<<(OStream &os, const std::pair<T, U> &pa) { return os << '(' << pa.first << ',' << pa.second << ')'; } template <class OStream, class TK, class TV> OStream &operator<<(OStream &os, const std::map<TK, TV> &mp) { os << '{'; for (auto v : mp) os << v.first << "=>" << v.second << ','; os << '}'; return os; } template <class OStream, class TK, class TV, class TH> OStream &operator<<(OStream &os, const std::unordered_map<TK, TV, TH> &mp) { os << '{'; for (auto v : mp) os << v.first << "=>" << v.second << ','; os << '}'; return os; } #ifdef HITONANODE_LOCAL const string COLOR_RESET = "\033[0m", BRIGHT_GREEN = "\033[1;32m", BRIGHT_RED = "\033[1;31m", BRIGHT_CYAN = "\033[1;36m", NORMAL_CROSSED = "\033[0;9;37m", RED_BACKGROUND = "\033[1;41m", NORMAL_FAINT = "\033[0;2m"; #define dbg(x) std::cerr << BRIGHT_CYAN << #x << COLOR_RESET << " = " << (x) << NORMAL_FAINT << " (L" << __LINE__ << ") " << __FILE__ << COLOR_RESET << std::endl #define dbgif(cond, x) ((cond) ? std::cerr << BRIGHT_CYAN << #x << COLOR_RESET << " = " << (x) << NORMAL_FAINT << " (L" << __LINE__ << ") " << __FILE__ << COLOR_RESET << std::endl : std::cerr) #else #define dbg(x) ((void)0) #define dbgif(cond, x) ((void)0) #endif template <int md> struct ModInt { using lint = long long; constexpr static int mod() { return md; } static int get_primitive_root() { static int primitive_root = 0; if (!primitive_root) { primitive_root = [&]() { std::set<int> fac; int v = md - 1; for (lint i = 2; i * i <= v; i++) while (v % i == 0) fac.insert(i), v /= i; if (v > 1) fac.insert(v); for (int g = 1; g < md; g++) { bool ok = true; for (auto i : fac) if (ModInt(g).pow((md - 1) / i) == 1) { ok = false; break; } if (ok) return g; } return -1; }(); } return primitive_root; } int val_; int val() const noexcept { return val_; } constexpr ModInt() : val_(0) {} constexpr ModInt &_setval(lint v) { return val_ = (v >= md ? v - md : v), *this; } constexpr ModInt(lint v) { _setval(v % md + md); } constexpr explicit operator bool() const { return val_ != 0; } constexpr ModInt operator+(const ModInt &x) const { return ModInt()._setval((lint)val_ + x.val_); } constexpr ModInt operator-(const ModInt &x) const { return ModInt()._setval((lint)val_ - x.val_ + md); } constexpr ModInt operator*(const ModInt &x) const { return ModInt()._setval((lint)val_ * x.val_ % md); } constexpr ModInt operator/(const ModInt &x) const { return ModInt()._setval((lint)val_ * x.inv().val() % md); } constexpr ModInt operator-() const { return ModInt()._setval(md - val_); } constexpr ModInt &operator+=(const ModInt &x) { return *this = *this + x; } constexpr ModInt &operator-=(const ModInt &x) { return *this = *this - x; } constexpr ModInt &operator*=(const ModInt &x) { return *this = *this * x; } constexpr ModInt &operator/=(const ModInt &x) { return *this = *this / x; } friend constexpr ModInt operator+(lint a, const ModInt &x) { return ModInt(a) + x; } friend constexpr ModInt operator-(lint a, const ModInt &x) { return ModInt(a) - x; } friend constexpr ModInt operator*(lint a, const ModInt &x) { return ModInt(a) * x; } friend constexpr ModInt operator/(lint a, const ModInt &x) { return ModInt(a) / x; } constexpr bool operator==(const ModInt &x) const { return val_ == x.val_; } constexpr bool operator!=(const ModInt &x) const { return val_ != x.val_; } constexpr bool operator<(const ModInt &x) const { return val_ < x.val_; } // To use std::map<ModInt, T> friend std::istream &operator>>(std::istream &is, ModInt &x) { lint t; return is >> t, x = ModInt(t), is; } constexpr friend std::ostream &operator<<(std::ostream &os, const ModInt &x) { return os << x.val_; } constexpr ModInt pow(lint n) const { ModInt ans = 1, tmp = *this; while (n) { if (n & 1) ans *= tmp; tmp *= tmp, n >>= 1; } return ans; } static constexpr int cache_limit = std::min(md, 1 << 21); static std::vector<ModInt> facs, facinvs, invs; constexpr static void _precalculation(int N) { const int l0 = facs.size(); if (N > md) N = md; if (N <= l0) return; facs.resize(N), facinvs.resize(N), invs.resize(N); for (int i = l0; i < N; i++) facs[i] = facs[i - 1] * i; facinvs[N - 1] = facs.back().pow(md - 2); for (int i = N - 2; i >= l0; i--) facinvs[i] = facinvs[i + 1] * (i + 1); for (int i = N - 1; i >= l0; i--) invs[i] = facinvs[i] * facs[i - 1]; } constexpr ModInt inv() const { if (this->val_ < cache_limit) { if (facs.empty()) facs = {1}, facinvs = {1}, invs = {0}; while (this->val_ >= int(facs.size())) _precalculation(facs.size() * 2); return invs[this->val_]; } else { return this->pow(md - 2); } } constexpr ModInt fac() const { while (this->val_ >= int(facs.size())) _precalculation(facs.size() * 2); return facs[this->val_]; } constexpr ModInt facinv() const { while (this->val_ >= int(facs.size())) _precalculation(facs.size() * 2); return facinvs[this->val_]; } constexpr ModInt doublefac() const { lint k = (this->val_ + 1) / 2; return (this->val_ & 1) ? ModInt(k * 2).fac() / (ModInt(2).pow(k) * ModInt(k).fac()) : ModInt(k).fac() * ModInt(2).pow(k); } constexpr ModInt nCr(int r) const { if (r < 0 or this->val_ < r) return ModInt(0); return this->fac() * (*this - r).facinv() * ModInt(r).facinv(); } constexpr ModInt nPr(int r) const { if (r < 0 or this->val_ < r) return ModInt(0); return this->fac() * (*this - r).facinv(); } static ModInt binom(int n, int r) { static long long bruteforce_times = 0; if (r < 0 or n < r) return ModInt(0); if (n <= bruteforce_times or n < (int)facs.size()) return ModInt(n).nCr(r); r = std::min(r, n - r); ModInt ret = ModInt(r).facinv(); for (int i = 0; i < r; ++i) ret *= n - i; bruteforce_times += r; return ret; } // Multinomial coefficient, (k_1 + k_2 + ... + k_m)! / (k_1! k_2! ... k_m!) // Complexity: O(sum(ks)) template <class Vec> static ModInt multinomial(const Vec &ks) { ModInt ret{1}; int sum = 0; for (int k : ks) { assert(k >= 0); ret *= ModInt(k).facinv(), sum += k; } return ret * ModInt(sum).fac(); } // Catalan number, C_n = binom(2n, n) / (n + 1) // C_0 = 1, C_1 = 1, C_2 = 2, C_3 = 5, C_4 = 14, ... // https://oeis.org/A000108 // Complexity: O(n) static ModInt catalan(int n) { if (n < 0) return ModInt(0); return ModInt(n * 2).fac() * ModInt(n + 1).facinv() * ModInt(n).facinv(); } ModInt sqrt() const { if (val_ == 0) return 0; if (md == 2) return val_; if (pow((md - 1) / 2) != 1) return 0; ModInt b = 1; while (b.pow((md - 1) / 2) == 1) b += 1; int e = 0, m = md - 1; while (m % 2 == 0) m >>= 1, e++; ModInt x = pow((m - 1) / 2), y = (*this) * x * x; x *= (*this); ModInt z = b.pow(m); while (y != 1) { int j = 0; ModInt t = y; while (t != 1) j++, t *= t; z = z.pow(1LL << (e - j - 1)); x *= z, z *= z, y *= z; e = j; } return ModInt(std::min(x.val_, md - x.val_)); } }; template <int md> std::vector<ModInt<md>> ModInt<md>::facs = {1}; template <int md> std::vector<ModInt<md>> ModInt<md>::facinvs = {1}; template <int md> std::vector<ModInt<md>> ModInt<md>::invs = {0}; using mint = ModInt<998244353>; #include <algorithm> #include <array> #include <cassert> #include <tuple> #include <vector> // Integer convolution for arbitrary mod // with NTT (and Garner's algorithm) for ModInt / ModIntRuntime class. // We skip Garner's algorithm if `skip_garner` is true or mod is in `nttprimes`. // input: a (size: n), b (size: m) // return: vector (size: n + m - 1) template <typename MODINT> std::vector<MODINT> nttconv(std::vector<MODINT> a, std::vector<MODINT> b, bool skip_garner); constexpr int nttprimes[3] = {998244353, 167772161, 469762049}; // Integer FFT (Fast Fourier Transform) for ModInt class // (Also known as Number Theoretic Transform, NTT) // is_inverse: inverse transform // ** Input size must be 2^n ** template <typename MODINT> void ntt(std::vector<MODINT> &a, bool is_inverse = false) { int n = a.size(); if (n == 1) return; static const int mod = MODINT::mod(); static const MODINT root = MODINT::get_primitive_root(); assert(__builtin_popcount(n) == 1 and (mod - 1) % n == 0); static std::vector<MODINT> w{1}, iw{1}; for (int m = w.size(); m < n / 2; m *= 2) { MODINT dw = root.pow((mod - 1) / (4 * m)), dwinv = 1 / dw; w.resize(m * 2), iw.resize(m * 2); for (int i = 0; i < m; i++) w[m + i] = w[i] * dw, iw[m + i] = iw[i] * dwinv; } if (!is_inverse) { for (int m = n; m >>= 1;) { for (int s = 0, k = 0; s < n; s += 2 * m, k++) { for (int i = s; i < s + m; i++) { MODINT x = a[i], y = a[i + m] * w[k]; a[i] = x + y, a[i + m] = x - y; } } } } else { for (int m = 1; m < n; m *= 2) { for (int s = 0, k = 0; s < n; s += 2 * m, k++) { for (int i = s; i < s + m; i++) { MODINT x = a[i], y = a[i + m]; a[i] = x + y, a[i + m] = (x - y) * iw[k]; } } } int n_inv = MODINT(n).inv().val(); for (auto &v : a) v *= n_inv; } } template <int MOD> std::vector<ModInt<MOD>> nttconv_(const std::vector<int> &a, const std::vector<int> &b) { int sz = a.size(); assert(a.size() == b.size() and __builtin_popcount(sz) == 1); std::vector<ModInt<MOD>> ap(sz), bp(sz); for (int i = 0; i < sz; i++) ap[i] = a[i], bp[i] = b[i]; ntt(ap, false); if (a == b) bp = ap; else ntt(bp, false); for (int i = 0; i < sz; i++) ap[i] *= bp[i]; ntt(ap, true); return ap; } long long garner_ntt_(int r0, int r1, int r2, int mod) { using mint2 = ModInt<nttprimes[2]>; static const long long m01 = 1LL * nttprimes[0] * nttprimes[1]; static const long long m0_inv_m1 = ModInt<nttprimes[1]>(nttprimes[0]).inv().val(); static const long long m01_inv_m2 = mint2(m01).inv().val(); int v1 = (m0_inv_m1 * (r1 + nttprimes[1] - r0)) % nttprimes[1]; auto v2 = (mint2(r2) - r0 - mint2(nttprimes[0]) * v1) * m01_inv_m2; return (r0 + 1LL * nttprimes[0] * v1 + m01 % mod * v2.val()) % mod; } template <typename MODINT> std::vector<MODINT> nttconv(std::vector<MODINT> a, std::vector<MODINT> b, bool skip_garner) { if (a.empty() or b.empty()) return {}; int sz = 1, n = a.size(), m = b.size(); while (sz < n + m) sz <<= 1; if (sz <= 16) { std::vector<MODINT> ret(n + m - 1); for (int i = 0; i < n; i++) { for (int j = 0; j < m; j++) ret[i + j] += a[i] * b[j]; } return ret; } int mod = MODINT::mod(); if (skip_garner or std::find(std::begin(nttprimes), std::end(nttprimes), mod) != std::end(nttprimes)) { a.resize(sz), b.resize(sz); if (a == b) { ntt(a, false); b = a; } else { ntt(a, false), ntt(b, false); } for (int i = 0; i < sz; i++) a[i] *= b[i]; ntt(a, true); a.resize(n + m - 1); } else { std::vector<int> ai(sz), bi(sz); for (int i = 0; i < n; i++) ai[i] = a[i].val(); for (int i = 0; i < m; i++) bi[i] = b[i].val(); auto ntt0 = nttconv_<nttprimes[0]>(ai, bi); auto ntt1 = nttconv_<nttprimes[1]>(ai, bi); auto ntt2 = nttconv_<nttprimes[2]>(ai, bi); a.resize(n + m - 1); for (int i = 0; i < n + m - 1; i++) a[i] = garner_ntt_(ntt0[i].val(), ntt1[i].val(), ntt2[i].val(), mod); } return a; } template <typename MODINT> std::vector<MODINT> nttconv(const std::vector<MODINT> &a, const std::vector<MODINT> &b) { return nttconv<MODINT>(a, b, false); } #include <algorithm> #include <cassert> #include <optional> #include <vector> namespace fps_nttmod { // Calculate the inverse of f(x) mod x^d // f(x) * g(x) = 1 mod x^d // If d = -1, d is set to f.size() // Complexity: O(d log d) template <class NTTModInt> std::vector<NTTModInt> inv(const std::vector<NTTModInt> &f, int d = -1) { assert(d >= -1); const int n = f.size(); if (d == -1) d = n; if (d == 0) return {}; assert(f.front() != NTTModInt(0)); using F = std::vector<NTTModInt>; F res{f.front().inv()}; // f(x) g_m(x) = 1 mod x^m for (int m = 1; m < d; m *= 2) { // g_2m = (2g_m - f g_m^2) mod x^2m F g_m{res.cbegin(), res.cbegin() + m}; g_m.resize(2 * m); ntt(g_m, false); F f_{f.cbegin(), f.cbegin() + std::min(n, 2 * m)}; f_.resize(2 * m); ntt(f_, false); for (int i = 0; i < 2 * m; ++i) f_.at(i) *= g_m.at(i); ntt(f_, true); std::rotate(f_.begin(), f_.begin() + m, f_.end()); for (int i = m; i < 2 * m; ++i) f_.at(i) = 0; ntt(f_, false); for (int i = 0; i < 2 * m; ++i) f_.at(i) *= g_m.at(i); ntt(f_, true); for (int i = 0; i < m; ++i) f_.at(i) = -f_.at(i); res.insert(res.end(), f_.begin(), f_.begin() + m); } res.resize(d); return res; } // Calculate the integral of f(x) // Complexity: O(len(f)) template <class NTTModInt> void integ_inplace(std::vector<NTTModInt> &f) { if (f.empty()) return; for (int i = (int)f.size() - 1; i > 0; --i) f.at(i) = f.at(i - 1) * NTTModInt(i).inv(); f.front() = NTTModInt(0); } // Calculate the derivative of f(x) // Complexity: O(len(f)) template <class NTTModInt> void deriv_inplace(std::vector<NTTModInt> &f) { if (f.empty()) return; for (int i = 1; i < (int)f.size(); ++i) f.at(i - 1) = f.at(i) * i; f.back() = NTTModInt(0); } // Calculate log f(x) mod x^d // Require f(0) = 1 mod x^d // Complexity: O(d log d) template <class NTTModInt> std::vector<NTTModInt> log(const std::vector<NTTModInt> &f, int d = -1) { assert(d >= -1); const int n = f.size(); if (d < 0) d = n; if (d == 0) return {}; assert(f.front() == NTTModInt(1)); std::vector<NTTModInt> inv_f = inv(f, d), df{f.cbegin(), f.cbegin() + std::min(d, n)}; deriv_inplace(df); auto ret = nttconv(inv_f, df); ret.resize(d); integ_inplace(ret); return ret; } template <class NTTModInt> std::vector<NTTModInt> exp(const std::vector<NTTModInt> &h, int d = -1) { assert(d >= -1); const int n = h.size(); if (d < 0) d = n; if (d == 0) return {}; assert(h.empty() or h.front() == NTTModInt(0)); using F = std::vector<NTTModInt>; F g{1}, g_fft; std::vector<NTTModInt> ret(d); ret.front() = 1; auto h_deriv = h; h_deriv.resize(d); deriv_inplace(h_deriv); for (int m = 1; m < d; m *= 2) { F f_fft = ret; f_fft.resize(m * 2); ntt(f_fft, false); // 2a if (m > 1) { F tmp{f_fft.cbegin(), f_fft.cbegin() + m}; for (int i = 0; i < m; ++i) tmp.at(i) *= g_fft.at(i); ntt(tmp, true); tmp.erase(tmp.begin(), tmp.begin() + m / 2); tmp.resize(m); ntt(tmp, false); for (int i = 0; i < m; ++i) tmp.at(i) *= -g_fft.at(i); ntt(tmp, true); tmp.resize(m / 2); g.insert(g.end(), tmp.cbegin(), tmp.cbegin() + m / 2); } // F t{ret.cbegin(), ret.cbegin() + m}; deriv_inplace(t); { F r{h_deriv.cbegin(), h_deriv.cbegin() + m - 1}; r.resize(m); ntt(r, false); for (int i = 0; i < m; ++i) r.at(i) *= f_fft.at(i); ntt(r, true); for (int i = 0; i < m; ++i) t.at(i) -= r.at(i); std::rotate(t.begin(), t.end() - 1, t.end()); } // t.resize(2 * m); ntt(t, false); g_fft = g; g_fft.resize(2 * m); ntt(g_fft, false); for (int i = 0; i < 2 * m; ++i) t.at(i) *= g_fft.at(i); ntt(t, true); t.resize(m); // F v{h.begin() + std::min(m, n), h.begin() + std::min({d, 2 * m, n})}; v.resize(m); t.insert(t.begin(), m - 1, 0); t.push_back(0); integ_inplace(t); for (int i = 0; i < m; ++i) v.at(i) -= t.at(m + i); // v.resize(2 * m); ntt(v, false); for (int i = 0; i < 2 * m; ++i) v.at(i) *= f_fft.at(i); ntt(v, true); v.resize(m); for (int i = 0; i < std::min(d - m, m); ++i) ret.at(m + i) = v.at(i); } return ret; } // Calculate f(x)^k mod x^d // assume 0^0 = 1 template <class NTTModInt> std::vector<NTTModInt> pow(const std::vector<NTTModInt> &A, long long k, int d = -1) { assert(d >= -1); const int n = A.size(); if (d < 0) d = n; if (k == 0) { std::vector<NTTModInt> ret{NTTModInt(1)}; // assume 0^0 = 1 ret.resize(d); return ret; } int l = 0; long long shift = 0; while (l < (int)A.size() and A.at(l) == NTTModInt(0) and shift < d) { ++l; shift += k; } if (l == (int)A.size() or shift >= d) return std::vector<NTTModInt>(d, 0); const NTTModInt cpow = A.at(l).pow(k), cinv = A.at(l).inv(); std::vector<NTTModInt> tmp{A.cbegin() + l, A.cbegin() + std::min<int>(n, d - l * k + l)}; for (auto &x : tmp) x *= cinv; tmp = log(tmp, d - l * k); for (auto &x : tmp) x *= k; tmp = exp(tmp, d - l * k); for (auto &x : tmp) x *= cpow; tmp.insert(tmp.begin(), l * k, NTTModInt(0)); tmp.resize(d); return tmp; } } // namespace fps_nttmod int main() { mint P1, P2, Q1, Q2; int T; cin >> P1 >> P2 >> Q1 >> Q2 >> T; const mint p = P1 / P2, q = Q1 / Q2; // g[t] : t step 後の生存確率 // f[t] : t step 後に増える確率 // g / (1 - f) vector<mint> g(T + 2), f(T + 2); mint alive = 1; g.at(0) = 1; for (int t = 1; t <= T + 1; ++t) { g.at(t) = alive; f.at(t) = alive * p; alive *= q.pow(t); } g.erase(g.begin()); auto den = f; for (auto &x : den) x = -x; den.at(0) += 1; den = fps_nttmod::inv(den, T + 1); auto ret = nttconv(g, den); cout << ret.at(T) << '\n'; }