#![allow(
clippy::cast_possible_truncation,
clippy::cast_sign_loss,
clippy::cast_possible_wrap,
clippy::many_single_char_names,
clippy::too_many_lines,
clippy::similar_names,
clippy::module_name_repetitions
)]
use crate::core::error::{IgraphError, IgraphResult};
use crate::core::graph::{Graph, VertexId};
#[derive(Debug, Clone, PartialEq, Eq)]
pub struct McsResult {
pub alpha: Vec<u32>,
pub alpham1: Vec<u32>,
}
#[derive(Debug, Clone, PartialEq, Eq)]
pub struct ChordalResult {
pub chordal: bool,
pub fill_in: Vec<(u32, u32)>,
}
pub fn maximum_cardinality_search(graph: &Graph) -> IgraphResult<McsResult> {
let n = graph.vcount() as usize;
if n == 0 {
return Ok(McsResult {
alpha: vec![],
alpham1: vec![],
});
}
let adj = build_simple_adj(graph)?;
let mut alpha = vec![0u32; n];
let mut alpham1 = vec![0u32; n];
let mut size = vec![0i32; n];
let mut head = vec![0u32; n];
let mut next = vec![0u32; n];
let mut prev = vec![0u32; n];
head[0] = 1; for v in 0..n {
next[v] = u32::try_from(v + 2).unwrap_or(0);
prev[v] = v as u32;
}
next[n - 1] = 0;
let mut i = n; let mut j: i32 = 0;
while i >= 1 {
let v = (head[j as usize])
.checked_sub(1)
.ok_or_else(|| IgraphError::InvalidArgument("MCS: empty bucket unexpectedly".into()))?
as usize;
let x = next[v];
head[j as usize] = x;
if x != 0 {
prev[(x - 1) as usize] = 0;
}
i -= 1;
alpha[v] = i as u32;
alpham1[i] = v as u32;
size[v] = -1;
for &w in &adj[v] {
let wu = w as usize;
let ws = size[wu];
if ws < 0 {
continue;
}
let nw = next[wu];
let pw = prev[wu];
if nw != 0 {
prev[(nw - 1) as usize] = pw;
}
if pw != 0 {
next[(pw - 1) as usize] = nw;
} else {
head[ws as usize] = nw;
}
size[wu] += 1;
let ws_new = size[wu] as usize;
let nw_new = head[ws_new];
next[wu] = nw_new;
prev[wu] = 0;
if nw_new != 0 {
prev[(nw_new - 1) as usize] = w + 1;
}
head[ws_new] = w + 1;
}
j += 1;
if (j as usize) < n {
while j >= 0 && head[j as usize] == 0 {
j -= 1;
}
}
}
Ok(McsResult { alpha, alpham1 })
}
pub fn is_chordal(graph: &Graph, alpha: Option<&[u32]>) -> IgraphResult<ChordalResult> {
let n = graph.vcount() as usize;
if n == 0 {
return Ok(ChordalResult {
chordal: true,
fill_in: vec![],
});
}
if let Some(a) = alpha {
if a.len() != n {
return Err(IgraphError::InvalidArgument(format!(
"is_chordal: alpha length {} != vcount {}",
a.len(),
n
)));
}
let mut inv = vec![0u32; n];
for (v, &rank) in a.iter().enumerate() {
let ri = rank as usize;
if ri >= n {
return Err(IgraphError::InvalidArgument(format!(
"is_chordal: alpha[{v}] = {rank} out of range [0, {n})"
)));
}
inv[ri] = v as u32;
}
is_chordal_impl(graph, n, a, &inv)
} else {
let mcs = maximum_cardinality_search(graph)?;
is_chordal_impl(graph, n, &mcs.alpha, &mcs.alpham1)
}
}
fn is_chordal_impl(
graph: &Graph,
n: usize,
alpha: &[u32],
alpham1: &[u32],
) -> IgraphResult<ChordalResult> {
let adj = build_simple_adj(graph)?;
let mut f = vec![0u32; n]; let mut index = vec![0u32; n];
let mut mark = vec![0u32; n]; let mut fill_in: Vec<(u32, u32)> = Vec::new();
let mut chordal = true;
for (i, &alpham1_i) in alpham1.iter().enumerate() {
let w = alpham1_i as usize;
f[w] = w as u32;
index[w] = i as u32;
for &nei in &adj[w] {
mark[nei as usize] = (w + 1) as u32;
}
for &nei in &adj[w] {
let v = nei as usize;
if alpha[v] as usize >= i {
continue;
}
let mut x = v;
while (index[x] as usize) < i {
index[x] = i as u32;
if mark[x] != (w + 1) as u32 {
chordal = false;
fill_in.push((x as u32, w as u32));
}
x = f[x] as usize;
}
if f[x] == x as u32 {
f[x] = w as u32;
}
}
}
Ok(ChordalResult { chordal, fill_in })
}
fn build_simple_adj(graph: &Graph) -> IgraphResult<Vec<Vec<VertexId>>> {
let n = graph.vcount() as usize;
let mut adj: Vec<Vec<VertexId>> = vec![Vec::new(); n];
let ec = graph.ecount();
for e in 0..ec {
let eid = u32::try_from(e)
.map_err(|_| IgraphError::InvalidArgument("edge count exceeds u32".into()))?;
let (u, v) = graph.edge(eid)?;
if u == v {
continue;
}
adj[u as usize].push(v);
adj[v as usize].push(u);
}
for list in &mut adj {
list.sort_unstable();
list.dedup();
}
Ok(adj)
}
#[cfg(test)]
mod tests {
use super::*;
use crate::create;
fn make(n: u32, edges: &[(u32, u32)]) -> Graph {
create(edges, n, false).expect("make")
}
#[test]
fn test_mcs_empty() {
let g = make(0, &[]);
let r = maximum_cardinality_search(&g).expect("ok");
assert!(r.alpha.is_empty());
assert!(r.alpham1.is_empty());
}
#[test]
fn test_mcs_singleton() {
let g = make(1, &[]);
let r = maximum_cardinality_search(&g).expect("ok");
assert_eq!(r.alpha, vec![0]);
assert_eq!(r.alpham1, vec![0]);
}
#[test]
fn test_mcs_triangle() {
let g = make(3, &[(0, 1), (1, 2), (2, 0)]);
let r = maximum_cardinality_search(&g).expect("ok");
assert_eq!(r.alpha.len(), 3);
assert_eq!(r.alpham1.len(), 3);
for (v, &rank) in r.alpha.iter().enumerate() {
assert_eq!(r.alpham1[rank as usize], v as u32);
}
}
#[test]
fn test_mcs_path() {
let g = make(4, &[(0, 1), (1, 2), (2, 3)]);
let r = maximum_cardinality_search(&g).expect("ok");
for (v, &rank) in r.alpha.iter().enumerate() {
assert_eq!(r.alpham1[rank as usize], v as u32);
}
}
#[test]
fn test_mcs_disconnected() {
let g = make(4, &[(0, 1)]);
let r = maximum_cardinality_search(&g).expect("ok");
assert_eq!(r.alpha.len(), 4);
for (v, &rank) in r.alpha.iter().enumerate() {
assert_eq!(r.alpham1[rank as usize], v as u32);
}
}
#[test]
fn test_mcs_self_loops() {
let g = make(3, &[(0, 0), (0, 1), (1, 2), (2, 2)]);
let r = maximum_cardinality_search(&g).expect("ok");
assert_eq!(r.alpha.len(), 3);
for (v, &rank) in r.alpha.iter().enumerate() {
assert_eq!(r.alpham1[rank as usize], v as u32);
}
}
#[test]
fn test_chordal_empty() {
let g = make(0, &[]);
let r = is_chordal(&g, None).expect("ok");
assert!(r.chordal);
assert!(r.fill_in.is_empty());
}
#[test]
fn test_chordal_singleton() {
let g = make(1, &[]);
let r = is_chordal(&g, None).expect("ok");
assert!(r.chordal);
assert!(r.fill_in.is_empty());
}
#[test]
fn test_chordal_triangle() {
let g = make(3, &[(0, 1), (1, 2), (2, 0)]);
let r = is_chordal(&g, None).expect("ok");
assert!(r.chordal);
assert!(r.fill_in.is_empty());
}
#[test]
fn test_chordal_complete_k4() {
let g = make(4, &[(0, 1), (0, 2), (0, 3), (1, 2), (1, 3), (2, 3)]);
let r = is_chordal(&g, None).expect("ok");
assert!(r.chordal);
assert!(r.fill_in.is_empty());
}
#[test]
fn test_not_chordal_4cycle() {
let g = make(4, &[(0, 1), (1, 2), (2, 3), (3, 0)]);
let r = is_chordal(&g, None).expect("ok");
assert!(!r.chordal);
assert!(!r.fill_in.is_empty());
}
#[test]
fn test_not_chordal_5cycle() {
let g = make(5, &[(0, 1), (1, 2), (2, 3), (3, 4), (4, 0)]);
let r = is_chordal(&g, None).expect("ok");
assert!(!r.chordal);
}
#[test]
fn test_chordal_tree() {
let g = make(5, &[(0, 1), (0, 2), (1, 3), (1, 4)]);
let r = is_chordal(&g, None).expect("ok");
assert!(r.chordal);
assert!(r.fill_in.is_empty());
}
#[test]
fn test_chordal_isolated_vertices() {
let g = make(5, &[]);
let r = is_chordal(&g, None).expect("ok");
assert!(r.chordal);
}
#[test]
fn test_chordal_with_alpha() {
let g = make(4, &[(0, 1), (1, 2), (2, 3), (3, 0)]);
let mcs = maximum_cardinality_search(&g).expect("mcs");
let r = is_chordal(&g, Some(&mcs.alpha)).expect("ok");
assert!(!r.chordal);
}
#[test]
fn test_chordal_wrong_alpha_size() {
let g = make(3, &[(0, 1)]);
assert!(is_chordal(&g, Some(&[0, 1])).is_err());
}
#[test]
fn test_chordal_igraph_lmu_graph() {
let g = make(
6,
&[
(0, 1),
(0, 2),
(1, 1),
(1, 3),
(2, 0),
(2, 0),
(2, 3),
(3, 4),
(3, 4),
],
);
let r = is_chordal(&g, None).expect("ok");
assert!(!r.chordal);
assert!(!r.fill_in.is_empty());
}
#[test]
fn test_chordal_directed_ignored() {
let g = create(&[(0, 1), (1, 2), (2, 0)], 3, true).expect("directed");
let r = is_chordal(&g, None).expect("ok");
assert!(r.chordal);
}
#[test]
fn test_chordal_fan_graph() {
let g = make(5, &[(0, 1), (0, 2), (0, 3), (0, 4), (1, 2), (2, 3), (3, 4)]);
let r = is_chordal(&g, None).expect("ok");
assert!(r.chordal);
}
#[test]
fn test_fill_in_makes_chordal() {
let g = make(4, &[(0, 1), (1, 2), (2, 3), (3, 0)]);
let r = is_chordal(&g, None).expect("ok");
assert!(!r.chordal);
let mut all_edges: Vec<(u32, u32)> = vec![(0, 1), (1, 2), (2, 3), (3, 0)];
all_edges.extend(r.fill_in.iter());
let g2 = make(4, &all_edges);
let r2 = is_chordal(&g2, None).expect("ok");
assert!(r2.chordal);
}
}
#[cfg(test)]
#[cfg(all(test, feature = "proptest-harness"))]
mod proptests {
use super::*;
use crate::create;
use proptest::prelude::*;
proptest! {
#[test]
fn prop_mcs_inverse_property(n in 1u32..30, edge_count in 0usize..60) {
let mut edges = Vec::new();
let mut rng_state = 0x1234_5678u64;
for _ in 0..edge_count {
rng_state = rng_state.wrapping_mul(6_364_136_223_846_793_005).wrapping_add(1);
let u = (rng_state >> 32) as u32 % n;
rng_state = rng_state.wrapping_mul(6_364_136_223_846_793_005).wrapping_add(1);
let v = (rng_state >> 32) as u32 % n;
if u != v {
edges.push((u, v));
}
}
let g = create(&edges, n, false).expect("graph");
let r = maximum_cardinality_search(&g).expect("mcs");
prop_assert_eq!(r.alpha.len(), n as usize);
prop_assert_eq!(r.alpham1.len(), n as usize);
for (v, &rank) in r.alpha.iter().enumerate() {
prop_assert!((rank as usize) < n as usize, "rank out of range");
prop_assert_eq!(r.alpham1[rank as usize], v as u32, "inverse mismatch");
}
}
#[test]
fn prop_mcs_ranks_are_permutation(n in 1u32..30, edge_count in 0usize..60) {
let mut edges = Vec::new();
let mut rng_state = 0xABCD_EF01u64;
for _ in 0..edge_count {
rng_state = rng_state.wrapping_mul(6_364_136_223_846_793_005).wrapping_add(1);
let u = (rng_state >> 32) as u32 % n;
rng_state = rng_state.wrapping_mul(6_364_136_223_846_793_005).wrapping_add(1);
let v = (rng_state >> 32) as u32 % n;
if u != v {
edges.push((u, v));
}
}
let g = create(&edges, n, false).expect("graph");
let r = maximum_cardinality_search(&g).expect("mcs");
let mut sorted_ranks: Vec<u32> = r.alpha.clone();
sorted_ranks.sort_unstable();
let expected: Vec<u32> = (0..n).collect();
prop_assert_eq!(sorted_ranks, expected, "ranks should be a permutation of 0..n");
}
#[test]
fn prop_fill_in_makes_chordal(n in 3u32..20, edge_count in 3usize..40) {
let mut edges = Vec::new();
let mut rng_state = 0xDEAD_BEEFu64;
for _ in 0..edge_count {
rng_state = rng_state.wrapping_mul(6_364_136_223_846_793_005).wrapping_add(1);
let u = (rng_state >> 32) as u32 % n;
rng_state = rng_state.wrapping_mul(6_364_136_223_846_793_005).wrapping_add(1);
let v = (rng_state >> 32) as u32 % n;
if u != v {
edges.push((u, v));
}
}
let g = create(&edges, n, false).expect("graph");
let r = is_chordal(&g, None).expect("chordal");
if !r.chordal {
let mut all_edges = edges.clone();
all_edges.extend(r.fill_in.iter());
let g2 = create(&all_edges, n, false).expect("graph2");
let r2 = is_chordal(&g2, None).expect("chordal2");
prop_assert!(r2.chordal, "adding fill-in should make graph chordal");
}
}
#[test]
fn prop_tree_is_chordal(n in 2u32..30) {
let mut edges = Vec::new();
let mut rng_state = 0xCAFE_BABEu64;
for v in 1..n {
rng_state = rng_state.wrapping_mul(6_364_136_223_846_793_005).wrapping_add(1);
let parent = (rng_state >> 32) as u32 % v;
edges.push((parent, v));
}
let g = create(&edges, n, false).expect("tree");
let r = is_chordal(&g, None).expect("chordal");
prop_assert!(r.chordal, "every tree is chordal");
prop_assert!(r.fill_in.is_empty(), "tree needs no fill-in");
}
#[test]
fn prop_complete_graph_is_chordal(n in 2u32..15) {
let mut edges = Vec::new();
for u in 0..n {
for v in (u + 1)..n {
edges.push((u, v));
}
}
let g = create(&edges, n, false).expect("complete");
let r = is_chordal(&g, None).expect("chordal");
prop_assert!(r.chordal, "complete graph is chordal");
}
}
}