Skip to main content

faer/sparse/linalg/
mod.rs

1use crate::internal_prelude_sp::*;
2const CHOLESKY_SUPERNODAL_RATIO_FACTOR: f64 = 40.0;
3const QR_SUPERNODAL_RATIO_FACTOR: f64 = 40.0;
4const LU_SUPERNODAL_RATIO_FACTOR: f64 = 40.0;
5/// tuning parameters for the supernodal factorizations
6#[derive(Copy, Clone, Debug)]
7pub struct SymbolicSupernodalParams<'a> {
8	/// supernode relaxation thresholds
9	///
10	/// let `n` be the total number of columns in two adjacent supernodes.
11	/// let `z` be the fraction of zero entries in the two supernodes if they
12	/// are merged (z includes zero entries from prior amalgamations). the
13	/// two supernodes are merged if:
14	///
15	/// `(n <= relax[0].0 && z < relax[0].1) || (n <= relax[1].0 && z <
16	/// relax[1].1) || ...`
17	pub relax: Option<&'a [(usize, f64)]>,
18}
19const DEFAULT_RELAX: &'static [(usize, f64)] =
20	&[(4, 1.0), (16, 0.8), (48, 0.1), (usize::MAX, 0.05)];
21impl Default for SymbolicSupernodalParams<'_> {
22	#[inline]
23	fn default() -> Self {
24		Self {
25			relax: Some(DEFAULT_RELAX),
26		}
27	}
28}
29/// nonnegative threshold controlling when the supernodal factorization is used
30///
31/// increasing it makes it more likely for the simplicial factorization to be
32/// used, while decreasing it makes it more likely for the supernodal
33/// factorization to be used
34///
35/// `1.0` is the default value
36#[derive(Copy, Clone, Debug, PartialEq)]
37pub struct SupernodalThreshold(pub f64);
38impl Default for SupernodalThreshold {
39	#[inline]
40	fn default() -> Self {
41		Self(1.0)
42	}
43}
44impl SupernodalThreshold {
45	/// determine automatically which variant to select
46	pub const AUTO: Self = Self(1.0);
47	/// simplicial factorization is always selected
48	pub const FORCE_SIMPLICIAL: Self = Self(f64::INFINITY);
49	/// supernodal factorization is always selected
50	pub const FORCE_SUPERNODAL: Self = Self(0.0);
51}
52/// sparse $ll^\top$ error
53#[derive(Copy, Clone, Debug)]
54pub enum LltError {
55	/// numerical error
56	Numeric(linalg::cholesky::llt::factor::LltError),
57	/// non algorithmic error
58	Generic(FaerError),
59}
60impl core::fmt::Display for LltError {
61	#[inline]
62	fn fmt(&self, f: &mut core::fmt::Formatter<'_>) -> core::fmt::Result {
63		core::fmt::Debug::fmt(self, f)
64	}
65}
66impl core::error::Error for LltError {}
67impl From<linalg::cholesky::llt::factor::LltError> for LltError {
68	fn from(value: linalg::cholesky::llt::factor::LltError) -> Self {
69		Self::Numeric(value)
70	}
71}
72/// sparse $lu$ error.
73#[derive(Copy, Clone, Debug)]
74pub enum LuError {
75	/// rank deficient symbolic structure
76	SymbolicSingular {
77		/// iteration at which a pivot could not be found
78		index: usize,
79	},
80	/// rank deficient numeric values: every candidate pivot of a column is
81	/// exactly zero
82	NumericSingular {
83		/// iteration at which the pivot is zero
84		index: usize,
85	},
86	/// non algorithmic error
87	Generic(FaerError),
88}
89impl core::fmt::Display for LuError {
90	#[inline]
91	fn fmt(&self, f: &mut core::fmt::Formatter<'_>) -> core::fmt::Result {
92		core::fmt::Debug::fmt(self, f)
93	}
94}
95impl core::error::Error for LuError {}
96impl<T: Into<FaerError>> From<T> for LltError {
97	fn from(value: T) -> Self {
98		Self::Generic(value.into())
99	}
100}
101impl<T: Into<FaerError>> From<T> for LuError {
102	fn from(value: T) -> Self {
103		Self::Generic(value.into())
104	}
105}
106pub mod amd;
107pub mod cholesky;
108pub mod colamd;
109pub mod lu;
110/// sparse matrix multiplication
111pub mod matmul;
112pub mod qr;
113/// high-level sparse matrix solvers
114#[path = "../solvers.rs"]
115pub mod solvers;
116/// sparse matrix triangular solve
117pub mod triangular_solve;
118mod ghost {
119	use crate::Index;
120	pub use crate::utils::bound::*;
121	pub const NONE_BYTE: u8 = u8::MAX;
122	#[inline]
123	pub fn fill_zero<'n, 'a, I: Index>(
124		slice: &'a mut [I],
125		size: Dim<'n>,
126	) -> &'a mut [Idx<'n, I>] {
127		let len = slice.len();
128		if len > 0 {
129			assert!(*size > 0);
130		}
131		unsafe {
132			core::ptr::write_bytes(slice.as_mut_ptr(), 0u8, len);
133			&mut *(slice as *mut _ as *mut _)
134		}
135	}
136	#[inline]
137	pub fn fill_none<'n, 'a, I: Index>(
138		slice: &'a mut [I::Signed],
139		size: Dim<'n>,
140	) -> &'a mut [MaybeIdx<'n, I>] {
141		let _ = size;
142		let len = slice.len();
143		unsafe { core::ptr::write_bytes(slice.as_mut_ptr(), NONE_BYTE, len) };
144		unsafe { &mut *(slice as *mut _ as *mut _) }
145	}
146	#[inline]
147	pub fn copy_slice<'n, 'a, I: Index>(
148		dst: &'a mut [I],
149		src: &[Idx<'n, I>],
150	) -> &'a mut [Idx<'n, I>] {
151		let dst: &mut [Idx<'_, I>] = unsafe { &mut *(dst as *mut _ as *mut _) };
152		dst.copy_from_slice(src);
153		dst
154	}
155}