m1une's library

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

View on GitHub

:heavy_check_mark: Replacement Paths
(graph/replacement_paths.hpp)

Overview

Replacement paths answer every single-failure shortest-path query along one fixed shortest path $P=(p_0=s,p_1,\ldots,p_k=t)$. One call computes either the shortest $s$-$t$ distance after deleting each logical path edge $(p_i,p_{i+1})$, one edge at a time, or the shortest distance after deleting each path vertex $p_i$, one vertex at a time.

The graph must be undirected and built with Graph<T>::add_edge. Every active edge must have strictly positive cost. Inactive edges are ignored. The input graph and its alive flags are never changed.

The automatic overloads select one shortest path. The GraphPath overloads preserve the caller’s exact shortest path, including its logical edge ids. Edge ids are necessary because parallel edges may join the same pair of vertices.

Path and Results

GraphPath has the following members:

Member Type Meaning
vertices std::vector<int> The ordered vertices $p_0,\ldots,p_k$.
edges std::vector<int> edges[i] is the active logical edge joining vertices[i] and vertices[i + 1].

It satisfies vertices.size() == edges.size() + 1. A supplied path must be nonempty and simple, and its stored total cost must equal the shortest distance between its endpoints.

EdgeReplacementPathsResult<T> contains:

Member / method Type / Signature Meaning
path GraphPath The selected or supplied fixed path.
replacement_dist std::vector<T> Entry i is the shortest distance after deleting logical edge path.edges[i]. Its size is path.edges.size().
inf T The unreachable-distance sentinel.
reachable bool reachable(int path_edge_index) const Tests whether the corresponding replacement distance is finite.

VertexReplacementPathsResult<T> contains:

Member / method Type / Signature Meaning
path GraphPath The selected or supplied fixed path.
replacement_dist std::vector<T> Entry i is the shortest distance after deleting vertex path.vertices[i]. Its size is path.vertices.size().
inf T The unreachable-distance sentinel.
reachable bool reachable(int path_vertex_index) const Tests whether the corresponding replacement distance is finite.

Deleting the source or target makes the route unavailable, so the first and last vertex-replacement entries are inf. If s == t, the path has one vertex and no edges, the edge result is empty, and the vertex result contains one inf.

Functions

Let $N$ be the number of vertices, $M$ the number of logical edges, and $K$ the number of edges in the fixed path.

Function Signature Description Complexity
edge_replacement_paths template <class T> EdgeReplacementPathsResult<T> edge_replacement_paths(const Graph<T>& g, int s, int t, T inf = std::numeric_limits<T>::max() / T(4)) Selects one shortest path and processes every path-edge failure. $O((N+M)\log N)$ time, $O(N+M)$ memory.
edge_replacement_paths template <class T> EdgeReplacementPathsResult<T> edge_replacement_paths(const Graph<T>& g, const GraphPath& path, T inf = std::numeric_limits<T>::max() / T(4)) Processes every edge failure on the exact supplied shortest path. $O((N+M)\log N)$ time, $O(N+M)$ memory.
vertex_replacement_paths template <class T> VertexReplacementPathsResult<T> vertex_replacement_paths(const Graph<T>& g, int s, int t, T inf = std::numeric_limits<T>::max() / T(4)) Selects one shortest path and processes every path-vertex failure. $O(K(N+M)\log N)$ time, $O(N+M)$ memory.
vertex_replacement_paths template <class T> VertexReplacementPathsResult<T> vertex_replacement_paths(const Graph<T>& g, const GraphPath& path, T inf = std::numeric_limits<T>::max() / T(4)) Processes every vertex failure on the exact supplied shortest path. $O(K(N+M)\log N)$ time, $O(N+M)$ memory.

All additions saturate at inf. Choose inf larger than every finite answer that must be represented.

Edge Algorithm and Correctness

The edge solver runs Dijkstra from both endpoints. It constructs a shortest-path tree rooted at $s$ that contains the exact fixed path: fixed-path parent edges are forced, and every other reachable vertex receives a tight predecessor. Strictly positive costs make predecessor distances decrease strictly toward the root, so the forced tree is acyclic even when shortest paths tie.

Conceptually remove the fixed-path edges from this tree. Each remaining tree component contains exactly one path vertex. Define block[v] = i when v lies in the component containing $p_i$.

Consider an active non-path logical edge joining blocks $a<b$, oriented from u in block $a$ to v in block $b$. The tree shortest path from $s$ to u stays on the source side of every fixed-path edge with index in $[a,b)`. Thus

dist_s[u] + cost(u, v) + dist_t[v]

is a candidate for every such edge failure. If the chosen shortest suffix from v crosses the failed cut, replacing the offending portion by the appropriate tree/path suffix gives an available route with no greater length. Conversely, every replacement path has a first edge that crosses from the source-side block to the target-side block; its prefix and suffix are no shorter than the two Dijkstra distances. Taking minima therefore gives the exact answer. Range chmin on $[a,b)$ processes all failures together.

The exact logical ids in path.edges are excluded as detour edges: a failed edge cannot replace itself. A parallel edge has a different id and remains a valid detour, even when its endpoints and cost are equal.

For a vertex $p_i$, an edge between blocks $a<b$ can similarly participate in a bypass only for internal indices $a<i<b$, the integer range $[a+1,b)`, when its prefix and suffix avoid $p_i$. Unlike edge failure, these one-edge ranges are not complete: an optimal vertex-avoiding route can enter the failed vertex’s tree block through one detour edge and leave through another without visiting the failed vertex. The current vertex solver therefore uses a proven Dijkstra run that ignores each path vertex rather than making an incorrect near-linear claim.

