weavatrix_graph/algo/
cycles.rs

1use super::strongly_connected_components_filtered;
2use crate::IndexGraphView;
3use crate::Vec;
4
5#[derive(Debug, Clone, PartialEq, Eq)]
6pub struct CycleEnumeration<Node> {
7    paths: Vec<Vec<Node>>,
8    truncated: bool,
9}
10
11impl<Node> CycleEnumeration<Node> {
12    #[must_use]
13    pub fn paths(&self) -> &[Vec<Node>] {
14        &self.paths
15    }
16
17    #[must_use]
18    pub const fn truncated(&self) -> bool {
19        self.truncated
20    }
21}
22
23/// Enumerates elementary directed circuits with Johnson's blocked-set algorithm.
24///
25/// Every returned circuit repeats its start node at the end. `max_cycles`
26/// bounds exponential output.
27pub fn johnson_cycles<G>(graph: &G, max_cycles: usize) -> CycleEnumeration<G::Node>
28where
29    G: IndexGraphView,
30{
31    if max_cycles == 0 {
32        return CycleEnumeration {
33            paths: Vec::new(),
34            truncated: false,
35        };
36    }
37    let (by_slot, adjacency) = indexed_adjacency(graph);
38    let mut state = JohnsonState::new(graph.node_bound(), max_cycles);
39    let mut lower = 0;
40    while let Some(component) = least_cyclic_component(graph, &adjacency, lower) {
41        let Some(start) = component.iter().min().copied() else {
42            break;
43        };
44        let mut allowed = vec![false; graph.node_bound()];
45        for &node in &component {
46            allowed[node] = true;
47        }
48        state.reset();
49        circuit(start, start, &adjacency, &allowed, &by_slot, &mut state);
50        if state.truncated {
51            break;
52        }
53        lower = start + 1;
54    }
55    CycleEnumeration {
56        paths: state.results,
57        truncated: state.truncated,
58    }
59}
60
61fn indexed_adjacency<G: IndexGraphView>(graph: &G) -> (Vec<Option<G::Node>>, Vec<Vec<usize>>) {
62    let mut by_slot = vec![None; graph.node_bound()];
63    let mut adjacency = vec![Vec::new(); graph.node_bound()];
64    for node in graph.node_indices() {
65        let source = G::node_slot(node);
66        by_slot[source] = Some(node);
67        adjacency[source] = graph
68            .outgoing_edges(node)
69            .filter_map(|edge| graph.edge_endpoints(edge))
70            .map(|endpoints| G::node_slot(endpoints.target()))
71            .collect();
72        adjacency[source].sort_unstable();
73        adjacency[source].dedup();
74    }
75    (by_slot, adjacency)
76}
77
78fn least_cyclic_component<G>(
79    graph: &G,
80    adjacency: &[Vec<usize>],
81    lower: usize,
82) -> Option<Vec<usize>>
83where
84    G: IndexGraphView,
85{
86    strongly_connected_components_filtered(graph, |edge| {
87        graph.edge_endpoints(edge).is_some_and(|endpoints| {
88            G::node_slot(endpoints.source()) >= lower && G::node_slot(endpoints.target()) >= lower
89        })
90    })
91    .into_iter()
92    .map(|component| component.into_iter().map(G::node_slot).collect::<Vec<_>>())
93    .filter(|component| {
94        component.iter().all(|node| *node >= lower)
95            && (component.len() > 1
96                || component
97                    .first()
98                    .is_some_and(|node| adjacency[*node].binary_search(node).is_ok()))
99    })
100    .min_by_key(|component| component.iter().min().copied())
101}
102
103struct JohnsonState<Node> {
104    blocked: Vec<bool>,
105    blocked_by: Vec<Vec<usize>>,
106    stack: Vec<usize>,
107    results: Vec<Vec<Node>>,
108    limit: usize,
109    truncated: bool,
110}
111
112impl<Node> JohnsonState<Node> {
113    fn new(bound: usize, limit: usize) -> Self {
114        Self {
115            blocked: vec![false; bound],
116            blocked_by: vec![Vec::new(); bound],
117            stack: Vec::new(),
118            results: Vec::new(),
119            limit,
120            truncated: false,
121        }
122    }
123
124    fn reset(&mut self) {
125        self.blocked.fill(false);
126        for dependencies in &mut self.blocked_by {
127            dependencies.clear();
128        }
129        self.stack.clear();
130    }
131}
132
133fn circuit<Node: Copy>(
134    node: usize,
135    start: usize,
136    adjacency: &[Vec<usize>],
137    allowed: &[bool],
138    by_slot: &[Option<Node>],
139    state: &mut JohnsonState<Node>,
140) -> bool {
141    let Some(start_node) = by_slot[start] else {
142        return false;
143    };
144    let mut found = false;
145    state.stack.push(node);
146    state.blocked[node] = true;
147    for &neighbor in &adjacency[node] {
148        if !allowed[neighbor] {
149            continue;
150        }
151        if neighbor == start {
152            if state.results.len() == state.limit {
153                state.truncated = true;
154                break;
155            }
156            let mut cycle = state
157                .stack
158                .iter()
159                .filter_map(|slot| by_slot[*slot])
160                .collect::<Vec<_>>();
161            cycle.push(start_node);
162            state.results.push(cycle);
163            found = true;
164        } else if !state.blocked[neighbor]
165            && circuit(neighbor, start, adjacency, allowed, by_slot, state)
166        {
167            found = true;
168        }
169        if state.truncated {
170            break;
171        }
172    }
173    if found {
174        unblock(node, state);
175    } else {
176        for &neighbor in &adjacency[node] {
177            if allowed[neighbor] && !state.blocked_by[neighbor].contains(&node) {
178                state.blocked_by[neighbor].push(node);
179            }
180        }
181    }
182    state.stack.pop();
183    found
184}
185
186fn unblock<Node>(node: usize, state: &mut JohnsonState<Node>) {
187    state.blocked[node] = false;
188    let dependencies = core::mem::take(&mut state.blocked_by[node]);
189    for dependency in dependencies {
190        if state.blocked[dependency] {
191            unblock(dependency, state);
192        }
193    }
194}
weavatrix_graph/algo/cycles.rs

weavatrix_graph/algo/
cycles.rs
