diffsol_la/linear_solver/faer/
mod.rs

1pub mod lu;
2pub mod sparse_lu;
3
4#[cfg(test)]
5mod tests {
6    use crate::linear_solver::tests::diagonal_op_n;
7    use crate::{
8        Context, DenseMatrix, FaerContext, FaerMat, FaerSparseMat, FaerVec, LinearSolver, Matrix,
9        Vector,
10    };
11    use faer::Par;
12
13    /// Runs every faer kernel that can take a `Par`: dense/sparse `gemv`, LU factorisation and LU
14    /// solve.
15    fn exercise<M, LS>(ctx: FaerContext)
16    where
17        M: Matrix<T = f64, V = FaerVec<f64>, C = FaerContext>,
18        LS: LinearSolver<M>,
19    {
20        let n = 64;
21        let op = diagonal_op_n::<M>(n, 2.0, ctx);
22        let mut s = LS::default();
23        s.set_sparsity(&op);
24        s.set_linearisation(&op);
25        let b = ctx.vector_from_element::<FaerVec<f64>>(n, 1.0);
26        let x = s.solve(&b).unwrap();
27        let a = M::from_diagonal(&b);
28        let mut y = ctx.vector_zeros::<FaerVec<f64>>(n);
29        a.gemv(1.0, &x, 0.0, &mut y);
30        y.assert_eq_st(&ctx.vector_from_element::<FaerVec<f64>>(n, 0.5), 1e-10);
31    }
32
33    /// `gemv_cols` is a second macro with its own `Par` argument, and dense only.
34    fn exercise_dense(ctx: FaerContext) {
35        let n = 64;
36        let a = FaerMat::<f64>::from_diagonal(&ctx.vector_from_element::<FaerVec<f64>>(n, 2.0));
37        let mut y = ctx.vector_zeros::<FaerVec<f64>>(n);
38        let nc = crate::matrix::MAX_SMALL_COLS;
39        a.gemv_cols(0, nc, 1.0, &vec![1.0; nc], 0.0, &mut y);
40        let expected = (0..n).map(|i| if i < nc { 2.0 } else { 0.0 }).collect();
41        y.assert_eq_st(&ctx.vector_from_vec::<FaerVec<f64>>(expected), 1e-10);
42    }
43
44    /// No faer call site may fall back on faer's process-global parallelism setting: they must all
45    /// read [`FaerContext::par`].  `get_global_parallelism` panics once the global is disabled, so
46    /// a call site that still queries it fails here.
47    ///
48    /// The global is process-wide, so this briefly affects any test running concurrently -- which
49    /// is only safe because no diffsol code path reads it any more.  Restore it before returning.
50    #[test]
51    fn faer_kernels_never_query_the_global_parallelism() {
52        for par in par_settings() {
53            let mut ctx = FaerContext::default();
54            ctx.par = par;
55            faer::disable_global_parallelism();
56            exercise::<FaerMat<f64>, super::lu::LU<f64>>(ctx);
57            exercise::<FaerSparseMat<f64>, super::sparse_lu::FaerSparseLU<f64>>(ctx);
58            exercise_dense(ctx);
59            faer::set_global_parallelism(Par::Seq);
60        }
61    }
62
63    #[cfg(not(target_arch = "wasm32"))]
64    fn par_settings() -> [Par; 2] {
65        [Par::Seq, Par::rayon(2)]
66    }
67
68    #[cfg(target_arch = "wasm32")]
69    fn par_settings() -> [Par; 1] {
70        [Par::Seq]
71    }
72}
diffsol_la/linear_solver/faer/mod.rs

diffsol_la/linear_solver/faer/
mod.rs