Preconditions

Assertions check that:

Directed graphs, graphs containing active zero-cost edges, and replacement-path reconstruction are not supported.

Examples

Automatic path selection:

#include "graph/replacement_paths.hpp"

int main() {
    int n = 5;
    m1une::graph::Graph<long long> g(n);
    g.add_edge(0, 1, 2);
    g.add_edge(1, 4, 2);
    g.add_edge(0, 2, 3);
    g.add_edge(2, 4, 3);

    auto edge_result = m1une::graph::edge_replacement_paths(g, 0, 4);
    for (long long distance : edge_result.replacement_dist) {
        (void)distance;
    }
}

An exact externally supplied shortest path:

#include "graph/replacement_paths.hpp"

int main() {
    m1une::graph::Graph<long long> g(4);
    int e01 = g.add_edge(0, 1, 2);
    int e13 = g.add_edge(1, 3, 2);
    g.add_edge(0, 2, 3);
    g.add_edge(2, 3, 3);

    m1une::graph::GraphPath fixed_path{
        .vertices = {0, 1, 3},
        .edges = {e01, e13},
    };
    auto vertex_result = m1une::graph::vertex_replacement_paths(g, fixed_path);
    (void)vertex_result;
}

Depends on

Required by

Verified with

Code

#ifndef M1UNE_GRAPH_REPLACEMENT_PATHS_HPP
#define M1UNE_GRAPH_REPLACEMENT_PATHS_HPP 1

#include <algorithm>
#include <cassert>
#include <functional>
#include <limits>
#include <queue>
#include <utility>
#include <vector>

#include "dijkstra.hpp"
#include "graph.hpp"

namespace m1une {
namespace graph {

struct GraphPath {
    std::vector<int> vertices;
    std::vector<int> edges;
};

template <class T>
struct EdgeReplacementPathsResult {
    GraphPath path;
    std::vector<T> replacement_dist;
    T inf;

    bool reachable(int path_edge_index) const {
        assert(0 <= path_edge_index && path_edge_index < int(replacement_dist.size()));
        return replacement_dist[path_edge_index] != inf;
    }
};

template <class T>
struct VertexReplacementPathsResult {
    GraphPath path;
    std::vector<T> replacement_dist;
    T inf;

