diff --git a/Cargo.lock b/Cargo.lock index c9c28c997..0463c3b2b 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -163,6 +163,12 @@ version = "0.2.137" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "fc7fcc620a3bff7cdd7a365be3376c97191aeaccc2a603e600951e452615bf89" +[[package]] +name = "libm" +version = "0.2.6" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "348108ab3fba42ec82ff6e9564fc4ca0247bdccdc68dd8af9764bbc79c3c8ffb" + [[package]] name = "lock_api" version = "0.4.9" @@ -248,6 +254,7 @@ source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "578ede34cf02f8924ab9447f50c28075b4d3e5b269972345e7e0372b38c6cdcd" dependencies = [ "autocfg", + "libm", ] [[package]] @@ -451,6 +458,16 @@ dependencies = [ "getrandom", ] +[[package]] +name = "rand_distr" +version = "0.4.3" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "32cb0b9bc82b0a0876c2dd994a7e7a2683d3e7390ca40e6886785ef0c7e3ee31" +dependencies = [ + "num-traits", + "rand", +] + [[package]] name = "rand_pcg" version = "0.3.1" @@ -533,6 +550,7 @@ dependencies = [ "petgraph", "priority-queue", "rand", + "rand_distr", "rayon", ] diff --git a/releasenotes/notes/configuration-model-a0f0a7677c88f8ed.yaml b/releasenotes/notes/configuration-model-a0f0a7677c88f8ed.yaml new file mode 100644 index 000000000..20e55ece2 --- /dev/null +++ b/releasenotes/notes/configuration-model-a0f0a7677c88f8ed.yaml @@ -0,0 +1,14 @@ +--- +features: + - | + Added a new generator function, :func:`~rustworkx.undirected_configuration_model` + that will generate a random undirected graph from a given degree sequence. + For example: + + .. jupyter-execute:: + + import rustworkx + from rustworkx.visualization import mpl_draw + + graph = rustworkx.undirected_configuration_model([1, 2, 3, 4]) + mpl_draw(graph) diff --git a/rustworkx-core/Cargo.toml b/rustworkx-core/Cargo.toml index 6c9a00306..a8125ed86 100644 --- a/rustworkx-core/Cargo.toml +++ b/rustworkx-core/Cargo.toml @@ -17,6 +17,7 @@ petgraph = "0.6.2" rayon = "1.6" num-traits = "0.2" priority-queue = "1.2" +rand = "0.8" [dependencies.hashbrown] version = "0.12" @@ -28,3 +29,4 @@ features = ["rayon"] [dev-dependencies] rand = "0.8" +rand_distr = "0.4.3" diff --git a/rustworkx-core/src/generators/configuration_model.rs b/rustworkx-core/src/generators/configuration_model.rs new file mode 100644 index 000000000..6351f91d6 --- /dev/null +++ b/rustworkx-core/src/generators/configuration_model.rs @@ -0,0 +1,344 @@ +// 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 crate::generators::InvalidInputError; +use num_traits::{Num, ToPrimitive}; +use petgraph::data::{Build, Create}; +use petgraph::visit::{Data, GraphProp, NodeIndexable}; +use petgraph::Undirected; +use rand::distributions::Distribution; +use rand::prelude::SliceRandom; +use rand::Rng; +use std::iter::repeat; + +/// A degree sequence is a vector of integers that sum up to an even number. +/// +/// A graph with a given degree sequence can be generated with +/// [configuration_model](fn.configuration_model.html). +#[derive(Debug, Clone, PartialEq, Eq)] +pub struct DegreeSequence { + values: Vec, +} + +impl DegreeSequence { + /// Returns a new degree sequence. + /// + /// Arguments: + /// * `vec`: A vector of integers that sum up to an even number. + pub fn new(vec: Vec) -> Result { + let sum: usize = vec.iter().sum(); + if sum % 2 != 0 { + return Err(DegreeSequenceError::OddSum); + } + Ok(Self { values: vec }) + } + + /// Generates a random degree sequence with the given length and values + /// sampled from the given distribution. + /// + /// Arguments: + /// * `rng`: A random number generator. + /// * `length`: The length of the degree sequence. + /// * `distribution`: A distribution that samples values that are convertible to `usize`. + /// + /// NOTE: If the sampled values are floating point numbers, they will be truncated. + /// + /// # Example + /// ```rust + /// use rand::SeedableRng; + /// use rustworkx_core::petgraph::graph::UnGraph; + /// use rustworkx_core::generators::DegreeSequence; + /// use rand_distr::Poisson; + /// + /// let mut rng = rand::rngs::StdRng::seed_from_u64(42); + /// let poisson = Poisson::new(2.0).unwrap(); + /// let degree_sequence = DegreeSequence::from_distribution(&mut rng, 10, &poisson).unwrap(); + /// let expected_sequence = vec![3, 0, 4, 0, 3, 1, 0, 1, 4, 2]; + /// assert_eq!(*degree_sequence, expected_sequence); + /// ``` + pub fn from_distribution( + rng: &mut R, + length: usize, + distribution: &impl Distribution, + ) -> Result { + let mut degree_sequence: Vec = Vec::with_capacity(length); + let mut sum: usize = 0; + for value in distribution.sample_iter(&mut *rng).take(length) { + if value < I::zero() { + return Err(DegreeSequenceError::NegativeDegree); + } + let value = match value.to_usize() { + Some(value) => value, + None => return Err(DegreeSequenceError::Conversion), + }; + degree_sequence.push(value); + sum += value; + } + if sum % 2 == 1 { + let index = rng.gen_range(0..length); + degree_sequence[index] += 1; + } + Ok(Self { + values: degree_sequence, + }) + } +} + +impl std::ops::Deref for DegreeSequence { + type Target = [usize]; + + fn deref(&self) -> &Self::Target { + &self.values + } +} + +/// An error that can occur when generating a degree sequence. +#[derive(Clone, Copy, Debug, PartialEq, Eq)] +pub enum DegreeSequenceError { + NegativeDegree, + OddSum, + Conversion, +} + +impl std::fmt::Display for DegreeSequenceError { + fn fmt(&self, f: &mut std::fmt::Formatter) -> std::fmt::Result { + f.write_str(match self { + DegreeSequenceError::NegativeDegree => "Sampled value was negative", + DegreeSequenceError::OddSum => "Degree sequence has an odd sum", + DegreeSequenceError::Conversion => "Sampled value could not be converted to usize", + }) + } +} + +impl std::error::Error for DegreeSequenceError {} + +/// Generates a random undirected graph with the given degree sequence. +/// +/// Arguments: +/// +/// * `rng`: A random number generator. +/// * `degree_sequence`: A degree sequence. +/// * `weights` - A `Vec` of node weight objects. If specified, the length must be +/// equal to the length of the `degree_sequence`. +/// * `default_node_weight` - A callable that will return the weight to use +/// for newly created nodes. This is ignored if `weights` is specified, +/// as the weights from that argument will be used instead. +/// * `default_edge_weight` - A callable that will return the weight object +/// to use for newly created edges. +/// +/// The algorithm is based on "M.E.J. Newman, "The structure and function of complex networks", +/// SIAM REVIEW 45-2, pp 167-256, 2003.". +/// +/// The graph construction process might attempt to insert parallel edges and self-loops. +/// If the graph type does not support either of them, the resulting graph might not have the +/// exact degree sequence. However, the probability of parallel edges and self-loops tends to +/// converge towards zero as the number of nodes increases. +/// +/// Time complexity: **O(n + m)** +/// +/// # Example +/// +/// Approximating Erdős-Rényi G(n, p) model by using a Poisson degree distribution: +/// ```rust +/// use rustworkx_core::generators::{undirected_configuration_model, DegreeSequence}; +/// use petgraph::graph::UnGraph; +/// use rand_distr::Poisson; +/// +/// let mut rng = rand::thread_rng(); +/// let n = 200; +/// let p = 0.01; +/// let lambda = n as f64 * p; +/// let poisson = Poisson::new(lambda).unwrap(); +/// let degree_sequence = DegreeSequence::from_distribution(&mut rng, n, &poisson).unwrap(); +/// let graph: UnGraph<(), ()> = undirected_configuration_model(&mut rng, °ree_sequence, None, || (), || ()).unwrap(); +/// assert_eq!(graph.node_count(), 200); +/// ``` +pub fn undirected_configuration_model( + rng: &mut R, + degree_sequence: &DegreeSequence, + weights: Option>, + mut default_node_weight: F, + mut default_edge_weight: H, +) -> Result +where + G: GraphProp + + Build + + Create + + Data + + NodeIndexable, + F: FnMut() -> T, + H: FnMut() -> M, + R: Rng + ?Sized, +{ + let num_nodes = degree_sequence.len(); + + let mut graph = G::with_capacity(num_nodes, num_nodes); + let mut nodes: Vec = Vec::with_capacity(num_nodes); + + match weights { + Some(weights) => { + if weights.len() != num_nodes { + return Err(InvalidInputError {}); + } + for weight in weights { + nodes.push(graph.add_node(weight)); + } + } + None => { + for _ in 0..num_nodes { + nodes.push(graph.add_node(default_node_weight())); + } + } + }; + + let mut stubs: Vec = nodes + .into_iter() + .zip(degree_sequence.iter()) + .flat_map(|(node, degree)| repeat(node).take(*degree)) + .collect(); + stubs.shuffle(rng); + + let (sources, targets) = stubs.split_at(stubs.len() / 2); + for (&source, &target) in sources.iter().zip(targets) { + graph.add_edge(source, target, default_edge_weight()); + } + Ok(graph) +} + +#[cfg(test)] +mod tests { + use super::*; + use petgraph::graph::UnGraph; + use petgraph::visit::EdgeRef; + use rand::rngs::StdRng; + use rand::SeedableRng; + + #[test] + fn test_degree_sequence() { + let degree_sequence = DegreeSequence::new(vec![1, 2, 3, 4]).unwrap(); + assert_eq!(*degree_sequence, vec![1, 2, 3, 4]); + } + + #[test] + fn test_degree_sequence_empty() { + let degree_sequence = DegreeSequence::new(vec![]).unwrap(); + assert_eq!(*degree_sequence, vec![]); + } + + #[test] + fn test_degree_sequence_odd_sum() { + let result = DegreeSequence::new(vec![1, 2, 3, 5]); + assert_eq!(result, Err(DegreeSequenceError::OddSum)); + } + + #[test] + fn test_degree_sequence_from_discrete_distribution() { + let mut rng = StdRng::seed_from_u64(42); + let poisson = rand_distr::Poisson::new(2.0).unwrap(); + let degree_sequence = DegreeSequence::from_distribution(&mut rng, 10, &poisson).unwrap(); + assert_eq!(*degree_sequence, vec![3, 0, 4, 0, 3, 1, 0, 1, 4, 2]); + } + + #[test] + fn test_degree_sequence_from_continuous_distribution() { + let mut rng = StdRng::seed_from_u64(42); + let exp = rand_distr::Exp::new(0.5).unwrap(); + let degree_sequence = DegreeSequence::from_distribution(&mut rng, 10, &exp).unwrap(); + assert_eq!(*degree_sequence, vec![1, 1, 0, 0, 0, 1, 0, 4, 0, 1]); + } + + #[test] + fn test_degree_sequence_from_distribution_empty() { + let mut rng = StdRng::seed_from_u64(42); + let uniform = rand_distr::Uniform::new(0, 10); + let degree_sequence = DegreeSequence::from_distribution(&mut rng, 0, &uniform).unwrap(); + assert_eq!(*degree_sequence, vec![]); + } + + #[test] + fn test_degree_sequence_from_distribution_conversion_error() { + let mut rng = StdRng::seed_from_u64(42); + let exp = rand_distr::Exp::new(0.0).unwrap(); // samples always yield infinity + let result = DegreeSequence::from_distribution(&mut rng, 1, &exp); + assert_eq!(result, Err(DegreeSequenceError::Conversion)); + } + + #[test] + fn test_degree_sequence_from_distribution_negative_degree() { + let mut rng = StdRng::seed_from_u64(42); + let poisson = rand_distr::Uniform::new(-1, 0); + let result = DegreeSequence::from_distribution(&mut rng, 10, &poisson); + assert_eq!(result, Err(DegreeSequenceError::NegativeDegree)); + } + + #[test] + fn test_configuration_model() { + let mut rng = StdRng::seed_from_u64(42); + let degree_sequence = DegreeSequence::new(vec![1, 2, 3, 4]).unwrap(); + let graph: UnGraph<(), ()> = + undirected_configuration_model(&mut rng, °ree_sequence, None, || (), || ()).unwrap(); + assert_eq!(graph.node_count(), 4); + assert_eq!(graph.edge_count(), 5); + assert_eq!( + vec![(3, 3), (2, 3), (3, 2), (2, 1), (0, 1)], + graph + .edge_references() + .map(|edge| (edge.source().index(), edge.target().index())) + .collect::>(), + ); + } + + #[test] + fn test_configuration_model_empty() { + let mut rng = StdRng::seed_from_u64(42); + let degree_sequence = DegreeSequence::new(vec![]).unwrap(); + let graph: UnGraph<(), ()> = + undirected_configuration_model(&mut rng, °ree_sequence, None, || (), || ()).unwrap(); + assert_eq!(graph.node_count(), 0); + assert_eq!(graph.edge_count(), 0); + } + + #[test] + fn test_configuration_model_with_weights() { + let mut rng = StdRng::seed_from_u64(42); + let degree_sequence = DegreeSequence::new(vec![1, 2, 3, 4]).unwrap(); + let graph: UnGraph = undirected_configuration_model( + &mut rng, + °ree_sequence, + Some(vec![0, 1, 2, 3]), + || 1, + || (), + ) + .unwrap(); + assert_eq!( + vec![0, 1, 2, 3], + graph.node_weights().copied().collect::>(), + ); + } + + #[test] + fn test_configuration_model_degree_sequence_and_weights_have_different_lengths() { + let mut rng = StdRng::seed_from_u64(42); + let degree_sequence = DegreeSequence::new(vec![1, 2, 3, 4]).unwrap(); + let result: Result, _> = undirected_configuration_model( + &mut rng, + °ree_sequence, + Some(vec![0, 1, 2]), + || 1, + || (), + ); + match result { + Ok(_) => panic!("Expected error"), + Err(e) => assert_eq!(e, InvalidInputError), + } + } +} diff --git a/rustworkx-core/src/generators/mod.rs b/rustworkx-core/src/generators/mod.rs index 40237685d..6f1390452 100644 --- a/rustworkx-core/src/generators/mod.rs +++ b/rustworkx-core/src/generators/mod.rs @@ -15,6 +15,7 @@ mod barbell_graph; mod binomial_tree_graph; mod complete_graph; +mod configuration_model; mod cycle_graph; mod full_rary_tree_graph; mod grid_graph; @@ -46,6 +47,7 @@ impl fmt::Display for InvalidInputError { pub use barbell_graph::barbell_graph; pub use binomial_tree_graph::binomial_tree_graph; pub use complete_graph::complete_graph; +pub use configuration_model::{undirected_configuration_model, DegreeSequence}; pub use cycle_graph::cycle_graph; pub use full_rary_tree_graph::full_rary_tree_graph; pub use grid_graph::grid_graph; diff --git a/src/lib.rs b/src/lib.rs index a9ad4f4dd..b70c0c559 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -425,6 +425,7 @@ fn rustworkx(py: Python<'_>, m: &PyModule) -> PyResult<()> { m.add_wrapped(wrap_pyfunction!(undirected_gnp_random_graph))?; m.add_wrapped(wrap_pyfunction!(directed_gnm_random_graph))?; m.add_wrapped(wrap_pyfunction!(undirected_gnm_random_graph))?; + m.add_wrapped(wrap_pyfunction!(undirected_configuration_model))?; m.add_wrapped(wrap_pyfunction!(random_geometric_graph))?; m.add_wrapped(wrap_pyfunction!(cycle_basis))?; m.add_wrapped(wrap_pyfunction!(simple_cycles))?; diff --git a/src/random_graph.rs b/src/random_graph.rs index dd97f33de..0e03d5962 100644 --- a/src/random_graph.rs +++ b/src/random_graph.rs @@ -26,6 +26,7 @@ use petgraph::prelude::*; use rand::distributions::{Distribution, Uniform}; use rand::prelude::*; use rand_pcg::Pcg64; +use rustworkx_core::generators as core_generators; /// Return a :math:`G_{np}` directed random graph, also known as an /// Erdős-Rényi graph or a binomial graph. @@ -493,3 +494,67 @@ pub fn random_geometric_graph( }; Ok(graph) } + +/// Generate an undirected configuration model graph. +/// +/// The configuration model is a random graph with a given degree sequence. +/// +/// The algorithm is based on "M.E.J. Newman, "The structure and function of complex networks", +/// SIAM REVIEW 45-2, pp 167-256, 2003.". +/// +/// The graph construction process may insert parallel edges and self-loops. +/// However, the probability of parallel edges and self-loops tends to +/// converge towards zero as the number of nodes increases. +/// +/// This algorithm has a time complexity of :math:`O(n + m)` +/// +/// :param list degree_sequence: A list of positive integers specifying the degree sequence +/// :param list weights: An optional list of node weights. If specified, the length of the list +/// must be equal to the length of the degree sequence +/// :param int seed: An optional seed to use for the random number generator +/// +/// :return: A PyGraph object +/// :rtype: PyGraph +#[pyfunction] +#[pyo3(text_signature = "(degree_sequence, /, weights=None, seed=None)")] +pub fn undirected_configuration_model( + py: Python, + degree_sequence: Vec, + weights: Option>, + seed: Option, +) -> PyResult { + let default_fn = || py.None(); + let mut rng: Pcg64 = match seed { + Some(seed) => Pcg64::seed_from_u64(seed), + None => Pcg64::from_entropy(), + }; + let degree_sequence = match core_generators::DegreeSequence::new(degree_sequence) { + Ok(degree_sequence) => degree_sequence, + Err(_) => { + return Err(PyValueError::new_err( + "degree_sequence must sum up to an even number", + )); + } + }; + + let graph: StablePyGraph = match core_generators::undirected_configuration_model( + &mut rng, + °ree_sequence, + weights, + default_fn, + default_fn, + ) { + Ok(graph) => graph, + Err(_) => { + return Err(PyValueError::new_err( + "degree_sequence and weights have different lengths", + )); + } + }; + Ok(graph::PyGraph { + graph, + node_removed: false, + multigraph: true, + attrs: py.None(), + }) +} diff --git a/tests/rustworkx_tests/generators/test_configuration_model.py b/tests/rustworkx_tests/generators/test_configuration_model.py new file mode 100644 index 000000000..fd6addd7e --- /dev/null +++ b/tests/rustworkx_tests/generators/test_configuration_model.py @@ -0,0 +1,43 @@ +# 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. + +import unittest + +import rustworkx + + +class TestConfigurationModel(unittest.TestCase): + def test_undirected_configuration_model(self): + graph = rustworkx.undirected_configuration_model([0, 1, 1, 1, 1, 0], seed=42) + self.assertEqual(len(graph), 6) + self.assertEqual(len(graph.edges()), 2) + + def test_undirected_configuration_model_empty(self): + graph = rustworkx.undirected_configuration_model([], seed=42) + self.assertEqual(len(graph), 0) + self.assertEqual(len(graph.edges()), 0) + + def test_undirected_configuration_model_weights(self): + graph = rustworkx.undirected_configuration_model( + [1, 2, 3, 4], weights=list(range(4)), seed=42 + ) + self.assertEqual(len(graph), 4) + self.assertEqual([x for x in range(4)], graph.nodes()) + self.assertEqual(len(graph.edges()), 5) + + def test_undirected_configuration_model_length_mismatch(self): + with self.assertRaises(ValueError): + rustworkx.undirected_configuration_model([1, 2, 3, 4], weights=list(range(3))) + + def test_undirected_configuration_model_odd_sum(self): + with self.assertRaises(ValueError): + rustworkx.undirected_configuration_model([1, 2, 3, 5])