Skip to content

tree_diameter.hpp

SECTIONGraph INCLUDEnoya/tree_diameter.hpp

求带权或无权树的直径端点、长度与路径;适合树上最远距离和中心类问题。

Complexity: Time: O(n) for a tree diameter; O(1) diameter-monoid merge. Space: O(n).

AC 记录:tree_diameter

跳到代码 · GitHub ↗

Implementation

当前头文件,省略 include guard;依赖见 #include

/// @complexity Time: O(n) for a tree diameter; O(1) diameter-monoid merge.
/// Space: O(n).

#include "noya/lowest_common_ancestor.hpp"
#include "noya/shortest_path.hpp"

#include <algorithm>
#include <array>
#include <cassert>
#include <cstdint>
#include <tuple>
#include <utility>
#include <vector>

namespace noya {

/// @brief Unweighted tree diameter via double BFS. @return (diameter, u, v, eccentricity[]).
inline std::tuple<int, int, int, std::vector<int>>
tree_diam(std::vector<std::vector<int>> &g) {
  if (g.empty())
    return {-1, -1, -1, {}};
  auto d0 = bfs_unweighted(g, 0).first;
  int p = std::max_element(d0.begin(), d0.end()) - d0.begin();
  auto dp = bfs_unweighted(g, p).first;
  int q = std::max_element(dp.begin(), dp.end()) - dp.begin();
  auto dq = bfs_unweighted(g, q).first;
  int n = int(g.size());
  std::vector<int> ecc(n);
  for (int i = 0; i < n; i++)
    ecc[i] = std::max(dp[i], dq[i]);
  return {dp[q], p, q, ecc};
}

/// @brief Weighted tree diameter and one realizing vertex path.
/// @return (distance, path from one endpoint to the other).
template <class Weight>
std::pair<Weight, std::vector<int>> weighted_tree_diameter(
    const std::vector<std::vector<std::pair<int, Weight>>> &G) {
  if (G.empty()) {
    return {Weight{}, {}};
  }
  auto dfs = [&](int s) {
    std::vector<Weight> dis(G.size());
    std::vector<int> fa(G.size(), -1);
    std::vector<int> stk = {s};
    fa[s] = s;
    for (int idx = 0; idx < int(stk.size()); idx++) {
      int u = stk[idx];
      for (auto [nxt, w] : G[u]) {
        if (fa[nxt] != -1) {
          continue;
        }
        fa[nxt] = u;
        dis[nxt] = dis[u] + w;
        stk.push_back(nxt);
      }
    }
    int far = s;
    for (int u = 0; u < int(G.size()); u++) {
      if (dis[far] < dis[u]) {
        far = u;
      }
    }
    return std::tuple{far, std::move(dis), std::move(fa)};
  };

  auto [a1, id1, ip] = dfs(0);
  auto [b1, dis, fa] = dfs(a1);
  std::vector<int> pth;
  for (int u = b1;; u = fa[u]) {
    pth.push_back(u);
    if (u == a1) {
      break;
    }
  }
  std::reverse(pth.begin(), pth.end());
  return {dis[b1], pth};
}

/// @brief Diameter monoid for segment tree. Merge two vertex sets and track the farthest pair.
struct diameter_monoid {
  using value_type = std::pair<int64_t, std::array<int, 2>>;
  using S = value_type;
  using X = value_type;

  static constexpr value_type identity = {-1, {-1, -1}};
  static constexpr bool commute = true;

  diameter_monoid() = default;
  explicit diameter_monoid(const fastlca &l) { set_lca(l); }

  static const fastlca *&get_lca() {
    static const fastlca *lca = nullptr;
    return lca;
  }

  static void set_lca(const fastlca &l) { get_lca() = &l; }

  static value_type unit() { return identity; }
  static value_type e() { return unit(); }

  static value_type make(int v) { return {0, {v, v}}; }
  static value_type from_vertex(int v) { return make(v); }

  static value_type op(value_type a, value_type b) {
    if (a == unit()) return b;
    if (b == unit()) return a;
    const fastlca *lca = get_lca();
    assert(lca != nullptr);
    value_type c = std::max(a, b);
    for (auto x : a.second)
      for (auto y : b.second) {
        int64_t d = lca->distance(x, y);
        if (d > c.first)
          c = {d, {x, y}};
      }
    return c;
  }

