HLIBpro  2.6
Algebra

Classes

struct  fac_options_t
 options for matrix factorisations More...
 
class  TLU
 Computes LU factorisation \( A = LU \). More...
 
class  TLDU
 Computes LDU factorisation \( A = LDU \). More...
 
class  TLL
 computes Cholesky factorisation \( A = LL^T \) or \( A=LL^H \) More...
 
class  TLDL
 computes LDL factorisation \( A = LDL^T \) or \( A = LDL^H \) More...
 
struct  inv_options_t
 options for matrix inversion More...
 
struct  solve_option_t
 options for how to solve with given matrix More...
 
struct  eval_option_t
 options for how to evaluate given matrix More...
 
class  TLowRankApx
 base class for all low rank approximation techniques More...
 
class  TZeroLRApx
 Approximate all low-rank blocks by zero, e.g. for nearfield only. More...
 
class  TDenseLRApx< T >
 Computes dense matrix block without approximation. More...
 
class  TSVDLRApx< T >
 Uses exact SVD to compute low rank approximation (WARNING: O(n³) complexity) More...
 
class  TRandSVDLRApx< T >
 Uses randomized SVD to compute low rank approximation (WARNING: O(n²) complexity) More...
 
class  TRRQRLRApx< T >
 Uses rank-revealing QR to compute low rank approximation. More...
 
class  TACA< T >
 Defines interface for all ACA algorithms and implements classical ACA. More...
 
class  TACAPlus< T >
 Implements ACA+, which corrects some of the deficits of the original ACA algorithm. More...
 
class  TACAFull< T >
 ACA with full pivot search (complexity: O(n²)) More...
 
class  THCA< T >
 uses hybrid cross approximation (HCA) for computing low rank approximation More...
 
class  TPermHCAGeneratorFn< T_val >
 base class for HCA generator functions using row/column permutations More...
 

Functions

void eval_diag (const TMatrix *A, TVector *v, const matop_t op, const eval_option_t &options)
 evaluate D·x=y with (block) diagonal D More...
 
void eval (const TMatrix *A, TVector *v, const matop_t op, const eval_option_t &options)
 evaluate L·U·x=y with lower triangular L and upper triangular U More...
 
void eval_lower (const TMatrix *A, TVector *v, const matop_t op, const eval_option_t &options)
 evaluate A·x=y with lower triangular A More...
 
void eval_upper (const TMatrix *A, TVector *v, const matop_t op, const eval_option_t &options)
 evaluate A·x = y with upper triangular A More...
 
void invert (TMatrix *A, const TTruncAcc &acc, const inv_options_t &opts=inv_options_t())
 compute \( A^{-1} \) More...
 
void invert_ll (TMatrix *A, const TTruncAcc &acc, const inv_options_t &opts=inv_options_t())
 compute \( A^{-1} \) More...
 
void invert_ll (TMatrix *L, TMatrix *X, const TTruncAcc &acc, const inv_options_t &opts=inv_options_t())
 compute \( X = L^{-1} \) More...
 
void invert_ur (TMatrix *A, const TTruncAcc &acc, const inv_options_t &opts=inv_options_t())
 compute \( A^{-1} \) More...
 
void invert_diag (TMatrix *A, const TTruncAcc &acc, const inv_options_t &opts=inv_options_t())
 compute \( A^{-1} \) More...
 
TVectorinverse_diag (TMatrix *A, const TTruncAcc &acc, const inv_options_t &opts=inv_options_t())
 Compute diagonal of \( A^{-1} \). More...
 
void gauss_elim (TMatrix *A, TMatrix *C, const TTruncAcc &acc, const inv_options_t &opts=inv_options_t())
 Compute \(A^{-1} \) via Gaussian elimination. More...
 
void gauss_elim (TMatrix *A, const TTruncAcc &acc, const inv_options_t &opts=inv_options_t())
 Compute \(A^{-1} \) via Gaussian elimination. More...
 
void invert (TDenseMatrix *A)
 compute \( A^{-1} \) More...
 
std::unique_ptr< TDenseMatrixinverse (const TDenseMatrix *A)
 compute and return \( A^{-1} \) More...
 
