library

This documentation is automatically generated by online-judge-tools/verification-helper

View the Project on GitHub maspypy/library

:warning: graph/ds/tree_wavelet_matrix.hpp

Depends on

Code

#include "ds/wavelet_matrix/wavelet_matrix.hpp"
#include "graph/tree.hpp"

// https://atcoder.jp/contests/pakencamp-2022-day1/tasks/pakencamp_2022_day1_j
// https://atcoder.jp/contests/utpc2011/tasks/utpc2011_12
template <typename TREE, bool edge, typename T, bool COMPRESS,
          typename Monoid = Monoid_Add<T>>
struct Tree_Wavelet_Matrix {
  TREE& tree;
  int N;
  using WM = Wavelet_Matrix<T, COMPRESS, Monoid_Add<T>>;
  using X = typename Monoid::value_type;
  WM wm;

  Tree_Wavelet_Matrix(TREE& tree, vc<T> A, vc<X> SUM_data = {}, int log = -1)
      : tree(tree), N(tree.N) {
    vc<X>& S = SUM_data;
    vc<T> A1;
    vc<X> S1;
    A1.resize(N);
    if (!S.empty()) S1.resize(N);
    if (!edge) {
      assert(len(A) == N && (len(S) == 0 || len(S) == N));
      FOR(v, N) A1[tree.LID[v]] = A[v];
      if (len(S) == N) { FOR(v, N) S1[tree.LID[v]] = S[v]; }
      wm.build(A1, S1, log);
    } else {
      assert(len(A) == N - 1 && (len(S) == 0 || len(S) == N - 1));
      if (!S.empty()) {
        FOR(e, N - 1) { S1[tree.LID[tree.e_to_v(e)]] = S[e]; }
      }
      FOR(e, N - 1) { A1[tree.LID[tree.e_to_v(e)]] = A[e]; }
      wm.build(A1, S1, log);
    }
  }

  // xor した結果で [a, b) に収まるものを数える
  int count_path(int s, int t, T a, T b, T xor_val = 0) {
    return wm.count(get_segments(s, t), a, b, xor_val);
  }

  // xor した結果で [a, b) に収まるものを数える
  int count_subtree(int u, T a, T b, T xor_val = 0) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.count(l + edge, r, a, b, xor_val);
  }

  // xor した結果で、[L, R) の中で k>=0 番目と prefix sum
  pair<T, X> kth_value_and_sum_path(int s, int t, int k, T xor_val = 0) {
    return wm.kth_value_and_sum(get_segments(s, t), k, xor_val);
  }

  // xor した結果で、[L, R) の中で k>=0 番目と prefix sum
  pair<T, X> kth_value_and_sum_subtree(int u, int k, T xor_val = 0) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.kth_value_and_sum(l + edge, r, k, xor_val);
  }

  // xor した結果で、[L, R) の中で k>=0 番目
  T kth_path(int s, int t, int k, T xor_val = 0) {
    return wm.kth(get_segments(s, t), k, xor_val);
  }

  // xor した結果で、[L, R) の中で k>=0 番目
  T kth_subtree(int u, int k, T xor_val = 0) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.kth(l + edge, r, k, xor_val);
  }

  // xor した結果で、[L, R) の中で中央値。
  // LOWER = true:下側中央値、false:上側中央値
  T median_path(bool UPPER, int s, int t, T xor_val = 0) {
    return wm.median(UPPER, get_segments(s, t), xor_val);
  }

  T median_subtree(bool UPPER, int u, T xor_val = 0) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.median(UPPER, l + edge, r, xor_val);
  }

  // xor した結果で [k1, k2) 番目であるところの SUM_data の和
  X sum_path(int s, int t, int k1, int k2, T xor_val = 0) {
    return wm.sum(get_segments(s, t), k1, k2, xor_val);
  }

  // xor した結果で [k1, k2) 番目であるところの SUM_data の和
  X sum_subtree(int u, int k1, int k2, T xor_val = 0) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.sum(l + edge, r, k1, k2, xor_val);
  }

  X sum_all_path(int s, int t) { return wm.sum_all(get_segments(s, t)); }

  X sum_all_subtree(int u) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.sum_all(l + edge, r);
  }

private:
  vc<pair<int, int>> get_segments(int s, int t) {
    vc<pair<int, int>> segments = tree.get_path_decomposition(s, t, edge);
    for (auto&& [a, b]: segments) {
      if (a >= b) swap(a, b);
      ++b;
    }
    return segments;
  }
};
#line 1 "graph/ds/tree_wavelet_matrix.hpp"

#line 1 "ds/bit_vector.hpp"
struct Bit_Vector {
  int n;
  bool prepared = 0;
  vc<pair<u64, u32>> dat;
  Bit_Vector(int n = 0) : n(n) { dat.assign((n + 127) >> 6, {0, 0}); }
  void set(int i) {
    assert(!prepared && (0 <= i && i < n));
    dat[i >> 6].fi |= u64(1) << (i & 63);
  }
  void reset() {
    fill(all(dat), pair<u64, u32>{0, 0});
    prepared = 0;
  }
  void build() {
    prepared = 1;
    FOR(i, len(dat) - 1) dat[i + 1].se = dat[i].se + popcnt(dat[i].fi);
  }
  bool operator[](int i) const { return dat[i >> 6].fi >> (i & 63) & 1; }
  // [0, k) 内の 1 の個数
  int count_prefix(int k, bool f = true) const {
    assert(prepared);
    auto [a, b] = dat[k >> 6];
    int ret = b + popcnt(a & ((u64(1) << (k & 63)) - 1));
    return (f ? ret : k - ret);
  }
  int count(int L, int R, bool f = true) const {
    return count_prefix(R, f) - count_prefix(L, f);
  }
  string to_string() const {
    string ans;
    FOR(i, n) ans += '0' + (dat[i / 64].fi >> (i % 64) & 1);
    return ans;
  }
};
#line 1 "alg/monoid/dummy.hpp"
struct Monoid_Dummy {
  using value_type = char;
  static constexpr bool commute = true;
  static value_type op(value_type, value_type) { return 0; }
  static value_type unit() { return 0; }
};
#line 2 "ds/dummy_data_structure.hpp"

