rivrs_sparse/ordering/
perm.rs

1//! Permutation construction utilities for sparse matrix orderings.
2//!
3//! Provides [`perm_from_forward`] to bridge ordering algorithm output (forward
4//! permutation arrays) to faer's [`Perm<usize>`](faer::perm::Perm) type by
5//! computing the inverse mapping automatically.
6//!
7//! Reference: Hogg, Duff & Lopez (2020), "A New Sparse LDL^T Solver Using
8//! A Posteriori Threshold Pivoting", SIAM J. Sci. Comput. 42(4) — permutation
9//! handling in the context of fill-reducing orderings (Section 2.1).
10
11use faer::perm::Perm;
12
13use crate::error::SparseError;
14use crate::validate;
15
16/// Construct a [`Perm<usize>`] from a forward permutation array, computing
17/// the inverse automatically.
18///
19/// Takes ownership of the forward array to avoid an extra allocation when
20/// converting to `Box<[usize]>` for faer's `Perm::new_checked`.
21///
22/// # Errors
23///
24/// Returns [`SparseError::InvalidInput`] if the array contains duplicate
25/// indices or out-of-bounds values.
26///
27/// # Examples
28///
29/// ```
30/// use rivrs_sparse::ordering::perm_from_forward;
31///
32/// let fwd = vec![2, 0, 1];
33/// let perm = perm_from_forward(fwd).unwrap();
34/// let (fwd_arr, inv_arr) = perm.as_ref().arrays();
35/// assert_eq!(fwd_arr, &[2, 0, 1]);
36/// assert_eq!(inv_arr, &[1, 2, 0]);
37/// ```
38pub fn perm_from_forward(fwd: Vec<usize>) -> Result<Perm<usize>, SparseError> {
39    let n = fwd.len();
40    validate::validate_permutation(&fwd, n)?;
41
42    // Compute inverse: inv[fwd[i]] = i
43    let mut inv = vec![0usize; n];
44    for (i, &f) in fwd.iter().enumerate() {
45        inv[f] = i;
46    }
47
48    let fwd_box = fwd.into_boxed_slice();
49    let inv_box = inv.into_boxed_slice();
50    Ok(Perm::new_checked(fwd_box, inv_box, n))
51}
52
53#[cfg(test)]
54mod tests {
55    use super::*;
56
57    // ---- perm_from_forward unit tests ----
58
59    #[test]
60    fn identity_permutation() {
61        let fwd = vec![0, 1, 2, 3, 4];
62        let perm = perm_from_forward(fwd).unwrap();
63        let (fwd_arr, inv_arr) = perm.as_ref().arrays();
64        assert_eq!(fwd_arr, &[0, 1, 2, 3, 4]);
65        assert_eq!(inv_arr, &[0, 1, 2, 3, 4]);
66    }
67
68    #[test]
69    fn nontrivial_permutation_round_trip() {
70        let fwd = vec![3, 1, 4, 0, 2];
71        let perm = perm_from_forward(fwd).unwrap();
72        let (fwd_arr, inv_arr) = perm.as_ref().arrays();
73
74        // Verify round-trip: inv[fwd[i]] == i for all i
75        for i in 0..5 {
76            assert_eq!(inv_arr[fwd_arr[i]], i);
77        }
78        // And fwd[inv[i]] == i for all i
79        for i in 0..5 {
80            assert_eq!(fwd_arr[inv_arr[i]], i);
81        }
82    }
83
84    #[test]
85    fn composition_of_two_permutations() {
86        let p1 = perm_from_forward(vec![1, 0, 2]).unwrap();
87        let p2 = perm_from_forward(vec![2, 1, 0]).unwrap();
88        let (p1_fwd, _) = p1.as_ref().arrays();
89        let (p2_fwd, _) = p2.as_ref().arrays();
90
91        // Compose: apply p1 then p2 → composed_fwd[i] = p2_fwd[p1_fwd[i]]
92        let composed: Vec<usize> = (0..3).map(|i| p2_fwd[p1_fwd[i]]).collect();
93        let composed_perm = perm_from_forward(composed).unwrap();
94        let (comp_fwd, comp_inv) = composed_perm.as_ref().arrays();
95
96        // Verify: composition of swap(0,1) then reverse → should be [2, 0, 1]
97        // p1: [1, 0, 2] — swap 0 and 1
98        // p2: [2, 1, 0] — reverse
99        // composed[0] = p2[p1[0]] = p2[1] = 1
100        // composed[1] = p2[p1[1]] = p2[0] = 2
101        // composed[2] = p2[p1[2]] = p2[2] = 0
102        assert_eq!(comp_fwd, &[1, 2, 0]);
103
104        // Verify forward/inverse relationship
105        for i in 0..3 {
106            assert_eq!(comp_inv[comp_fwd[i]], i);
107        }
108    }
109
110    #[test]
111    fn empty_array_dimension_0() {
112        let perm = perm_from_forward(vec![]).unwrap();
113        let (fwd_arr, inv_arr) = perm.as_ref().arrays();
114        assert!(fwd_arr.is_empty());
115        assert!(inv_arr.is_empty());
116    }
117
118    #[test]
119    fn single_element() {
120        let perm = perm_from_forward(vec![0]).unwrap();
121        let (fwd_arr, inv_arr) = perm.as_ref().arrays();
122        assert_eq!(fwd_arr, &[0]);
123        assert_eq!(inv_arr, &[0]);
124    }
125
126    #[test]
127    fn error_on_duplicate_indices() {
128        let result = perm_from_forward(vec![0, 1, 1]);
129        assert!(result.is_err());
130    }
131
132    #[test]
133    fn error_on_out_of_bounds() {
134        let result = perm_from_forward(vec![0, 5, 2]);
135        assert!(result.is_err());
136    }
137
138    #[test]
139    fn forward_inverse_relationship() {
140        let fwd = vec![4, 2, 0, 3, 1];
141        let perm = perm_from_forward(fwd).unwrap();
142        let (fwd_arr, inv_arr) = perm.as_ref().arrays();
143
144        for i in 0..5 {
145            assert_eq!(inv_arr[fwd_arr[i]], i, "inv[fwd[{}]] != {}", i, i);
146            assert_eq!(fwd_arr[inv_arr[i]], i, "fwd[inv[{}]] != {}", i, i);
147        }
148    }
149
150    #[test]
151    fn scale_test_n_10000() {
152        // perm_from_forward at n=10,000, verify round-trip identity
153        let n = 10_000;
154        // Create a permutation: rotate by 1
155        let fwd: Vec<usize> = (0..n).map(|i| (i + 1) % n).collect();
156        let perm = perm_from_forward(fwd).unwrap();
157        let (fwd_arr, inv_arr) = perm.as_ref().arrays();
158
159        // Verify round-trip
160        for i in 0..n {
161            assert_eq!(inv_arr[fwd_arr[i]], i);
162            assert_eq!(fwd_arr[inv_arr[i]], i);
163        }
164    }
165}
rivrs_sparse/ordering/perm.rs

rivrs_sparse/ordering/
perm.rs
