tree_diameter.hpp¶
Unweighted tree diameter via double BFS.
Verified by tree_diameter.
求带权或无权树的直径端点、长度与路径；适合树上最远距离和中心类问题。
Implementation¶
#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>>> &graph) {
if (graph.empty()) {
return {Weight{}, {}};
}
auto traverse = [&](int start) {
std::vector<Weight> distance(graph.size());
std::vector<int> parent(graph.size(), -1);
std::vector<int> stack = {start};
parent[start] = start;
for (int index = 0; index < int(stack.size()); index++) {
int vertex = stack[index];
for (auto [next, weight] : graph[vertex]) {
if (parent[next] != -1) {
continue;
}
parent[next] = vertex;
distance[next] = distance[vertex] + weight;
stack.push_back(next);
}
}
int farthest = start;
for (int vertex = 0; vertex < int(graph.size()); vertex++) {
if (distance[farthest] < distance[vertex]) {
farthest = vertex;
}
}
return std::tuple{farthest, std::move(distance), std::move(parent)};
};
auto [first, ignored_distance, ignored_parent] = traverse(0);
auto [second, distance, parent] = traverse(first);
std::vector<int> path;
for (int vertex = second;; vertex = parent[vertex]) {
path.push_back(vertex);
if (vertex == first) {
break;
}
}
std::reverse(path.begin(), path.end());
return {distance[second], path};
}
/// @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 Operation, auto Identity,
bool Legacy = !std::is_same_v<decltype(Operation), std::nullptr_t>>
struct sparse_table_algebra;
template <class Semilattice, auto Operation, auto Identity>
struct sparse_table_algebra<Semilattice, Operation, Identity, false> {
static_assert(std::is_same_v<decltype(Identity), 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 &left, const value_type &right) {
return Semilattice::op(left, right);
}
};
template <class Value, auto Operation, auto Identity>
struct sparse_table_algebra<Value, Operation, Identity, true> {
using value_type = Value;
static value_type unit() {
if constexpr (std::is_same_v<decltype(Identity), std::nullptr_t>) {
return value_type{};
} else {
return value_type(Identity);
}
}
static value_type op(const value_type &left, const value_type &right) {
return Operation(left, right);
}
};
} // 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 Operation = nullptr, auto Identity = nullptr>
class sparse_table {
using algebra_type =
internal::sparse_table_algebra<Algebra, Operation, Identity>;
public:
using value_type = typename algebra_type::value_type;
sparse_table() = default;
explicit sparse_table(const std::vector<value_type> &values) {
build(values);
}
/// @brief Rebuild from a static array.
void build(const std::vector<value_type> &values) {
size_ = int(values.size());
int levels = int(std::bit_width(static_cast<unsigned>(size_)));
table_.assign(levels, {});
if (size_ == 0) {
return;
}
table_[0] = values;
for (int level = 1; level < levels; level++) {
int width = 1 << level;
int half = width >> 1;
table_[level].resize(size_ - width + 1);
for (int left = 0; left + width <= size_; left++) {
table_[level][left] = algebra_type::op(table_[level - 1][left],
table_[level - 1][left + half]);
}
}
}
int size() const { return size_; }
bool empty() const { return size_ == 0; }
/// @brief Return the idempotent product over [left, right).
value_type prod(int left, int right) const {
assert(0 <= left && left <= right && right <= size_);
if (left == right) {
return algebra_type::unit();
}
int level = int(std::bit_width(static_cast<unsigned>(right - left))) - 1;
int width = 1 << level;
return algebra_type::op(table_[level][left], table_[level][right - width]);
}
private:
int size_ = 0;
std::vector<std::vector<value_type>> table_;
};
} // 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 weighted;
std::vector<int> fa;
std::vector<std::vector<int>> ancestors;
template <class T>
fastlca(const std::vector<T> &g = {}, const bool &_weighted = false,
const int &root = 0) {
weighted = _weighted;
if (!g.empty())
build(g, root);
};
template <class T>
void build(const std::vector<std::vector<T>> &g, const int &root = 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, root);
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, root);
}
void build(const std::vector<std::vector<int>> &g = {}, const int &root = 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 parent) -> void {
fa[u] = parent;
siz[u] = 1;
dfn[u] = idx++;
a.push_back(parent == -1 ? std::pair{-1, -1}
: std::pair{dfn[parent], parent});
for (int v : g[u]) {
if (v == parent) {
continue;
}
d[v] = d[u] + 1;
self(self, v, u);
siz[u] += siz[v];
}
};
dfs(dfs, root, -1);
ST.build(a);
int levels = 1;
while ((1 << levels) <= std::max(1, n)) {
levels++;
}
ancestors.assign(levels, std::vector<int>(n, -1));
ancestors[0] = fa;
for (int level = 1; level < levels; level++) {
for (int vertex = 0; vertex < n; vertex++) {
int middle = ancestors[level - 1][vertex];
if (middle != -1) {
ancestors[level][vertex] = ancestors[level - 1][middle];
}
}
}
}
/// @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 vertex, int k) const {
assert(0 <= vertex && vertex < n);
if (k < 0 || k > d[vertex]) {
return -1;
}
for (int level = 0; k > 0; level++, k >>= 1) {
if (k & 1) {
vertex = ancestors[level][vertex];
}
}
return vertex;
}
/// @brief Return the k-th vertex on path first->second, zero-indexed.
int kth_on_path(int first, int second, int k) const {
int ancestor = lca(first, second);
int first_length = d[first] - d[ancestor];
int second_length = d[second] - d[ancestor];
if (k < 0 || k > first_length + second_length) {
return -1;
}
if (k <= first_length) {
return kth_ancestor(first, k);
}
return kth_ancestor(second, first_length + second_length - 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 (!weighted) {
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 start) {
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[start] = 0;
pre[start] = start;
std::vector<int> que{start};
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 start) {
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[start] = 0, start);
pre[start] = start;
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 start) {
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[start] = 0;
pre[start] = start;
que.push_front(start);
while (!que.empty()) {
int u = que.front();
que.pop_front();
for (const auto &[v, weight] : g[u]) {
assert(weight == 0 || weight == 1);
if (dis[v] <= dis[u] + weight) {
continue;
}
dis[v] = dis[u] + weight;
pre[v] = u;
if (weight == 0) {
que.push_front(v);
} else {
que.push_back(v);
}
}
}
return {dis, pre};
}
struct bellman_ford_result {
std::vector<int64_t> distance;
std::vector<int> predecessor;
std::vector<bool> negative_infinite;
};
/// @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>> &edges) {
assert(n >= 0);
std::vector<__int128> distance(n);
std::vector<int> predecessor_edge(n, -1);
int changed_vertex = -1;
for (int iteration = 0; iteration < n; iteration++) {
changed_vertex = -1;
for (int id = 0; id < int(edges.size()); id++) {
auto [from, to, weight] = edges[id];
assert(0 <= from && from < n);
assert(0 <= to && to < n);
__int128 candidate = distance[from] + __int128(weight);
if (candidate < distance[to]) {
distance[to] = candidate;
predecessor_edge[to] = id;
changed_vertex = to;
}
}
}
if (changed_vertex == -1) {
return {};
}
for (int step = 0; step < n; step++) {
int id = predecessor_edge[changed_vertex];
assert(id != -1);
changed_vertex = std::get<0>(edges[id]);
}
int start = changed_vertex;
std::vector<int> cycle;
do {
int id = predecessor_edge[changed_vertex];
assert(id != -1);
cycle.push_back(id);
changed_vertex = std::get<0>(edges[id]);
} while (changed_vertex != start);
std::reverse(cycle.begin(), cycle.end());
return cycle;
}
/// @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 start) {
const int N = int(g.size());
const int64_t INF = std::numeric_limits<int64_t>::max();
const int64_t NEG_INF = std::numeric_limits<int64_t>::lowest();
const __int128 WIDE_INF = __int128(1) << 120;
std::vector<__int128> wide_distance(N, WIDE_INF);
std::vector<int> pre(N, -1);
std::vector<bool> negative_infinite(N);
wide_distance[start] = 0;
pre[start] = start;
for (int iteration = 0; iteration < N; iteration++) {
bool updated = false;
for (int u = 0; u < N; u++) {
if (wide_distance[u] == WIDE_INF) {
continue;
}
for (const auto &[v, weight] : g[u]) {
__int128 candidate = wide_distance[u] + static_cast<__int128>(weight);
if (candidate >= wide_distance[v]) {
continue;
}
wide_distance[v] = candidate;
pre[v] = u;
updated = true;
if (iteration == N - 1) {
negative_infinite[v] = true;
}
}
}
if (!updated) {
break;
}
}
std::vector<int> que;
for (int u = 0; u < N; u++) {
if (negative_infinite[u]) {
que.push_back(u);
}
}
for (int i = 0; i < int(que.size()); i++) {
int u = que[i];
for (const auto &[v, weight] : g[u]) {
(void)weight;
if (!negative_infinite[v]) {
negative_infinite[v] = true;
que.push_back(v);
}
}
}
std::vector<int64_t> dis(N, INF);
for (int u = 0; u < N; u++) {
if (negative_infinite[u] || wide_distance[u] < __int128(NEG_INF)) {
dis[u] = NEG_INF;
} else if (wide_distance[u] < WIDE_INF) {
dis[u] = wide_distance[u] > __int128(INF)
? INF
: static_cast<int64_t>(wide_distance[u]);
}
}
return {dis, pre, negative_infinite};
}
inline std::vector<int> find_path(std::vector<int> &pre, int s, int t) {
std::vector<int> path;
int cur = t;
while (cur != s) {
assert(cur >= 0);
path.push_back(cur);
cur = pre[cur];
}
path.push_back(s);
std::reverse(path.begin(), path.end());
return path;
}
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>>> &graph) {
if (graph.empty()) {
return {Weight{}, {}};
}
auto traverse = [&](int start) {
std::vector<Weight> distance(graph.size());
std::vector<int> parent(graph.size(), -1);
std::vector<int> stack = {start};
parent[start] = start;
for (int index = 0; index < int(stack.size()); index++) {
int vertex = stack[index];
for (auto [next, weight] : graph[vertex]) {
if (parent[next] != -1) {
continue;
}
parent[next] = vertex;
distance[next] = distance[vertex] + weight;
stack.push_back(next);
}
}
int farthest = start;
for (int vertex = 0; vertex < int(graph.size()); vertex++) {
if (distance[farthest] < distance[vertex]) {
farthest = vertex;
}
}
return std::tuple{farthest, std::move(distance), std::move(parent)};
};
auto [first, ignored_distance, ignored_parent] = traverse(0);
auto [second, distance, parent] = traverse(first);
std::vector<int> path;
for (int vertex = second;; vertex = parent[vertex]) {
path.push_back(vertex);
if (vertex == first) {
break;
}
}
std::reverse(path.begin(), path.end());
return {distance[second], path};
}
/// @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