struct Dummy_Data_Structure {
  using MX = Monoid_Dummy;
  using T = typename MX::value_type;
  void build(const vc<T>& A) {}
};
#line 3 "ds/wavelet_matrix/wavelet_matrix.hpp"

template <typename Y, typename SEGTREE>
struct Uncompressed_Wavelet_Matrix {
  using Mono = typename SEGTREE::MX;
  using T = typename Mono::value_type;
  static_assert(Mono::commute);
  static_assert(is_same_v<Y, int> || is_same_v<Y, ll>);
  int n = 0, log = 0;
  vc<int> mid;
  vc<Bit_Vector> bv;
  vc<SEGTREE> seg;
  Y limit;

  Uncompressed_Wavelet_Matrix() = default;

  // f(i) = {A[i], dat[i]}

  template <typename F>
  Uncompressed_Wavelet_Matrix(int n, F f, int log = -1) {
    build(n, f, log);
  }
  Uncompressed_Wavelet_Matrix(const vc<Y>& A, int log = -1) {
    static_assert(is_same_v<SEGTREE, Dummy_Data_Structure>);
    build(
        len(A), [&](int i) -> pair<Y, T> { return {A[i], Mono::unit()}; }, log);
  }

  template <typename F>
  void build(int n, F f, int log = -1) {
    this->n = n;
    vc<Y> A(n);
    vc<T> S(n);
    FOR(i, n) tie(A[i], S[i]) = f(i);
    if (log == -1) {
      log = (n == 0 ? 0 : topbit(MAX(A)) + 1);
    } else {
      for (auto& x : A) assert(0 <= x && topbit(x) < log);
    }
    this->log = log;
    limit = Y(1) << log;
    if constexpr (is_same_v<Y, int>) assert(0 <= log && log <= 30);
    if constexpr (is_same_v<Y, ll>) assert(0 <= log && log <= 62);
    mid.resize(log), bv.assign(log, Bit_Vector(n));
    vc<Y> A0(n), A1(n);
    vc<T> S0(n), S1(n);
    seg.resize(log + 1);
    seg[log].build(S);
    for (int d = log - 1; d >= 0; --d) {
      int p0 = 0, p1 = 0;
      for (int i = 0; i < n; ++i) {
        if (A[i] >> d & 1) {
          bv[d].set(i), A1[p1] = A[i], S1[p1] = S[i], p1++;
        } else {
          A0[p0] = A[i], S0[p0] = S[i], p0++;
        }
      }
      swap(A, A0), swap(S, S0);
      move(A1.begin(), A1.begin() + p1, A.begin() + p0);
      move(S1.begin(), S1.begin() + p1, S.begin() + p0);
      mid[d] = p0, bv[d].build(), seg[d].build(S);
    }
  }

  tuple<int, int, int, int> get_subtree(int d, int L, int R) const {
    assert(1 <= d && d <= log);
    int a = bv[d - 1].count_prefix(L), b = bv[d - 1].count_prefix(R);
    return {L - a, R - b, mid[d - 1] + a, mid[d - 1] + b};
  }

  template <typename F>
  void work_point(F f, int i) {
    assert(0 <= i && i < n);
    f(log, i);
    FOR_R(d, log) {
      int a = bv[d].count_prefix(i);
      if (bv[d][i]) {
        i = mid[d] + a;
      } else {
        i = i - a;
      }
      f(d, i);
    }
  }

  template <typename F>
  void work_prefix(F f, int L, int R, Y y) const {
    chmin(y, limit);
    if (y == 0) return;
    if (y == limit) {
      f(log, L, R);
      return;
    }
    FOR_R(d, log) {
      auto [L0, R0, L1, R1] = get_subtree(d + 1, L, R);
      if (y >> d & 1) {
        f(d, L0, R0);
        L = L1, R = R1;
      } else {
        L = L0, R = R0;
      }
    }
  }

  template <typename F>
  void work_range(F f, int L, int R, Y y1, Y y2) const {
    chmin(y2, limit);
    if (y1 >= y2) return;
    assert(0 <= y1 && y1 <= y2 && y2 <= limit);
    if (y1 == 0) return work_prefix(f, L, R, y2);
    auto dfs = [&](auto& dfs, int d, int L, int R, Y y1, Y y2) -> void {
      if (y1 == y2) return;
      if (y1 == 0 && y2 == Y(1) << d) {
        f(d, L, R);
        return;
      }
      assert(d > 0);
      auto [L0, R0, L1, R1] = get_subtree(d, L, R);
      Y m = (Y(1) << (d - 1));

      if (y2 <= m) {
        dfs(dfs, d - 1, L0, R0, y1, y2);
      } else if (y1 >= m) {
        dfs(dfs, d - 1, L1, R1, y1 - m, y2 - m);
      } else {
        dfs(dfs, d - 1, L0, R0, y1, m);
        dfs(dfs, d - 1, L1, R1, 0, y2 - m);
      }
    };
    dfs(dfs, log, L, R, y1, y2);
  }

  // [L,R) x [0,y)