template<typename value_t >
void multiply_diag (const value_t alpha, const matop_t op_A, const TMatrix *A, const matop_t op_D, const TMatrix *D, const matop_t op_B, const TMatrix *B, const value_t beta, TMatrix *C, const TTruncAcc &acc, TProgressBar *progress=NULL)
 compute C ≔ β·C + α·op(A)·op(D)·op(B) with (block) diagonal matrix D
 
template<typename value_t >
void multiply_diag_accu (const value_t alpha, const matop_t op_A, const TMatrix *A, const matop_t op_D, const TMatrix *D, const matop_t op_B, const TMatrix *B, TMatrix *C, const TTruncAcc &acc)
 compute C ≔ C + α·op(A)·op(D)·op(B) with (block) diagonal matrix D using accumulators
 
template<typename value_t >
std::unique_ptr< TMatrixmultiply_diag (const value_t alpha, const matop_t op_A, const TMatrix *A, const matop_t op_D, const TMatrix *D, const matop_t op_B, const TMatrix *B)
 compute C ≔ α·op(A)·op(D)·op(B) with (block) diagonal matrix D More...
 
size_t multiply_diag_steps (const matop_t op_A, const TMatrix *A, const matop_t op_D, const TMatrix *D, const matop_t op_B, const TMatrix *B, const TMatrix *C)
 return number of steps for computing C ≔ C + op(A)·op(D)·op(B)
 
void mul_diag_left (const TScalarVector &v, TMatrix *A)
 compute B = diag(v)·A and overwrite A
 
void mul_diag_right (TMatrix *A, const TScalarVector &v)
 compute B = A·diag(v) and overwrite A
 

Fourier Transformation

Functions for the forward and backward Fourier transformation of vectors.

void fft (TVector *v)
 perform FFT for vector v (inplace)
 
void ifft (TVector *v)
 perform inverse FFT for vector v (inplace)
 

Matrix Addition

Functions for matrix addition

template<typename T_value >
void add (const T_value alpha, const TMatrix *A, const T_value beta, TMatrix *C, const TTruncAcc &acc)
 C ≔ α·A + β·C. More...
 
template<typename T_value >
void add_identity (TMatrix *A, const T_value lambda)
 compute A ≔ A + λ·I
 

Matrix Factorisation

Functions related to matrix factorisation, e.g. LU or Cholesky factorisation.