    bool reachable(int path_vertex_index) const {
        assert(0 <= path_vertex_index && path_vertex_index < int(replacement_dist.size()));
        return replacement_dist[path_vertex_index] != inf;
    }
};

namespace internal {

template <class T>
T replacement_paths_safe_add(T a, T b, T inf) {
    if (a >= inf || b >= inf) return inf;
    if (a > inf - b) return inf;
    return a + b;
}

template <class T>
DijkstraResult<T> replacement_paths_dijkstra(const Graph<T>& g, int s, T inf) {
    int n = g.size();
    assert(0 <= s && s < n);
    DijkstraResult<T> result;
    result.dist.assign(n, inf);
    result.reached.assign(n, false);
    result.parent.assign(n, -1);
    result.parent_edge.assign(n, -1);
    result.inf = inf;

    using P = std::pair<T, int>;
    std::priority_queue<P, std::vector<P>, std::greater<P>> que;
    result.dist[s] = T(0);
    result.reached[s] = true;
    que.emplace(T(0), s);
    while (!que.empty()) {
        auto [d, v] = que.top();
        que.pop();
        if (result.dist[v] != d) continue;
        for (const auto& e : g[v]) {
            if (!e.alive) continue;
            T nd = replacement_paths_safe_add(d, e.cost, inf);
            if (result.dist[e.to] <= nd) continue;
            result.reached[e.to] = true;
            result.dist[e.to] = nd;
            result.parent[e.to] = v;
            result.parent_edge[e.to] = e.id;
            que.emplace(nd, e.to);
        }
    }
    return result;
}

template <class T>
std::vector<Edge<T>> replacement_paths_validate_graph(const Graph<T>& g, T inf) {
    assert(T(0) < inf);
    std::vector<int> occurrence(g.edge_count(), 0);
    std::vector<Edge<T>> edge_by_id(g.edge_count());
    for (int v = 0; v < g.size(); v++) {
        for (const auto& e : g[v]) {
            assert(e.from == v);
            assert(0 <= e.to && e.to < g.size());
            assert(0 <= e.id && e.id < g.edge_count());
            if (e.alive) assert(T(0) < e.cost);
            if (occurrence[e.id] == 0) {
                edge_by_id[e.id] = e;
            } else {
                assert(occurrence[e.id] == 1);
                const auto& other = edge_by_id[e.id];
                assert(e.from == other.to && e.to == other.from);
                assert(e.cost == other.cost && e.alive == other.alive);
            }
            occurrence[e.id]++;
        }
    }

    for (int id = 0; id < g.edge_count(); id++) {
        // add_edge creates exactly two mutually reversed arcs with one logical id.
        assert(occurrence[id] == 2);
    }
    return edge_by_id;
}

template <class T>
void replacement_paths_validate_path(const Graph<T>& g, const GraphPath& path,
                                     const std::vector<Edge<T>>& edge_by_id,
                                     const DijkstraResult<T>& from_s, T inf) {
    assert(!path.vertices.empty());
    assert(path.edges.size() + 1 == path.vertices.size());
    std::vector<char> used_vertex(g.size(), false);
    for (int v : path.vertices) {
        assert(0 <= v && v < g.size());
        assert(!used_vertex[v]);
        used_vertex[v] = true;
    }

    T path_cost = T(0);
    for (int i = 0; i < int(path.edges.size()); i++) {
        int id = path.edges[i];
        assert(0 <= id && id < g.edge_count());
        assert(g.is_edge_alive(id));
        const auto& e = edge_by_id[id];
        int u = path.vertices[i];
        int v = path.vertices[i + 1];
        assert((e.from == u && e.to == v) || (e.from == v && e.to == u));
        assert(T(0) < e.cost);
        path_cost = replacement_paths_safe_add(path_cost, e.cost, inf);
    }
    assert(from_s.reachable(path.vertices.back()));
    assert(path_cost == from_s.dist[path.vertices.back()]);
}

template <class T>
GraphPath replacement_paths_restore_path(const DijkstraResult<T>& result, int s, int t) {
    assert(result.reachable(t));
    GraphPath path;
    for (int v = t; v != s; v = result.parent[v]) {
        assert(v != -1 && result.parent[v] != -1 && result.parent_edge[v] != -1);
        path.vertices.push_back(v);
        path.edges.push_back(result.parent_edge[v]);
    }
    path.vertices.push_back(s);
    std::reverse(path.vertices.begin(), path.vertices.end());
    std::reverse(path.edges.begin(), path.edges.end());
    return path;
}

template <class T>
struct ReplacementPathsData {
    GraphPath path;
    std::vector<T> dist_s;
    std::vector<T> dist_t;
    std::vector<int> block;
    std::vector<char> is_path_edge;
    std::vector<Edge<T>> edge_by_id;
    T inf;
};

template <class T>
ReplacementPathsData<T> replacement_paths_prepare(const Graph<T>& g, const GraphPath& path,
                                                   T inf, const DijkstraResult<T>* known_from_s) {
    auto edge_by_id = replacement_paths_validate_graph(g, inf);
    int s = path.vertices.front();
    int t = path.vertices.back();
    auto computed_from_s = known_from_s == nullptr
                               ? replacement_paths_dijkstra(g, s, inf)
                               : DijkstraResult<T>();
    const auto& from_s = known_from_s == nullptr ? computed_from_s : *known_from_s;
    replacement_paths_validate_path(g, path, edge_by_id, from_s, inf);
    auto from_t = replacement_paths_dijkstra(g, t, inf);

    int n = g.size();
    std::vector<int> path_position(n, -1);
    std::vector<char> is_path_edge(g.edge_count(), false);
    for (int i = 0; i < int(path.vertices.size()); i++) path_position[path.vertices[i]] = i;
    for (int id : path.edges) is_path_edge[id] = true;

    std::vector<int> parent(n, -1);
    for (int i = 0; i < int(path.edges.size()); i++) {
        int v = path.vertices[i + 1];
        parent[v] = path.vertices[i];
        const auto& e = edge_by_id[path.edges[i]];
        assert(replacement_paths_safe_add(from_s.dist[parent[v]], e.cost, inf) == from_s.dist[v]);
    }
    for (int v = 0; v < n; v++) {
        if (!from_s.reachable(v) || v == s || path_position[v] != -1) continue;
        for (const auto& e : g[v]) {
            if (!e.alive || !from_s.reachable(e.to)) continue;
            if (replacement_paths_safe_add(from_s.dist[e.to], e.cost, inf) != from_s.dist[v]) {
                continue;
            }
            parent[v] = e.to;
            break;
        }
        assert(parent[v] != -1);
        assert(from_s.dist[parent[v]] < from_s.dist[v]);
    }

    std::vector<std::vector<int>> children(n);
    for (int v = 0; v < n; v++) {
        if (parent[v] != -1) children[parent[v]].push_back(v);
    }
    std::vector<int> block(n, -1);
    block[s] = 0;
    std::vector<int> stack = {s};
    while (!stack.empty()) {
        int v = stack.back();
        stack.pop_back();
        for (int to : children[v]) {
            block[to] = path_position[to] == -1 ? block[v] : path_position[to];
            stack.push_back(to);
        }
    }
    for (int v = 0; v < n; v++) assert(!from_s.reachable(v) || block[v] != -1);

    return {path, from_s.dist, from_t.dist, block, is_path_edge, edge_by_id, inf};
}

template <class T>
class ReplacementPathsRangeChmin {
   private:
    int _size;
    std::vector<T> _lazy;

   public:
    ReplacementPathsRangeChmin(int n, T inf) : _size(1) {
        while (_size < n) _size <<= 1;
        _lazy.assign(2 * _size, inf);
    }

    void apply(int l, int r, T value) {
        assert(0 <= l && l <= r && r <= _size);
        for (l += _size, r += _size; l < r; l >>= 1, r >>= 1) {
            if (l & 1) {
                _lazy[l] = std::min(_lazy[l], value);
                l++;
            }
            if (r & 1) {
                --r;
                _lazy[r] = std::min(_lazy[r], value);
            }
        }
    }