  int prefix_count(int L, int R, Y y) const {
    int cnt = 0;
    work_prefix([&](int d, int a, int b) { cnt += b - a; }, L, R, y);
    return cnt;
  }

  // [L,R) x [y1,y2)

  int count(int L, int R, Y y1, Y y2) const {
    return prefix_count(L, R, y2) - prefix_count(L, R, y1);
  }

  // [L,R) x [0,y)

  T prefix_prod(int L, int R, Y y) const {
    T ans = Mono::unit();
    work_prefix(
        [&](int d, int a, int b) { ans = Mono::op(ans, seg[d].prod(a, b)); }, L,
        R, y);
    return ans;
  }
  // [L,R) x [y1,y2)

  T prod(int L, int R, Y y1, Y y2) const {
    T ans = Mono::unit();
    work_range(
        [&](int d, int a, int b) { ans = Mono::op(ans, seg[d].prod(a, b)); }, L,
        R, y1, y2);
    return ans;
  }
  T prod_all(int L, int R) const { return seg[log].prod(L, R); }

  // [L,R) x [0,y)

  pair<int, T> prefix_count_and_prod(int L, int R, Y y) const {
    pair<int, T> ans = {0, Mono::unit()};
    work_prefix(
        [&](int d, int a, int b) {
          ans.fi += b - a;
          ans.se = Mono::op(ans.se, seg[d].prod(a, b));
        },
        L, R, y);
    return ans;
  }
  // [L,R) x [y1,y2)

  pair<int, T> count_and_prod(int L, int R, Y y1, Y y2) const {
    pair<int, T> ans = {0, Mono::unit()};
    work_range(
        [&](int d, int a, int b) {
          ans.fi += b - a;
          ans.se = Mono::op(ans.se, seg[d].prod(a, b));
        },
        L, R, y1, y2);
    return ans;
  }

  Y kth(int L, int R, int k) const {
    assert(0 <= k && k < R - L);
    Y ans = 0;
    for (int d = log - 1; d >= 0; --d) {
      auto [L0, R0, L1, R1] = get_subtree(d + 1, L, R);
      if (k < R0 - L0) {
        L = L0, R = R0;
      } else {
        ans |= Y(1) << d;
        k -= R0 - L0, L = L1, R = R1;
      }
    }
    return ans;
  }

  template <bool upper>
  Y median(int L, int R) const {
    assert(0 <= L && L < R && R <= n);
    int k = (upper ? (R - L) / 2 : (R - L - 1) / 2);
    return kth(L, R, k);
  }

  void set(int i, T t) {
    assert(0 <= i && i < n);
    work_point([&](int d, int i) { seg[d].set(i, t); }, i);
  }
  void multiply(int i, T t) {
    assert(0 <= i && i < n);
    work_point([&](int d, int i) { seg[d].multiply(i, t); }, i);
  }
  void add(int i, T t) {
    assert(0 <= i && i < n);
    work_point([&](int d, int i) { seg[d].add(i, t); }, i);
  }

  // [L,R) x [0,y) での check(y, cnt, prod) が true となる最大の (Y,cnt,prod)

  template <typename F>
  tuple<Y, int, T> max_right(F check, int L, int R) const {
    assert(limit < infty<Y>);
    int cnt = 0;
    Y y = 0;
    T t = Mono::unit();
    T t_all = seg[log].prod(L, R);
    assert(check(0, 0, Mono::unit()));
    if (check(limit, R - L, t_all)) {
      y = binary_search([&](Y y) -> bool { return check(y, R - L, t_all); },
                        limit, infty<Y> + 1);
      return {y, R - L, t_all};
    }
    for (int d = log - 1; d >= 0; --d) {
      auto [L0, R0, L1, R1] = get_subtree(d + 1, L, R);
      Y y1 = y | Y(1) << d;
      int cnt1 = cnt + R0 - L0;
      T t1 = Mono::op(t, seg[d].prod(L0, R0));
      if (check(y1, cnt1, t1)) {
        y = y1, cnt = cnt1, t = t1, L = L1, R = R1;
      } else {
        L = L0, R = R0;
      }
    }
    return {y, cnt, t};
  }

  // [L,R) x [0,y) での check(y, cnt, prod) が true となる最大の (Y,cnt,prod)

  template <typename F>
  tuple<Y, int, T> max_right_many(F check, vc<pair<int, int>> LR) const {
    assert(limit < infty<Y>);
    int cnt = 0;
    Y y = 0;
    T t = Mono::unit();
    T t_all = Mono::unit();
    int cnt_all = 0;
    for (auto& [l, r] : LR)
      t_all = Mono::op(t_all, prod_all(l, r)), cnt_all += r - l;
    assert(check(0, 0, Mono::unit()));
    if (check(limit, cnt_all, t_all)) {
      y = binary_search([&](Y y) -> bool { return check(y, cnt_all, t_all); },
                        limit, infty<Y> + 1);
      return {y, cnt_all, t_all};
    }
    for (int d = log - 1; d >= 0; --d) {
      Y y1 = Y(1) << d;
      T t1 = t;
      int cnt1 = 0;
      for (auto& [L, R] : LR) {
        auto [L0, R0, L1, R1] = get_subtree(d + 1, L, R);
        cnt1 += R0 - L0;
        t1 = Mono::op(t1, seg[d].prod(L0, R0));
      }
      if (check(y1, cnt1, t1)) {
        y = y1, cnt = cnt1, t = t1;
        for (auto& [L, R] : LR) {
          auto [L0, R0, L1, R1] = get_subtree(d + 1, L, R);
          L = L1, R = R1;
        }
      } else {
        for (auto& [L, R] : LR) {
          auto [L0, R0, L1, R1] = get_subtree(d + 1, L, R);
          L = L0, R = R0;
        }
      }
    }
    return {y, cnt, t};
  }

