#include <bits/stdc++.h>
using i32 = int;
using u32 = unsigned int;
using i64 = long long;
using u64 = unsigned long long;
using i128 = __int128_t;
using u128 = __uint128_t;
using f64 = double;
using f80 = long double;
using f128 = __float128;
constexpr i32 operator"" _i32(u64 v) { return v; }
constexpr u32 operator"" _u32(u64 v) { return v; }
constexpr i64 operator"" _i64(u64 v) { return v; }
constexpr u64 operator"" _u64(u64 v) { return v; }
constexpr f64 operator"" _f64(f80 v) { return v; }
constexpr f80 operator"" _f80(f80 v) { return v; }
using Istream = std::istream;
using Ostream = std::ostream;
using Str = std::string;
template<typename T> using Lt = std::less<T>;
template<typename T> using Gt = std::greater<T>;
template<int n> using BSet = std::bitset<n>;
template<typename T1, typename T2> using Pair = std::pair<T1, T2>;
template<typename... Ts> using Tup = std::tuple<Ts...>;
template<typename T, int N> using Arr = std::array<T, N>;
template<typename... Ts> using Deq = std::deque<Ts...>;
template<typename... Ts> using Set = std::set<Ts...>;
template<typename... Ts> using MSet = std::multiset<Ts...>;
template<typename... Ts> using USet = std::unordered_set<Ts...>;
template<typename... Ts> using UMSet = std::unordered_multiset<Ts...>;
template<typename... Ts> using Map = std::map<Ts...>;
template<typename... Ts> using MMap = std::multimap<Ts...>;
template<typename... Ts> using UMap = std::unordered_map<Ts...>;
template<typename... Ts> using UMMap = std::unordered_multimap<Ts...>;
template<typename... Ts> using Vec = std::vector<Ts...>;
template<typename... Ts> using Stack = std::stack<Ts...>;
template<typename... Ts> using Queue = std::queue<Ts...>;
template<typename T> using MaxHeap = std::priority_queue<T>;
template<typename T> using MinHeap = std::priority_queue<T, Vec<T>, Gt<T>>;
template<typename T> using Opt = std::optional<T>;
template<typename... Ts> using Span = std::span<Ts...>;
constexpr bool LOCAL = false;
template<typename T> static constexpr T OjLocal(T oj, T local) { return LOCAL ? local : oj; }
template<typename T> constexpr T LIMMIN = std::numeric_limits<T>::min();
template<typename T> constexpr T LIMMAX = std::numeric_limits<T>::max();
template<typename T> constexpr T INF = (LIMMAX<T> - 1) / 2;
template<typename T = i64> constexpr T TEN(int N) { return N == 0 ? T{1} : TEN<T>(N - 1) * T{10}; }
constexpr auto ABS(auto x) { return (x >= 0 ? x : -x); }
constexpr auto makePair(const auto& x1, const auto& x2) { return std::make_pair(x1, x2); }
constexpr auto makeTup(const auto&... xs) { return std::make_tuple(xs...); }
template<typename T> constexpr bool chmin(T& x, const T& y, auto comp) { return (comp(y, x) ? (x = y, true) : false); }
template<typename T> constexpr bool chmin(T& x, const T& y) { return chmin(x, y, Lt<T>{}); }
template<typename T> constexpr bool chmax(T& x, const T& y, auto comp) { return (comp(x, y) ? (x = y, true) : false); }
template<typename T> constexpr bool chmax(T& x, const T& y) { return chmax(x, y, Lt<T>{}); }
constexpr i64 floorDiv(i64 x, i64 y)
{
assert(y != 0);
if (y < 0) { x = -x, y = -y; }
return x >= 0 ? x / y : (x - y + 1) / y;
}
constexpr i64 ceilDiv(i64 x, i64 y)
{
assert(y != 0);
if (y < 0) { x = -x, y = -y; }
return x >= 0 ? (x + y - 1) / y : x / y;
}
template<typename T> constexpr T powerSemiGroup(const T& x, i64 N, auto mul)
{
assert(N > 0);
if (N == 1) { return x; }
return (N % 2 == 1 ? mul(x, powerSemiGroup(x, N - 1, mul)) : powerSemiGroup(mul(x, x), N / 2, mul));
}
template<typename T> constexpr auto powerSemiGroup(const auto& v, i64 N) { return powerSemiGroup(v, N, std::multiplies<T>{}); }
template<typename T> constexpr T powerMonoid(const T& x, i64 N, const T& e, auto mul)
{
assert(N >= 0);
return (N == 0 ? e : powerSemiGroup(x, N, mul));
}
template<typename T> constexpr T powerMonoid(T x, i64 N, const T& e) { return powerMonoid(x, N, e, std::multiplies<T>{}); }
template<typename T> constexpr T powerInt(T x, i64 N) { return powerMonoid(x, N, T{1}); }
constexpr u64 powerMod(u64 x, i64 N, u64 mod)
{
assert(0 < mod);
return powerMonoid(x, N, u64{1}, [&](u64 x, u64 y) {
if (mod <= (u64)LIMMAX<u32>) {
return x * y % mod;
} else {
return (u64)((u128)x * y % mod);
}
});
}
constexpr auto sumAll(const auto& xs) { return std::accumulate(std::begin(xs), std::end(xs), decltype(xs[0]){}); }
constexpr int lbInd(const auto& xs, const auto& x) { return std::ranges::lower_bound(xs, x) - std::begin(xs); }
constexpr int ubInd(const auto& xs, const auto& x) { return std::ranges::upper_bound(xs, x) - std::begin(xs); }
constexpr int find(const auto& xs, const auto& x)
{
const int i = lbInd(xs, x);
if (i < std::ssize(xs) and xs[i] == x) { return i; }
return std::ssize(xs);
}
constexpr void concat(auto& xs1, const auto& xs2) { std::ranges::copy(xs2, std::back_inserter(xs1)); }
constexpr auto concatCopy(const auto& xs1, const auto& xs2)
{
auto Ans = xs1;
return concat(Ans, xs2), Ans;
}
template<typename Ts, typename T> constexpr void fillAll(Ts& arr, const T& x)
{
if constexpr (std::is_convertible<T, Ts>::value) {
arr = x;
} else {
for (auto& subarr : arr) { fillAll(subarr, x); }
}
}
template<typename T> constexpr Vec<T> genVec(int N, auto gen)
{
Vec<T> ans;
std::ranges::generate_n(std::back_inserter(ans), N, gen);
return ans;
}
constexpr auto minAll(const auto& xs, auto... args) { return std::ranges::min(xs, args...); }
constexpr auto maxAll(const auto& xs, auto... args) { return std::ranges::max(xs, args...); }
constexpr int minInd(const auto& xs, auto... args) { return std::ranges::min_element(xs, args...) - std::begin(xs); }
constexpr int maxInd(const auto& xs, auto... args) { return std::ranges::max_element(xs, args...) - std::begin(xs); }
template<typename T = int> constexpr Vec<T> iotaVec(int N, T offset = 0)
{
Vec<T> ans(N);
std::iota(std::begin(ans), std::end(ans), offset);
return ans;
}
constexpr void plusAll(auto& xs, const auto& x)
{
std::ranges::for_each(xs, [&](auto& e) { e += x; });
}
constexpr void reverseAll(auto& xs) { std::ranges::reverse(xs); }
constexpr void sortAll(auto& xs, auto... args) { std::ranges::sort(xs, args...); }
constexpr auto runLengthEncode(const auto& xs)
{
using T = typename std::decay<decltype(xs[0])>::type;
Vec<Pair<T, int>> Ans;
auto [l, px] = makePair(0, T{});
for (const T& x : xs) {
if (l == 0 or x != px) {
if (l > 0) { Ans.push_back({px, l}); }
l = 1, px = x;
} else {
l++;
}
}
if (l > 0) { Ans.push_back({px, l}); }
return Ans;
}
inline Ostream& operator<<(Ostream& os, u128 v)
{
Str ans;
if (v == 0) { ans = "0"; }
while (v) { ans.push_back('0' + v % 10), v /= 10; }
std::reverse(ans.begin(), ans.end());
return os << ans;
}
inline Ostream& operator<<(Ostream& os, i128 v)
{
bool minus = false;
if (v < 0) { minus = true, v = -v; }
return os << (minus ? "-" : "") << (u128)v;
}
constexpr bool isBitOn(u64 x, int i) { return assert(0 <= i and i < 64), ((x >> i) & 1_u64); }
constexpr bool isBitOff(u64 x, int i) { return assert(0 <= i and i < 64), (not isBitOn(x, i)); }
constexpr u64 bitMask(int w) { return assert(0 <= w and w <= 64), (w == 64 ? ~0_u64 : (1_u64 << w) - 1); }
constexpr u64 bitMask(int s, int e) { return assert(0 <= s and s <= e and e <= 64), (bitMask(e - s) << s); }
constexpr int floorLog2(u64 x) { return 63 - std::countl_zero(x); }
constexpr int ceilLog2(u64 x) { return x == 0 ? -1 : std::bit_width(x - 1); }
constexpr int order2(u64 x) { return std::countr_zero(x); }
template<typename F> struct Fix : F
{
constexpr Fix(F&& f) : F{std::forward<F>(f)} {}
template<typename... Args> constexpr auto operator()(Args&&... args) const
{
return F::operator()(*this, std::forward<Args>(args)...);
}
};
class irange
{
private:
struct itr
{
constexpr itr(i64 start, i64 end, i64 step) : m_cnt{start}, m_step{step}, m_end{end} {}
constexpr bool operator!=(const itr&) const { return (m_step > 0 ? m_cnt < m_end : m_end < m_cnt); }
constexpr i64 operator*() { return m_cnt; }
constexpr itr& operator++() { return m_cnt += m_step, *this; }
i64 m_cnt, m_step, m_end;
};
i64 m_start, m_end, m_step;
public:
constexpr irange(i64 start, i64 end, i64 step = 1) : m_start{start}, m_end{end}, m_step{step} { assert(step != 0); }
constexpr itr begin() const { return itr{m_start, m_end, m_step}; }
constexpr itr end() const { return itr{m_start, m_end, m_step}; }
};
constexpr irange rep(i64 end) { return irange(0, end, 1); }
constexpr irange per(i64 rend) { return irange(rend - 1, -1, -1); }
template<typename Engine> class RNG
{
public:
using result_type = typename Engine::result_type;
using U = result_type;
static constexpr U min() { return Engine::min(); }
static constexpr U max() { return Engine::max(); }
RNG() : RNG(std::random_device{}()) {}
RNG(U seed) : m_rng(seed) {}
U operator()() { return m_rng(); }
template<typename T>
requires std::is_integral_v<T>
T val(T min, T max)
{
return std::uniform_int_distribution<T>(min, max)(m_rng);
}
template<typename T> Vec<T> vec(int N, T min, T max)
{
return genVec<T>(N, [&]() { return val<T>(min, max); });
}
private:
Engine m_rng;
};
inline RNG<std::mt19937> rng;
inline RNG<std::mt19937_64> rng64;
constexpr i64 binSearch(i64 ng, i64 ok, auto check)
{
while (ABS(ok - ng) > 1) {
const i64 mid = (ok + ng) / 2;
(check(mid) ? ok : ng) = mid;
}
return ok;
}
constexpr f80 binSearch(f80 ng, f80 ok, auto check, int times)
{
for (auto _temp_name_0 [[maybe_unused]] : rep(times)) {
const f80 mid = (ok + ng) / 2;
(check(mid) ? ok : ng) = mid;
}
return ok;
}
constexpr u64 intNthRoot(u64 A, int K)
{
assert(K > 0);
if (A == 0) { return 0; }
if (K == 1) { return A; }
if (K > 64) { return 1; }
return binSearch(1_i64 << 32, 1_i64, [&](i64 a) {
u64 x = (u64)1;
for (auto _temp_name_1 [[maybe_unused]] : rep(K)) {
if (x > A / a) { return false; }
x *= a;
}
return true;
});
}
constexpr u64 intSqrt(u64 A) { return intNthRoot(A, 2); }
constexpr Pair<i64, i64> extgcd(i64 a, i64 b)
{
assert(a != 0 or b != 0);
const i64 A = ABS(a), B = ABS(b);
auto [x, y, g] = Fix([&](auto self, i64 a, i64 b) -> Tup<i64, i64, i64> {
assert(0 <= a and a < b);
if (a == 0) { return {0, 1, b}; }
const auto [px, py, pg] = self(b % a, a);
return {py - (b / a) * px, px, pg};
})(std::ranges::min(A, B), std::ranges::max(A, B));
if (A > B) { std::swap(x, y); }
if (a < 0) { x = -x; }
if (b < 0) { y = -y; }
if (x < 0) { x += B / g, y -= (b > 0 ? a / g : -a / g); }
return {x, y};
}
constexpr i64 inverseMod(i64 a, i64 mod)
{
assert(mod > 0 and a % mod != 0);
return extgcd(a % mod, mod).first;
}
template<u64 Mod, bool dynamic = false>
requires(dynamic or (0 < Mod and Mod <= (u64)LIMMAX<i64>))
class modint
{
public:
static constexpr bool isDynamic() { return dynamic; }
static constexpr u64 mod()
{
if constexpr (dynamic) {
return modRef();
} else {
return Mod;
}
}
static constexpr void setMod(u64 m)
requires dynamic
{
assert(0 < m and m <= LIMMAX<i64>), modRef() = m;
}
constexpr modint() : m_val{0} {}
constexpr modint(i64 v) : m_val{normll(v)} {}
constexpr friend modint operator-(const modint& x) { return modint{0} - x; }
constexpr friend modint& operator+=(modint& x1, const modint& x2) { return x1.m_val = norm(x1.m_val + x2.m_val), x1; }
constexpr friend modint& operator-=(modint& x1, const modint& x2) { return x1.m_val = norm(x1.m_val + mod() - x2.m_val), x1; }
constexpr friend modint& operator*=(modint& x1, const modint& x2)
{
if constexpr (dynamic) {
if (mod() <= (u64)LIMMAX<u32>) {
return (x1.m_val *= x2.m_val) %= mod(), x1;
} else {
return x1.m_val = (u64)((u128)x1.m_val * (u128)x2.m_val % (u128)mod()), x1;
}
} else {
if constexpr (Mod <= (u64)LIMMAX<u32>) {
return (x1.m_val *= x2.m_val) %= mod(), x1;
} else {
return x1.m_val = (u64)((u128)x1.m_val * (u128)x2.m_val % (u128)mod()), x1;
}
}
}
constexpr friend modint& operator/=(modint& x1, const modint& x2) { return x1 *= x2.inv(); }
constexpr friend modint operator+(const modint& x1, const modint& x2)
{
auto ans = x1;
return ans += x2;
}
constexpr friend modint operator-(const modint& x1, const modint& x2)
{
auto ans = x1;
return ans -= x2;
}
constexpr friend modint operator*(const modint& x1, const modint& x2)
{
auto ans = x1;
return ans *= x2;
}
constexpr friend modint operator/(const modint& x1, const modint& x2)
{
auto ans = x1;
return ans /= x2;
}
constexpr friend bool operator==(const modint& x1, const modint& x2) { return x1.m_val == x2.m_val; }
friend Istream& operator>>(Istream& is, modint& x)
{
i64 v;
return is >> v, x = v, is;
}
friend Ostream& operator<<(Ostream& os, const modint& x) { return os << x.m_val; }
constexpr u64 val() const { return m_val; }
constexpr modint pow(i64 n) const { return powerInt(*this, n); }
constexpr modint inv() const { return inverseMod(m_val, mod()); }
static modint sinv(int n)
{
assert(1 <= n);
auto& is = sinvRef();
for (int i : irange((int)is.size(), n + 1)) { is.push_back(-is[mod() % i] * (mod() / i)); }
return is[n];
}
static modint fact(int n)
{
assert(0 <= n);
auto& fs = factRef();
for (int i : irange((int)fs.size(), n + 1)) { fs.push_back(fs.back() * i); }
return fs[n];
}
static modint ifact(int n)
{
auto& ifs = ifactRef();
for (int i : irange((int)ifs.size(), n + 1)) { ifs.push_back(ifs.back() * sinv(i)); }
return ifs[n];
}
static modint perm(int n, int k) { return k > n or k < 0 ? modint{0} : fact(n) * ifact(n - k); }
static modint comb(int n, int k) { return k > n or k < 0 ? modint{0} : fact(n) * ifact(n - k) * ifact(k); }
private:
static u64& modRef()
requires dynamic
{
static u64 mod_ = 0;
return mod_;
}
static Vec<modint>& sinvRef()
{
static Vec<modint> is{1, 1};
return is;
}
static Vec<modint>& factRef()
{
static Vec<modint> fs{1, 1};
return fs;
}
static Vec<modint>& ifactRef()
{
static Vec<modint> ifs{1, 1};
return ifs;
}
static constexpr u64 norm(u64 x) { return x < mod() ? x : x - mod(); }
static constexpr u64 normll(i64 x)
{
x %= (i64)mod();
return norm(u64(x < 0 ? x + (i64)mod() : x));
}
u64 m_val;
};
using modint_1000000007 = modint<1000000007, false>;
using modint_998244353 = modint<998244353, false>;
template<u64 id>
requires(id < (u64)LIMMAX<i64>)
using modint_dynamic = modint<id, true>;
template<u64 id>
requires(id < (u64)LIMMAX<i64>)
using modint_dynamic_reserved = modint<id | (1_u64 << 63), true>;
constexpr bool millerRabin(u64 X, const Vec<u64>& as)
{
using mint = modint_dynamic_reserved<81165>;
mint::setMod(X);
const u64 d = (X - 1) >> order2(X - 1);
for (u64 a : as) {
if (X <= a) { break; }
u64 s = d;
mint x = mint(a).pow(s);
while (x != 1 and x != X - 1 and s != X - 1) { s *= 2, x *= x; }
if (x != X - 1 and s % 2 == 0) { return false; }
}
return true;
}
constexpr bool isPrime(u64 X)
{
if (X == 1) { return false; }
if (X % 2 == 0) { return X == 2; }
if (X < 4759123141_u64) {
return millerRabin(X, {2_u64, 7_u64, 61_u64});
} else {
return millerRabin(X, {2_u64, 325_u64, 9375_u64, 28178_u64, 450775_u64, 9780504_u64, 1795265022_u64});
}
}
constexpr u64 pollardRho(u64 X)
{
assert(1 <= X and X <= (u64)LIMMAX<i64>);
if (X % 2 == 0) { return 2; }
if (X == 1 or isPrime(X)) { return X; }
using mint = modint_dynamic_reserved<77726>;
mint::setMod(X);
auto f = [&](mint x, mint c) { return x * x + c; };
const u64 gcdBlock = intNthRoot(X, 8);
while (true) {
const u64 a = rng64.val<u64>(0, X - 1);
const mint c = rng64.val(2_u64, X - 1);
mint x = a, y = a, sx = x, sy = y;
mint p = 1;
u64 g = 1;
while (g == 1) {
sx = x, sy = y;
for (auto _temp_name_2 [[maybe_unused]] : rep(gcdBlock)) { x = f(x, c), y = f(f(y, c), c), p *= (x - y); }
g = std::gcd(X, p.val());
}
if (g == X) {
x = sx, y = sy, g = 1;
while (g == 1) { x = f(x, c), y = f(f(y, c), c), g = std::gcd(X, (x - y).val()); }
}
if (g != X) { return g; }
}
return X;
}
constexpr Vec<Pair<u64, int>> primeFactors(u64 X)
{
Vec<u64> Ans;
Fix([&](auto dfs, u64 x) -> void {
while (x % 2 == 0) { x /= 2, Ans.push_back(2); }
if (x == 1) { return; }
const u64 d = pollardRho(x);
if (d == x) { return Ans.push_back(d), void(); }
dfs(d), dfs(x / d);
})(X);
sortAll(Ans);
return runLengthEncode(Ans);
}
constexpr Vec<u64> divisors(const Vec<Pair<u64, int>>& factors)
{
Vec<u64> Ans{1};
for (const auto& [p, e] : factors) {
const int dn = (int)Ans.size();
u64 pe = p;
for (auto _temp_name_3 [[maybe_unused]] : rep(e)) {
for (int j : rep(dn)) { Ans.push_back(Ans[j] * pe); }
pe *= p;
}
}
return sortAll(Ans), Ans;
}
constexpr u64 primitiveRoot(u64 P)
{
assert(P >= 2);
Vec<u64> ps;
for (const auto& e : primeFactors(P - 1)) { ps.push_back(e.first); }
for (u64 r = 1; r < P; r++) {
bool ok = true;
for (u64 p : ps) {
const u64 Q = powerMod(r, (P - 1) / p, P);
if (Q == 1) {
ok = false;
break;
}
}
if (ok) { return r; }
}
return 0;
}
template<typename mint> class NumberTheoriticTransform
{
private:
static constexpr u64 mod() { return mint::mod(); }
static constexpr u64 root() { return primitiveRoot(mint::mod()); }
static constexpr int order() { return order2(mint::mod() - 1); }
public:
static constexpr Vec<mint> convolute(Vec<mint> F, Vec<mint> G)
{
const int A = F.size(), B = G.size();
const int C = A + B - 1;
const int N = (int)std::bit_ceil((u64)C);
F.resize(N, 0), G.resize(N, 0);
transform(F, false), transform(G, false);
for (int i : rep(N)) { F[i] *= G[i]; }
transform(F, true), F.resize(C);
return F;
}
static constexpr void transform(Vec<mint>& F, bool rev)
{
const int N = F.size();
assert((N & (N - 1)) == 0);
assert(N <= (1 << order()));
if (N == 1) { return; }
const int L = floorLog2(N);
const auto l_range = (rev ? irange(1, L + 1, 1) : irange(L, 0, -1));
for (int l : l_range) {
const int H = (1 << l), B = N / H;
for (int b : rep(B)) {
const mint W = zeta(l, rev);
for (mint W_h = 1; int h : rep(H / 2)) {
const int y1 = H * b + h, y2 = y1 + H / 2;
const mint f1 = F[y1], f2 = F[y2];
const mint nf1 = (rev ? f1 + f2 * W_h : f1 + f2), nf2 = (rev ? f1 - f2 * W_h : (f1 - f2) * W_h);
F[y1] = nf1, F[y2] = nf2;
W_h *= W;
}
}
}
if (rev) {
const mint iN = mint{1} / N;
for (auto& a : F) { a *= iN; }
}
}
private:
static mint zeta(int i, bool rev)
{
static Vec<mint> zs, izs;
if (zs.empty()) {
const auto MOD = mod(), ROOT = root();
const auto LMAX = order();
zs.resize(LMAX + 1, 1), izs.resize(LMAX + 1, 1);
zs[LMAX] = mint(ROOT).pow((MOD - 1) / (1 << LMAX)), izs[LMAX] = zs[LMAX].inv();
for (int l : per(LMAX)) { zs[l] = zs[l + 1] * zs[l + 1], izs[l] = izs[l + 1] * izs[l + 1]; }
}
return (rev ? izs[i] : zs[i]);
}
};
class Garner final
{
public:
explicit Garner() = delete;
~Garner() = delete;
template<typename mint, typename mint1, typename mint2> static constexpr mint restoreMod(const mint1& x1, const mint2& x2)
{
const u64 m1 = mint1::mod();
const auto [y0, y1] = coeff(x1, x2);
return mint(y0) + mint(y1) * m1;
}
template<typename mint, typename mint1, typename mint2, typename mint3>
static constexpr mint restoreMod(const mint1& x1, const mint2& x2, const mint3& x3)
{
const u64 m1 = mint1::mod(), m2 = mint2::mod();
const auto [y0, y1, y2] = coeff(x1, x2, x3);
return mint(y0) + mint(y1) * m1 + mint(y2) * m1 * m2;
}
template<typename mint1, typename mint2> static constexpr i64 restoreInt(const mint1& x1, const mint2& x2)
{
const u64 m1 = mint1::mod(), m2 = mint2::mod();
const auto [y0, y1] = coeff(x1, x2);
const i128 M = (i128)m1 * m2;
const i128 X = i128(y0) + i128(y1) * m1;
return (i64)(M / 2 >= X ? X : X - M);
}
template<typename mint1, typename mint2, typename mint3>
static constexpr i64 restoreInt(const mint1& x1, const mint2& x2, const mint3& x3)
{
const u64 m1 = mint1::mod(), m2 = mint2::mod(), m3 = mint3::mod();
const auto [y0, y1, y2] = coeff(x1, x2, x3);
const i128 M = (i128)m1 * m2 * m3;
const i128 X = i128(y0) + i128(y1) * m1 + i128(y2) * m1 * m2;
return (i64)(M / 2 >= X ? X : X - M);
}
private:
template<typename mint1, typename mint2> static constexpr Pair<i64, i64> coeff(const mint1& x1, const mint2& x2)
{
const u64 m1 = mint1::mod();
const mint2 m1_inv = mint2(m1).inv();
const i64 y0 = x1.val();
const i64 y1 = ((x2 - mint2(y0)) * m1_inv).val();
return {y0, y1};
}
template<typename mint1, typename mint2, typename mint3>
static constexpr Tup<i64, i64, i64> coeff(const mint1& x1, const mint2& x2, const mint3& x3)
{
const u64 m1 = mint1::mod(), m2 = mint2::mod();
const mint2 m1_inv = mint2(m1).inv();
const mint3 m1m2_inv = (mint3(m1) * mint3(m2)).inv();
const i64 y0 = x1.val();
const i64 y1 = ((x2 - mint2(y0)) * m1_inv).val();
const i64 y2 = ((x3 - mint3(y0) - mint3(y1) * m1) * m1m2_inv).val();
return {y0, y1, y2};
}
};
template<typename mint> constexpr Vec<mint> convoluteMod(const Vec<mint>& F, const Vec<mint>& G)
{
const int LMAX = order2(mint::mod() - 1), NMAX = (1 << LMAX);
const int A = F.size(), B = G.size();
if (A == 0 or B == 0) { return {}; }
const int N = A + B - 1;
if (std::min(A, B) <= 100) {
Vec<mint> ans(N, 0);
for (int i : rep(A)) {
for (int j : rep(B)) { ans[i + j] += F[i] * G[j]; }
}
return ans;
}
if (N <= NMAX) {
return NumberTheoriticTransform<mint>::convolute(F, G);
} else {
assert(N <= (1 << 24));
using submint1 = modint<469762049>;
using submint2 = modint<167772161>;
using submint3 = modint<754974721>;
Vec<submint1> as1(A), bs1(B);
Vec<submint2> as2(A), bs2(B);
Vec<submint3> as3(A), bs3(B);
for (int i : rep(A)) { as1[i] = F[i].val(), as2[i] = F[i].val(), as3[i] = F[i].val(); }
for (int i : rep(B)) { bs1[i] = G[i].val(), bs2[i] = G[i].val(), bs3[i] = G[i].val(); }
const auto cs1 = NumberTheoriticTransform<submint1>::convolute(as1, bs1);
const auto cs2 = NumberTheoriticTransform<submint2>::convolute(as2, bs2);
const auto cs3 = NumberTheoriticTransform<submint3>::convolute(as3, bs3);
Vec<mint> cs(N);
for (int i : rep(N)) { cs[i] = Garner::restoreMod<mint>(cs1[i], cs2[i], cs3[i]); }
return cs;
}
}
constexpr Vec<i64> convoluteInt(const Vec<i64>& F, const Vec<i64>& G)
{
const int A = F.size(), B = G.size();
if (A == 0 or B == 0) { return {}; }
const int N = A + B - 1;
if (std::min(A, B) <= 100) {
Vec<i64> ans(N, 0);
for (int i : rep(A)) {
for (int j : rep(B)) { ans[i + j] += F[i] * G[j]; }
}
return ans;
}
assert(N <= (1 << 24));
using submint1 = modint<469762049>;
using submint2 = modint<167772161>;
using submint3 = modint<754974721>;
Vec<submint1> as1(A), bs1(B);
Vec<submint2> as2(A), bs2(B);
Vec<submint3> as3(A), bs3(B);
for (int i : rep(A)) { as1[i] = F[i], as2[i] = F[i], as3[i] = F[i]; }
for (int i : rep(B)) { bs1[i] = G[i], bs2[i] = G[i], bs3[i] = G[i]; }
const auto cs1 = NumberTheoriticTransform<submint1>::convolute(as1, bs1);
const auto cs2 = NumberTheoriticTransform<submint2>::convolute(as2, bs2);
const auto cs3 = NumberTheoriticTransform<submint3>::convolute(as3, bs3);
Vec<i64> cs(N);
for (int i : rep(N)) { cs[i] = Garner::restoreInt(cs1[i], cs2[i], cs3[i]); }
return cs;
}
template<typename mint> constexpr Vec<mint> convoluteModReverse(Vec<mint> F, const Vec<mint>& G)
{
const int A = (int)F.size(), B = (int)G.size();
reverseAll(F);
const auto cs = convoluteMod(F, G);
Vec<mint> ans(B);
for (int i : rep(B)) { ans[i] = cs[i + A - 1]; }
return ans;
}
constexpr Vec<i64> convoluteIntReverse(Vec<i64> F, const Vec<i64>& G)
{
const int A = (int)F.size(), B = (int)G.size();
reverseAll(F);
const auto cs = convoluteInt(F, G);
Vec<i64> ans(B);
for (int i : rep(B)) { ans[i] = cs[i + A - 1]; }
return ans;
}
template<typename mint> class Polynomial : public Vec<mint>
{
using Vec<mint>::resize;
using Vec<mint>::push_back;
using Vec<mint>::pop_back;
using Vec<mint>::back;
public:
using Vec<mint>::size;
explicit constexpr Polynomial(const auto&... args) : Vec<mint>{args...} {}
constexpr const Vec<mint>& asVec() const { return static_cast<const Vec<mint>&>(*this); }
constexpr int deg() const { return (int)size() - 1; }
constexpr void shrink(i64 n)
{
if (n >= (i64)size()) { return; }
resize(n), normalize();
}
constexpr Polynomial low(i64 n) const
{
const int sz = (int)std::min(n, (i64)size());
return Polynomial{this->begin(), this->begin() + (int)sz};
}
constexpr mint& operator[](int n)
{
if (n >= (int)size()) { resize(n + 1, mint{}); }
return Vec<mint>::operator[](n);
}
constexpr mint at(i64 n) const { return (n < (i64)size() ? (*this)[n] : mint{0}); }
constexpr friend Polynomial operator-(const Polynomial& F)
{
Polynomial ans = F;
for (auto& v : ans) { v = -v; }
return ans;
}
constexpr friend Polynomial& operator+=(Polynomial& F, const Polynomial& G)
{
for (int i : rep(G.size())) { F[i] += G[i]; }
return F;
}
constexpr friend Polynomial& operator-=(Polynomial& F, const Polynomial& G)
{
for (int i : rep(G.size())) { F[i] -= G[i]; }
return F;
}
constexpr friend Polynomial& operator*=(Polynomial& F, const Polynomial& G) { return F = F * G; }
constexpr friend Polynomial& operator/=(Polynomial& F, const Polynomial& G) { return F = F / G; }
constexpr friend Polynomial& operator%=(Polynomial& F, const Polynomial& G) { return F = F % G; }
constexpr friend Polynomial& operator<<=(Polynomial& F, int s) { return F = (F << s); }
constexpr friend Polynomial& operator>>=(Polynomial& F, int s) { return F = (F >> s); }
constexpr friend Polynomial operator+(const Polynomial& F, const Polynomial& G)
{
auto ans{F};
return ans += G;
}
constexpr friend Polynomial operator-(const Polynomial& F, const Polynomial& G)
{
auto ans{F};
return ans -= G;
}
constexpr friend Polynomial operator*(const Polynomial& F, const Polynomial& G) { return Polynomial{convoluteMod(F, G)}; }
constexpr friend Polynomial operator/(const Polynomial& F, const Polynomial& G) { return quot(F, G); }
constexpr friend Polynomial operator%(const Polynomial& F, const Polynomial& G) { return F - G * (F / G); }
constexpr friend Polynomial operator<<(const Polynomial& F, int s)
{
assert(s >= 0);
Polynomial ans;
for (int i : irange(s, (int)F.size())) { ans[i - s] = F[i]; }
return ans;
}
constexpr friend Polynomial operator>>(const Polynomial& F, int s)
{
assert(s >= 0);
Polynomial ans;
for (int i : rep(F.size())) { ans[i + s] = F[i]; }
return ans;
}
constexpr int order() const
{
for (int i : rep(size())) {
if ((*this)[i] != 0) { return i; }
}
return size();
}
constexpr bool isZero() const { return (size() == 1) and ((*this)[0] == 0); }
constexpr Polynomial derivative() const
{
Polynomial ans;
for (int i : irange(1, (int)size())) { ans[i - 1] = (*this)[i] * i; }
return ans;
}
constexpr Polynomial integral() const
{
Polynomial ans;
for (int i : irange(1, (int)size() + 1)) { ans[i] = (*this)[i - 1] / i; }
return ans;
}
constexpr Polynomial inv(int n) const
{
assert(n >= 1);
assert((*this)[0].val() != 0);
Polynomial G{(*this)[0].inv()};
for (int m = 1; m < n; m *= 2) {
const auto F = low(m * 2);
G = ((-(F * G).low(m * 2) + Polynomial{2}) * G).low(2 * m);
}
return G.low(n);
}
constexpr Polynomial log(int n) const
{
assert(n >= 1);
if (n == 1) { return Polynomial{0}; }
assert((*this)[0] == 1);
return (low(n).derivative() * inv(n)).low(n - 1).integral();
}
constexpr Polynomial exp(int n) const
{
assert(n >= 1);
assert((*this)[0] == 0);
Polynomial G{1};
for (int m = 1; m < n; m *= 2) {
const auto F = low(m * 2);
G = (G * (F - G.log(m * 2) + Polynomial{1})).low(m * 2);
}
return G.low(n);
}
constexpr Polynomial pow(i64 m, int n) const
{
assert(n >= 1);
if (m == 0) { return Polynomial{1}; }
if (isZero()) { return Polynomial{0}; }
const int k = order();
if ((i64)k >= ceilDiv((i64)n, m)) { return Polynomial{}; }
n -= k * m;
auto f = ((*this) << k).low(n);
const mint c = f[0];
f *= Polynomial{c.inv()};
return ((f.log(n) * Polynomial{m}).exp(n) * Polynomial{c.pow(m)}) >> (k * m);
}
constexpr Polynomial taylorShift(const mint& c) const
{
const int N = (int)size();
Vec<mint> as(N), bs(N);
for (mint ce = 1; int i : rep(N)) {
as[i] = ce * mint::ifact(i), bs[i] = (*this)[i] * mint::fact(i);
ce *= c;
}
auto cs = convoluteModReverse(as, bs);
for (int i : rep(N)) { cs[i] *= mint::ifact(i); }
return Polynomial{cs};
}
constexpr mint eval(const mint& x) const
{
mint ans = 0;
for (mint xe = 1; int i : rep(size())) { ans += xe * (*this)[i], xe *= x; }
return ans;
}
private:
const mint& operator[](const int n) const { return assert(n < (int)size()), Vec<mint>::operator[](n); }
static constexpr Polynomial quot(const Polynomial& F, const Polynomial& G)
{
const int A = (int)F.deg(), B = (int)G.deg();
if (A < B) { return Polynomial{0}; }
Polynomial Frev = F, Grev = G;
reverseAll(Frev), reverseAll(Grev);
Polynomial qrev = (Frev.low(A - B + 1) * Grev.inv(A - B + 1)).low(A - B + 1);
return reverseAll(qrev), qrev;
}
constexpr void normalize()
{
while (size() > 0 and back() == 0) { pop_back(); }
if (size() == 0) { push_back(0); }
}
};
template<typename T> constexpr Vec<T> berlekampMassey(const Vec<T>& As)
{
const int N = (int)As.size();
Vec<T> C{1};
Tup<Vec<T>, int, T> lastFailure = {{}, -1, 0};
auto eval = [&](const Vec<T>& C, int i) {
assert(not C.empty());
const int D = (int)C.size() - 1;
if (i < D) { return T{}; }
T s = 0;
for (int j : rep(D + 1)) { s += As[i - j] * C[j]; }
return s;
};
for (int i : rep(N)) {
const T d = eval(C, i);
if (d == 0) { continue; }
const Tup<Vec<T>, int, T> nlastFailure = {C, i, d};
const auto& [pC, pi, pd] = lastFailure;
const int nl = i - (int)C.size(), pl = pi - (int)pC.size();
if (pi != -1) {
const bool updateFailure = (pl < nl);
const auto [c, s] = makePair(d / pd, i - pi);
if (C.size() < pC.size() + s) { C.resize(pC.size() + s); }
for (int j : rep(pC.size())) { C[j + s] -= pC[j] * c; }
if (updateFailure) { lastFailure = nlastFailure; }
} else {
C = Vec<T>(i + 2), C[0] = 1;
lastFailure = nlastFailure;
}
}
return C;
}
template<typename mint> constexpr mint bostanMori(Polynomial<mint> F, Polynomial<mint> G, i64 N)
{
assert(not G.isZero());
const int ford = F.order(), gord = G.order();
assert(gord <= ford);
if (N < (i64)(ford - gord)) { return 0; }
F <<= ford, G <<= gord, N -= (ford - gord);
while (N > 0) {
F.shrink(N + 1), G.shrink(N + 1);
Polynomial<mint> mG;
for (int i : rep(G.size())) { mG[i] = G[i] * (i % 2 == 0 ? 1 : -1); }
const auto nG = G * mG;
const auto nF = F * mG;
for (int i : rep(ceilDiv(nF.size(), 2))) { F[i] = nF.at(2 * i + (N % 2)); }
for (int i : rep(ceilDiv(nG.size(), 2))) { G[i] = nG.at(2 * i); }
N /= 2;
}
return F[0] / G[0];
}
template<typename mint> constexpr mint guessNthTerm(const Vec<mint>& As, i64 N)
{
const Polynomial<mint> G{berlekampMassey(As)};
const int L = G.size();
if (L == 1) { return As[0]; }
const auto F = (Polynomial<mint>{As} * G).low(L - 1);
return bostanMori(F, G, N);
}
class Printer
{
public:
Printer(Ostream& os = std::cout) : m_os{os} { m_os << std::fixed << std::setprecision(15); }
int operator()(const auto&... args) { return dump(args...), 0; }
int ln(const auto&... args) { return dump(args...), m_os << '\n', 0; }
int el(const auto&... args) { return dump(args...), m_os << std::endl, 0; }
private:
void dump(const auto& v) { m_os << v; }
int dump(const auto& v, const auto&... args) { return dump(v), m_os << ' ', dump(args...), 0; }
template<typename... Args> void dump(const Vec<Args...>& vs)
{
for (Str delim = ""; const auto& v : vs) { m_os << std::exchange(delim, " "), dump(v); }
}
Ostream& m_os;
};
inline Printer out;
class Scanner
{
public:
Scanner(Istream& is = std::cin) : m_is{is} { m_is.tie(nullptr)->sync_with_stdio(false); }
template<typename T> T val()
{
T v;
return m_is >> v, v;
}
template<typename T> T val(T offset) { return val<T>() - offset; }
template<typename T> Vec<T> vec(int N)
{
return genVec<T>(N, [&]() { return val<T>(); });
}
template<typename T> Vec<T> vec(int N, T offset)
{
return genVec<T>(N, [&]() { return val<T>(offset); });
}
template<typename T> Vec<Vec<T>> vvec(int n, int m)
{
return genVec<Vec<T>>(n, [&]() { return vec<T>(m); });
}
template<typename T> Vec<Vec<T>> vvec(int n, int m, const T offset)
{
return genVec<Vec<T>>(n, [&]() { return vec<T>(m, offset); });
}
template<typename... Args> auto tup() { return Tup<Args...>{val<Args>()...}; }
template<typename... Args> auto tup(Args... offsets) { return Tup<Args...>{val<Args>(offsets)...}; }
private:
Istream& m_is;
};
inline Scanner in;
int main()
{
using mint = modint_998244353;
const auto [D, K] = in.tup<int, i64>();
const auto as = in.vec<mint>(D);
const auto cs = concatCopy(Vec<mint>{1}, in.vec<mint>(D));
Polynomial<mint> f{as};
Polynomial<mint> g{cs};
for (int i : irange(1, D + 1)) {
g[i] *= (-1);
}
(f *= g).shrink(D);
const auto ans = bostanMori(f, g, K);
out.ln(ans.val());
return 0;
}