    std::vector<T> values(int n) {
        for (int v = 1; v < _size; v++) {
            _lazy[2 * v] = std::min(_lazy[2 * v], _lazy[v]);
            _lazy[2 * v + 1] = std::min(_lazy[2 * v + 1], _lazy[v]);
        }
        return std::vector<T>(_lazy.begin() + _size, _lazy.begin() + _size + n);
    }
};

template <class T>
std::vector<T> replacement_paths_solve_edges(const ReplacementPathsData<T>& data) {
    int answer_size = int(data.path.edges.size());
    ReplacementPathsRangeChmin<T> range_chmin(answer_size, data.inf);
    for (const auto& e : data.edge_by_id) {
        if (!e.alive || data.is_path_edge[e.id]) continue;
        int u = e.from;
        int v = e.to;
        if (data.block[u] == -1 || data.block[v] == -1 || data.block[u] == data.block[v]) continue;
        if (data.block[u] > data.block[v]) std::swap(u, v);
        int a = data.block[u];
        int b = data.block[v];
        T candidate = replacement_paths_safe_add(data.dist_s[u], e.cost, data.inf);
        candidate = replacement_paths_safe_add(candidate, data.dist_t[v], data.inf);
        if (candidate == data.inf) continue;
        range_chmin.apply(a, b, candidate);
    }
    return range_chmin.values(answer_size);
}

template <class T>
T replacement_paths_without_vertex(const Graph<T>& g, int s, int t, int removed, T inf) {
    if (s == removed || t == removed) return inf;
    std::vector<T> dist(g.size(), inf);
    using P = std::pair<T, int>;
    std::priority_queue<P, std::vector<P>, std::greater<P>> que;
    dist[s] = T(0);
    que.emplace(T(0), s);
    while (!que.empty()) {
        auto [d, v] = que.top();
        que.pop();
        if (dist[v] != d) continue;
        for (const auto& e : g[v]) {
            if (!e.alive || e.to == removed) continue;
            T nd = replacement_paths_safe_add(d, e.cost, inf);
            if (dist[e.to] <= nd) continue;
            dist[e.to] = nd;
            que.emplace(nd, e.to);
        }
    }
    return dist[t];
}

template <class T>
std::vector<T> replacement_paths_solve_vertices(const Graph<T>& g,
                                                const ReplacementPathsData<T>& data) {
    // One edge can cross an edge cut, but a vertex-avoiding path may enter and
    // leave the failed vertex's tree block through two different detour edges.
    int path_size = int(data.path.vertices.size());
    std::vector<T> answer(path_size, data.inf);
    int s = data.path.vertices.front();
    int t = data.path.vertices.back();
    for (int i = 1; i + 1 < path_size; i++) {
        answer[i] = replacement_paths_without_vertex(g, s, t, data.path.vertices[i], data.inf);
    }
    return answer;
}

}  // namespace internal

template <class T>
EdgeReplacementPathsResult<T> edge_replacement_paths(
    const Graph<T>& g, const GraphPath& path, T inf = std::numeric_limits<T>::max() / T(4)) {
    assert(!path.vertices.empty());
    auto data = internal::replacement_paths_prepare(
        g, path, inf, static_cast<const DijkstraResult<T>*>(nullptr));
    auto replacement_dist = internal::replacement_paths_solve_edges(data);
    return {path, replacement_dist, inf};
}

template <class T>
EdgeReplacementPathsResult<T> edge_replacement_paths(
    const Graph<T>& g, int s, int t, T inf = std::numeric_limits<T>::max() / T(4)) {
    assert(0 <= s && s < g.size());
    assert(0 <= t && t < g.size());
    auto from_s = internal::replacement_paths_dijkstra(g, s, inf);
    assert(from_s.reachable(t));
    auto path = internal::replacement_paths_restore_path(from_s, s, t);
    auto data = internal::replacement_paths_prepare(g, path, inf, &from_s);
    auto replacement_dist = internal::replacement_paths_solve_edges(data);
    return {path, replacement_dist, inf};
}

template <class T>
VertexReplacementPathsResult<T> vertex_replacement_paths(
    const Graph<T>& g, const GraphPath& path, T inf = std::numeric_limits<T>::max() / T(4)) {
    assert(!path.vertices.empty());
    auto data = internal::replacement_paths_prepare(
        g, path, inf, static_cast<const DijkstraResult<T>*>(nullptr));
    auto replacement_dist = internal::replacement_paths_solve_vertices(g, data);
    return {path, replacement_dist, inf};
}

template <class T>
VertexReplacementPathsResult<T> vertex_replacement_paths(
    const Graph<T>& g, int s, int t, T inf = std::numeric_limits<T>::max() / T(4)) {
    assert(0 <= s && s < g.size());
    assert(0 <= t && t < g.size());
    auto from_s = internal::replacement_paths_dijkstra(g, s, inf);
    assert(from_s.reachable(t));
    auto path = internal::replacement_paths_restore_path(from_s, s, t);
    auto data = internal::replacement_paths_prepare(g, path, inf, &from_s);
    auto replacement_dist = internal::replacement_paths_solve_vertices(g, data);
    return {path, replacement_dist, inf};
}

}  // namespace graph
}  // namespace m1une

#endif  // M1UNE_GRAPH_REPLACEMENT_PATHS_HPP
#line 1 "graph/replacement_paths.hpp"



#include <algorithm>
#include <cassert>
#include <functional>
#include <limits>
#include <queue>
#include <utility>
#include <vector>

#line 1 "graph/dijkstra.hpp"



#line 8 "graph/dijkstra.hpp"

#line 1 "graph/graph.hpp"



#include <array>
#line 8 "graph/graph.hpp"

namespace m1une {
namespace graph {

template <class T = int>
struct Edge {
    using cost_type = T;

    int from;
    int to;
    T cost;
    int id;
    bool alive;

    Edge() : from(-1), to(-1), cost(T()), id(-1), alive(true) {}
    Edge(int from_, int to_, T cost_ = T(1), int id_ = -1, bool alive_ = true)
        : from(from_), to(to_), cost(cost_), id(id_), alive(alive_) {}

    int other(int v) const {
        assert(v == from || v == to);
        return from ^ to ^ v;
    }
};

template <class T = int>
struct Graph {
    using edge_type = Edge<T>;
    using cost_type = T;