  // [L,R) x [y, inf) での check(y, cnt, prod) が true となる最小の (y,cnt,prod)

  // cnt==0 だと true であることは仮定する

  // https://qoj.ac/contest/1047/problem/5094

  template <typename F>
  tuple<Y, int, T> min_left_many(F check, vc<pair<int, int>> LR) const {
    assert(check(limit, 0, Mono::unit()));
    int cnt = 0;
    Y y = limit;
    T t = Mono::unit();
    T t_all = Mono::unit();
    int cnt_all = 0;
    for (auto& [l, r] : LR)
      t_all = Mono::op(t_all, prod_all(l, r)), cnt_all += r - l;
    if (check(0, cnt_all, t_all)) {
      return {0, cnt_all, t_all};
    }
    for (int d = log - 1; d >= 0; --d) {
      Y y1 = y - (Y(1) << d);
      T t1 = t;
      int cnt1 = cnt;
      for (auto& [L, R] : LR) {
        auto [L0, R0, L1, R1] = get_subtree(d + 1, L, R);
        cnt1 += R1 - L1;
        t1 = Mono::op(t1, seg[d].prod(L1, R1));
      }
      if (check(y1, cnt1, t1)) {
        y = y1, cnt = cnt1, t = t1;
        SHOW(y);
        for (auto& [L, R] : LR) {
          auto [L0, R0, L1, R1] = get_subtree(d + 1, L, R);
          L = L0, R = R0;
        }
      } else {
        for (auto& [L, R] : LR) {
          auto [L0, R0, L1, R1] = get_subtree(d + 1, L, R);
          L = L1, R = R1;
        }
      }
    }
    SHOW(y, cnt, t);
    return {y, cnt, t};
  }
};

template <typename Y, typename SEGTREE>
struct Compressed_Wavelet_Matrix {
  using Mono = typename SEGTREE::MX;
  using T = typename Mono::value_type;

  int n = 0;
  vc<Y> key;
  Uncompressed_Wavelet_Matrix<int, SEGTREE> wm;

  Compressed_Wavelet_Matrix() = default;

  // f(i) = {A[i], dat[i]}

  template <typename F>
  Compressed_Wavelet_Matrix(int n, F f) {
    build(n, f);
  }

  Compressed_Wavelet_Matrix(const vc<Y>& A) {
    static_assert(is_same_v<SEGTREE, Dummy_Data_Structure>);
    build(A);
  }

  template <typename F>
  void build(int n, F f) {
    this->n = n;
    vc<Y> A(n);
    vc<T> S(n);
    FOR(i, n) tie(A[i], S[i]) = f(i);

    key = A;
    UNIQUE(key);

    wm.build(n, [&](int i) -> pair<int, T> {
      int k = LB(key, A[i]);
      return {k, S[i]};
    });
  }

  void build(const vc<Y>& A) {
    static_assert(is_same_v<SEGTREE, Dummy_Data_Structure>);
    n = len(A);
    key = A;
    UNIQUE(key);

    wm.build(n, [&](int i) -> pair<int, T> {
      int k = LB(key, A[i]);
      return {k, Mono::unit()};
    });
  }

  Y kth(int L, int R, int k) const { return key[wm.kth(L, R, k)]; }

  template <bool upper>
  Y median(int L, int R) const {
    return key[wm.template median<upper>(L, R)];
  }

  // [L,R) x [-inf,y)

  int prefix_count(int L, int R, Y y) const {
    return wm.prefix_count(L, R, LB(key, y));
  }

  // [L,R) x [y1,y2)

  int count(int L, int R, Y y1, Y y2) const {
    return wm.count(L, R, LB(key, y1), LB(key, y2));
  }

  // [L,R) x [-inf,y)

  T prefix_prod(int L, int R, Y y) const {
    return wm.prefix_prod(L, R, LB(key, y));
  }

  // [L,R) x [y1,y2)

  T prod(int L, int R, Y y1, Y y2) const {
    return wm.prod(L, R, LB(key, y1), LB(key, y2));
  }

  T prod_all(int L, int R) const { return wm.prod_all(L, R); }

  // [L,R) x [-inf,y)

  pair<int, T> prefix_count_and_prod(int L, int R, Y y) const {
    return wm.prefix_count_and_prod(L, R, LB(key, y));
  }

  // [L,R) x [y1,y2)

  pair<int, T> count_and_prod(int L, int R, Y y1, Y y2) const {
    return wm.count_and_prod(L, R, LB(key, y1), LB(key, y2));
  }

  void set(int i, T t) { wm.set(i, t); }

  void multiply(int i, T t) { wm.multiply(i, t); }

  void add(int i, T t) { wm.add(i, t); }
};

template <typename Y, bool compress, typename SEGTREE = Dummy_Data_Structure>
using Wavelet_Matrix =
    conditional_t<compress, Compressed_Wavelet_Matrix<Y, SEGTREE>,
                  Uncompressed_Wavelet_Matrix<Y, SEGTREE>>;
#line 2 "graph/tree.hpp"

#line 2 "ds/hashmap.hpp"

// u64 -> Val

template <typename Val>
struct HashMap {
  // n は入れたいものの個数で ok

  HashMap(u32 n = 0) { build(n); }
  void build(u32 n) {
    u32 k = 8;
    while (k < n * 2) k *= 2;
    cap = k / 2, mask = k - 1;
    key.resize(k), val.resize(k), used.assign(k, 0);
  }