  static value_type merge(value_type a, value_type b) { return op(a, b); }
};

} // namespace noya
#ifndef NOYA_TREE_DIAMETER_HPP
#define NOYA_TREE_DIAMETER_HPP 1

/// @complexity Time: O(n) for a tree diameter; O(1) diameter-monoid merge.
/// Space: O(n).

#include "noya/lowest_common_ancestor.hpp"
#include "noya/shortest_path.hpp"

#include <algorithm>
#include <array>
#include <cassert>
#include <cstdint>
#include <tuple>
#include <utility>
#include <vector>

namespace noya {

/// @brief Unweighted tree diameter via double BFS. @return (diameter, u, v, eccentricity[]).
inline std::tuple<int, int, int, std::vector<int>>
tree_diam(std::vector<std::vector<int>> &g) {
  if (g.empty())
    return {-1, -1, -1, {}};
  auto d0 = bfs_unweighted(g, 0).first;
  int p = std::max_element(d0.begin(), d0.end()) - d0.begin();
  auto dp = bfs_unweighted(g, p).first;
  int q = std::max_element(dp.begin(), dp.end()) - dp.begin();
  auto dq = bfs_unweighted(g, q).first;
  int n = int(g.size());
  std::vector<int> ecc(n);
  for (int i = 0; i < n; i++)
    ecc[i] = std::max(dp[i], dq[i]);
  return {dp[q], p, q, ecc};
}

/// @brief Weighted tree diameter and one realizing vertex path.
/// @return (distance, path from one endpoint to the other).
template <class Weight>
std::pair<Weight, std::vector<int>> weighted_tree_diameter(
    const std::vector<std::vector<std::pair<int, Weight>>> &G) {
  if (G.empty()) {
    return {Weight{}, {}};
  }
  auto dfs = [&](int s) {
    std::vector<Weight> dis(G.size());
    std::vector<int> fa(G.size(), -1);
    std::vector<int> stk = {s};
    fa[s] = s;
    for (int idx = 0; idx < int(stk.size()); idx++) {
      int u = stk[idx];
      for (auto [nxt, w] : G[u]) {
        if (fa[nxt] != -1) {
          continue;
        }
        fa[nxt] = u;
        dis[nxt] = dis[u] + w;
        stk.push_back(nxt);
      }
    }
    int far = s;
    for (int u = 0; u < int(G.size()); u++) {
      if (dis[far] < dis[u]) {
        far = u;
      }
    }
    return std::tuple{far, std::move(dis), std::move(fa)};
  };

  auto [a1, id1, ip] = dfs(0);
  auto [b1, dis, fa] = dfs(a1);
  std::vector<int> pth;
  for (int u = b1;; u = fa[u]) {
    pth.push_back(u);
    if (u == a1) {
      break;
    }
  }
  std::reverse(pth.begin(), pth.end());
  return {dis[b1], pth};
}

/// @brief Diameter monoid for segment tree. Merge two vertex sets and track the farthest pair.
struct diameter_monoid {
  using value_type = std::pair<int64_t, std::array<int, 2>>;
  using S = value_type;
  using X = value_type;

  static constexpr value_type identity = {-1, {-1, -1}};
  static constexpr bool commute = true;

  diameter_monoid() = default;
  explicit diameter_monoid(const fastlca &l) { set_lca(l); }

  static const fastlca *&get_lca() {
    static const fastlca *lca = nullptr;
    return lca;
  }

  static void set_lca(const fastlca &l) { get_lca() = &l; }

  static value_type unit() { return identity; }
  static value_type e() { return unit(); }

  static value_type make(int v) { return {0, {v, v}}; }
  static value_type from_vertex(int v) { return make(v); }

  static value_type op(value_type a, value_type b) {
    if (a == unit()) return b;
    if (b == unit()) return a;
    const fastlca *lca = get_lca();
    assert(lca != nullptr);
    value_type c = std::max(a, b);
    for (auto x : a.second)
      for (auto y : b.second) {
        int64_t d = lca->distance(x, y);
        if (d > c.first)
          c = {d, {x, y}};
      }
    return c;
  }

