strumpack::HODLR::HODLRMatrix< scalar_t > Class Template Reference

Hierarchically low-rank matrix representation. More...

#include <HODLRMatrix.hpp>

Inheritance diagram for strumpack::HODLR::HODLRMatrix< scalar_t >:
Collaboration diagram for strumpack::HODLR::HODLRMatrix< scalar_t >:

Public Types

using real_t = typename RealType< scalar_t >::value_type
 
using mult_t = typename std::function< void(Trans, const DenseM_t &, DenseM_t &)>
 
using delem_blocks_t = typename std::function< void(VecVec_t &I, VecVec_t &J, std::vector< DistMW_t > &B, ExtractionMeta &)>
 
using elem_blocks_t = typename std::function< void(VecVec_t &I, VecVec_t &J, std::vector< DenseMW_t > &B, ExtractionMeta &)>
 

Public Member Functions

 HODLRMatrix ()
 
 HODLRMatrix (const MPIComm &c, kernel::Kernel< real_t > &K, const opts_t &opts)
 
 HODLRMatrix (const MPIComm &c, const structured::ClusterTree &tree, const std::function< scalar_t(int i, int j)> &Aelem, const opts_t &opts)
 
 HODLRMatrix (const MPIComm &c, const structured::ClusterTree &tree, const elem_blocks_t &Aelem, const opts_t &opts)
 
 HODLRMatrix (const MPIComm &c, const structured::ClusterTree &tree, const std::function< void(Trans op, const DenseM_t &R, DenseM_t &S)> &Amult, const opts_t &opts)
 
 HODLRMatrix (const MPIComm &c, const structured::ClusterTree &tree, const opts_t &opts)
 
template<typename integer_t >
 HODLRMatrix (const MPIComm &c, const structured::ClusterTree &tree, const CSRGraph< integer_t > &graph, const opts_t &opts)
 
 HODLRMatrix (const HODLRMatrix< scalar_t > &h)=delete
 
 HODLRMatrix (HODLRMatrix< scalar_t > &&h)
 
virtual ~HODLRMatrix ()
 
HODLRMatrix< scalar_t > & operator= (const HODLRMatrix< scalar_t > &h)=delete
 
HODLRMatrix< scalar_t > & operator= (HODLRMatrix< scalar_t > &&h)
 
std::size_t rows () const override
 
std::size_t cols () const override
 
std::size_t lrows () const
 
std::size_t local_rows () const override
 
std::size_t begin_row () const override
 
std::size_t end_row () const override
 
const std::vector< int > & dist () const override
 
const MPICommComm () const
 
std::size_t memory () const override
 
std::size_t nonzeros () const override
 
std::size_t total_memory () const
 
std::size_t factor_memory () const
 
std::size_t total_factor_memory () const
 
std::size_t rank () const override
 
std::size_t max_rank () const
 
double get_stat (const std::string &name) const
 
void print_stats ()
 
void compress (const std::function< void(Trans op, const DenseM_t &R, DenseM_t &S)> &Amult)
 
void compress (const std::function< void(Trans op, const DenseM_t &R, DenseM_t &S)> &Amult, int rank_guess)
 
void compress (const delem_blocks_t &Aelem)
 
void compress (const elem_blocks_t &Aelem)
 
void mult (Trans op, const DenseM_t &X, DenseM_t &Y) const override
 
void mult (Trans op, const DistM_t &X, DistM_t &Y) const override
 
void factor () override
 
long long int solve (const DenseM_t &B, DenseM_t &X) const
 
void solve (DenseM_t &B) const override
 
long long int solve (const DistM_t &B, DistM_t &X) const
 
void solve (DistM_t &B) const override
 
long long int inv_mult (Trans op, const DenseM_t &B, DenseM_t &X) const
 
void extract_elements (const VecVec_t &I, const VecVec_t &J, std::vector< DistM_t > &B)
 
void extract_elements (const VecVec_t &I, const VecVec_t &J, std::vector< DenseM_t > &B)
 
void extract_elements (const Vec_t &I, const Vec_t &J, DenseM_t &B)
 