  // size を保ったまま. size=0 にするときは build すること.

  void clear() {
    used.assign(len(used), 0);
    cap = (mask + 1) / 2;
  }
  int size() { return len(used) / 2 - cap; }

  int index(const u64& k) {
    int i = 0;
    for (i = hash(k); used[i] && key[i] != k; i = (i + 1) & mask) {}
    return i;
  }

  Val& operator[](const u64& k) {
    if (cap == 0) extend();
    int i = index(k);
    if (!used[i]) { used[i] = 1, key[i] = k, val[i] = Val{}, --cap; }
    return val[i];
  }

  Val get(const u64& k, Val default_value) {
    int i = index(k);
    return (used[i] ? val[i] : default_value);
  }

  bool count(const u64& k) {
    int i = index(k);
    return used[i] && key[i] == k;
  }

  // f(key, val)

  template <typename F>
  void enumerate_all(F f) {
    FOR(i, len(used)) if (used[i]) f(key[i], val[i]);
  }

private:
  u32 cap, mask;
  vc<u64> key;
  vc<Val> val;
  vc<bool> used;

  u64 hash(u64 x) {
    static const u64 FIXED_RANDOM = std::chrono::steady_clock::now().time_since_epoch().count();
    x += FIXED_RANDOM;
    x = (x ^ (x >> 30)) * 0xbf58476d1ce4e5b9;
    x = (x ^ (x >> 27)) * 0x94d049bb133111eb;
    return (x ^ (x >> 31)) & mask;
  }

  void extend() {
    vc<pair<u64, Val>> dat;
    dat.reserve(len(used) / 2 - cap);
    FOR(i, len(used)) {
      if (used[i]) dat.eb(key[i], val[i]);
    }
    build(2 * len(dat));
    for (auto& [a, b]: dat) (*this)[a] = b;
  }
};
#line 3 "graph/base.hpp"

template <typename T>
struct Edge {
  int frm, to;
  T cost;
  int id;
};

template <typename T = int, bool directed = false>
struct Graph {
  static constexpr bool is_directed = directed;
  int N, M;
  using cost_type = T;
  using edge_type = Edge<T>;
  vector<edge_type> edges;
  vector<int> indptr;
  vector<edge_type> csr_edges;
  vc<int> vc_deg, vc_indeg, vc_outdeg;
  bool prepared;

  class OutgoingEdges {
  public:
    OutgoingEdges(const Graph* G, int l, int r) : G(G), l(l), r(r) {}

    const edge_type* begin() const {
      if (l == r) { return 0; }
      return &G->csr_edges[l];
    }

    const edge_type* end() const {
      if (l == r) { return 0; }
      return &G->csr_edges[r];
    }

  private:
    const Graph* G;
    int l, r;
  };

  bool is_prepared() { return prepared; }

  Graph() : N(0), M(0), prepared(0) {}
  Graph(int N) : N(N), M(0), prepared(0) {}

  void build(int n) {
    N = n, M = 0;
    prepared = 0;
    edges.clear();
    indptr.clear();
    csr_edges.clear();
    vc_deg.clear();
    vc_indeg.clear();
    vc_outdeg.clear();
  }

  void add(int frm, int to, T cost = 1, int i = -1) {
    assert(!prepared);
    assert(0 <= frm && 0 <= to && to < N);
    if (i == -1) i = M;
    auto e = edge_type({frm, to, cost, i});
    edges.eb(e);
    ++M;
  }

#ifdef FASTIO
  // wt, off
  void read_tree(bool wt = false, int off = 1) { read_graph(N - 1, wt, off); }

  void read_graph(int M, bool wt = false, int off = 1) {
    for (int m = 0; m < M; ++m) {
      INT(a, b);
      a -= off, b -= off;
      if (!wt) {
        add(a, b);
      } else {
        T c;
        read(c);
        add(a, b, c);
      }
    }
    build();
  }
#endif

  void build() {
    assert(!prepared);
    prepared = true;
    indptr.assign(N + 1, 0);
    for (auto&& e: edges) {
      indptr[e.frm + 1]++;
      if (!directed) indptr[e.to + 1]++;
    }
    for (int v = 0; v < N; ++v) { indptr[v + 1] += indptr[v]; }
    auto counter = indptr;
    csr_edges.resize(indptr.back() + 1);
    for (auto&& e: edges) {
      csr_edges[counter[e.frm]++] = e;
      if (!directed) csr_edges[counter[e.to]++] = edge_type({e.to, e.frm, e.cost, e.id});
    }
  }

  OutgoingEdges operator[](int v) const {
    assert(prepared);
    return {this, indptr[v], indptr[v + 1]};
  }

  vc<int> deg_array() {
    if (vc_deg.empty()) calc_deg();
    return vc_deg;
  }

  pair<vc<int>, vc<int>> deg_array_inout() {
    if (vc_indeg.empty()) calc_deg_inout();
    return {vc_indeg, vc_outdeg};
  }

  int deg(int v) {
    if (vc_deg.empty()) calc_deg();
    return vc_deg[v];
  }

  int in_deg(int v) {
    if (vc_indeg.empty()) calc_deg_inout();
    return vc_indeg[v];
  }

  int out_deg(int v) {
    if (vc_outdeg.empty()) calc_deg_inout();
    return vc_outdeg[v];
  }

#ifdef FASTIO
  void debug() {
#ifdef LOCAL
    print("Graph");
    if (!prepared) {
      print("frm to cost id");
      for (auto&& e: edges) print(e.frm, e.to, e.cost, e.id);
    } else {
      print("indptr", indptr);
      print("frm to cost id");
      FOR(v, N) for (auto&& e: (*this)[v]) print(e.frm, e.to, e.cost, e.id);
    }
#endif
  }
#endif

