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
51#[cfg(feature = "exact")]
52mod exact;
53mod ldlt;
54mod lu;
55mod matrix;
56mod vector;
57
58use core::fmt;
59
60/// Default absolute threshold used for singularity/degeneracy detection.
61///
62/// This is intentionally conservative for geometric predicates and small systems.
63///
64/// Conceptually, this is an absolute bound for deciding when a scalar should be treated
65/// as "numerically zero" (e.g. LU pivots, LDLT diagonal entries).
66pub const DEFAULT_SINGULAR_TOL: f64 = 1e-12;
67
68/// Default absolute pivot magnitude threshold used for LU pivot selection / singularity detection.
69///
70/// This name is kept for backwards compatibility; prefer [`DEFAULT_SINGULAR_TOL`] when the
71/// tolerance is not specifically about pivot selection.
72pub const DEFAULT_PIVOT_TOL: f64 = DEFAULT_SINGULAR_TOL;
73
74/// Linear algebra errors.
75#[derive(Clone, Copy, Debug, PartialEq, Eq)]
76pub enum LaError {
77    /// The matrix is (numerically) singular.
78    Singular {
79        /// The factorization column/step where a suitable pivot/diagonal could not be found.
80        pivot_col: usize,
81    },
82    /// A non-finite value (NaN/∞) was encountered.
83    NonFinite {
84        /// The column where a non-finite value was detected.
85        col: usize,
86    },
87}
88
89impl fmt::Display for LaError {
90    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
91        match *self {
92            Self::Singular { pivot_col } => {
93                write!(f, "singular matrix at pivot column {pivot_col}")
94            }
95            Self::NonFinite { col } => {
96                write!(f, "non-finite value encountered at column {col}")
97            }
98        }
99    }
100}
101
102impl std::error::Error for LaError {}
103
104pub use ldlt::Ldlt;
105pub use lu::Lu;
106pub use matrix::Matrix;
107pub use vector::Vector;
108
109/// Common imports for ergonomic usage.
110///
111/// This prelude re-exports the primary types and constants: [`Matrix`], [`Vector`], [`Lu`],
112/// [`Ldlt`], [`LaError`], [`DEFAULT_PIVOT_TOL`], and [`DEFAULT_SINGULAR_TOL`].
113pub mod prelude {
114    pub use crate::{DEFAULT_PIVOT_TOL, DEFAULT_SINGULAR_TOL, LaError, Ldlt, Lu, Matrix, Vector};
115}
116
117#[cfg(test)]
118mod tests {
119    use super::*;
120
121    use approx::assert_abs_diff_eq;
122
123    #[test]
124    fn default_singular_tol_is_expected() {
125        assert_abs_diff_eq!(DEFAULT_SINGULAR_TOL, 1e-12, epsilon = 0.0);
126        assert_abs_diff_eq!(DEFAULT_PIVOT_TOL, DEFAULT_SINGULAR_TOL, epsilon = 0.0);
127    }
128
129    #[test]
130    fn laerror_display_formats_singular() {
131        let err = LaError::Singular { pivot_col: 3 };
132        assert_eq!(err.to_string(), "singular matrix at pivot column 3");
133    }
134
135    #[test]
136    fn laerror_display_formats_nonfinite() {
137        let err = LaError::NonFinite { col: 2 };
138        assert_eq!(err.to_string(), "non-finite value encountered at column 2");
139    }
140
141    #[test]
142    fn laerror_is_std_error_with_no_source() {
143        let err = LaError::Singular { pivot_col: 0 };
144        let e: &dyn std::error::Error = &err;
145        assert!(e.source().is_none());
146    }
147
148    #[test]
149    fn prelude_reexports_compile_and_work() {
150        use crate::prelude::*;
151
152        // Use the items so we know they are in scope and usable.
153        let m = Matrix::<2>::identity();
154        let v = Vector::<2>::new([1.0, 2.0]);
155        let _ = m.lu(DEFAULT_PIVOT_TOL).unwrap().solve_vec(v).unwrap();
156        let _ = m.ldlt(DEFAULT_SINGULAR_TOL).unwrap().solve_vec(v).unwrap();
157    }
158}
la_stack/lib.rs

la_stack/
lib.rs
