la_stack/
lib.rs

1#![forbid(unsafe_code)]
2#![warn(missing_docs)]
3#![doc = include_str!("../README.md")]
4
5#[cfg(doc)]
6mod readme_doctests {
7    //! Executable versions of README examples.
8    /// ```rust
9    /// use la_stack::prelude::*;
10    ///
11    /// // This system requires pivoting (a[0][0] = 0), so it's a good LU demo.
12    /// let a = Matrix::<5>::from_rows([
13    ///     [0.0, 1.0, 1.0, 1.0, 1.0],
14    ///     [1.0, 0.0, 1.0, 1.0, 1.0],
15    ///     [1.0, 1.0, 0.0, 1.0, 1.0],
16    ///     [1.0, 1.0, 1.0, 0.0, 1.0],
17    ///     [1.0, 1.0, 1.0, 1.0, 0.0],
18    /// ]);
19    ///
20    /// let b = Vector::<5>::new([14.0, 13.0, 12.0, 11.0, 10.0]);
21    ///
22    /// let lu = a.lu(DEFAULT_PIVOT_TOL).unwrap();
23    /// let x = lu.solve_vec(b).unwrap().into_array();
24    ///
25    /// // Floating-point rounding is expected; compare with a tolerance.
26    /// let expected = [1.0, 2.0, 3.0, 4.0, 5.0];
27    /// for (x_i, e_i) in x.iter().zip(expected.iter()) {
28    ///     assert!((*x_i - *e_i).abs() <= 1e-12);
29    /// }
30    /// ```
31    fn solve_5x5_example() {}
32
33    /// ```rust
34    /// use la_stack::prelude::*;
35    ///
36    /// // This matrix is symmetric positive-definite (A = L*L^T) so LDLT works without pivoting.
37    /// let a = Matrix::<5>::from_rows([
38    ///     [1.0, 1.0, 0.0, 0.0, 0.0],
39    ///     [1.0, 2.0, 1.0, 0.0, 0.0],
40    ///     [0.0, 1.0, 2.0, 1.0, 0.0],
41    ///     [0.0, 0.0, 1.0, 2.0, 1.0],
42    ///     [0.0, 0.0, 0.0, 1.0, 2.0],
43    /// ]);
44    ///
45    /// let det = a.ldlt(DEFAULT_SINGULAR_TOL).unwrap().det();
46    /// assert!((det - 1.0).abs() <= 1e-12);
47    /// ```
48    fn det_5x5_ldlt_example() {}
49}
50
51mod ldlt;
52mod lu;
53mod matrix;
54mod vector;
55
56use core::fmt;
57
58/// Default absolute threshold used for singularity/degeneracy detection.
59///
60/// This is intentionally conservative for geometric predicates and small systems.
61///
62/// Conceptually, this is an absolute bound for deciding when a scalar should be treated
63/// as "numerically zero" (e.g. LU pivots, LDLT diagonal entries).
64pub const DEFAULT_SINGULAR_TOL: f64 = 1e-12;
65
66/// Default absolute pivot magnitude threshold used for LU pivot selection / singularity detection.
67///
68/// This name is kept for backwards compatibility; prefer [`DEFAULT_SINGULAR_TOL`] when the
69/// tolerance is not specifically about pivot selection.
70pub const DEFAULT_PIVOT_TOL: f64 = DEFAULT_SINGULAR_TOL;
71
72/// Linear algebra errors.
73#[derive(Clone, Copy, Debug, PartialEq, Eq)]
74pub enum LaError {
75    /// The matrix is (numerically) singular.
76    Singular {
77        /// The factorization column/step where a suitable pivot/diagonal could not be found.
78        pivot_col: usize,
79    },
80    /// A non-finite value (NaN/∞) was encountered.
81    NonFinite {
82        /// The column being processed when non-finite values were detected.
83        pivot_col: usize,
84    },
85}
86
87impl fmt::Display for LaError {
88    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
89        match *self {
90            Self::Singular { pivot_col } => {
91                write!(f, "singular matrix at pivot column {pivot_col}")
92            }
93            Self::NonFinite { pivot_col } => {
94                write!(
95                    f,
96                    "non-finite value encountered at pivot column {pivot_col}"
97                )
98            }
99        }
100    }
101}
102
103impl std::error::Error for LaError {}
104
105pub use ldlt::Ldlt;
106pub use lu::Lu;
107pub use matrix::Matrix;
108pub use vector::Vector;
109
110/// Common imports for ergonomic usage.
111///
112/// This prelude re-exports the primary types and constants: [`Matrix`], [`Vector`], [`Lu`],
113/// [`Ldlt`], [`LaError`], [`DEFAULT_PIVOT_TOL`], and [`DEFAULT_SINGULAR_TOL`].
114pub mod prelude {
115    pub use crate::{DEFAULT_PIVOT_TOL, DEFAULT_SINGULAR_TOL, LaError, Ldlt, Lu, Matrix, Vector};
116}
117
118#[cfg(test)]
119mod tests {
120    use super::*;
121
122    use approx::assert_abs_diff_eq;
123
124    #[test]
125    fn default_singular_tol_is_expected() {
126        assert_abs_diff_eq!(DEFAULT_SINGULAR_TOL, 1e-12, epsilon = 0.0);
127        assert_abs_diff_eq!(DEFAULT_PIVOT_TOL, DEFAULT_SINGULAR_TOL, epsilon = 0.0);
128    }
129
130    #[test]
131    fn laerror_display_formats_singular() {
132        let err = LaError::Singular { pivot_col: 3 };
133        assert_eq!(err.to_string(), "singular matrix at pivot column 3");
134    }
135
136    #[test]
137    fn laerror_display_formats_nonfinite() {
138        let err = LaError::NonFinite { pivot_col: 2 };
139        assert_eq!(
140            err.to_string(),
141            "non-finite value encountered at pivot column 2"
142        );
143    }
144
145    #[test]
146    fn laerror_is_std_error_with_no_source() {
147        let err = LaError::Singular { pivot_col: 0 };
148        let e: &dyn std::error::Error = &err;
149        assert!(e.source().is_none());
150    }
151
152    #[test]
153    fn prelude_reexports_compile_and_work() {
154        use crate::prelude::*;
155
156        // Use the items so we know they are in scope and usable.
157        let m = Matrix::<2>::identity();
158        let v = Vector::<2>::new([1.0, 2.0]);
159        let _ = m.lu(DEFAULT_PIVOT_TOL).unwrap().solve_vec(v).unwrap();
160        let _ = m.ldlt(DEFAULT_SINGULAR_TOL).unwrap().solve_vec(v).unwrap();
161    }
162}
la_stack/lib.rs

la_stack/
lib.rs