  static value_type merge(value_type a, value_type b) { return op(a, b); }
};

} // namespace noya

#endif // NOYA_TREE_DIAMETER_HPP
#include <algorithm>
#include <array>
#include <bit>
#include <cassert>
#include <cstddef>
#include <cstdint>
#include <deque>
#include <functional>
#include <limits>
#include <queue>
#include <tuple>
#include <type_traits>
#include <utility>
#include <vector>

/// @complexity Time: O(n) for a tree diameter; O(1) diameter-monoid merge.
/// Space: O(n).

/// @complexity Time: O(n log n) build and O(1) LCA/ancestor queries.
/// Space: O(n log n).

/// @complexity Time: O(n log n) build and O(1) idempotent range query.
/// Space: O(n log n).

namespace noya {

namespace internal {

template <class Algebra, auto fn, auto e,
          bool old = !std::is_same_v<decltype(fn), std::nullptr_t>>
struct sparse_table_algebra;

template <class Semilattice, auto fn, auto e>
struct sparse_table_algebra<Semilattice, fn, e, false> {
  static_assert(std::is_same_v<decltype(e), std::nullptr_t>,
                "a policy sparse_table gets its identity from unit()");

  using value_type = typename Semilattice::value_type;

  static value_type unit() { return Semilattice::unit(); }
  static value_type op(const value_type &l, const value_type &r) {
    return Semilattice::op(l, r);
  }
};

template <class Value, auto fn, auto e>
struct sparse_table_algebra<Value, fn, e, true> {
  using value_type = Value;

  static value_type unit() {
    if constexpr (std::is_same_v<decltype(e), std::nullptr_t>) {
      return value_type{};
    } else {
      return value_type(e);
    }
  }

  static value_type op(const value_type &l, const value_type &r) {
    return fn(l, r);
  }
};

} // namespace internal

/// @brief Static O(1) range product for an idempotent semilattice.
/// The preferred interface is sparse_table<Semilattice>, where Semilattice
/// provides value_type, unit(), and an associative, idempotent op(). The legacy
/// sparse_table<Value, operation, identity> spelling remains supported.
/// @details For k = floor(log2(r-l)), [l, r) is covered by the length-2^k
/// blocks beginning at l and ending at r. They may overlap, so combining them
/// is correct precisely because op(x, x) = x.
template <class Algebra, auto fn = nullptr, auto e = nullptr>
class sparse_table {
  using algebra_type = internal::sparse_table_algebra<Algebra, fn, e>;

public:
  using value_type = typename algebra_type::value_type;

  sparse_table() = default;

  explicit sparse_table(const std::vector<value_type> &a) { build(a); }

  /// @brief Rebuild from a static array.
  void build(const std::vector<value_type> &a) {
    n = int(a.size());
    int lg = int(std::bit_width(static_cast<unsigned>(n)));
    st.assign(lg, {});
    if (n == 0) {
      return;
    }

    st[0] = a;
    for (int dep = 1; dep < lg; dep++) {
      int w = 1 << dep;
      int hf = w >> 1;
      st[dep].resize(n - w + 1);
      for (int l = 0; l + w <= n; l++) {
        st[dep][l] = algebra_type::op(st[dep - 1][l], st[dep - 1][l + hf]);
      }
    }
  }

  int size() const { return n; }
  bool empty() const { return n == 0; }

  /// @brief Return the idempotent product over [l, r).
  value_type prod(int l, int r) const {
    assert(0 <= l && l <= r && r <= n);
    if (l == r) {
      return algebra_type::unit();
    }
    int dep = int(std::bit_width(static_cast<unsigned>(r - l))) - 1;
    int w = 1 << dep;
    return algebra_type::op(st[dep][l], st[dep][r - w]);
  }

private:
  int n = 0;
  std::vector<std::vector<value_type>> st;
};

} // namespace noya

namespace noya {
/// @brief Sparse-table-based LCA with O(n log n) build and O(1) query.
struct fastlca {
  struct minimum_pair_semilattice {
    using value_type = std::pair<int, int>;
    static value_type unit() {
      return {std::numeric_limits<int>::max(), -1};
    }
    static value_type op(value_type a, value_type b) { return std::min(a, b); }
  };
  int n;
  std::vector<int> dfn;
  std::vector<int> d;
  std::vector<int64_t> len;
  std::vector<int> siz;
  sparse_table<minimum_pair_semilattice> ST;
  bool wt;

