#if defined(__x86_64__) && defined(__GNUC__) && !defined(__clang__) #include #pragma GCC target("sse3,ssse3,sse4.1,sse4.2,popcnt,cx16,sahf,avx,avx2,bmi,bmi2,fma,f16c,lzcnt,movbe,pclmul,vpclmulqdq") #endif // https://yukicoder.me/problems/no/650 // // 提出は次を実装する。 // struct Solver { // Solver(int n, const vector> &edges); // 辺 i は edges[i]。行列は単位行列で始まる // void set(int e, const array &x); // 辺 e の行列を (x0 x1; x2 x3) にする // array product(int u, int v); // u (先祖) から v への道の辺の行列の積。根側が左 // }; // // 剰余は 1e9+7 で、入出力は素の整数で渡す。積は可換でないので、根側を左にして // 掛ける。根は 0。辺に載った値を扱うので、HLD は子側の頂点に持たせ、Link-Cut 木は // 辺を頂点にして間に挟む。 // // 構築を計測区間に入れてある。HLD は木を列に潰し、Link-Cut 木は 2(N-1) 回の link // を積むので、構築そのものが比較の対象になる。 // すべての問題のハーネスと提出が共通で使うもの。 // bits/stdc++.h は Apple clang に無いので、必要なものを名指しで include する。 #include #include #include #include #include #include #include using namespace std; // scanf の %lld と揃えたいので int64_t ではなく long long を使う。 using i64 = long long; using u64 = unsigned long long; using u32 = unsigned int; // 読み取れた項目数が想定と違ったら落とす。黙って進むと WA の原因が見えない。 inline void must_scan(int got, int want) { if (got != want) { fprintf(stderr, "input format error\n"); exit(1); } } // 空白で区切られた次のトークンを読む。長さが分からない入力 (文字列) 用。 // scanf と同じ stdin を読むので、混ぜて使ってよい。 inline string read_token() { string s; int c = getchar(); while (c != EOF && isspace(c)) c = getchar(); if (c == EOF) { fprintf(stderr, "input format error\n"); exit(1); } while (c != EOF && !isspace(c)) { s.push_back((char)c); c = getchar(); } return s; } // 整数を n 個読む。 inline vector read_ints(int n) { vector a(n); for (auto &x : a) must_scan(scanf("%lld", &x), 1); return a; } // 数の列を書く。既定は 1 行 1 個で、sep を渡すと区切りを変えられる。 // 末尾は sep に関わらず改行にする。空なら何も書かない。 // // 整形は計測区間の外に置きたいので、いったん文字列に溜めてから 1 回で書く。 template inline void print_all(const vector &v, char sep = '\n') { if (v.empty()) return; string out; out.reserve(v.size() * 12); for (size_t i = 0; i + 1 < v.size(); ++i) { out += to_string(v[i]); out += sep; } out += to_string(v.back()); out += '\n'; fwrite(out.data(), 1, out.size(), stdout); } // 計測区間の時間を stderr に 1 行の JSON で出す。実行側はこの接頭辞の行だけを拾う。 // メモリはハーネスでは測らない。Linux では pj が差し込む rss_preload が終了時の VmHWM を // 同じ形の行で報告し、macOS では pj が ru_maxrss を読む (pj/execute.py の peak_rss_kb)。 inline void report_metrics(long long algo_time_ns) { fprintf(stderr, "PJ_METRICS {\"algo_time_ns\":%lld}\n", algo_time_ns); } // ライブラリを使う提出が共有するもの。比べたいのは木の持ち方なので、載せる // モノイド (2x2 行列の積) は 1 つに揃える。 #include #include #include #include template constexpr inline Uint mod_inv(Uint a, Uint mod) { std::make_signed_t x= 1, y= 0, z= 0; for (Uint q= 0, b= mod, c= 0; b;) z= x, x= y, y= z - y * (q= a / b), c= a, a= b, b= c - b * q; return assert(a == 1), x < 0 ? mod - (-x) % mod : x % mod; } namespace math_internal { using namespace std; using u8= unsigned char; using u32= unsigned; using i64= long long; using u64= unsigned long long; using u128= __uint128_t; struct MP_Na { // mod < 2^32 u32 mod; constexpr MP_Na(): mod(0) {} constexpr MP_Na(u32 m): mod(m) {} constexpr inline u32 mul(u32 l, u32 r) const { return u64(l) * r % mod; } constexpr inline u32 set(u32 n) const { return n; } constexpr inline u32 get(u32 n) const { return n; } constexpr inline u32 norm(u32 n) const { return n; } constexpr inline u32 plus(u64 l, u32 r) const { return l+= r, l < mod ? l : l - mod; } constexpr inline u32 diff(u64 l, u32 r) const { return l-= r, l >> 63 ? l + mod : l; } }; template struct MP_Mo { // mod < 2^32, mod < 2^62 u_t mod; constexpr MP_Mo(): mod(0), iv(0), r2(0) {} constexpr MP_Mo(u_t m): mod(m), iv(inv(m)), r2(-du_t(mod) % mod) {} constexpr inline u_t mul(u_t l, u_t r) const { return reduce(du_t(l) * r); } constexpr inline u_t set(u_t n) const { return mul(n, r2); } constexpr inline u_t get(u_t n) const { return n= reduce(n), n >= mod ? n - mod : n; } constexpr inline u_t norm(u_t n) const { return n >= mod ? n - mod : n; } constexpr inline u_t plus(u_t l, u_t r) const { return l+= r, l < (mod << 1) ? l : l - (mod << 1); } constexpr inline u_t diff(u_t l, u_t r) const { return l-= r, l >> (B - 1) ? l + (mod << 1) : l; } private: u_t iv, r2; static constexpr u_t inv(u_t n, int e= 6, u_t x= 1) { return e ? inv(n, e - 1, x * (2 - x * n)) : x; } constexpr inline u_t reduce(const du_t &w) const { return u_t(w >> B) + mod - ((du_t(u_t(w) * iv) * mod) >> B); } }; using MP_Mo32= MP_Mo; using MP_Mo64= MP_Mo; struct MP_Br { // 2^20 < mod <= 2^41 u64 mod; constexpr MP_Br(): mod(0), x(0) {} constexpr MP_Br(u64 m): mod(m), x((u128(1) << 84) / m) {} constexpr inline u64 mul(u64 l, u64 r) const { return rem(u128(l) * r); } static constexpr inline u64 set(u64 n) { return n; } constexpr inline u64 get(u64 n) const { return n >= mod ? n - mod : n; } constexpr inline u64 norm(u64 n) const { return n >= mod ? n - mod : n; } constexpr inline u64 plus(u64 l, u64 r) const { return l+= r, l < (mod << 1) ? l : l - (mod << 1); } constexpr inline u64 diff(u64 l, u64 r) const { return l-= r, l >> 63 ? l + (mod << 1) : l; } private: u64 x; constexpr inline u128 quo(const u128 &n) const { return (n * x) >> 84; } constexpr inline u64 rem(const u128 &n) const { return n - quo(n) * mod; } }; template struct MP_D2B1 { // mod < 2^63, mod < 2^64 u64 mod; constexpr MP_D2B1(): mod(0), s(0), d(0), v(0) {} constexpr MP_D2B1(u64 m): mod(m), s(__builtin_clzll(m)), d(m << s), v(u128(-1) / d) {} constexpr inline u64 mul(u64 l, u64 r) const { return rem((u128(l) * r) << s) >> s; } constexpr inline u64 set(u64 n) const { return n; } constexpr inline u64 get(u64 n) const { return n; } constexpr inline u64 norm(u64 n) const { return n; } constexpr inline u64 plus(du_t l, u64 r) const { return l+= r, l < mod ? l : l - mod; } constexpr inline u64 diff(du_t l, u64 r) const { return l-= r, l >> B ? l + mod : l; } private: u8 s; u64 d, v; constexpr inline u64 rem(const u128 &u) const { u128 q= (u >> 64) * v + u; u64 r= u64(u) - (q >> 64) * d - d; if (r > u64(q)) r+= d; if (r >= d) r-= d; return r; } }; using MP_D2B1_1= MP_D2B1; using MP_D2B1_2= MP_D2B1; template constexpr u_t pow(u_t x, u64 k, const MP &md) { for (u_t ret= md.set(1);; x= md.mul(x, x)) if (k & 1 ? ret= md.mul(ret, x) : 0; !(k>>= 1)) return ret; } } #include namespace math_internal { struct m_b {}; struct s_b: m_b {}; } template constexpr bool is_modint_v= std::is_base_of_v; template constexpr bool is_staticmodint_v= std::is_base_of_v; namespace math_internal { template struct SB: s_b { protected: static constexpr MP md= MP(MOD); }; template struct MInt: public B { using Uint= U; static constexpr inline auto mod() { return B::md.mod; } constexpr MInt(): x(0) {} template && !is_same_v>> constexpr MInt(T v): x(B::md.set(v.val() % B::md.mod)) {} constexpr MInt(__int128_t n): x(B::md.set((n < 0 ? ((n= (-n) % B::md.mod) ? B::md.mod - n : n) : n % B::md.mod))) {} constexpr MInt operator-() const { return MInt() - *this; } #define FUNC(name, op) \ constexpr MInt name const { \ MInt ret; \ return ret.x= op, ret; \ } FUNC(operator+(const MInt & r), B::md.plus(x, r.x)) FUNC(operator-(const MInt & r), B::md.diff(x, r.x)) FUNC(operator*(const MInt & r), B::md.mul(x, r.x)) FUNC(pow(u64 k), math_internal::pow(x, k, B::md)) #undef FUNC constexpr MInt operator/(const MInt& r) const { return *this * r.inv(); } constexpr MInt& operator+=(const MInt& r) { return *this= *this + r; } constexpr MInt& operator-=(const MInt& r) { return *this= *this - r; } constexpr MInt& operator*=(const MInt& r) { return *this= *this * r; } constexpr MInt& operator/=(const MInt& r) { return *this= *this / r; } constexpr bool operator==(const MInt& r) const { return B::md.norm(x) == B::md.norm(r.x); } constexpr bool operator!=(const MInt& r) const { return !(*this == r); } constexpr bool operator<(const MInt& r) const { return B::md.norm(x) < B::md.norm(r.x); } constexpr inline MInt inv() const { return mod_inv(val(), B::md.mod); } constexpr inline Uint val() const { return B::md.get(x); } friend ostream& operator<<(ostream& os, const MInt& r) { return os << r.val(); } friend istream& operator>>(istream& is, MInt& r) { i64 v; return is >> v, r= MInt(v), is; } private: Uint x; }; template using MP_B= conditional_t < (MOD < (1 << 30)) & MOD, MP_Mo32, conditional_t < MOD < (1ull << 32), MP_Na, conditional_t<(MOD < (1ull << 62)) & MOD, MP_Mo64, conditional_t>>>>; template using ModInt= MInt < conditional_t, SB, MOD>>; } using math_internal::ModInt; #include #include #include #include #include #include #include namespace _la_internal { using namespace std; template struct Vector { valarray dat; Vector()= default; Vector(size_t n): dat(n) {} Vector(size_t n, const R &v): dat(v, n) {} Vector(const initializer_list &v): dat(v) {} R &operator[](int i) { return dat[i]; } const R &operator[](int i) const { return dat[i]; } bool operator==(const Vector &r) const { if (dat.size() != r.dat.size()) return false; for (int i= dat.size(); i--;) if (dat[i] != r.dat[i]) return false; return true; } bool operator!=(const Vector &r) const { return !(*this == r); } explicit operator bool() const { return dat.size(); } Vector operator-() const { return Vector(dat.size())-= *this; } Vector &operator+=(const Vector &r) { return dat+= r.dat, *this; } Vector &operator-=(const Vector &r) { return dat-= r.dat, *this; } Vector &operator*=(const R &r) { return dat*= r, *this; } Vector operator+(const Vector &r) const { return Vector(*this)+= r; } Vector operator-(const Vector &r) const { return Vector(*this)-= r; } Vector operator*(const R &r) const { return Vector(*this)*= r; } size_t size() const { return dat.size(); } friend R dot(const Vector &a, const Vector &b) { return assert(a.size() == b.size()), (a.dat * b.dat).sum(); } }; using u128= __uint128_t; using u64= uint64_t; using u8= uint8_t; class Ref { u128 *ref; u8 i; public: Ref(u128 *ref, u8 i): ref(ref), i(i) {} Ref &operator=(const Ref &r) { return *this= bool(r); } Ref &operator=(bool b) { return *ref&= ~(u128(1) << i), *ref|= u128(b) << i, *this; } Ref &operator|=(bool b) { return *ref|= u128(b) << i, *this; } Ref &operator&=(bool b) { return *ref&= ~(u128(!b) << i), *this; } Ref &operator^=(bool b) { return *ref^= u128(b) << i, *this; } operator bool() const { return (*ref >> i) & 1; } }; template <> class Vector { size_t n; public: valarray dat; Vector(): n(0) {} Vector(size_t n): n(n), dat((n + 127) >> 7) {} Vector(size_t n, bool b): n(n), dat(-u128(b), (n + 127) >> 7) { if (int k= n & 127; k) dat[dat.size() - 1]&= (u128(1) << k) - 1; } Vector(const initializer_list &v): n(v.size()), dat((n + 127) >> 7) { int i= 0; for (bool b: v) dat[i >> 7]|= u128(b) << (i & 127), ++i; } Ref operator[](int i) { return {begin(dat) + (i >> 7), u8(i & 127)}; } bool operator[](int i) const { return (dat[i >> 7] >> (i & 127)) & 1; } bool operator==(const Vector &r) const { if (dat.size() != r.dat.size()) return false; for (int i= dat.size(); i--;) if (dat[i] != r.dat[i]) return false; return true; } bool operator!=(const Vector &r) const { return !(*this == r); } explicit operator bool() const { return n; } Vector operator-() const { return Vector(*this); } Vector &operator+=(const Vector &r) { return dat^= r.dat, *this; } Vector &operator-=(const Vector &r) { return dat^= r.dat, *this; } Vector &operator*=(bool b) { return dat*= b, *this; } Vector operator+(const Vector &r) const { return Vector(*this)+= r; } Vector operator-(const Vector &r) const { return Vector(*this)-= r; } Vector operator*(bool b) const { return Vector(*this)*= b; } size_t size() const { return n; } friend bool dot(const Vector &a, const Vector &b) { assert(a.size() == b.size()); u128 v= 0; for (int i= a.dat.size(); i--;) v^= a.dat[i] & b.dat[i]; return __builtin_parityll(v >> 64) ^ __builtin_parityll(u64(v)); } }; template Vector operator*(const R &r, const Vector &v) { return v * r; } template ostream &operator<<(ostream &os, const Vector &v) { os << '['; for (int _= 0, __= v.size(); _ < __; ++_) os << (_ ? ", " : "") << v[_]; return os << ']'; } } using _la_internal::Vector; namespace _la_internal { template struct Mat { Mat(): W(0) {} Mat(size_t h, size_t w): W(w), dat(h * w) {} Mat(size_t h, size_t w, R v): W(w), dat(v, h * w) {} Mat(initializer_list> v): W(v.size() ? v.begin()->size() : 0), dat(v.size() * W) { auto it= begin(dat); for(const auto& r: v) { assert(r.size() == W); for(R x: r) *it++= x; } } size_t width() const { return W; } size_t height() const { return W ? dat.size() / W : 0; } auto operator[](int i) { return begin(dat) + i * W; } auto operator[](int i) const { return begin(dat) + i * W; } protected: size_t W; valarray dat; void add(const Mat& r) { assert(dat.size() == r.dat.size()), assert(W == r.W), dat+= r.dat; } D mul(const Mat& r) const { const size_t h= height(), w= r.W, l= W; assert(l == r.height()); D ret(h, w); auto a= begin(dat); auto c= begin(ret.dat); for(int i= h; i--; c+= w) { auto b= begin(r.dat); for(int k= l; k--; ++a) { auto d= c; auto v= *a; for(int j= w; j--; ++b, ++d) *d+= v * *b; } } return ret; } Vector mul(const Vector& r) const { assert(W == r.size()); const size_t h= height(); Vector ret(h); auto a= begin(dat); for(size_t i= 0; i < h; ++i) for(size_t k= 0; k < W; ++k, ++a) ret[i]+= *a * r[k]; return ret; } }; template struct Mat { struct Array { u128* bg; Array(u128* it): bg(it) {} Ref operator[](int i) { return Ref{bg + (i >> 7), u8(i & 127)}; } bool operator[](int i) const { return (bg[i >> 7] >> (i & 127)) & 1; } }; struct ConstArray { const u128* bg; ConstArray(const u128* it): bg(it) {} bool operator[](int i) const { return (bg[i >> 7] >> (i & 127)) & 1; } }; Mat(): H(0), W(0), m(0) {} Mat(size_t h, size_t w): H(h), W(w), m((w + 127) >> 7), dat(h * m) {} Mat(size_t h, size_t w, bool b): H(h), W(w), m((w + 127) >> 7), dat(-u128(b), h * m) { if(size_t i= h, k= w & 127; k) for(u128 s= (u128(1) << k) - 1; i--;) dat[i * m]&= s; } Mat(const initializer_list>& v): H(v.size()), W(H ? v.begin()->size() : 0), m((W + 127) >> 7), dat(H * m) { auto it= begin(dat); for(const auto& r: v) { assert(r.size() == W); int i= 0; for(bool b: r) it[i >> 7]|= u128(b) << (i & 127), ++i; it+= m; } } size_t width() const { return W; } size_t height() const { return H; } Array operator[](int i) { return {begin(dat) + i * m}; } ConstArray operator[](int i) const { return {begin(dat) + i * m}; } ConstArray get(int i) const { return {begin(dat) + i * m}; } protected: size_t H, W, m; valarray dat; void add(const Mat& r) { assert(H == r.H), assert(W == r.W), dat^= r.dat; } D mul(const Mat& r) const { assert(W == r.H); D ret(H, r.W); valarray tmp(r.m << 8); auto y= begin(r.dat); for(size_t l= 0; l < W; l+= 8) { auto t= begin(tmp) + r.m; for(int i= 0, n= min(8, W - l); i < n; ++i, y+= r.m) { auto u= begin(tmp); for(int s= 1 << i; s--;) { auto z= y; for(int j= r.m; j--; ++u, ++t, ++z) *t= *u ^ *z; } } auto a= begin(dat) + (l >> 7); auto c= begin(ret.dat); for(int i= H; i--; a+= m) { auto u= begin(tmp) + ((*a >> (l & 127)) & 255) * r.m; for(int j= r.m; j--; ++c, ++u) *c^= *u; } } return ret; } Vector mul(const Vector& r) const { assert(W == r.size()); Vector ret(H); auto a= begin(dat); for(size_t i= 0; i < H; ++i) { u128 v= 0; for(size_t j= 0; j < m; ++j, ++a) v^= *a & r.dat[j]; ret[i]= __builtin_parityll(v >> 64) ^ __builtin_parityll(u64(v)); } return ret; } }; template struct Matrix: public Mat> { using Mat>::Mat; explicit operator bool() const { return this->W; } static Matrix identity(int n) { Matrix ret(n, n); for(; n--;) ret[n][n]= R(true); return ret; } Matrix submatrix(const vector& rows, const vector& cols) const { Matrix ret(rows.size(), cols.size()); for(int i= rows.size(); i--;) for(int j= cols.size(); j--;) ret[i][j]= (*this)[rows[i]][cols[j]]; return ret; } Matrix submatrix_rm(vector rows, vector cols) const { sort(begin(rows), end(rows)), sort(begin(cols), end(cols)), rows.erase(unique(begin(rows), end(rows)), end(rows)), cols.erase(unique(begin(cols), end(cols)), end(cols)); const int H= this->height(), W= this->width(), n= rows.size(), m= cols.size(); vector rs(H - n), cs(W - m); for(int i= 0, j= 0, k= 0; i < H; ++i) if(j < n && rows[j] == i) ++j; else rs[k++]= i; for(int i= 0, j= 0, k= 0; i < W; ++i) if(j < m && cols[j] == i) ++j; else cs[k++]= i; return submatrix(rs, cs); } bool operator==(const Matrix& r) const { if(this->width() != r.width() || this->height() != r.height()) return false; for(int i= this->dat.size(); i--;) if(this->dat[i] != r.dat[i]) return false; return true; } bool operator!=(const Matrix& r) const { return !(*this == r); } Matrix& operator*=(const Matrix& r) { return *this= this->mul(r); } Matrix operator*(const Matrix& r) const { return this->mul(r); } Matrix& operator*=(R r) { return this->dat*= r, *this; } template Matrix operator*(T r) const { static_assert(is_convertible_v); return Matrix(*this)*= r; } Matrix& operator+=(const Matrix& r) { return this->add(r), *this; } Matrix operator+(const Matrix& r) const { return Matrix(*this)+= r; } Vector operator*(const Vector& r) const { return this->mul(r); } Vector operator()(const Vector& r) const { return this->mul(r); } Matrix pow(uint64_t k) const { size_t W= this->width(); assert(W == this->height()); for(Matrix ret= identity(W), b= *this;; b*= b) if(k& 1 ? ret*= b, !(k>>= 1) : !(k>>= 1)) return ret; } }; template Matrix operator*(const T& r, const Matrix& m) { return m * r; } template ostream& operator<<(ostream& os, const Matrix& m) { os << "\n["; for(int i= 0, h= m.height(); i < h; os << ']', ++i) { if(i) os << "\n "; os << '['; for(int j= 0, w= m.width(); j < w; ++j) os << (j ? ", " : "") << m[i][j]; } return os << ']'; } template static bool is_zero(K x) { if constexpr(is_floating_point_v) return abs(x) < 1e-8; else return x == K(); } } using _la_internal::Matrix; using Mint = ModInt<1000000007>; using Mat = Matrix; inline const Mat IDENTITY = Mat::identity(2); inline Mat to_mat(const array &x) { Mat m(2, 2); m[0][0] = Mint(x[0]), m[0][1] = Mint(x[1]); m[1][0] = Mint(x[2]), m[1][1] = Mint(x[3]); return m; } inline array from_mat(const Mat &m) { return {(i64)m[0][0].val(), (i64)m[0][1].val(), (i64)m[1][0].val(), (i64)m[1][1].val()}; } // 根側を左にして掛ける。可換でないので、Link-Cut 木の側は反転したときの積も持つ。 struct MatProd { using T = Mat; static T ti() { return IDENTITY; } static T op(const T &l, const T &r) { return l * r; } }; #include #include #include #include #define _LR(name, IT, CT) \ template struct name { \ using Iterator= typename std::vector::IT; \ Iterator bg, ed; \ Iterator begin() const { return bg; } \ Iterator end() const { return ed; } \ size_t size() const { return std::distance(bg, ed); } \ CT &operator[](int i) const { return bg[i]; } \ } _LR(ListRange, iterator, T); _LR(ConstListRange, const_iterator, const T); #undef _LR template struct CSRArray { std::vector dat; std::vector p; size_t size() const { return p.size() - 1; } ListRange operator[](int i) { return {dat.begin() + p[i], dat.begin() + p[i + 1]}; } ConstListRange operator[](int i) const { return {dat.cbegin() + p[i], dat.cbegin() + p[i + 1]}; } }; template