   private:
    struct EdgePositions {
        std::array<std::pair<int, int>, 2> value{};
        int size = 0;

        void push_back(std::pair<int, int> position) {
            assert(size < 2);
            value[size++] = position;
        }
    };

    int _n;
    int _edge_count;
    std::vector<std::vector<edge_type>> _g;
    std::vector<EdgePositions> _edge_positions;

   public:
    Graph() : _n(0), _edge_count(0) {}
    explicit Graph(int n) : _n(n), _edge_count(0), _g(n) {
        assert(0 <= n);
    }

    int size() const {
        return _n;
    }

    bool empty() const {
        return _n == 0;
    }

    int edge_count() const {
        return _edge_count;
    }

    int add_vertex() {
        _g.emplace_back();
        return _n++;
    }

    int add_directed_edge(int from, int to, T cost = T(1)) {
        assert(0 <= from && from < _n);
        assert(0 <= to && to < _n);
        int id = _edge_count++;
        int idx = int(_g[from].size());
        _g[from].push_back(edge_type(from, to, cost, id));
        _edge_positions.emplace_back();
        _edge_positions.back().push_back({from, idx});
        return id;
    }

    int add_edge(int u, int v, T cost = T(1)) {
        assert(0 <= u && u < _n);
        assert(0 <= v && v < _n);
        int id = _edge_count++;
        int u_idx = int(_g[u].size());
        _g[u].push_back(edge_type(u, v, cost, id));
        int v_idx = int(_g[v].size());
        _g[v].push_back(edge_type(v, u, cost, id));
        _edge_positions.emplace_back();
        _edge_positions.back().push_back({u, u_idx});
        _edge_positions.back().push_back({v, v_idx});
        return id;
    }

    void set_edge_alive(int id, bool alive) {
        assert(0 <= id && id < _edge_count);
        for (int i = 0; i < _edge_positions[id].size; ++i) {
            auto [v, idx] = _edge_positions[id].value[i];
            _g[v][idx].alive = alive;
        }
    }

    void erase_edge(int id) {
        set_edge_alive(id, false);
    }

    void revive_edge(int id) {
        set_edge_alive(id, true);
    }

    bool is_edge_alive(int id) const {
        assert(0 <= id && id < _edge_count);
        assert(_edge_positions[id].size != 0);
        auto [v, idx] = _edge_positions[id].value[0];
        return _g[v][idx].alive;
    }

    const std::vector<edge_type>& operator[](int v) const {
        assert(0 <= v && v < _n);
        return _g[v];
    }

    std::vector<edge_type>& operator[](int v) {
        assert(0 <= v && v < _n);
        return _g[v];
    }

    const std::vector<std::vector<edge_type>>& adjacency() const {
        return _g;
    }

    std::vector<std::vector<edge_type>>& adjacency() {
        return _g;
    }

    std::vector<edge_type> edges(bool include_inactive = false) const {
        std::vector<edge_type> result;
        result.reserve(_edge_count);
        std::vector<char> used(_edge_count, false);
        for (int v = 0; v < _n; v++) {
            for (const auto& e : _g[v]) {
                if (!include_inactive && !e.alive) continue;
                if (0 <= e.id && e.id < _edge_count) {
                    if (used[e.id]) continue;
                    used[e.id] = true;
                }
                result.push_back(e);
            }
        }
        return result;
    }

    Graph reversed() const {
        Graph result(_n);
        result._edge_count = _edge_count;
        result._edge_positions.assign(_edge_count, {});
        for (int v = 0; v < _n; v++) {
            for (const auto& e : _g[v]) {
                int idx = int(result._g[e.to].size());
                result._g[e.to].push_back(edge_type(e.to, e.from, e.cost, e.id, e.alive));
                if (0 <= e.id && e.id < _edge_count) result._edge_positions[e.id].push_back({e.to, idx});
            }
        }
        return result;
    }
};

}  // namespace graph
}  // namespace m1une


#line 10 "graph/dijkstra.hpp"

namespace m1une {
namespace graph {

template <class T>
struct DijkstraResult {
    std::vector<T> dist;
    std::vector<char> reached;
    std::vector<int> parent;
    std::vector<int> parent_edge;
    T inf = T();

    bool reachable(int v) const {
        assert(0 <= v && v < int(dist.size()));
        return reached[v];
    }

    std::vector<int> path(int t) const {
        assert(reachable(t));
        std::vector<int> result;
        for (int v = t; v != -1; v = parent[v]) result.push_back(v);
        std::reverse(result.begin(), result.end());
        return result;
    }
};

namespace internal {

template <class T>
class DijkstraHeap {
   private:
    const std::vector<T>& dist_;
    std::vector<int> heap_;
    std::vector<int> position_;

    bool less(int first, int second) const {
        return dist_[heap_[first]] < dist_[heap_[second]];
    }

    void swap_nodes(int first, int second) {
        std::swap(heap_[first], heap_[second]);
        position_[heap_[first]] = first;
        position_[heap_[second]] = second;
    }

    void sift_up(int index) {
        while (index != 0) {
            const int parent = (index - 1) / 2;
            if (!less(index, parent)) break;
            swap_nodes(index, parent);
            index = parent;
        }
    }

    void sift_down(int index) {
        while (2 * index + 1 < int(heap_.size())) {
            int child = 2 * index + 1;
            if (child + 1 < int(heap_.size()) && less(child + 1, child)) {
                ++child;
            }
            if (!less(child, index)) break;
            swap_nodes(index, child);
            index = child;
        }
    }