  std::vector<int> fa;
  std::vector<std::vector<int>> up;

  template <class T>
  fastlca(const std::vector<T> &g = {}, const bool &wt1 = false,
          const int &rt = 0) {
    wt = wt1;
    if (!g.empty())
      build(g, rt);
  };

  template <class T>
  void build(const std::vector<std::vector<T>> &g, const int &rt = 0) {
    n = int(g.size());
    std::vector<std::vector<int>> g2(n);

    for (int u = 0; u < n; u++) {
      for (auto &[v, w] : g[u]) {
        g2[u].push_back(v);
      }
    }
    build(g2, rt);

    len.assign(n, 0);
    auto dfs = [&](auto &self, int u) -> void {
      for (auto &[v, w] : g[u]) {
        if (fa[v] == u) {
          len[v] = len[u] + w;
          self(self, v);
        }
      }
    };
    dfs(dfs, rt);
  }

  void build(const std::vector<std::vector<int>> &g = {}, const int &rt = 0) {
    n = int(g.size());
    d.assign(n, 0);
    dfn.assign(n, -1);
    siz.assign(n, 0);
    fa.assign(n, -1);

    int idx = 0;
    std::vector<std::pair<int, int>> a;
    a.reserve(n);
    auto dfs = [&](auto &self, int u, int par) -> void {
      fa[u] = par;
      siz[u] = 1;
      dfn[u] = idx++;
      a.push_back(par == -1 ? std::pair{-1, -1}
                               : std::pair{dfn[par], par});
      for (int v : g[u]) {
        if (v == par) {
          continue;
        }
        d[v] = d[u] + 1;
        self(self, v, u);
        siz[u] += siz[v];
      }
    };
    dfs(dfs, rt, -1);
    ST.build(a);
    int dep = 1;
    while ((1 << dep) <= std::max(1, n)) {
      dep++;
    }
    up.assign(dep, std::vector<int>(n, -1));
    up[0] = fa;
    for (int de1 = 1; de1 < dep; de1++) {
      for (int u1 = 0; u1 < n; u1++) {
        int mid = up[de1 - 1][u1];
        if (mid != -1) {
          up[de1][u1] = up[de1 - 1][mid];
        }
      }
    }
  }

  /// @brief Check if node a is in the subtree of node b.
  bool is_subtree(int a, int b) {
    if (dfn[b] <= dfn[a] && dfn[a] < dfn[b] + siz[b]) {
      return true;
    } else {
      return false;
    }
  }

  /// @brief Return the lowest common ancestor of nodes u and v.
  int lca(int u, int v) const {
    assert(0 <= u && u < n);
    assert(0 <= v && v < n);

    if (u == v) {
      return u;
    } else {
      int a = dfn[u];
      int b = dfn[v];
      if (a > b) {
        std::swap(a, b);
      }
      return ST.prod(a + 1, b + 1).second;
    }
  }

  /// @brief Return the LCA when the tree is re-rooted at c.
  int rooted_lca(int a, int b, int c) const {
    return lca(a, b) ^ lca(a, c) ^ lca(b, c);
  }

  /// @brief Return the k-th ancestor of vertex, or -1 when it does not exist.
  int kth_ancestor(int u1, int k) const {
    assert(0 <= u1 && u1 < n);
    if (k < 0 || k > d[u1]) {
      return -1;
    }
    for (int de1 = 0; k > 0; de1++, k >>= 1) {
      if (k & 1) {
        u1 = up[de1][u1];
      }
    }
    return u1;
  }

  /// @brief Return the k-th vertex on path first->second, zero-indexed.
  int kth_on_path(int lhs, int rhs, int k) const {
    int anc = lca(lhs, rhs);
    int al = d[lhs] - d[anc];
    int bl = d[rhs] - d[anc];
    if (k < 0 || k > al + bl) {
      return -1;
    }
    if (k <= al) {
      return kth_ancestor(lhs, k);
    }
    return kth_ancestor(rhs, al + bl - k);
  }