DistM_t dense (const BLACSGrid *g) const
 
void set_sampling_parameter (double sample_param)
 
void set_BACA_block (int bsize)
 
DenseM_t redistribute_2D_to_1D (const DistM_t &R) const
 
void redistribute_2D_to_1D (const DistM_t &R2D, DenseM_t &R1D) const
 
void redistribute_1D_to_2D (const DenseM_t &S1D, DistM_t &S2D) const
 
DenseM_t gather_from_1D (const DenseM_t &A) const
 
DenseM_t all_gather_from_1D (const DenseM_t &A) const
 
DenseM_t scatter_to_1D (const DenseM_t &A) const
 
const std::vector< int > & perm () const
 
const std::vector< int > & iperm () const
 
- Public Member Functions inherited from strumpack::structured::StructuredMatrix< scalar_t >
virtual ~StructuredMatrix ()=default
 
virtual const std::vector< int > & rdist () const
 
virtual const std::vector< int > & cdist () const
 
void mult (Trans op, int m, const scalar_t *x, int ldx, scalar_t *y, int ldy) const
 
virtual void solve (int nrhs, scalar_t *b, int ldb) const
 
virtual void shift (scalar_t s)
 

Friends

template<typename S >
class ButterflyMatrix
 

Detailed Description

template<typename scalar_t>
class strumpack::HODLR::HODLRMatrix< scalar_t >

Hierarchically low-rank matrix representation.

This requires MPI support.

There are 3 different ways to create an HODLRMatrix

  • By specifying a matrix-(multiple)vector multiplication routine.
  • By specifying an element extraction routine.
  • By specifying a strumpack::kernel::Kernel matrix, defined by a collection of points and a kernel function.
Template Parameters
scalar_tCan be double, or std::complex<double>.
See also
HSS::HSSMatrix

Constructor & Destructor Documentation

◆ HODLRMatrix() [1/9]

template<typename scalar_t >
strumpack::HODLR::HODLRMatrix< scalar_t >::HODLRMatrix ( )
inline

Default constructor, makes an empty 0 x 0 matrix.

◆ HODLRMatrix() [2/9]

template<typename scalar_t >
strumpack::HODLR::HODLRMatrix< scalar_t >::HODLRMatrix ( const MPIComm c,
kernel::Kernel< real_t > &  K,
const opts_t opts 
)

Construct an HODLR approximation for the kernel matrix K.

Parameters
cMPI communicator, this communicator is copied internally.
KKernel matrix object. The data associated with this kernel will be permuted according to the clustering algorithm selected by the HODLROptions objects. The permutation will be stored in the kernel object.
optsobject containing a number of HODLR options

◆ HODLRMatrix() [3/9]

template<typename scalar_t >
strumpack::HODLR::HODLRMatrix< scalar_t >::HODLRMatrix ( const MPIComm c,
const structured::ClusterTree tree,
const std::function< scalar_t(int i, int j)> &  Aelem,
const opts_t opts 
)