   public:
    DijkstraHeap(const std::vector<T>& dist, int size)
        : dist_(dist), position_(size, -1) {
        heap_.reserve(size);
    }

    bool empty() const {
        return heap_.empty();
    }

    void push_or_decrease(int vertex) {
        int& position = position_[vertex];
        if (position == -1) {
            position = int(heap_.size());
            heap_.push_back(vertex);
        }
        sift_up(position);
    }

    int pop_min() {
        const int result = heap_.front();
        position_[result] = -1;
        if (heap_.size() == 1) {
            heap_.pop_back();
            return result;
        }
        heap_.front() = heap_.back();
        position_[heap_.front()] = 0;
        heap_.pop_back();
        sift_down(0);
        return result;
    }
};

}  // namespace internal

template <class T>
DijkstraResult<T> dijkstra(const Graph<T>& g,
                           const std::vector<int>& sources) {
    int n = g.size();
    DijkstraResult<T> result;
    result.dist.resize(n);
    result.reached.assign(n, false);
    result.parent.assign(n, -1);
    result.parent_edge.assign(n, -1);

    internal::DijkstraHeap<T> que(result.dist, n);
    for (int s : sources) {
        assert(0 <= s && s < n);
        if (result.reached[s]) continue;
        result.reached[s] = true;
        result.dist[s] = T();
        que.push_or_decrease(s);
    }

    while (!que.empty()) {
        const int current = que.pop_min();
        for (const auto& e : g[current]) {
            if (!e.alive) continue;
            T nd = result.dist[current] + e.cost;
            if (result.reached[e.to] && !(nd < result.dist[e.to])) continue;
            result.reached[e.to] = true;
            result.dist[e.to] = std::move(nd);
            result.parent[e.to] = current;
            result.parent_edge[e.to] = e.id;
            que.push_or_decrease(e.to);
        }
    }

    return result;
}

template <class T>
DijkstraResult<T> dijkstra(const Graph<T>& g, int s) {
    return dijkstra(g, std::vector<int>{s});
}

// Compatibility overload: unreachable distances are replaced by inf after the
// search. Reachability itself never depends on this sentinel.
template <class T>
DijkstraResult<T> dijkstra(const Graph<T>& g,
                           const std::vector<int>& sources, const T& inf) {
    DijkstraResult<T> result = dijkstra(g, sources);
    result.inf = inf;
    for (int v = 0; v < int(result.dist.size()); v++) {
        if (!result.reachable(v)) result.dist[v] = inf;
    }
    return result;
}

template <class T>
DijkstraResult<T> dijkstra(const Graph<T>& g, int s, const T& inf) {
    return dijkstra(g, std::vector<int>{s}, inf);
}

}  // namespace graph
}  // namespace m1une


#line 14 "graph/replacement_paths.hpp"

namespace m1une {
namespace graph {

struct GraphPath {
    std::vector<int> vertices;
    std::vector<int> edges;
};

template <class T>
struct EdgeReplacementPathsResult {
    GraphPath path;
    std::vector<T> replacement_dist;
    T inf;

    bool reachable(int path_edge_index) const {
        assert(0 <= path_edge_index && path_edge_index < int(replacement_dist.size()));
        return replacement_dist[path_edge_index] != inf;
    }
};

template <class T>
struct VertexReplacementPathsResult {
    GraphPath path;
    std::vector<T> replacement_dist;
    T inf;