  /// @brief Return the (weighted or unweighted) distance between nodes a and b.
  int64_t distance(int a, int b) const {
    int c = lca(a, b);
    if (!wt) {
      return d[a] + d[b] - d[c] * 2;
    } else {
      return len[a] + len[b] - len[c] * 2;
    }
  }

  /// @brief Compute the intersection of paths (a,b) and (c,d) as a pair of endpoints.
  std::pair<int, int> intersection(int a, int b, int c, int d) const {
    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;
    int y = ab ^ ad ^ bd;
    if (x != y) {
      return {x, y};
    }
    int z = ac ^ ad ^ cd;
    if (x != z) {
      x = -1;
    }
    return {x, x};
  }

  std::pair<int, int> intersection(std::pair<int, int> a,
                                   std::pair<int, int> b) const {
    return intersection(a.first, a.second, b.first, b.second);
  }

  /// @brief Return the list of vertices on the path from a to b.
  std::vector<int> path(int a, int b) {
    int c = lca(a, b);
    std::vector<int> ac;
    while (a != c) {
      ac.push_back(a);
      a = fa[a];
    }
    std::vector<int> bc;
    while (b != c) {
      bc.push_back(b);
      b = fa[b];
    }
    std::vector<int> res = std::move(ac);
    res.push_back(c);
    res.insert(res.end(), bc.rbegin(), bc.rend());
    return res;
  }
};
} // namespace noya

/// @complexity Time: O(V + E) BFS/0-1 BFS, O((V + E) log V) Dijkstra, O(VE) signed-cycle routines.
/// Space: O(V + E).

