dune-istl 2.10
Loading...
Searching...
No Matches

Functions

template<class M >
void blockDILUDecomposition (M &A, std::vector< typename M::block_type > &Dinv_)
 
template<class M , class X , class Y >
void blockDILUBacksolve (const M &A, const std::vector< typename M::block_type > Dinv_, X &v, const Y &d)
 

Function Documentation

◆ blockDILUBacksolve()

template<class M , class X , class Y >
void Dune::DILU::blockDILUBacksolve ( const M & A,
const std::vector< typename M::block_type > Dinv_,
X & v,
const Y & d )

DILU backsolve

M = (D + L_A) D^-1 (D + U_A) (a LU decomposition of M) where L_A and U_A are the strictly lower and upper parts of A.

Working with residual d = b - Ax and update v = x_{n+1} - x_n solving the product M^-1(d) using upper and lower triangular solve:

v = M^{-1} d = (D + U_A)^{-1} D (D + L_A)^{-1} d

define y = (D + L_A)^{-1} d

lower triangular solve: (D + L_A) y = d upper triangular solve: (D + U_A) v = D y

◆ blockDILUDecomposition()

template<class M >
void Dune::DILU::blockDILUDecomposition ( M & A,
std::vector< typename M::block_type > & Dinv_ )

compute DILU decomposition of A

The preconditioner matrix M has the property

diag(A) = diag(M) = diag((D + L_A) D^{-1} (D + U_A)) = diag(D + L_A D^{-1} U_A)

such that the diagonal matrix D can be constructed:

D_11 = A_11 D_22 = A22 - L_A_{21} D_{11}^{-1} U_A_{12} and etc

we store the inverse of D to be used when applying the preconditioner.

For more details, see: R. Barrett et al., "Templates for the Solution of Linear Systems: Building Blocks for Iterative Methods", 1994. Available from: https://www.netlib.org/templates/templates.pdf