Add chain_decomposition function. (Qiskit#444)

* Add  function.

This commit adds a new function that finds
a chain decomposition of an undirected .
It's defined in https://doi.org/10.1016/j.ipl.2013.01.016
and can be used to compute all bridges and cut
vertices of the input graph.

One current limitation is that it's a recursive implementation
since it's based on the .

* switch to `usize` + remove `drain` + a note about multigraphs

* use petgraph traits

* move to retworkx-core + add docs

* use iterators to collect the final chains

Co-authored-by: Matthew Treinish <[email protected]>

* update docs

* add custom return type

* run cargo fmt

* add tests for new custom return type

* Update citation name in python docstring

docs/source/api.rst

@@ -130,6 +130,7 @@ Connectivity and Cycles
    retworkx.is_weakly_connected
    retworkx.cycle_basis
    retworkx.digraph_find_cycle
+   retworkx.chain_decomposition
 
 .. _other-algorithms:
 
@@ -328,3 +329,4 @@ Custom Return Types
    retworkx.AllPairsPathMapping
    retworkx.AllPairsPathLengthMapping
    retworkx.CentralityMapping
+   retworkx.Chains

releasenotes/notes/chain-decomposition-3fdf4b283b5b9ad1.yaml

@@ -0,0 +1,40 @@
+---
+features:
+  - |
+    Added a new function :func:`~retworkx.chain_decomposition` that finds
+    a chain decomposition of an undirected :class:`~retworkx.PyGraph`.
+    A chain decomposition is a set of cycles or paths derived from the
+    set of fundamental cycles of a depth-first tree. It's defined in
+    https://doi.org/10.1016/j.ipl.2013.01.016
+    For example:
+
+    .. jupyter-execute::
+
+      import retworkx
+      from retworkx.visualization import mpl_draw
+
+      graph = retworkx.PyGraph()
+      graph.extend_from_edge_list([
+          (0, 1), (0, 2), (1, 2), (2, 3),
+          (3, 4), (3, 5), (4, 5),
+      ])
+      chains = retworkx.chain_decomposition(graph)
+
+      def color_edges(graph, chains): 
+          COLORS = ['blue', 'red']
+
+          edges_in_chain = {}
+          for idx, chain in enumerate(chains):
+              for edge in chain:
+                  edge = tuple(sorted(edge))
+                  edges_in_chain[edge] = COLORS[idx]
+              
+          edge_colors = []
+          for edge in graph.edge_list():
+              edge = tuple(sorted(edge))
+              edge_colors += [edges_in_chain.get(edge, 'black')]
+
+          return edge_colors
+
+      mpl_draw(graph, node_color='black', node_size=150,
+              edge_color=color_edges(graph, chains))

retworkx-core/src/connectivity/chain.rs

@@ -0,0 +1,172 @@
+// Licensed under the Apache License, Version 2.0 (the "License"); you may
+// not use this file except in compliance with the License. You may obtain
+// a copy of the License at
+//
+//     http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS, WITHOUT
+// WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the
+// License for the specific language governing permissions and limitations
+// under the License.
+
+use std::cmp::Eq;
+use std::hash::Hash;
+
+use hashbrown::HashMap;
+
+use petgraph::visit::{
+    depth_first_search, DfsEvent, GraphProp, IntoNeighbors,
+    IntoNodeIdentifiers, NodeCount, NodeIndexable, VisitMap, Visitable,
+};
+use petgraph::Undirected;
+
+fn _build_chain<G, VM: VisitMap<G::NodeId>>(
+    graph: G,
+    parent: &[usize],
+    mut u_id: G::NodeId,
+    mut v_id: G::NodeId,
+    visited: &mut VM,
+) -> Vec<(G::NodeId, G::NodeId)>
+where
+    G: Visitable + NodeIndexable,
+{
+    let mut chain = Vec::new();
+    while visited.visit(v_id) {
+        chain.push((u_id, v_id));
+        u_id = v_id;
+        let u = graph.to_index(u_id);
+        let v = parent[u];
+        v_id = graph.from_index(v);
+    }
+    chain.push((u_id, v_id));
+
+    chain
+}
+
+/// Returns the chain decomposition of a graph.
+///
+/// The *chain decomposition* of a graph with respect to a depth-first
+/// search tree is a set of cycles or paths derived from the set of
+/// fundamental cycles of the tree in the following manner. Consider
+/// each fundamental cycle with respect to the given tree, represented
+/// as a list of edges beginning with the nontree edge oriented away
+/// from the root of the tree. For each fundamental cycle, if it
+/// overlaps with any previous fundamental cycle, just take the initial
+/// non-overlapping segment, which is a path instead of a cycle. Each
+/// cycle or path is called a *chain*. For more information,
+/// see [`Schmidt`](https://doi.org/10.1016/j.ipl.2013.01.016).
+///
+/// The graph should be undirected. If `source` is specified only the chain
+/// decomposition for the connected component containing this node will be returned.
+/// This node indicates the root of the depth-first search tree. If it's not
+/// specified, a source will be chosen arbitrarly and repeated until all components
+/// of the graph are searched.
+///
+/// Returns a list of list of edges where each inner list is a chain.
+///
+/// # Note
+/// The function implicitly assumes that there are no parallel edges
+/// or self loops. It may produce incorrect/unexpected results if the
+/// input graph has self loops or parallel edges.
+///
+/// # Example
+/// ```rust
+/// use retworkx_core::connectivity::chain_decomposition;
+/// use retworkx_core::petgraph::graph::{NodeIndex, UnGraph};
+///
+/// let mut graph : UnGraph<(), ()> = UnGraph::new_undirected();
+/// let a = graph.add_node(()); // node with no weight
+/// let b = graph.add_node(());
+/// let c = graph.add_node(());
+/// let d = graph.add_node(());
+/// let e = graph.add_node(());
+/// let f = graph.add_node(());
+/// let g = graph.add_node(());
+/// let h = graph.add_node(());
+///
+/// graph.extend_with_edges(&[
+///     (a, b),
+///     (b, c),
+///     (c, d),
+///     (d, a),
+///     (e, f),
+///     (b, e),
+///     (f, g),
+///     (g, h),
+///     (h, e)
+/// ]);
+/// // a ---- b ---- e ---- f
+/// // |      |      |      |
+/// // d ---- c      h ---- g
+///
+/// let chains = chain_decomposition(&graph, None);
+/// assert_eq!(
+///     chains,
+///     vec![
+///         vec![(a, d), (d, c), (c, b), (b, a)],
+///         vec![(e, h), (h, g), (g, f), (f, e)]
+///     ]
+/// );
+/// ```
+pub fn chain_decomposition<G>(
+    graph: G,
+    source: Option<G::NodeId>,
+) -> Vec<Vec<(G::NodeId, G::NodeId)>>
+where
+    G: IntoNodeIdentifiers
+        + IntoNeighbors
+        + Visitable
+        + NodeIndexable
+        + NodeCount
+        + GraphProp<EdgeType = Undirected>,
+    G::NodeId: Eq + Hash,
+{
+    let roots = match source {
+        Some(node) => vec![node],
+        None => graph.node_identifiers().collect(),
+    };
+
+    let mut parent = vec![std::usize::MAX; graph.node_bound()];
+    let mut back_edges: HashMap<G::NodeId, Vec<G::NodeId>> = HashMap::new();
+
+    // depth-first-index (DFI) ordered nodes.
+    let mut nodes = Vec::with_capacity(graph.node_count());
+    depth_first_search(graph, roots, |event| match event {
+        DfsEvent::Discover(u, _) => {
+            nodes.push(u);
+        }
+        DfsEvent::TreeEdge(u, v) => {
+            let u = graph.to_index(u);
+            let v = graph.to_index(v);
+            parent[v] = u;
+        }
+        DfsEvent::BackEdge(u_id, v_id) => {
+            let u = graph.to_index(u_id);
+            let v = graph.to_index(v_id);
+
+            // do *not* consider ``(u, v)`` as a back edge if ``(v, u)`` is a tree edge.
+            if parent[u] != v {
+                back_edges
+                    .entry(v_id)
+                    .and_modify(|v_edges| v_edges.push(u_id))
+                    .or_insert(vec![u_id]);
+            }
+        }
+        _ => {}
+    });
+
+    let visited = &mut graph.visit_map();
+    nodes
+        .into_iter()
+        .filter_map(|u| {
+            visited.visit(u);
+            back_edges.get(&u).map(|vs| {
+                vs.iter()
+                    .map(|v| _build_chain(graph, &parent, u, *v, visited))
+                    .collect::<Vec<Vec<(G::NodeId, G::NodeId)>>>()
+            })
+        })
+        .flatten()
+        .collect()
+}

retworkx-core/src/connectivity/mod.rs

@@ -0,0 +1,17 @@
+// Licensed under the Apache License, Version 2.0 (the "License"); you may
+// not use this file except in compliance with the License. You may obtain
+// a copy of the License at
+//
+//     http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS, WITHOUT
+// WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the
+// License for the specific language governing permissions and limitations
+// under the License.
+
+//! Module for connectivity and cut algorithms.
+
+mod chain;
+
+pub use chain::chain_decomposition;

retworkx-core/src/lib.rs

@@ -47,6 +47,7 @@
 //! The crate is organized into
 //!
 //! * [`centrality`](./centrality/index.html)
+//! * [`connectivity`](./connectivity/index.html)
 //! * [`max_weight_matching`](./max_weight_matching/index.html)
 //! * [`shortest_path`](./shortest_path/index.html)
 //!
@@ -68,6 +69,7 @@ pub type Result<T, E = Infallible> = core::result::Result<T, E>;
 
 /// Module for centrality algorithms
 pub mod centrality;
+pub mod connectivity;
 /// Module for depth first search edge methods
 pub mod dfs_edges;
 /// Module for maximum weight matching algorithmss

src/connectivity/mod.rs

@@ -33,7 +33,8 @@ use petgraph::visit::{EdgeRef, IntoEdgeReferences, NodeCount, NodeIndexable};
 use ndarray::prelude::*;
 use numpy::IntoPyArray;
 
-use crate::iterators::EdgeList;
+use crate::iterators::{Chains, EdgeList};
+use retworkx_core::connectivity;
 
 /// Return a list of cycles which form a basis for cycles of a given PyGraph
 ///
@@ -662,3 +663,57 @@ pub fn digraph_core_number(
 ) -> PyResult<PyObject> {
     core_number::core_number(py, &graph.graph)
 }
+
+/// Returns the chain decomposition of a graph.
+///
+/// The *chain decomposition* of a graph with respect to a depth-first
+/// search tree is a set of cycles or paths derived from the set of
+/// fundamental cycles of the tree in the following manner. Consider
+/// each fundamental cycle with respect to the given tree, represented
+/// as a list of edges beginning with the nontree edge oriented away
+/// from the root of the tree. For each fundamental cycle, if it
+/// overlaps with any previous fundamental cycle, just take the initial
+/// non-overlapping segment, which is a path instead of a cycle. Each
+/// cycle or path is called a *chain*. For more information, see [Schmidt]_.
+///
+/// .. note::
+///
+///     The function implicitly assumes that there are no parallel edges
+///     or self loops. It may produce incorrect/unexpected results if the
+///     input graph has self loops or parallel edges. It's also a recursive
+///     implementation and might run out of memory in large graphs.
+///
+/// :param PyGraph: The undirected graph to be used
+/// :param int source: An optional node index in the graph. If specified,
+///     only the chain decomposition for the connected component containing
+///     this node will be returned. This node indicates the root of the depth-first
+///     search tree. If this is not specified then a source will be chosen
+///     arbitrarly and repeated until all components of the graph are searched.
+/// :returns: A list of list of edges where each inner list is a chain.
+/// :rtype: list of EdgeList
+///
+/// .. [Schmidt] Jens M. Schmidt (2013). "A simple test on 2-vertex-
+///       and 2-edge-connectivity." *Information Processing Letters*,
+///       113, 241–244. Elsevier. <https://doi.org/10.1016/j.ipl.2013.01.016>
+#[pyfunction]
+#[pyo3(text_signature = "(graph, /, source=None)")]
+pub fn chain_decomposition(
+    graph: graph::PyGraph,
+    source: Option<usize>,
+) -> Chains {
+    let chains = connectivity::chain_decomposition(
+        &graph.graph,
+        source.map(NodeIndex::new),
+    );
+    Chains {
+        chains: chains
+            .into_iter()
+            .map(|chain| EdgeList {
+                edges: chain
+                    .into_iter()
+                    .map(|(a, b)| (a.index(), b.index()))
+                    .collect(),
+            })
+            .collect(),
+    }
+}

src/iterators.rs

@@ -698,6 +698,56 @@ custom_vec_iter_impl!(
 );
 default_pygc_protocol_impl!(EdgeIndices);
 
+impl PyHash for EdgeList {
+    fn hash<H: Hasher>(&self, py: Python, state: &mut H) -> PyResult<()> {
+        PyHash::hash(&self.edges, py, state)?;
+        Ok(())
+    }
+}
+
+impl PyEq<PyAny> for EdgeList {
+    #[inline]
+    fn eq(&self, other: &PyAny, py: Python) -> PyResult<bool> {
+        PyEq::eq(&self.edges, other.downcast::<PySequence>()?, py)
+    }
+}
+
+impl PyDisplay for EdgeList {
+    fn str(&self, py: Python) -> PyResult<String> {
+        Ok(format!("EdgeList{}", self.edges.str(py)?))
+    }
+}
+
+custom_vec_iter_impl!(
+    Chains,
+    chains,
+    EdgeList,
+    "A custom class for the return of a list of list of edges.
+
+    This class is a container class for the results of functions that
+    return a list of list of edges. It implements the Python sequence
+    protocol. So you can treat the return as a read-only sequence/list
+    that is integer indexed. If you want to use it as an iterator you
+    can by wrapping it in an ``iter()`` that will yield the results in
+    order.
+
+    For example::
+
+        import retworkx
+
+        graph = retworkx.generators.hexagonal_lattice_graph(2, 2)
+        chains = retworkx.chain_decomposition(graph)
+        # Index based access
+        third_chain = chains[2]
+        # Use as iterator
+        chains_iter = iter(chains)
+        first_chain = next(chains_iter)
+        second_chain = next(chains_iter)
+
+    "
+);
+default_pygc_protocol_impl!(Chains);
+
 macro_rules! py_object_protocol_impl {
     ($name:ident, $data:ident) => {
         #[pyproto]

src/lib.rs

@@ -301,9 +301,11 @@ fn retworkx(py: Python<'_>, m: &PyModule) -> PyResult<()> {
     ))?;
     m.add_wrapped(wrap_pyfunction!(metric_closure))?;
     m.add_wrapped(wrap_pyfunction!(steiner_tree))?;
+    m.add_wrapped(wrap_pyfunction!(chain_decomposition))?;
     m.add_class::<digraph::PyDiGraph>()?;
     m.add_class::<graph::PyGraph>()?;
     m.add_class::<iterators::BFSSuccessors>()?;
+    m.add_class::<iterators::Chains>()?;
     m.add_class::<iterators::NodeIndices>()?;
     m.add_class::<iterators::EdgeIndices>()?;
     m.add_class::<iterators::EdgeList>()?;