faer-precond 0.2.0

Numerical preconditioners (Jacobi, block-Jacobi, ILU(0), IC(0)) for iterative solvers, built on faer.
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37

//! Triangular solve dispatch for ILU(0).

use faer::sparse::linalg::triangular_solve;
use faer::{Conj, MatMut, Par};
use faer_traits::{ComplexField, Index};

use super::numeric::Ilu0;

/// Apply `M^{-1} = U^{-1} L^{-1}` (forward then back) to `rhs` in place.
pub(crate) fn solve_in_place<I, T>(ilu: &Ilu0<I, T>, conj: Conj, rhs: MatMut<'_, T>, par: Par)
where
    I: Index,
    T: ComplexField,
{
    let l = ilu.l_view();
    let u = ilu.u_view();
    let mut rhs = rhs;
    triangular_solve::solve_unit_lower_triangular_in_place(l, conj, rhs.as_mut(), par);
    triangular_solve::solve_upper_triangular_in_place(u, conj, rhs.as_mut(), par);
}

/// Apply `M^{-T} = L^{-T} U^{-T}` (so first solve U^T then L^T) to `rhs` in place.
pub(crate) fn solve_transpose_in_place<I, T>(
    ilu: &Ilu0<I, T>,
    conj: Conj,
    rhs: MatMut<'_, T>,
    par: Par,
) where
    I: Index,
    T: ComplexField,
{
    let l = ilu.l_view();
    let u = ilu.u_view();
    let mut rhs = rhs;
    triangular_solve::solve_upper_triangular_transpose_in_place(u, conj, rhs.as_mut(), par);
    triangular_solve::solve_unit_lower_triangular_transpose_in_place(l, conj, rhs.as_mut(), par);
}