Construct an HODLR approximation using a routine to evaluate individual matrix elements. This will construct an approximation of a permuted matrix, with the permutation available through this->perm() (and it's inverse this->iperm()). This permutation is applied symmetrically to rows and columns.

Parameters
cMPI communicator, this communicator is copied internally.
ttree specifying the HODLR matrix partitioning
AelemRoutine, std::function, which can also be a lambda function or a functor (class object implementing the member "scalar_t operator()(int i, int j)"), that evaluates/returns the matrix element A(i,j)
optsobject containing a number of HODLR options

◆ HODLRMatrix() [4/9]

template<typename scalar_t >
strumpack::HODLR::HODLRMatrix< scalar_t >::HODLRMatrix ( const MPIComm c,
const structured::ClusterTree tree,
const elem_blocks_t &  Aelem,
const opts_t opts 
)

Construct an HODLR approximation using a routine to evaluate multiple sub-blocks of the matrix at once. This will construct an approximation of a permuted matrix, with the permutation available through this->perm() (and it's inverse this->iperm()). This permutation is applied symmetrically to rows and columns.

Parameters
cMPI communicator, this communicator is copied internally.
ttree specifying the HODLR matrix partitioning
AelemRoutine, std::function, which can also be a lambda function or a functor. This should have the signature: void(std::vector<std::vector<std::size_t>>& I, std::vector<std::vector<std::size_t>>& J, std::vector<DenseMatrixWrapper<scalar_t>>& B, ExtractionMeta&); The ExtractionMeta object can be ignored.
optsobject containing a number of HODLR options

◆ HODLRMatrix() [5/9]

template<typename scalar_t >
strumpack::HODLR::HODLRMatrix< scalar_t >::HODLRMatrix ( const MPIComm c,
const structured::ClusterTree tree,
const std::function< void(Trans op, const DenseM_t &R, DenseM_t &S)> &  Amult,
const opts_t opts 
)

Construct an HODLR matrix using a specified HODLR tree and matrix-vector multiplication routine.

Parameters
cMPI communicator, this communicator is copied internally.
ttree specifying the HODLR matrix partitioning
AmultRoutine for the matrix-vector product. Trans op argument will be N, T or C for none, transpose or complex conjugate. The const DenseM_t& argument is the random matrix R, and the final DenseM_t& argument S is what the user routine should compute as A*R, A^t*R or A^c*R. S will already be allocated.
optsobject containing a number of options for HODLR compression
See also
compress, HODLROptions

◆ HODLRMatrix() [6/9]

template<typename scalar_t >
strumpack::HODLR::HODLRMatrix< scalar_t >::HODLRMatrix ( const MPIComm c,
const structured::ClusterTree tree,
const opts_t opts 
)

Construct an HODLR matrix using a specified HODLR tree. After construction, the HODLR matrix will be empty, and can be filled by calling one of the compress member routines.

Parameters
cMPI communicator, this communicator is copied internally.
ttree specifying the HODLR matrix partitioning
optsobject containing a number of options for HODLR compression
See also
compress, HODLROptions

◆ HODLRMatrix() [7/9]

template<typename scalar_t >
template<typename integer_t >
strumpack::HODLR::HODLRMatrix< scalar_t >::HODLRMatrix ( const MPIComm c,
const structured::ClusterTree tree,
const CSRGraph< integer_t > &  graph,
const opts_t opts 
)

Construct an HODLR matrix using a specified HODLR tree. After construction, the HODLR matrix will be empty, and can be filled by calling one of the compress member routines.

Parameters
cMPI communicator, this communicator is copied internally.
ttree specifying the HODLR matrix partitioning
graphconnectivity info for the dofs
optsobject containing a number of options for HODLR compression
See also
compress, HODLROptions

◆ HODLRMatrix() [8/9]

template<typename scalar_t >
strumpack::HODLR::HODLRMatrix< scalar_t >::HODLRMatrix ( const HODLRMatrix< scalar_t > &  h)
delete

Copy constructor is not supported.

◆ HODLRMatrix() [9/9]

template<typename scalar_t >
strumpack::HODLR::HODLRMatrix< scalar_t >::HODLRMatrix ( HODLRMatrix< scalar_t > &&  h)
inline

Move constructor.

Parameters
hHODLRMatrix to move from, will be emptied.

◆ ~HODLRMatrix()

template<typename scalar_t >
virtual strumpack::HODLR::HODLRMatrix< scalar_t >::~HODLRMatrix ( )
virtual

Virtual destructor.

Member Function Documentation

◆ begin_row()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::begin_row ( ) const
inlineoverridevirtual

Return the first row of the local rows owned by this process.

Returns
Return first local row

Reimplemented from strumpack::structured::StructuredMatrix< scalar_t >.

◆ cols()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::cols ( ) const
inlineoverridevirtual

Return the number of columns in the matrix.

Returns
Global number of columns in the matrix.

Implements strumpack::structured::StructuredMatrix< scalar_t >.

◆ Comm()

template<typename scalar_t >
const MPIComm& strumpack::HODLR::HODLRMatrix< scalar_t >::Comm ( ) const
inline

Return MPI communicator wrapper object.

◆ compress() [1/4]

template<typename scalar_t >
void strumpack::HODLR::HODLRMatrix< scalar_t >::compress ( const delem_blocks_t &  Aelem)

Construct the compressed HODLR representation of the matrix, using only a element evaluation (multiple blocks at once).

Parameters
Aelemelement extraction routine, extracting multiple blocks at once.

◆ compress() [2/4]

template<typename scalar_t >
void strumpack::HODLR::HODLRMatrix< scalar_t >::compress ( const elem_blocks_t &  Aelem)

Construct the compressed HODLR representation of the matrix, using only a element evaluation (multiple blocks at once).

Parameters
Aelemelement extraction routine, extracting multiple blocks at once.

◆ compress() [3/4]

template<typename scalar_t >
void strumpack::HODLR::HODLRMatrix< scalar_t >::compress ( const std::function< void(Trans op, const DenseM_t &R, DenseM_t &S)> &  Amult)

Construct the compressed HODLR representation of the matrix, using only a matrix-(multiple)vector multiplication routine.

Parameters
AmultRoutine for the matrix-vector product. Trans op argument will be N, T or C for none, transpose or complex conjugate. The const DenseM_t& argument is the the random matrix R, and the final DenseM_t& argument S is what the user routine should compute as A*R, A^t*R or A^c*R. S will already be allocated.

◆ compress() [4/4]

template<typename scalar_t >
void strumpack::HODLR::HODLRMatrix< scalar_t >::compress ( const std::function< void(Trans op, const DenseM_t &R, DenseM_t &S)> &  Amult,
int  rank_guess 
)

Construct the compressed HODLR representation of the matrix, using only a matrix-(multiple)vector multiplication routine.

Parameters
AmultRoutine for the matrix-vector product. Trans op argument will be N, T or C for none, transpose or complex conjugate. The const DenseM_t& argument is the the random matrix R, and the final DenseM_t& argument S is what the user routine should compute as A*R, A^t*R or A^c*R. S will already be allocated.
rank_guessInitial guess for the rank

◆ dense()

template<typename scalar_t >
DistM_t strumpack::HODLR::HODLRMatrix< scalar_t >::dense ( const BLACSGrid g) const

Create a dense matrix from this HODLR compressed matrix. This is mainly for debugging.

Parameters
gBLACSFGrid to be used for the output
Returns
dense matrix representation of *this, in 2D bloc cyclic format.

◆ dist()

template<typename scalar_t >
const std::vector<int>& strumpack::HODLR::HODLRMatrix< scalar_t >::dist ( ) const
inlineoverridevirtual

Return vector describing the 1d block row distribution. dist()[rank]==begin_row() and dist()[rank+1]==end_row()

Returns
1D block row distribution

Reimplemented from strumpack::structured::StructuredMatrix< scalar_t >.

◆ end_row()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::end_row ( ) const
inlineoverridevirtual

Return last row (+1) of the local rows (begin_rows()+lrows())

Returns
Final local row (+1).

Reimplemented from strumpack::structured::StructuredMatrix< scalar_t >.

◆ extract_elements()

template<typename scalar_t >
void strumpack::HODLR::HODLRMatrix< scalar_t >::extract_elements ( const Vec_t &  I,
const Vec_t &  J,
DenseM_t B 
)

Extract a submatrix defined by index sets I (rows) and J (columns), and put the result in matrix B (only on the root!).

Parameters
Iset of row indices, needs to be specified on all ranks!
Jset of column indices, needs to be specified on all ranks!
Boutput, the extracted elements. B should have the correct size B.rows() == I.size() and B.cols() == J.size()

◆ factor()

template<typename scalar_t >
void strumpack::HODLR::HODLRMatrix< scalar_t >::factor ( )
overridevirtual

Compute the factorization of this HODLR matrix. The matrix can still be used for multiplication.

See also
solve, inv_mult

Reimplemented from strumpack::structured::StructuredMatrix< scalar_t >.

◆ factor_memory()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::factor_memory ( ) const
inline

Return the additional memory for the factorization of this HODLR matrix, in bytes. This memory is only for the additional storage for the factorization, not the HODLR matrix itself.

◆ get_stat()

template<typename scalar_t >
double strumpack::HODLR::HODLRMatrix< scalar_t >::get_stat ( const std::string &  name) const

Get certain statistics about the HODLR matrix. See the HODLR code at https://github.com/liuyangzhuan/ButterflyPACK for more info.

Parameters
namefi : Rank_max, Mem_Factor, Mem_Fill, Flop_Fill, Flop_Factor, Flop_C_Mult

◆ inv_mult()

template<typename scalar_t >
long long int strumpack::HODLR::HODLRMatrix< scalar_t >::inv_mult ( Trans  op,
const DenseM_t B,
DenseM_t X 
) const

Solve a system of linear equations op(A)*X=B, with possibly multiple right-hand sides, where op can be none, transpose or complex conjugate.

Parameters
BRight hand side. This is the local part of the distributed matrix B. Should be B.rows() == this.lrows().
XResult, should be X.cols() == B.cols(), X.rows() == this.lrows(). X should be allocated.
Returns
number of flops
See also
factor, solve, lrows, begin_row, end_row

◆ iperm()

template<typename scalar_t >
const std::vector<int>& strumpack::HODLR::HODLRMatrix< scalar_t >::iperm ( ) const
inline

The inverse permutation for the matrix, which is applied to both rows and columns. This is 1-based!.

◆ local_rows()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::local_rows ( ) const
inlineoverridevirtual

Return the number of local rows, owned by this process.

Returns
Number of local rows.

Reimplemented from strumpack::structured::StructuredMatrix< scalar_t >.

◆ lrows()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::lrows ( ) const
inline

Return the number of local rows, owned by this process.

Returns
Number of local rows.

◆ max_rank()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::max_rank ( ) const
inline

Return the maximal rank encountered in this HODLR matrix, taking the maximum over all processes. This is collective on Comm().

◆ memory()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::memory ( ) const
inlineoverridevirtual

Return the memory for this HODLR matrix, on this rank, in bytes.

Implements strumpack::structured::StructuredMatrix< scalar_t >.

◆ mult() [1/2]

template<typename scalar_t >
void strumpack::HODLR::HODLRMatrix< scalar_t >::mult ( Trans  op,
const DenseM_t X,
DenseM_t Y 
) const
overridevirtual

Multiply this HODLR matrix with a dense matrix: Y = op(this)*X, where op can be none, transpose or complex conjugate. X and Y are the local parts of block-row distributed matrices. The number of rows in X and Y should correspond to the distribution of this HODLR matrix.

Parameters
opTranspose, conjugate, or none.
XRight-hand side matrix. This is the local part of the distributed matrix X. Should be X.rows() == this.lrows().
YResult, should be Y.cols() == X.cols(), Y.rows() == this.lrows()
See also
lrows, begin_row, end_row, mult

Reimplemented from strumpack::structured::StructuredMatrix< scalar_t >.

◆ mult() [2/2]

template<typename scalar_t >
void strumpack::HODLR::HODLRMatrix< scalar_t >::mult ( Trans  op,
const DistM_t X,
DistM_t Y 
) const
overridevirtual

Multiply this HODLR matrix with a dense matrix: Y = op(this)*X, where op can be none, transpose or complex conjugate. X and Y are in 2D block cyclic distribution.

Parameters
opTranspose, conjugate, or none.
XRight-hand side matrix. Should be X.rows() == this.rows().
YResult, should be Y.cols() == X.cols(), Y.rows() == this.rows()
See also
mult

Reimplemented from strumpack::structured::StructuredMatrix< scalar_t >.

◆ nonzeros()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::nonzeros ( ) const
inlineoverridevirtual

Return the total number of nonzeros stored by this matrix.

Returns
Nonzeros in the matrix representation.
See also
memory

Implements strumpack::structured::StructuredMatrix< scalar_t >.

◆ operator=() [1/2]

template<typename scalar_t >
HODLRMatrix<scalar_t>& strumpack::HODLR::HODLRMatrix< scalar_t >::operator= ( const HODLRMatrix< scalar_t > &  h)
delete

Copy assignment operator is not supported.

◆ operator=() [2/2]

template<typename scalar_t >
HODLRMatrix<scalar_t>& strumpack::HODLR::HODLRMatrix< scalar_t >::operator= ( HODLRMatrix< scalar_t > &&  h)

Move assignment operator.

Parameters
hHODLRMatrix to move from, will be emptied.

◆ perm()

template<typename scalar_t >
const std::vector<int>& strumpack::HODLR::HODLRMatrix< scalar_t >::perm ( ) const
inline

The permutation for the matrix, which is applied to both rows and columns. This is 1-based!.

◆ rank()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::rank ( ) const
inlineoverridevirtual

Return the maximal rank encountered in this HODLR matrix.

Implements strumpack::structured::StructuredMatrix< scalar_t >.

◆ rows()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::rows ( ) const
inlineoverridevirtual

Return the number of rows in the matrix.

Returns
Global number of rows in the matrix.

Implements strumpack::structured::StructuredMatrix< scalar_t >.

◆ solve() [1/4]

template<typename scalar_t >
long long int strumpack::HODLR::HODLRMatrix< scalar_t >::solve ( const DenseM_t B,
DenseM_t X 
) const

Solve a system of linear equations A*X=B, with possibly multiple right-hand sides. X and B are distributed using 1D block row distribution.

Parameters
BRight hand side. This is the local part of the distributed matrix B. Should be B.rows() == this.lrows().
XResult, should be X.cols() == B.cols(), X.rows() == this.lrows(). X should be allocated.
Returns
number of flops
See also
factor, lrows, begin_row, end_row, inv_mult

◆ solve() [2/4]

template<typename scalar_t >
long long int strumpack::HODLR::HODLRMatrix< scalar_t >::solve ( const DistM_t B,
DistM_t X 
) const

Solve a system of linear equations A*X=B, with possibly multiple right-hand sides. X and B are in 2D block cyclic distribution.

Parameters
BRight hand side. This is the local part of the distributed matrix B. Should be B.rows() == this.rows().
XResult, should be X.cols() == B.cols(), X.rows() == this.rows(). X should be allocated.
Returns
number of flops
See also
factor, inv_mult

◆ solve() [3/4]

template<typename scalar_t >
void strumpack::HODLR::HODLRMatrix< scalar_t >::solve ( DenseM_t B) const
overridevirtual

Solve a system of linear equations A*X=B, with possibly multiple right-hand sides. X and B are distributed using 1D block row distribution. The solution X will overwrite the right-hand side vector B.

Parameters
BRight hand side. This is the local part of the distributed matrix B. Should be B.rows() == this.lrows(). Wil lbe overwritten with the solution.
See also
factor, lrows, begin_row, end_row, inv_mult

Reimplemented from strumpack::structured::StructuredMatrix< scalar_t >.

◆ solve() [4/4]

template<typename scalar_t >
void strumpack::HODLR::HODLRMatrix< scalar_t >::solve ( DistM_t B) const
overridevirtual

Solve a system of linear equations A*X=B, with possibly multiple right-hand sides. X and B are in 2D block cyclic distribution. The solution X will overwrite the righ-hand side B.

Parameters
BRight hand side. This is the local part of the distributed matrix B. Should be B.rows() == this.rows(). Will be overwritten by the solution.
See also
factor, inv_mult, solve

Reimplemented from strumpack::structured::StructuredMatrix< scalar_t >.

◆ total_factor_memory()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::total_factor_memory ( ) const
inline

Return the total additional memory for the factorization of this HODLR matrix, in bytes, summed over all processes in Comm(). This call is collective on Comm(). This memory is only for the additional storage for the factorization, not the HODLR matrix itself.

◆ total_memory()

template<typename scalar_t >
std::size_t strumpack::HODLR::HODLRMatrix< scalar_t >::total_memory ( ) const
inline

Return the total memory for this HODLR matrix, summed over all ranks, in bytes.


The documentation for this class was generated from the following file: