Skip to main content

constraint_theory_core/
percolation.rs

1//! Rigidity Percolation using Laman's Theorem
2
3use std::collections::HashMap;
4
5/// Result of rigidity percolation analysis
6///
7/// Provides comprehensive metrics about the rigidity properties of a graph
8/// based on Laman's theorem for structural rigidity.
9pub struct RigidityResult {
10    /// Whether the graph satisfies Laman's condition for minimal rigidity
11    pub is_rigid: bool,
12    /// The rank of the rigidity matrix (number of independent constraints)
13    pub rank: usize,
14    /// Deficiency from minimal rigidity (0 = minimally rigid)
15    pub deficiency: usize,
16    /// Number of connected clusters in the graph
17    pub n_clusters: usize,
18    /// Fraction of nodes that are part of rigid clusters
19    pub rigid_fraction: f32,
20}
21
22/// Fast union-find data structure for rigidity percolation
23///
24/// Uses path compression and union by rank for efficient
25/// connectivity analysis of large graphs.
26pub struct FastPercolation {
27    parent: Vec<usize>,
28    rank: Vec<usize>,
29    size: Vec<usize>,
30}
31
32impl FastPercolation {
33    /// Creates a new percolation structure for n nodes
34    pub fn new(n: usize) -> Self {
35        Self {
36            parent: (0..n).collect(),
37            rank: vec![0; n],
38            size: vec![1; n],
39        }
40    }
41
42    fn find(&mut self, x: usize) -> usize {
43        if self.parent[x] != x {
44            self.parent[x] = self.find(self.parent[x]);
45        }
46        self.parent[x]
47    }
48
49    fn union(&mut self, x: usize, y: usize) {
50        let xr = self.find(x);
51        let yr = self.find(y);
52
53        if xr == yr {
54            return;
55        }
56
57        if self.rank[xr] < self.rank[yr] {
58            self.parent[xr] = yr;
59            self.size[yr] += self.size[xr];
60        } else if self.rank[xr] > self.rank[yr] {
61            self.parent[yr] = xr;
62            self.size[xr] += self.size[yr];
63        } else {
64            self.parent[yr] = xr;
65            self.rank[xr] += 1;
66            self.size[xr] += self.size[yr];
67        }
68    }
69
70    /// Computes rigidity metrics for a graph
71    ///
72    /// Uses Laman's theorem to determine if a graph is minimally rigid.
73    /// A graph with V vertices is minimally rigid if it has exactly 2V-3 edges
74    /// and every subgraph with k vertices has at most 2k-3 edges.
75    ///
76    /// # Arguments
77    ///
78    /// * `edges` - Slice of (u, v) tuples representing edges
79    /// * `n_nodes` - Total number of nodes in the graph
80    ///
81    /// # Returns
82    ///
83    /// RigidityResult containing comprehensive rigidity metrics
84    pub fn compute_rigidity(&mut self, edges: &[(usize, usize)], n_nodes: usize) -> RigidityResult {
85        for &(u, v) in edges {
86            if u < n_nodes && v < n_nodes {
87                self.union(u, v);
88            }
89        }
90
91        let mut clusters: HashMap<usize, usize> = HashMap::new();
92        for i in 0..n_nodes {
93            let root = self.find(i);
94            *clusters.entry(root).or_insert(0) += 1;
95        }
96
97        let n_clusters = clusters.len();
98        let n_edges = edges.len();
99
100        // Laman's theorem (2V - 3) only applies for V >= 3
101        // For smaller graphs, use special cases
102        let (is_rigid, rank, deficiency) = if n_nodes < 3 {
103            // Graphs with < 3 vertices are trivially non-rigid
104            (false, n_edges, if n_edges > 0 { n_edges } else { 0 })
105        } else {
106            let expected_edges = 2 * n_nodes - 3;
107            let is_rigid = n_edges >= expected_edges;
108            let rank = n_edges.min(expected_edges);
109            let deficiency = if n_edges >= 2 * n_nodes - 2 {
110                n_edges - expected_edges
111            } else {
112                expected_edges - n_edges
113            };
114            (is_rigid, rank, deficiency)
115        };
116
117        let rigid_nodes: usize = clusters.values().filter(|&&s| s >= 3).sum();
118
119        let rigid_fraction = rigid_nodes as f32 / n_nodes as f32;
120
121        RigidityResult {
122            is_rigid,
123            rank,
124            deficiency,
125            n_clusters,
126            rigid_fraction,
127        }
128    }
129}
130
131#[cfg(test)]
132mod tests {
133    use super::*;
134
135    #[test]
136    fn test_percolation() {
137        let mut perc = FastPercolation::new(5);
138        let edges = [(0, 1), (1, 2), (2, 3), (3, 4)];
139        let result = perc.compute_rigidity(&edges, 5);
140
141        assert!(result.rigid_fraction > 0.0);
142    }
143}