namespace noya {
/// @brief BFS on unweighted graph. @return (distances, predecessors).
inline std::pair<std::vector<int>, std::vector<int>>
bfs_unweighted(std::vector<std::vector<int>> &g, int src) {
  int N = int(g.size());
  const int INF = std::numeric_limits<int>::max();
  std::vector<int> dis(N, INF);
  std::vector<int> pre(N, -1);
  dis[src] = 0;
  pre[src] = src;

  std::vector<int> que{src};
  for (int i = 0; i < int(que.size()); i++) {
    int u = que[i];
    for (auto v : g[u])
      if (dis[v] == INF) {
        dis[v] = dis[u] + 1;
        pre[v] = u;
        que.push_back(v);
      }
  }
  return {dis, pre};
}

/// @brief Dijkstra's algorithm. @return (distances, predecessors).
template <class T>
std::pair<std::vector<int64_t>, std::vector<int>>
dijkstra(std::vector<std::vector<T>> &g, int src) {
  int N = int(g.size());
  const int64_t INF = std::numeric_limits<int64_t>::max();
  std::vector<int64_t> dis(N, INF);
  std::vector<int> pre(N, -1);
  std::priority_queue<std::pair<int64_t, int>,
                      std::vector<std::pair<int64_t, int>>, std::greater<>>
      que;
  que.emplace(dis[src] = 0, src);
  pre[src] = src;
  while (!que.empty()) {
    auto [d, u] = que.top();
    que.pop();
    if (d > dis[u])
      continue;
    for (auto &[v, w] : g[u])
      if (dis[v] > dis[u] + w) {
        dis[v] = dis[u] + w;
        pre[v] = u;
        que.emplace(dis[v], v);
      }
  }
  return {dis, pre};
}

/// @brief Shortest paths in a graph whose edge weights are 0 or 1.
/// @return (distances, predecessors).
template <class T>
std::pair<std::vector<int64_t>, std::vector<int>>
bfs01(const std::vector<std::vector<T>> &g, int src) {
  const int N = int(g.size());
  const int64_t INF = std::numeric_limits<int64_t>::max();
  std::vector<int64_t> dis(N, INF);
  std::vector<int> pre(N, -1);
  std::deque<int> que;
  dis[src] = 0;
  pre[src] = src;
  que.push_front(src);
  while (!que.empty()) {
    int u = que.front();
    que.pop_front();
    for (const auto &[v, w1] : g[u]) {
      assert(w1 == 0 || w1 == 1);
      if (dis[v] <= dis[u] + w1) {
        continue;
      }
      dis[v] = dis[u] + w1;
      pre[v] = u;
      if (w1 == 0) {
        que.push_front(v);
      } else {
        que.push_back(v);
      }
    }
  }
  return {dis, pre};
}

struct bellman_ford_result {
  std::vector<int64_t> di1;
  std::vector<int> prv;
  std::vector<bool> neg;
};

/// @brief Return edge ids forming one negative directed cycle anywhere in the
/// graph, in traversal order, or an empty vector when no negative cycle exists.
template <class Weight>
std::vector<int>
find_negative_cycle(int n,
                    const std::vector<std::tuple<int, int, Weight>> &es) {
  assert(n >= 0);
  std::vector<__int128> di1(n);
  std::vector<int> pe(n, -1);
  int cv = -1;
  for (int ite = 0; ite < n; ite++) {
    cv = -1;
    for (int id = 0; id < int(es.size()); id++) {
      auto [u1, to, w1] = es[id];
      assert(0 <= u1 && u1 < n);
      assert(0 <= to && to < n);
      __int128 can = di1[u1] + __int128(w1);
      if (can < di1[to]) {
        di1[to] = can;
        pe[to] = id;
        cv = to;
      }
    }
  }
  if (cv == -1) {
    return {};
  }
  for (int stp = 0; stp < n; stp++) {
    int id = pe[cv];
    assert(id != -1);
    cv = std::get<0>(es[id]);
  }
  int src = cv;
  std::vector<int> cyc;
  do {
    int id = pe[cv];
    assert(id != -1);
    cyc.push_back(id);
    cv = std::get<0>(es[id]);
  } while (cv != src);
  std::reverse(cyc.begin(), cyc.end());
  return cyc;
}

/// @brief Bellman-Ford with propagation from source-reachable negative cycles.
/// A negative-infinite vertex has distance numeric_limits<int64_t>::lowest().
template <class T>
bellman_ford_result bellman_ford(const std::vector<std::vector<T>> &g,
                                 int src) {
  const int N = int(g.size());
  const int64_t INF = std::numeric_limits<int64_t>::max();
  const int64_t NIF = std::numeric_limits<int64_t>::lowest();
  const __int128 WIF = __int128(1) << 120;
  std::vector<__int128> wd(N, WIF);
  std::vector<int> pre(N, -1);
  std::vector<bool> neg(N);
  wd[src] = 0;
  pre[src] = src;

  for (int ite = 0; ite < N; ite++) {
    bool upd = false;
    for (int u = 0; u < N; u++) {
      if (wd[u] == WIF) {
        continue;
      }
      for (const auto &[v, w1] : g[u]) {
        __int128 can = wd[u] + static_cast<__int128>(w1);
        if (can >= wd[v]) {
          continue;
        }
        wd[v] = can;
        pre[v] = u;
        upd = true;
        if (ite == N - 1) {
          neg[v] = true;
        }
      }
    }
    if (!upd) {
      break;
    }
  }

  std::vector<int> que;
  for (int u = 0; u < N; u++) {
    if (neg[u]) {
      que.push_back(u);
    }
  }
  for (int i = 0; i < int(que.size()); i++) {
    int u = que[i];
    for (const auto &[v, w1] : g[u]) {
      (void)w1;
      if (!neg[v]) {
        neg[v] = true;
        que.push_back(v);
      }
    }
  }

  std::vector<int64_t> dis(N, INF);
  for (int u = 0; u < N; u++) {
    if (neg[u] || wd[u] < __int128(NIF)) {
      dis[u] = NIF;
    } else if (wd[u] < WIF) {
      dis[u] = wd[u] > __int128(INF)
                   ? INF
                   : static_cast<int64_t>(wd[u]);
    }
  }
  return {dis, pre, neg};
}

inline std::vector<int> find_path(std::vector<int> &pre, int s, int t) {
  std::vector<int> pth;
  int cur = t;
  while (cur != s) {
    assert(cur >= 0);
    pth.push_back(cur);
    cur = pre[cur];
  }
  pth.push_back(s);
  std::reverse(pth.begin(), pth.end());
  return pth;
}

template <class T>
std::vector<std::vector<int64_t>> floyd(std::vector<std::vector<T>> &g) {
  int N = int(g.size());
  const int64_t INF = std::numeric_limits<T>::max() / 2;
  std::vector<std::vector<int64_t>> f(N, std::vector<int64_t>(N, INF));
  for (int i = 0; i < N; i++)
    for (int j = 0; j < N; j++)
      f[i][j] = g[i][j];
  for (int i = 0; i < N; i++)
    f[i][i] = 0;
  for (int k = 0; k < N; k++)
    for (int i = 0; i < N; i++)
      for (int j = 0; j < N; j++)
        f[i][j] = std::min(f[i][j], f[i][k] + f[k][j]);
  return f;
}

} // namespace noya

