weavatrix_graph/algo/
shortest_extra.rs

1use super::{BellmanFord, WeightedPath};
2use crate::Vec;
3use crate::{GraphError, IndexGraphView, Result};
4use alloc::collections::{BinaryHeap, VecDeque};
5use core::cmp::Reverse;
6
7pub fn bidirectional_dijkstra<G, F>(
8    graph: &G,
9    source: G::Node,
10    target: G::Node,
11    edge_cost: F,
12) -> Option<WeightedPath<G::Node>>
13where
14    G: IndexGraphView,
15    F: Fn(G::Edge) -> u64,
16{
17    if !graph.contains_node(source) || !graph.contains_node(target) {
18        return None;
19    }
20    if source == target {
21        return Some(WeightedPath::from_parts(vec![source], 0));
22    }
23    let bound = graph.node_bound();
24    let mut by_slot = vec![None; bound];
25    for node in graph.node_indices() {
26        by_slot[G::node_slot(node)] = Some(node);
27    }
28    let mut forward = SearchSide::new(bound, G::node_slot(source));
29    let mut backward = SearchSide::new(bound, G::node_slot(target));
30    let mut best = None;
31    while !forward.queue.is_empty() && !backward.queue.is_empty() {
32        let forward_min = forward.queue.peek().map_or(u64::MAX, |entry| entry.0.0);
33        let backward_min = backward.queue.peek().map_or(u64::MAX, |entry| entry.0.0);
34        if best.is_some_and(|(cost, _)| forward_min.saturating_add(backward_min) >= cost) {
35            break;
36        }
37        if forward_min <= backward_min {
38            expand(
39                graph,
40                true,
41                &edge_cost,
42                &by_slot,
43                &mut forward,
44                &backward,
45                &mut best,
46            );
47        } else {
48            expand(
49                graph,
50                false,
51                &edge_cost,
52                &by_slot,
53                &mut backward,
54                &forward,
55                &mut best,
56            );
57        }
58    }
59    let (cost, meeting) = best?;
60    let nodes = reconstruct::<G>(
61        source,
62        target,
63        meeting,
64        &by_slot,
65        &forward.parent,
66        &backward.parent,
67    )?;
68    Some(WeightedPath::from_parts(nodes, cost))
69}
70
71/// Computes signed single-source paths with the queue-based SPFA algorithm.
72///
73/// # Errors
74///
75/// Returns an error for arithmetic overflow or a reachable negative cycle.
76pub fn spfa<G, F>(graph: &G, source: G::Node, edge_cost: F) -> Result<Option<BellmanFord<G::Node>>>
77where
78    G: IndexGraphView,
79    F: Fn(G::Edge) -> i64,
80{
81    spfa_filtered(graph, source, |edge| Some(edge_cost(edge)))
82}
83
84/// Computes SPFA while omitting edges whose cost is `None`.
85///
86/// # Errors
87///
88/// Returns an error for arithmetic overflow or a reachable negative cycle.
89pub fn spfa_filtered<G, F>(
90    graph: &G,
91    source: G::Node,
92    edge_cost: F,
93) -> Result<Option<BellmanFord<G::Node>>>
94where
95    G: IndexGraphView,
96    F: Fn(G::Edge) -> Option<i64>,
97{
98    if !graph.contains_node(source) {
99        return Ok(None);
100    }
101    let nodes = graph.node_indices().collect::<Vec<_>>();
102    let mut by_slot = vec![None; graph.node_bound()];
103    for &node in &nodes {
104        by_slot[G::node_slot(node)] = Some(node);
105    }
106    let mut distance = vec![0_i64; graph.node_bound()];
107    let mut reachable = vec![false; graph.node_bound()];
108    let mut parent = vec![None; graph.node_bound()];
109    let mut queued = vec![false; graph.node_bound()];
110    let mut relaxations = vec![0_usize; graph.node_bound()];
111    let source_slot = G::node_slot(source);
112    reachable[source_slot] = true;
113    queued[source_slot] = true;
114    let mut queue = VecDeque::from([source]);
115    while let Some(node) = queue.pop_front() {
116        let node_slot = G::node_slot(node);
117        queued[node_slot] = false;
118        for edge in graph.outgoing_edges(node) {
119            let Some(weight) = edge_cost(edge) else {
120                continue;
121            };
122            let Some(endpoints) = graph.edge_endpoints(edge) else {
123                continue;
124            };
125            let target = endpoints.target();
126            let target_slot = G::node_slot(target);
127            let candidate =
128                distance[node_slot]
129                    .checked_add(weight)
130                    .ok_or(GraphError::ArithmeticOverflow {
131                        operation: "SPFA relaxation",
132                    })?;
133            if !reachable[target_slot] || candidate < distance[target_slot] {
134                reachable[target_slot] = true;
135                distance[target_slot] = candidate;
136                parent[target_slot] = Some(node_slot);
137                relaxations[target_slot] += 1;
138                if relaxations[target_slot] >= nodes.len() {
139                    return Err(GraphError::NegativeCycle { algorithm: "SPFA" });
140                }
141                if !queued[target_slot] {
142                    queued[target_slot] = true;
143                    queue.push_back(target);
144                }
145            }
146        }
147    }
148    Ok(Some(BellmanFord::from_parts(
149        source,
150        nodes,
151        by_slot,
152        distance,
153        reachable,
154        parent,
155        G::node_slot,
156    )))
157}
158
159struct SearchSide {
160    distance: Vec<u64>,
161    parent: Vec<Option<usize>>,
162    settled: Vec<bool>,
163    queue: BinaryHeap<Reverse<(u64, usize)>>,
164}
165
166impl SearchSide {
167    fn new(bound: usize, start: usize) -> Self {
168        let mut distance = vec![u64::MAX; bound];
169        distance[start] = 0;
170        Self {
171            distance,
172            parent: vec![None; bound],
173            settled: vec![false; bound],
174            queue: BinaryHeap::from([Reverse((0, start))]),
175        }
176    }
177}
178
179#[allow(clippy::too_many_arguments)]
180fn expand<G, F>(
181    graph: &G,
182    outgoing: bool,
183    edge_cost: &F,
184    by_slot: &[Option<G::Node>],
185    side: &mut SearchSide,
186    other: &SearchSide,
187    best: &mut Option<(u64, usize)>,
188) where
189    G: IndexGraphView,
190    F: Fn(G::Edge) -> u64,
191{
192    let Some(Reverse((cost, slot))) = side.queue.pop() else {
193        return;
194    };
195    if cost != side.distance[slot] || side.settled[slot] {
196        return;
197    }
198    side.settled[slot] = true;
199    update_best(slot, side, other, best);
200    let Some(node) = by_slot[slot] else {
201        return;
202    };
203    let edges = if outgoing {
204        graph.outgoing_edges(node).collect::<Vec<_>>()
205    } else {
206        graph.incoming_edges(node).collect::<Vec<_>>()
207    };
208    for edge in edges {
209        let Some(endpoints) = graph.edge_endpoints(edge) else {
210            continue;
211        };
212        let neighbor = if outgoing {
213            endpoints.target()
214        } else {
215            endpoints.source()
216        };
217        let neighbor_slot = G::node_slot(neighbor);
218        let Some(candidate) = cost.checked_add(edge_cost(edge)) else {
219            continue;
220        };
221        if candidate < side.distance[neighbor_slot] {
222            side.distance[neighbor_slot] = candidate;
223            side.parent[neighbor_slot] = Some(slot);
224            side.queue.push(Reverse((candidate, neighbor_slot)));
225            update_best(neighbor_slot, side, other, best);
226        }
227    }
228}
229
230fn update_best(
231    slot: usize,
232    side: &SearchSide,
233    other: &SearchSide,
234    best: &mut Option<(u64, usize)>,
235) {
236    if side.distance[slot] == u64::MAX || other.distance[slot] == u64::MAX {
237        return;
238    }
239    let Some(cost) = side.distance[slot].checked_add(other.distance[slot]) else {
240        return;
241    };
242    if best.is_none_or(|current| cost < current.0 || (cost == current.0 && slot < current.1)) {
243        *best = Some((cost, slot));
244    }
245}
246
247fn reconstruct<G: IndexGraphView>(
248    source: G::Node,
249    target: G::Node,
250    meeting: usize,
251    by_slot: &[Option<G::Node>],
252    forward: &[Option<usize>],
253    backward: &[Option<usize>],
254) -> Option<Vec<G::Node>> {
255    let mut slots = vec![meeting];
256    while *slots.last()? != G::node_slot(source) {
257        slots.push(forward[*slots.last()?]?);
258    }
259    slots.reverse();
260    while *slots.last()? != G::node_slot(target) {
261        slots.push(backward[*slots.last()?]?);
262    }
263    slots
264        .into_iter()
265        .map(|slot| by_slot.get(slot).copied().flatten())
266        .collect::<Option<Vec<_>>>()
267}
weavatrix_graph/algo/shortest_extra.rs

weavatrix_graph/algo/
shortest_extra.rs