  vc<int> new_idx;
  vc<bool> used_e;

  // G における頂点 V[i] が、新しいグラフで i になるようにする
  // {G, es}
  // sum(deg(v)) の計算量になっていて、
  // 新しいグラフの n+m より大きい可能性があるので注意
  Graph<T, directed> rearrange(vc<int> V, bool keep_eid = 0) {
    if (len(new_idx) != N) new_idx.assign(N, -1);
    int n = len(V);
    FOR(i, n) new_idx[V[i]] = i;
    Graph<T, directed> G(n);
    vc<int> history;
    FOR(i, n) {
      for (auto&& e: (*this)[V[i]]) {
        if (len(used_e) <= e.id) used_e.resize(e.id + 1);
        if (used_e[e.id]) continue;
        int a = e.frm, b = e.to;
        if (new_idx[a] != -1 && new_idx[b] != -1) {
          history.eb(e.id);
          used_e[e.id] = 1;
          int eid = (keep_eid ? e.id : -1);
          G.add(new_idx[a], new_idx[b], e.cost, eid);
        }
      }
    }
    FOR(i, n) new_idx[V[i]] = -1;
    for (auto&& eid: history) used_e[eid] = 0;
    G.build();
    return G;
  }

  Graph<T, true> to_directed_tree(int root = -1) {
    if (root == -1) root = 0;
    assert(!is_directed && prepared && M == N - 1);
    Graph<T, true> G1(N);
    vc<int> par(N, -1);
    auto dfs = [&](auto& dfs, int v) -> void {
      for (auto& e: (*this)[v]) {
        if (e.to == par[v]) continue;
        par[e.to] = v, dfs(dfs, e.to);
      }
    };
    dfs(dfs, root);
    for (auto& e: edges) {
      int a = e.frm, b = e.to;
      if (par[a] == b) swap(a, b);
      assert(par[b] == a);
      G1.add(a, b, e.cost);
    }
    G1.build();
    return G1;
  }

  HashMap<int> MP_FOR_EID;

  int get_eid(u64 a, u64 b) {
    if (len(MP_FOR_EID) == 0) {
      MP_FOR_EID.build(N - 1);
      for (auto& e: edges) {
        u64 a = e.frm, b = e.to;
        u64 k = to_eid_key(a, b);
        MP_FOR_EID[k] = e.id;
      }
    }
    return MP_FOR_EID.get(to_eid_key(a, b), -1);
  }

  u64 to_eid_key(u64 a, u64 b) {
    if (!directed && a > b) swap(a, b);
    return N * a + b;
  }

private:
  void calc_deg() {
    assert(vc_deg.empty());
    vc_deg.resize(N);
    for (auto&& e: edges) vc_deg[e.frm]++, vc_deg[e.to]++;
  }

  void calc_deg_inout() {
    assert(vc_indeg.empty());
    vc_indeg.resize(N);
    vc_outdeg.resize(N);
    for (auto&& e: edges) { vc_indeg[e.to]++, vc_outdeg[e.frm]++; }
  }
};
#line 4 "graph/tree.hpp"

// HLD euler tour をとっていろいろ。

template <typename GT>
struct Tree {
  using Graph_type = GT;
  GT &G;
  using WT = typename GT::cost_type;
  int N;
  vector<int> LID, RID, head, V, parent, VtoE;
  vc<int> depth;
  vc<WT> depth_weighted;

  Tree(GT &G, int r = 0, bool hld = 1) : G(G) { build(r, hld); }

  void build(int r = 0, bool hld = 1) {
    if (r == -1) return; // build を遅延したいとき

    N = G.N;
    LID.assign(N, -1), RID.assign(N, -1), head.assign(N, r);
    V.assign(N, -1), parent.assign(N, -1), VtoE.assign(N, -1);
    depth.assign(N, -1), depth_weighted.assign(N, 0);
    assert(G.is_prepared());
    int t1 = 0;
    dfs_sz(r, -1, hld);
    dfs_hld(r, t1);
  }

  void dfs_sz(int v, int p, bool hld) {
    auto &sz = RID;
    parent[v] = p;
    depth[v] = (p == -1 ? 0 : depth[p] + 1);
    sz[v] = 1;
    int l = G.indptr[v], r = G.indptr[v + 1];
    auto &csr = G.csr_edges;
    // 使う辺があれば先頭にする

    for (int i = r - 2; i >= l; --i) {
      if (hld && depth[csr[i + 1].to] == -1) swap(csr[i], csr[i + 1]);
    }
    int hld_sz = 0;
    for (int i = l; i < r; ++i) {
      auto e = csr[i];
      if (depth[e.to] != -1) continue;
      depth_weighted[e.to] = depth_weighted[v] + e.cost;
      VtoE[e.to] = e.id;
      dfs_sz(e.to, v, hld);
      sz[v] += sz[e.to];
      if (hld && chmax(hld_sz, sz[e.to]) && l < i) { swap(csr[l], csr[i]); }
    }
  }

  void dfs_hld(int v, int &times) {
    LID[v] = times++;
    RID[v] += LID[v];
    V[LID[v]] = v;
    bool heavy = true;
    for (auto &&e: G[v]) {
      if (depth[e.to] <= depth[v]) continue;
      head[e.to] = (heavy ? head[v] : e.to);
      heavy = false;
      dfs_hld(e.to, times);
    }
  }