namespace noya {

/// @brief Unweighted tree diameter via double BFS. @return (diameter, u, v, eccentricity[]).
inline std::tuple<int, int, int, std::vector<int>>
tree_diam(std::vector<std::vector<int>> &g) {
  if (g.empty())
    return {-1, -1, -1, {}};
  auto d0 = bfs_unweighted(g, 0).first;
  int p = std::max_element(d0.begin(), d0.end()) - d0.begin();
  auto dp = bfs_unweighted(g, p).first;
  int q = std::max_element(dp.begin(), dp.end()) - dp.begin();
  auto dq = bfs_unweighted(g, q).first;
  int n = int(g.size());
  std::vector<int> ecc(n);
  for (int i = 0; i < n; i++)
    ecc[i] = std::max(dp[i], dq[i]);
  return {dp[q], p, q, ecc};
}

/// @brief Weighted tree diameter and one realizing vertex path.
/// @return (distance, path from one endpoint to the other).
template <class Weight>
std::pair<Weight, std::vector<int>> weighted_tree_diameter(
    const std::vector<std::vector<std::pair<int, Weight>>> &G) {
  if (G.empty()) {
    return {Weight{}, {}};
  }
  auto dfs = [&](int s) {
    std::vector<Weight> dis(G.size());
    std::vector<int> fa(G.size(), -1);
    std::vector<int> stk = {s};
    fa[s] = s;
    for (int idx = 0; idx < int(stk.size()); idx++) {
      int u = stk[idx];
      for (auto [nxt, w] : G[u]) {
        if (fa[nxt] != -1) {
          continue;
        }
        fa[nxt] = u;
        dis[nxt] = dis[u] + w;
        stk.push_back(nxt);
      }
    }
    int far = s;
    for (int u = 0; u < int(G.size()); u++) {
      if (dis[far] < dis[u]) {
        far = u;
      }
    }
    return std::tuple{far, std::move(dis), std::move(fa)};
  };

  auto [a1, id1, ip] = dfs(0);
  auto [b1, dis, fa] = dfs(a1);
  std::vector<int> pth;
  for (int u = b1;; u = fa[u]) {
    pth.push_back(u);
    if (u == a1) {
      break;
    }
  }
  std::reverse(pth.begin(), pth.end());
  return {dis[b1], pth};
}

/// @brief Diameter monoid for segment tree. Merge two vertex sets and track the farthest pair.
struct diameter_monoid {
  using value_type = std::pair<int64_t, std::array<int, 2>>;
  using S = value_type;
  using X = value_type;

  static constexpr value_type identity = {-1, {-1, -1}};
  static constexpr bool commute = true;

  diameter_monoid() = default;
  explicit diameter_monoid(const fastlca &l) { set_lca(l); }

  static const fastlca *&get_lca() {
    static const fastlca *lca = nullptr;
    return lca;
  }

  static void set_lca(const fastlca &l) { get_lca() = &l; }

  static value_type unit() { return identity; }
  static value_type e() { return unit(); }

  static value_type make(int v) { return {0, {v, v}}; }
  static value_type from_vertex(int v) { return make(v); }

  static value_type op(value_type a, value_type b) {
    if (a == unit()) return b;
    if (b == unit()) return a;
    const fastlca *lca = get_lca();
    assert(lca != nullptr);
    value_type c = std::max(a, b);
    for (auto x : a.second)
      for (auto y : b.second) {
        int64_t d = lca->distance(x, y);
        if (d > c.first)
          c = {d, {x, y}};
      }
    return c;
  }

  static value_type merge(value_type a, value_type b) { return op(a, b); }
};

} // namespace noya