    bool reachable(int path_vertex_index) const {
        assert(0 <= path_vertex_index && path_vertex_index < int(replacement_dist.size()));
        return replacement_dist[path_vertex_index] != inf;
    }
};

namespace internal {

template <class T>
T replacement_paths_safe_add(T a, T b, T inf) {
    if (a >= inf || b >= inf) return inf;
    if (a > inf - b) return inf;
    return a + b;
}

template <class T>
DijkstraResult<T> replacement_paths_dijkstra(const Graph<T>& g, int s, T inf) {
    int n = g.size();
    assert(0 <= s && s < n);
    DijkstraResult<T> result;
    result.dist.assign(n, inf);
    result.reached.assign(n, false);
    result.parent.assign(n, -1);
    result.parent_edge.assign(n, -1);
    result.inf = inf;

    using P = std::pair<T, int>;
    std::priority_queue<P, std::vector<P>, std::greater<P>> que;
    result.dist[s] = T(0);
    result.reached[s] = true;
    que.emplace(T(0), s);
    while (!que.empty()) {
        auto [d, v] = que.top();
        que.pop();
        if (result.dist[v] != d) continue;
        for (const auto& e : g[v]) {
            if (!e.alive) continue;
            T nd = replacement_paths_safe_add(d, e.cost, inf);
            if (result.dist[e.to] <= nd) continue;
            result.reached[e.to] = true;
            result.dist[e.to] = nd;
            result.parent[e.to] = v;
            result.parent_edge[e.to] = e.id;
            que.emplace(nd, e.to);
        }
    }
    return result;
}

template <class T>
std::vector<Edge<T>> replacement_paths_validate_graph(const Graph<T>& g, T inf) {
    assert(T(0) < inf);
    std::vector<int> occurrence(g.edge_count(), 0);
    std::vector<Edge<T>> edge_by_id(g.edge_count());
    for (int v = 0; v < g.size(); v++) {
        for (const auto& e : g[v]) {
            assert(e.from == v);
            assert(0 <= e.to && e.to < g.size());
            assert(0 <= e.id && e.id < g.edge_count());
            if (e.alive) assert(T(0) < e.cost);
            if (occurrence[e.id] == 0) {
                edge_by_id[e.id] = e;
            } else {
                assert(occurrence[e.id] == 1);
                const auto& other = edge_by_id[e.id];
                assert(e.from == other.to && e.to == other.from);
                assert(e.cost == other.cost && e.alive == other.alive);
            }
            occurrence[e.id]++;
        }
    }

    for (int id = 0; id < g.edge_count(); id++) {
        // add_edge creates exactly two mutually reversed arcs with one logical id.
        assert(occurrence[id] == 2);
    }
    return edge_by_id;
}

template <class T>
void replacement_paths_validate_path(const Graph<T>& g, const GraphPath& path,
                                     const std::vector<Edge<T>>& edge_by_id,
                                     const DijkstraResult<T>& from_s, T inf) {
    assert(!path.vertices.empty());
    assert(path.edges.size() + 1 == path.vertices.size());
    std::vector<char> used_vertex(g.size(), false);
    for (int v : path.vertices) {
        assert(0 <= v && v < g.size());
        assert(!used_vertex[v]);
        used_vertex[v] = true;
    }

    T path_cost = T(0);
    for (int i = 0; i < int(path.edges.size()); i++) {
        int id = path.edges[i];
        assert(0 <= id && id < g.edge_count());
        assert(g.is_edge_alive(id));
        const auto& e = edge_by_id[id];
        int u = path.vertices[i];
        int v = path.vertices[i + 1];
        assert((e.from == u && e.to == v) || (e.from == v && e.to == u));
        assert(T(0) < e.cost);
        path_cost = replacement_paths_safe_add(path_cost, e.cost, inf);
    }
    assert(from_s.reachable(path.vertices.back()));
    assert(path_cost == from_s.dist[path.vertices.back()]);
}

template <class T>
GraphPath replacement_paths_restore_path(const DijkstraResult<T>& result, int s, int t) {
    assert(result.reachable(t));
    GraphPath path;
    for (int v = t; v != s; v = result.parent[v]) {
        assert(v != -1 && result.parent[v] != -1 && result.parent_edge[v] != -1);
        path.vertices.push_back(v);
        path.edges.push_back(result.parent_edge[v]);
    }
    path.vertices.push_back(s);
    std::reverse(path.vertices.begin(), path.vertices.end());
    std::reverse(path.edges.begin(), path.edges.end());
    return path;
}

template <class T>
struct ReplacementPathsData {
    GraphPath path;
    std::vector<T> dist_s;
    std::vector<T> dist_t;
    std::vector<int> block;
    std::vector<char> is_path_edge;
    std::vector<Edge<T>> edge_by_id;
    T inf;
};

template <class T>
ReplacementPathsData<T> replacement_paths_prepare(const Graph<T>& g, const GraphPath& path,
                                                   T inf, const DijkstraResult<T>* known_from_s) {
    auto edge_by_id = replacement_paths_validate_graph(g, inf);
    int s = path.vertices.front();
    int t = path.vertices.back();
    auto computed_from_s = known_from_s == nullptr
                               ? replacement_paths_dijkstra(g, s, inf)
                               : DijkstraResult<T>();
    const auto& from_s = known_from_s == nullptr ? computed_from_s : *known_from_s;
    replacement_paths_validate_path(g, path, edge_by_id, from_s, inf);
    auto from_t = replacement_paths_dijkstra(g, t, inf);

    int n = g.size();
    std::vector<int> path_position(n, -1);
    std::vector<char> is_path_edge(g.edge_count(), false);
    for (int i = 0; i < int(path.vertices.size()); i++) path_position[path.vertices[i]] = i;
    for (int id : path.edges) is_path_edge[id] = true;

    std::vector<int> parent(n, -1);
    for (int i = 0; i < int(path.edges.size()); i++) {
        int v = path.vertices[i + 1];
        parent[v] = path.vertices[i];
        const auto& e = edge_by_id[path.edges[i]];
        assert(replacement_paths_safe_add(from_s.dist[parent[v]], e.cost, inf) == from_s.dist[v]);
    }
    for (int v = 0; v < n; v++) {
        if (!from_s.reachable(v) || v == s || path_position[v] != -1) continue;
        for (const auto& e : g[v]) {
            if (!e.alive || !from_s.reachable(e.to)) continue;
            if (replacement_paths_safe_add(from_s.dist[e.to], e.cost, inf) != from_s.dist[v]) {
                continue;
            }
            parent[v] = e.to;
            break;
        }
        assert(parent[v] != -1);
        assert(from_s.dist[parent[v]] < from_s.dist[v]);
    }

    std::vector<std::vector<int>> children(n);
    for (int v = 0; v < n; v++) {
        if (parent[v] != -1) children[parent[v]].push_back(v);
    }
    std::vector<int> block(n, -1);
    block[s] = 0;
    std::vector<int> stack = {s};
    while (!stack.empty()) {
        int v = stack.back();
        stack.pop_back();
        for (int to : children[v]) {
            block[to] = path_position[to] == -1 ? block[v] : path_position[to];
            stack.push_back(to);
        }
    }
    for (int v = 0; v < n; v++) assert(!from_s.reachable(v) || block[v] != -1);

    return {path, from_s.dist, from_t.dist, block, is_path_edge, edge_by_id, inf};
}

template <class T>
class ReplacementPathsRangeChmin {
   private:
    int _size;
    std::vector<T> _lazy;

   public:
    ReplacementPathsRangeChmin(int n, T inf) : _size(1) {
        while (_size < n) _size <<= 1;
        _lazy.assign(2 * _size, inf);
    }