  vc<int> heavy_path_at(int v) {
    vc<int> P = {v};
    while (1) {
      int a = P.back();
      for (auto &&e: G[a]) {
        if (e.to != parent[a] && head[e.to] == v) {
          P.eb(e.to);
          break;
        }
      }
      if (P.back() == a) break;
    }
    return P;
  }

  int heavy_child(int v) {
    int k = LID[v] + 1;
    if (k == N) return -1;
    int w = V[k];
    return (parent[w] == v ? w : -1);
  }

  vc<int> memo_tail;

  int tail(int v) {
    if (memo_tail.empty()) {
      memo_tail.assign(N, -1);
      FOR_R(i, N) {
        int v = V[i];
        int w = heavy_child(v);
        memo_tail[v] = (w == -1 ? v : memo_tail[w]);
      }
    }
    return memo_tail[v];
  }

  int e_to_v(int eid) {
    auto e = G.edges[eid];
    return (parent[e.frm] == e.to ? e.frm : e.to);
  }
  int v_to_e(int v) { return VtoE[v]; }
  int get_eid(int u, int v) {
    if (parent[u] != v) swap(u, v);
    assert(parent[u] == v);
    return VtoE[u];
  }

  int ELID(int v) { return 2 * LID[v] - depth[v]; }
  int ERID(int v) { return 2 * RID[v] - depth[v] - 1; }

  // 目標地点へ進む個数が k

  int LA(int v, int k) {
    assert(k <= depth[v]);
    while (1) {
      int u = head[v];
      if (LID[v] - k >= LID[u]) return V[LID[v] - k];
      k -= LID[v] - LID[u] + 1;
      v = parent[u];
    }
  }
  int la(int u, int v) { return LA(u, v); }

  int LCA(int u, int v) {
    for (;; v = parent[head[v]]) {
      if (LID[u] > LID[v]) swap(u, v);
      if (head[u] == head[v]) return u;
    }
  }

  int meet(int a, int b, int c) { return LCA(a, b) ^ LCA(a, c) ^ LCA(b, c); }
  int lca(int u, int v) { return LCA(u, v); }

  int subtree_size(int v, int root = -1) {
    if (root == -1) return RID[v] - LID[v];
    if (v == root) return N;
    int x = jump(v, root, 1);
    if (in_subtree(v, x)) return RID[v] - LID[v];
    return N - RID[x] + LID[x];
  }

  int dist(int a, int b) {
    int c = LCA(a, b);
    return depth[a] + depth[b] - 2 * depth[c];
  }

  WT dist_weighted(int a, int b) {
    int c = LCA(a, b);
    return depth_weighted[a] + depth_weighted[b] - WT(2) * depth_weighted[c];
  }

  // a is in b

  bool in_subtree(int a, int b) { return LID[b] <= LID[a] && LID[a] < RID[b]; }

  int jump(int a, int b, ll k) {
    if (k == 1) {
      if (a == b) return -1;
      return (in_subtree(b, a) ? LA(b, depth[b] - depth[a] - 1) : parent[a]);
    }
    int c = LCA(a, b);
    int d_ac = depth[a] - depth[c];
    int d_bc = depth[b] - depth[c];
    if (k > d_ac + d_bc) return -1;
    if (k <= d_ac) return LA(a, k);
    return LA(b, d_ac + d_bc - k);
  }

  vc<int> collect_child(int v) {
    vc<int> res;
    for (auto &&e: G[v])
      if (e.to != parent[v]) res.eb(e.to);
    return res;
  }

  vc<int> collect_subtree(int v) { return {V.begin() + LID[v], V.begin() + RID[v]}; }

  vc<int> collect_light(int v) {
    vc<int> res;
    bool skip = true;
    for (auto &&e: G[v])
      if (e.to != parent[v]) {
        if (!skip) res.eb(e.to);
        skip = false;
      }
    return res;
  }

  vc<pair<int, int>> get_path_decomposition(int u, int v, bool edge) {
    // [始点, 終点] の"閉"区間列。

    vc<pair<int, int>> up, down;
    while (1) {
      if (head[u] == head[v]) break;
      if (LID[u] < LID[v]) {
        down.eb(LID[head[v]], LID[v]);
        v = parent[head[v]];
      } else {
        up.eb(LID[u], LID[head[u]]);
        u = parent[head[u]];
      }
    }
    if (LID[u] < LID[v]) down.eb(LID[u] + edge, LID[v]);
    elif (LID[v] + edge <= LID[u]) up.eb(LID[u], LID[v] + edge);
    reverse(all(down));
    up.insert(up.end(), all(down));
    return up;
  }

  // 辺の列の情報 (frm,to,str)

  // str = "heavy_up", "heavy_down", "light_up", "light_down"

  vc<tuple<int, int, string>> get_path_decomposition_detail(int u, int v) {
    vc<tuple<int, int, string>> up, down;
    while (1) {
      if (head[u] == head[v]) break;
      if (LID[u] < LID[v]) {
        if (v != head[v]) down.eb(head[v], v, "heavy_down"), v = head[v];
        down.eb(parent[v], v, "light_down"), v = parent[v];
      } else {
        if (u != head[u]) up.eb(u, head[u], "heavy_up"), u = head[u];
        up.eb(u, parent[u], "light_up"), u = parent[u];
      }
    }
    if (LID[u] < LID[v]) down.eb(u, v, "heavy_down");
    elif (LID[v] < LID[u]) up.eb(u, v, "heavy_up");
    reverse(all(down));
    concat(up, down);
    return up;
  }

  vc<int> restore_path(int u, int v) {
    vc<int> P;
    for (auto &&[a, b]: get_path_decomposition(u, v, 0)) {
      if (a <= b) {
        FOR(i, a, b + 1) P.eb(V[i]);
      } else {
        FOR_R(i, b, a + 1) P.eb(V[i]);
      }
    }
    return P;
  }