std::unique_ptr< TFacMatrix > factorise (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 compute factorisation of A More...
 
std::unique_ptr< TFacInvMatrixfactorise_inv (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 compute factorisation of A and return inverse operator More...
 
void lu (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 compute LU factorisation
 
std::unique_ptr< TFacInvMatrixlu_inv (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 compute LU factorisation and return inverse operator
 
void ldu (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 compute LDU factorisation
 
std::unique_ptr< TFacInvMatrixldu_inv (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 compute LDU factorisation and return inverse operator
 
void ll (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 compute Cholesky factorisation
 
std::unique_ptr< TFacInvMatrixll_inv (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 compute Cholesky factorisation and return inverse operator
 
void ldl (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 compute LDL^H factorisation
 
std::unique_ptr< TFacInvMatrixldl_inv (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 compute LDL factorisation and return inverse operator
 
void test_convert_dense (TBlockMatrix *B, const uint i, const uint j)
 test, if given matrix B_ij is lowrank with too large rank and convert to dense
 
void chol (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 
std::unique_ptr< TFacInvMatrixchol_inv (TMatrix *A, const TTruncAcc &acc, const fac_options_t &options=fac_options_t())
 

Matrix Multiplication

Functions for matrix multiplication.

template<typename T_value >
void multiply (const T_value alpha, const matop_t op_A, const TMatrix *A, const matop_t op_B, const TMatrix *B, const T_value beta, TMatrix *C, const TTruncAcc &acc, TProgressBar *progress=NULL)
 compute C ≔ β·C + α·op(A)·op(B) More...
 
template<typename T_value >
void multiply (const T_value alpha, const TMatrix *A, const TMatrix *B, const T_value beta, TMatrix *C, const TTruncAcc &acc, TProgressBar *progress=NULL)
 compute C = β·C + α·A·B More...
 
template<typename value_t >
void multiply_accu (const value_t alpha, const matop_t op_A, const TMatrix *A, const matop_t op_B, const TMatrix *B, const value_t beta, TMatrix *C, const TTruncAcc &acc)
 compute C ≔ C + α·op(A)·op(B) using accumulators
 
template<typename T_value >
std::unique_ptr< TMatrixmultiply (const T_value alpha, const matop_t op_A, const TMatrix *A, const matop_t op_B, const TMatrix *B)
 compute C ≔ β·C + α·op(A)·op(B) More...
 
void multiply_ll_left (const real alpha, const TMatrix *L, TMatrix *A, const TTruncAcc &acc, const eval_option_t &opts)
 compute A ≔ α·L·A with lower left tridiagonal L
 
void multiply_ll_right (const real alpha, TMatrix *A, const TMatrix *L, const TTruncAcc &acc, const eval_option_t &opts)
 compute A ≔ α·A·L with lower left tridiagonal L
 
void multiply_ll_right (const real alpha, const TMatrix *A, const TMatrix *B, TMatrix *C, const TTruncAcc &acc, const eval_option_t &opts)
 compute C ≔ C + α·A·B with lower left tridiagonal B
 
void multiply_ur_left (const real alpha, const TMatrix *U, TMatrix *A, const TTruncAcc &acc, const eval_option_t &opts)
 compute A ≔ α·U·A with upper right tridiagonal U
 
void multiply_ur_right (const real alpha, TMatrix *A, const TMatrix *U, const TTruncAcc &acc, const eval_option_t &opts)
 compute A ≔ α·A·U with upper right tridiagonal U
 
void multiply_ur_ll (const TMatrix *U, const TMatrix *L, TMatrix *A, const TTruncAcc &acc)
 Compute A ≔ A + U·L. More...
 
void multiply_llt_d_ll (const TMatrix *L, const TMatrix *D, TMatrix *C, const matop_t op_L, const TTruncAcc &acc)
 Compute C ≔ C + L^T·D·L. More...
 
void multiply_llt_d_left (const real alpha, const TMatrix *L, const TMatrix *D, TMatrix *B, const matop_t op_L, const TTruncAcc &acc)
 Compute B ≔ op(L)·D·B inplace, overwriting B. More...
 
void multiply_llt_d_left (const real alpha, const TMatrix *L, const TMatrix *D, const TMatrix *B, TMatrix *C, const matop_t op_L, const TTruncAcc &acc)
 Compute C ≔ C + op(L)·D·B. More...
 
template<typename T_value >
void multiply (const T_value alpha, const TMatrixView &A, const TMatrix *B, const T_value beta, TMatrix *C, const TTruncAcc &acc, TProgressBar *progress=NULL)
 
template<typename T_value >
void multiply (const T_value alpha, const TMatrix *A, const TMatrixView &B, const T_value beta, TMatrix *C, const TTruncAcc &acc, TProgressBar *progress=NULL)
 
template<typename T_value >
void multiply (const T_value alpha, const TMatrixView &A, const TMatrixView &B, const T_value beta, TMatrix *C, const TTruncAcc &acc, TProgressBar *progress=NULL)
 
template<typename value_t >
void add_product (const value_t alpha, const matop_t op_A, const TMatrix *A, const matop_t op_B, const TMatrix *B, TMatrix *C, const TTruncAcc &acc, const bool lazy=CFG::Arith::lazy_eval)
 
template<typename T_value >
std::unique_ptr< TMatrixmultiply (const T_value alpha, const matop_t op_A, const TRkMatrix *A, const matop_t op_B, const TMatrix *B)
 
template<typename T_value >
std::unique_ptr< TMatrixmultiply (const T_value alpha, const matop_t op_A, const TMatrix *A, const matop_t op_B, const TRkMatrix *B)
 
template<typename T_value >
std::unique_ptr< TMatrixmultiply (const T_value alpha, const matop_t op_A, const TDenseMatrix *A, const matop_t op_B, const TMatrix *B)
 
template<typename T_value >
std::unique_ptr< TMatrixmultiply (const T_value alpha, const matop_t op_A, const TMatrix *A, const matop_t op_B, const TDenseMatrix *B)
 
size_t multiply_steps (const matop_t op_A, const TMatrix *A, const matop_t op_B, const TMatrix *B, const TMatrix *C)
 
size_t multiply_ll_left_steps (const TMatrix *A, const TMatrix *B)
 
size_t multiply_ll_right_steps (const TMatrix *A, const TMatrix *B)
 
size_t multiply_ur_left_steps (const TMatrix *A, const TMatrix *B)
 
size_t multiply_ur_right_steps (const TMatrix *A, const TMatrix *B)
 
size_t multiply_ur_ll_steps (const TMatrix *A)
 
size_t multiply_llt_d_ll_steps (const TMatrix *L)
 
size_t multiply_llt_d_left_steps (const TMatrix *L, const TMatrix *B)
 

Matrix-Vector Multiplication

Functions for (parallel) matrix-vector multiplication

void mul_vec (const TProcSet &ps, const real alpha, const TMatrix *A, const TVector *x, const real beta, TVector *y, const matop_t op)
 compute y ≔ α·A·x + β·y on (possibly) distributed matrix A on all processors in ps More...
 
void cmul_vec (const TProcSet &ps, const complex alpha, const TMatrix *A, const TVector *x, const complex beta, TVector *y, const matop_t op)
 same as mul_vec but with complex valued scalars ( More...
 
void mul_vec_diag (const real alpha, const matop_t op_A, const TMatrix *A, const matop_t op_D, const TMatrix *D, const TVector *x, const real beta, TVector *y)
 compute y ≔ α·A·D·x + β·y with diagonal matrix D More...
 

Detailed Description

This module provides most higher level algebra functions, e.g. matrix multiplication, inversion and factorisation. See also Basic Matrix Algebra and Matrix Factorisation for an introduction into H-arithmetic.

To include all algebra functions and classes add

#include <hlib-alg.hh>

to your source files.

Function Documentation

◆ add()

void HLIB::add ( const T_value  alpha,
const TMatrix A,
const T_value  beta,
TMatrix C,
const TTruncAcc acc 
)

The function computes the sum \( C := \alpha A + \beta C \) with a predefined accuracy acc. Thread parallel execution is supported.

Parameters
alphascaling factor for A
Aupdate matrix for C
betascaling factor for C
Cmatrix to update
accaccuracy of summation

◆ add_product()

void HLIB::add_product ( const value_t  alpha,
const matop_t  op_A,
const TMatrix A,
const matop_t  op_B,
const TMatrix B,
TMatrix C,
const TTruncAcc acc,
const bool  lazy = CFG::Arith::lazy_eval 
)

start accumulator based matrix multiplication C = C + α·A·B

◆ cmul_vec()

void HLIB::cmul_vec ( const TProcSet ps,
const complex  alpha,
const TMatrix A,
const TVector x,
const complex  beta,
TVector y,
const matop_t  op 
)
See also
mul_vec)

◆ eval()

void HLIB::eval ( const TMatrix A,
TVector v,
const matop_t  op,
const eval_option_t options 
)
  • L and U are both stored in A
  • on entry: v = x
  • on exit : v = y

◆ eval_diag()

void HLIB::eval_diag ( const TMatrix A,
TVector v,
const matop_t  op,
const eval_option_t options 
)
  • on entry: v = x
  • on exit : v = y

◆ eval_lower()

void HLIB::eval_lower ( const TMatrix A,
TVector v,
const matop_t  op,
const eval_option_t options 
)
  • on entry: v = x
  • on exit : v = y

◆ eval_upper()

void HLIB::eval_upper ( const TMatrix A,
TVector v,
const matop_t  op,
const eval_option_t options 
)
  • on entry: v = x
  • on exit : v = y

◆ factorise()

std::unique_ptr< TFacMatrix > HLIB::factorise ( TMatrix A,
const TTruncAcc acc,
const fac_options_t options = fac_options_t() 
)

Compute triangular factorisation of A while choosing appropriate factorisation method depending on the format of A, e.g. if unsymmetric, symmetric or hermitian.

The return value is an operator object representing the factorised form and suitable for evaluation, e.g. matrix-vector multiplication (see TFacMatrix). A will be overwritten with the actual factorisation data.

◆ factorise_inv()

std::unique_ptr< TFacInvMatrix > HLIB::factorise_inv ( TMatrix A,
const TTruncAcc acc,
const fac_options_t options = fac_options_t() 
)

Compute triangular factorisation of A as in factorise but instead of an operator for evaluation of the factorised A, an operator for evaluation of the inverse of A is returned (see TFacInvMatrix).

◆ gauss_elim() [1/2]

void HLIB::gauss_elim ( TMatrix A,
TMatrix C,
const TTruncAcc acc,
const inv_options_t opts = inv_options_t() 
)

A is overwritten by \(A^{-1} \). C is optional and may be used during computation as temporary space. If present, it should have same format as A.

◆ gauss_elim() [2/2]

void HLIB::gauss_elim ( TMatrix A,
const TTruncAcc acc,
const inv_options_t opts = inv_options_t() 
)

A is overwritten by \(A^{-1} \). The elimination is performed in-place, although local copy operations are needed. (ATTENTION: experimental)

◆ inverse()

std::unique_ptr< TDenseMatrix > HLIB::inverse ( const TDenseMatrix A)

Compute inverse of \( A \) (dense matrix version) but do not modify \( A \).

◆ inverse_diag()

TVector* HLIB::inverse_diag ( TMatrix A,
const TTruncAcc acc,
const inv_options_t opts = inv_options_t() 
)

Compute only the diagonal of A and return the resulting vector. A is modified during computation.

◆ invert() [1/2]

void HLIB::invert ( TMatrix A,
const TTruncAcc acc,
const inv_options_t opts = inv_options_t() 
)

Compute inverse of \( A \) with block-wise accuracy acc. A is overwritten by \( A^{-1} \) during inversion.

◆ invert() [2/2]

void HLIB::invert ( TDenseMatrix A)

Compute inverse of \( A \) (dense matrix version). A is overwritten by \( A^{-1} \) during inversion.

◆ invert_diag()

void HLIB::invert_diag ( TMatrix A,
const TTruncAcc acc,
const inv_options_t opts = inv_options_t() 
)

Compute inverse of diagonal matrix \( A \) with block-wise accuracy acc. A is overwritten by \( A^{-1} \) during inversion.

◆ invert_ll() [1/2]

void HLIB::invert_ll ( TMatrix A,
const TTruncAcc acc,
const inv_options_t opts = inv_options_t() 
)

Compute inverse of lower-left triangular matrix \( A \) with block-wise accuracy acc. A is overwritten by \( A^{-1} \) during inversion.

◆ invert_ll() [2/2]

void HLIB::invert_ll ( TMatrix L,
TMatrix X,
const TTruncAcc acc,
const inv_options_t opts = inv_options_t() 
)

Compute inverse of lower-left triangular matrix \( L \) with block-wise accuracy acc. L may be modified during inversion.

◆ invert_ur()

void HLIB::invert_ur ( TMatrix A,
const TTruncAcc acc,
const inv_options_t opts = inv_options_t() 
)

Compute inverse of upper-right triangular matrix \( A \) with block-wise accuracy acc. A is overwritten by \( A^{-1} \) during inversion.

◆ mul_vec()

void HLIB::mul_vec ( const TProcSet ps,
const real  alpha,
const TMatrix A,
const TVector x,
const real  beta,
TVector y,
const matop_t  op 
)
Parameters
psprocessor set defining processors involved in computation
alphascaling factor for multiplication
Amatrix to multiply with
xvector to multiply with
betascaling factor of updated vector
yvector to update with multiplication result
opdefines transformation of matrix, e.g. transposed, adjoint (
See also
matop_t)

◆ mul_vec_diag()

void HLIB::mul_vec_diag ( const real  alpha,
const matop_t  op_A,
const TMatrix A,
const matop_t  op_D,
const TMatrix D,
const TVector x,
const real  beta,
TVector y 
)
Parameters
alphascaling factor for update
op_Atransformation of matrix A, e.g. transposed, adjoint (
See also
matop_t)
Parameters
Aarbitrary matrix
op_Dtransformation of matrix D
Ddiagonal matrix
xsource vector
betascaling factor for destination
ydestination vector

◆ multiply() [1/3]

void HLIB::multiply ( const T_value  alpha,
const matop_t  op_A,
const TMatrix A,
const matop_t  op_B,
const TMatrix B,
const T_value  beta,
TMatrix C,
const TTruncAcc acc,
TProgressBar *  progress = NULL 
)

The function computes to matrix product \(C := \beta C + \alpha \tilde A \tilde B\) where \(\tilde A\) and \(\tilde B\) may be the non-modified, transposed or adjoint matrices \(A\) and \(B\) respectively.

The result of the multiplication is written to \(C\), whereby the block structure of C is not changed, e.g. the resulting block structure of the product is defined by C.

Thread-parallel execution is supported by the matrix multiplication.

This function is available in various versions without corresponding parameters, e.g. without op_A, op_B.

Parameters
alphascaling factor of product
op_Amatrix modifier for A
Afirst matrix factor
op_Bmatrix modifier for B
Bsecond matrix factor
betascaling factor for C
Cmatrix to update
accaccuracy of multiplication
progressoptional progress bar

◆ multiply() [2/3]

void HLIB::multiply ( const T_value  alpha,
const TMatrix A,
const TMatrix B,
const T_value  beta,
TMatrix C,
const TTruncAcc acc,
TProgressBar *  progress = NULL 
)
inline

matrices are provided as matrix view (or in combination with matrix)

◆ multiply() [3/3]

std::unique_ptr< TMatrix > HLIB::multiply ( const T_value  alpha,
const matop_t  op_A,
const TMatrix A,
const matop_t  op_B,
const TMatrix B 
)

The function computes to matrix product \(\alpha \tilde A \tilde B\) where \(\tilde A\) and \(\tilde B\) may be the non-modified, transposed or adjoint matrices \(A\) and \(B\) respectively.

The output format is based on the format of \(A\) or \(B\), which assumes that at least one of those is a non-blocked matrix. Otherwise, an exception is thrown!

Parameters
alphascaling factor of product
op_Amatrix modifier for A
Afirst matrix factor
op_Bmatrix modifier for B
Bsecond matrix factor

◆ multiply_diag()

std::unique_ptr< TMatrix > HLIB::multiply_diag ( const value_t  alpha,
const matop_t  op_A,
const TMatrix A,
const matop_t  op_D,
const TMatrix D,
const matop_t  op_B,
const TMatrix B 
)
  • not all matrices must be blocked!

◆ multiply_llt_d_left() [1/2]

void HLIB::multiply_llt_d_left ( const real  alpha,
const TMatrix L,
const TMatrix D,
TMatrix B,
const matop_t  op_L,
const TTruncAcc acc 
)

Compute B ≔ op(L)·D·B with lower left tridiagonal L and diagonal D.

◆ multiply_llt_d_left() [2/2]

void HLIB::multiply_llt_d_left ( const real  alpha,
const TMatrix L,
const TMatrix D,
const TMatrix B,
TMatrix C,
const matop_t  op_L,
const TTruncAcc acc 
)

Compute C ≔ C + op(L)·D·B with lower left tridiagonal L, diagonal D and general B

◆ multiply_llt_d_ll()

void HLIB::multiply_llt_d_ll ( const TMatrix L,
const TMatrix D,
TMatrix C,
const matop_t  op_L,
const TTruncAcc acc 
)

Compute C ≔ C + L^T·D·L with lower left tridiagonal L and diagonal D.

◆ multiply_ur_ll()

void HLIB::multiply_ur_ll ( const TMatrix U,
const TMatrix L,
TMatrix A,
const TTruncAcc acc 
)

Compute A ≔ A + U·L with lower left tridiagonal L and upper right tridiagonal U.