    void apply(int l, int r, T value) {
        assert(0 <= l && l <= r && r <= _size);
        for (l += _size, r += _size; l < r; l >>= 1, r >>= 1) {
            if (l & 1) {
                _lazy[l] = std::min(_lazy[l], value);
                l++;
            }
            if (r & 1) {
                --r;
                _lazy[r] = std::min(_lazy[r], value);
            }
        }
    }

    std::vector<T> values(int n) {
        for (int v = 1; v < _size; v++) {
            _lazy[2 * v] = std::min(_lazy[2 * v], _lazy[v]);
            _lazy[2 * v + 1] = std::min(_lazy[2 * v + 1], _lazy[v]);
        }
        return std::vector<T>(_lazy.begin() + _size, _lazy.begin() + _size + n);
    }
};

template <class T>
std::vector<T> replacement_paths_solve_edges(const ReplacementPathsData<T>& data) {
    int answer_size = int(data.path.edges.size());
    ReplacementPathsRangeChmin<T> range_chmin(answer_size, data.inf);
    for (const auto& e : data.edge_by_id) {
        if (!e.alive || data.is_path_edge[e.id]) continue;
        int u = e.from;
        int v = e.to;
        if (data.block[u] == -1 || data.block[v] == -1 || data.block[u] == data.block[v]) continue;
        if (data.block[u] > data.block[v]) std::swap(u, v);
        int a = data.block[u];
        int b = data.block[v];
        T candidate = replacement_paths_safe_add(data.dist_s[u], e.cost, data.inf);
        candidate = replacement_paths_safe_add(candidate, data.dist_t[v], data.inf);
        if (candidate == data.inf) continue;
        range_chmin.apply(a, b, candidate);
    }
    return range_chmin.values(answer_size);
}

template <class T>
T replacement_paths_without_vertex(const Graph<T>& g, int s, int t, int removed, T inf) {
    if (s == removed || t == removed) return inf;
    std::vector<T> dist(g.size(), inf);
    using P = std::pair<T, int>;
    std::priority_queue<P, std::vector<P>, std::greater<P>> que;
    dist[s] = T(0);
    que.emplace(T(0), s);
    while (!que.empty()) {
        auto [d, v] = que.top();
        que.pop();
        if (dist[v] != d) continue;
        for (const auto& e : g[v]) {
            if (!e.alive || e.to == removed) continue;
            T nd = replacement_paths_safe_add(d, e.cost, inf);
            if (dist[e.to] <= nd) continue;
            dist[e.to] = nd;
            que.emplace(nd, e.to);
        }
    }
    return dist[t];
}

template <class T>
std::vector<T> replacement_paths_solve_vertices(const Graph<T>& g,
                                                const ReplacementPathsData<T>& data) {
    // One edge can cross an edge cut, but a vertex-avoiding path may enter and
    // leave the failed vertex's tree block through two different detour edges.
    int path_size = int(data.path.vertices.size());
    std::vector<T> answer(path_size, data.inf);
    int s = data.path.vertices.front();
    int t = data.path.vertices.back();
    for (int i = 1; i + 1 < path_size; i++) {
        answer[i] = replacement_paths_without_vertex(g, s, t, data.path.vertices[i], data.inf);
    }
    return answer;
}

}  // namespace internal

template <class T>
EdgeReplacementPathsResult<T> edge_replacement_paths(
    const Graph<T>& g, const GraphPath& path, T inf = std::numeric_limits<T>::max() / T(4)) {
    assert(!path.vertices.empty());
    auto data = internal::replacement_paths_prepare(
        g, path, inf, static_cast<const DijkstraResult<T>*>(nullptr));
    auto replacement_dist = internal::replacement_paths_solve_edges(data);
    return {path, replacement_dist, inf};
}

template <class T>
EdgeReplacementPathsResult<T> edge_replacement_paths(
    const Graph<T>& g, int s, int t, T inf = std::numeric_limits<T>::max() / T(4)) {
    assert(0 <= s && s < g.size());
    assert(0 <= t && t < g.size());
    auto from_s = internal::replacement_paths_dijkstra(g, s, inf);
    assert(from_s.reachable(t));
    auto path = internal::replacement_paths_restore_path(from_s, s, t);
    auto data = internal::replacement_paths_prepare(g, path, inf, &from_s);
    auto replacement_dist = internal::replacement_paths_solve_edges(data);
    return {path, replacement_dist, inf};
}

template <class T>
VertexReplacementPathsResult<T> vertex_replacement_paths(
    const Graph<T>& g, const GraphPath& path, T inf = std::numeric_limits<T>::max() / T(4)) {
    assert(!path.vertices.empty());
    auto data = internal::replacement_paths_prepare(
        g, path, inf, static_cast<const DijkstraResult<T>*>(nullptr));
    auto replacement_dist = internal::replacement_paths_solve_vertices(g, data);
    return {path, replacement_dist, inf};
}

template <class T>
VertexReplacementPathsResult<T> vertex_replacement_paths(
    const Graph<T>& g, int s, int t, T inf = std::numeric_limits<T>::max() / T(4)) {
    assert(0 <= s && s < g.size());
    assert(0 <= t && t < g.size());
    auto from_s = internal::replacement_paths_dijkstra(g, s, inf);
    assert(from_s.reachable(t));
    auto path = internal::replacement_paths_restore_path(from_s, s, t);
    auto data = internal::replacement_paths_prepare(g, path, inf, &from_s);
    auto replacement_dist = internal::replacement_paths_solve_vertices(g, data);
    return {path, replacement_dist, inf};
}

}  // namespace graph
}  // namespace m1une
Back to top page