  // path [a,b] と [c,d] の交わり. 空ならば {-1,-1}.

  // https://codeforces.com/problemset/problem/500/G

  pair<int, int> path_intersection(int a, int b, int c, int d) {
    int ab = lca(a, b), ac = lca(a, c), ad = lca(a, d);
    int bc = lca(b, c), bd = lca(b, d), cd = lca(c, d);
    int x = ab ^ ac ^ bc, y = ab ^ ad ^ bd; // meet(a,b,c), meet(a,b,d)

    if (x != y) return {x, y};
    int z = ac ^ ad ^ cd;
    if (x != z) x = -1;
    return {x, x};
  }

  // uv path 上で check(v) を満たす最後の v

  // なければ (つまり check(v) が ng )-1

  template <class F>
  int max_path(F check, int u, int v) {
    if (!check(u)) return -1;
    auto pd = get_path_decomposition(u, v, false);
    for (auto [a, b]: pd) {
      if (!check(V[a])) return u;
      if (check(V[b])) {
        u = V[b];
        continue;
      }
      int c = binary_search([&](int c) -> bool { return check(V[c]); }, a, b, 0);
      return V[c];
    }
    return u;
  }
};
#line 4 "graph/ds/tree_wavelet_matrix.hpp"

// https://atcoder.jp/contests/pakencamp-2022-day1/tasks/pakencamp_2022_day1_j
// https://atcoder.jp/contests/utpc2011/tasks/utpc2011_12
template <typename TREE, bool edge, typename T, bool COMPRESS,
          typename Monoid = Monoid_Add<T>>
struct Tree_Wavelet_Matrix {
  TREE& tree;
  int N;
  using WM = Wavelet_Matrix<T, COMPRESS, Monoid_Add<T>>;
  using X = typename Monoid::value_type;
  WM wm;

  Tree_Wavelet_Matrix(TREE& tree, vc<T> A, vc<X> SUM_data = {}, int log = -1)
      : tree(tree), N(tree.N) {
    vc<X>& S = SUM_data;
    vc<T> A1;
    vc<X> S1;
    A1.resize(N);
    if (!S.empty()) S1.resize(N);
    if (!edge) {
      assert(len(A) == N && (len(S) == 0 || len(S) == N));
      FOR(v, N) A1[tree.LID[v]] = A[v];
      if (len(S) == N) { FOR(v, N) S1[tree.LID[v]] = S[v]; }
      wm.build(A1, S1, log);
    } else {
      assert(len(A) == N - 1 && (len(S) == 0 || len(S) == N - 1));
      if (!S.empty()) {
        FOR(e, N - 1) { S1[tree.LID[tree.e_to_v(e)]] = S[e]; }
      }
      FOR(e, N - 1) { A1[tree.LID[tree.e_to_v(e)]] = A[e]; }
      wm.build(A1, S1, log);
    }
  }

  // xor した結果で [a, b) に収まるものを数える
  int count_path(int s, int t, T a, T b, T xor_val = 0) {
    return wm.count(get_segments(s, t), a, b, xor_val);
  }

  // xor した結果で [a, b) に収まるものを数える
  int count_subtree(int u, T a, T b, T xor_val = 0) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.count(l + edge, r, a, b, xor_val);
  }

  // xor した結果で、[L, R) の中で k>=0 番目と prefix sum
  pair<T, X> kth_value_and_sum_path(int s, int t, int k, T xor_val = 0) {
    return wm.kth_value_and_sum(get_segments(s, t), k, xor_val);
  }

  // xor した結果で、[L, R) の中で k>=0 番目と prefix sum
  pair<T, X> kth_value_and_sum_subtree(int u, int k, T xor_val = 0) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.kth_value_and_sum(l + edge, r, k, xor_val);
  }

  // xor した結果で、[L, R) の中で k>=0 番目
  T kth_path(int s, int t, int k, T xor_val = 0) {
    return wm.kth(get_segments(s, t), k, xor_val);
  }

  // xor した結果で、[L, R) の中で k>=0 番目
  T kth_subtree(int u, int k, T xor_val = 0) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.kth(l + edge, r, k, xor_val);
  }

  // xor した結果で、[L, R) の中で中央値。
  // LOWER = true:下側中央値、false:上側中央値
  T median_path(bool UPPER, int s, int t, T xor_val = 0) {
    return wm.median(UPPER, get_segments(s, t), xor_val);
  }

  T median_subtree(bool UPPER, int u, T xor_val = 0) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.median(UPPER, l + edge, r, xor_val);
  }

  // xor した結果で [k1, k2) 番目であるところの SUM_data の和
  X sum_path(int s, int t, int k1, int k2, T xor_val = 0) {
    return wm.sum(get_segments(s, t), k1, k2, xor_val);
  }

  // xor した結果で [k1, k2) 番目であるところの SUM_data の和
  X sum_subtree(int u, int k1, int k2, T xor_val = 0) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.sum(l + edge, r, k1, k2, xor_val);
  }

  X sum_all_path(int s, int t) { return wm.sum_all(get_segments(s, t)); }

  X sum_all_subtree(int u) {
    int l = tree.LID[u], r = tree.RID[u];
    return wm.sum_all(l + edge, r);
  }

private:
  vc<pair<int, int>> get_segments(int s, int t) {
    vc<pair<int, int>> segments = tree.get_path_decomposition(s, t, edge);
    for (auto&& [a, b]: segments) {
      if (a >= b) swap(a, b);
      ++b;
    }
    return segments;
  }
};
Back to top page