MFEM  v4.6.0
Finite element discretization library
Public Member Functions | List of all members
mfem::HypreParMatrix Class Reference

Wrapper for hypre's ParCSR matrix class. More...

#include <hypre.hpp>

Inheritance diagram for mfem::HypreParMatrix:
[legend]
Collaboration diagram for mfem::HypreParMatrix:
[legend]

Public Member Functions

 HypreParMatrix ()
 An empty matrix to be used as a reference to an existing matrix. More...
 
void WrapHypreParCSRMatrix (hypre_ParCSRMatrix *a, bool owner=true)
 Converts hypre's format to HypreParMatrix. More...
 
 HypreParMatrix (hypre_ParCSRMatrix *a, bool owner=true)
 Converts hypre's format to HypreParMatrix. More...
 
 HypreParMatrix (MPI_Comm comm, HYPRE_BigInt glob_size, HYPRE_BigInt *row_starts, SparseMatrix *diag)
 Creates block-diagonal square parallel matrix. More...
 
 HypreParMatrix (MPI_Comm comm, HYPRE_BigInt global_num_rows, HYPRE_BigInt global_num_cols, HYPRE_BigInt *row_starts, HYPRE_BigInt *col_starts, SparseMatrix *diag)
 Creates block-diagonal rectangular parallel matrix. More...
 
 HypreParMatrix (MPI_Comm comm, HYPRE_BigInt global_num_rows, HYPRE_BigInt global_num_cols, HYPRE_BigInt *row_starts, HYPRE_BigInt *col_starts, SparseMatrix *diag, SparseMatrix *offd, HYPRE_BigInt *cmap, bool own_diag_offd=false)
 Creates general (rectangular) parallel matrix. More...
 
 HypreParMatrix (MPI_Comm comm, HYPRE_BigInt global_num_rows, HYPRE_BigInt global_num_cols, HYPRE_BigInt *row_starts, HYPRE_BigInt *col_starts, HYPRE_Int *diag_i, HYPRE_Int *diag_j, double *diag_data, HYPRE_Int *offd_i, HYPRE_Int *offd_j, double *offd_data, HYPRE_Int offd_num_cols, HYPRE_BigInt *offd_col_map, bool hypre_arrays=false)
 Creates general (rectangular) parallel matrix. More...
 
 HypreParMatrix (MPI_Comm comm, HYPRE_BigInt *row_starts, HYPRE_BigInt *col_starts, SparseMatrix *a)
 Creates a parallel matrix from SparseMatrix on processor 0. More...
 
 HypreParMatrix (MPI_Comm comm, HYPRE_BigInt global_num_rows, HYPRE_BigInt global_num_cols, HYPRE_BigInt *row_starts, HYPRE_BigInt *col_starts, Table *diag)
 Creates boolean block-diagonal rectangular parallel matrix. More...
 
 HypreParMatrix (MPI_Comm comm, int id, int np, HYPRE_BigInt *row, HYPRE_BigInt *col, HYPRE_Int *i_diag, HYPRE_Int *j_diag, HYPRE_Int *i_offd, HYPRE_Int *j_offd, HYPRE_BigInt *cmap, HYPRE_Int cmap_size)
 Creates boolean rectangular parallel matrix. More...
 
 HypreParMatrix (MPI_Comm comm, int nrows, HYPRE_BigInt glob_nrows, HYPRE_BigInt glob_ncols, int *I, HYPRE_BigInt *J, double *data, HYPRE_BigInt *rows, HYPRE_BigInt *cols)
 Creates a general parallel matrix from a local CSR matrix on each processor described by the I, J and data arrays. More...
 
 HypreParMatrix (const HypreParMatrix &P)
 Copy constructor for a ParCSR matrix which creates a deep copy of structure and data from P. More...
 
void MakeRef (const HypreParMatrix &master)
 Make this HypreParMatrix a reference to 'master'. More...
 
MPI_Comm GetComm () const
 MPI communicator. More...
 
 operator hypre_ParCSRMatrix * () const
 Typecasting to hypre's hypre_ParCSRMatrix*. More...
 
 operator HYPRE_ParCSRMatrix ()
 Typecasting to hypre's HYPRE_ParCSRMatrix, a.k.a. void *. More...
 
hypre_ParCSRMatrix * StealData ()
 Changes the ownership of the matrix. More...
 
void SetOwnerFlags (signed char diag, signed char offd, signed char colmap)
 Explicitly set the three ownership flags, see docs for diagOwner etc. More...
 
signed char OwnsDiag () const
 Get diag ownership flag. More...
 
signed char OwnsOffd () const
 Get offd ownership flag. More...
 
signed char OwnsColMap () const
 Get colmap ownership flag. More...
 
void CopyRowStarts ()
 
void CopyColStarts ()
 
HYPRE_BigInt NNZ () const
 Returns the global number of nonzeros. More...
 
HYPRE_BigIntRowPart ()
 Returns the row partitioning. More...
 
HYPRE_BigIntColPart ()
 Returns the column partitioning. More...
 
const HYPRE_BigIntRowPart () const
 Returns the row partitioning (const version) More...
 
const HYPRE_BigIntColPart () const
 Returns the column partitioning (const version) More...
 
HYPRE_BigInt M () const
 Returns the global number of rows. More...
 
HYPRE_BigInt N () const
 Returns the global number of columns. More...
 
void GetDiag (Vector &diag) const
 Get the local diagonal of the matrix. More...
 
void GetDiag (SparseMatrix &diag) const
 Get the local diagonal block. NOTE: 'diag' will not own any data. More...
 
void GetOffd (SparseMatrix &offd, HYPRE_BigInt *&cmap) const
 Get the local off-diagonal block. NOTE: 'offd' will not own any data. More...
 
void MergeDiagAndOffd (SparseMatrix &merged)
 Get a single SparseMatrix containing all rows from this processor, merged from the diagonal and off-diagonal blocks stored by the HypreParMatrix. More...
 
void AssembleDiagonal (Vector &diag) const override
 Return the diagonal of the matrix (Operator interface). More...
 
void GetBlocks (Array2D< HypreParMatrix *> &blocks, bool interleaved_rows=false, bool interleaved_cols=false) const
 
HypreParMatrixTranspose () const
 Returns the transpose of *this. More...
 
HypreParMatrixExtractSubmatrix (const Array< int > &indices, double threshold=0.0) const
 
int GetNumRows () const
 Returns the number of rows in the diagonal block of the ParCSRMatrix. More...
 
int GetNumCols () const
 Returns the number of columns in the diagonal block of the ParCSRMatrix. More...
 
HYPRE_BigInt GetGlobalNumRows () const
 Return the global number of rows. More...
 
HYPRE_BigInt GetGlobalNumCols () const
 Return the global number of columns. More...
 
HYPRE_BigIntGetRowStarts () const
 Return the parallel row partitioning array. More...
 
HYPRE_BigIntGetColStarts () const
 Return the parallel column partitioning array. More...
 
MemoryClass GetMemoryClass () const override
 Return the MemoryClass preferred by the Operator. More...
 
void EnsureMultTranspose () const
 Ensure the action of the transpose is performed fast. More...
 
void ResetTranspose () const
 Reset (destroy) the internal transpose matrix that is created by EnsureMultTranspose() and MultTranspose(). More...
 
HYPRE_Int Mult (HypreParVector &x, HypreParVector &y, double alpha=1.0, double beta=0.0) const
 Computes y = alpha * A * x + beta * y. More...
 
HYPRE_Int Mult (HYPRE_ParVector x, HYPRE_ParVector y, double alpha=1.0, double beta=0.0) const
 Computes y = alpha * A * x + beta * y. More...
 
HYPRE_Int MultTranspose (HypreParVector &x, HypreParVector &y, double alpha=1.0, double beta=0.0) const
 Computes y = alpha * A^t * x + beta * y. More...
 
void Mult (double a, const Vector &x, double b, Vector &y) const
 
void MultTranspose (double a, const Vector &x, double b, Vector &y) const
 Computes y = alpha * A^t * x + beta * y. More...
 
void Mult (const Vector &x, Vector &y) const override
 Operator application: y=A(x). More...
 
void MultTranspose (const Vector &x, Vector &y) const override
 Computes y = A^t * x. More...
 
void AddMult (const Vector &x, Vector &y, const double a=1.0) const override
 Operator application: y+=A(x) (default) or y+=a*A(x). More...
 
void AddMultTranspose (const Vector &x, Vector &y, const double a=1.0) const override
 Operator transpose application: y+=A^t(x) (default) or y+=a*A^t(x). More...
 
void AbsMult (double a, const Vector &x, double b, Vector &y) const
 Computes y = a * |A| * x + b * y, using entry-wise absolute values of the matrix A. More...
 
void AbsMultTranspose (double a, const Vector &x, double b, Vector &y) const
 Computes y = a * |At| * x + b * y, using entry-wise absolute values of the transpose of the matrix A. More...
 
void BooleanMult (int alpha, const int *x, int beta, int *y)
 The "Boolean" analog of y = alpha * A * x + beta * y, where elements in the sparsity pattern of the matrix are treated as "true". More...
 
void BooleanMultTranspose (int alpha, const int *x, int beta, int *y)
 The "Boolean" analog of y = alpha * A^T * x + beta * y, where elements in the sparsity pattern of the matrix are treated as "true". More...
 
HypreParMatrixoperator= (double value)
 Initialize all entries with value. More...
 
HypreParMatrixoperator+= (const HypreParMatrix &B)
 
HypreParMatrixAdd (const double beta, const HypreParMatrix &B)
 
HypreParMatrixLeftDiagMult (const SparseMatrix &D, HYPRE_BigInt *row_starts=NULL) const
 Multiply the HypreParMatrix on the left by a block-diagonal parallel matrix D and return the result as a new HypreParMatrix. More...
 
void ScaleRows (const Vector &s)
 Scale the local row i by s(i). More...
 
void InvScaleRows (const Vector &s)
 Scale the local row i by 1./s(i) More...
 
void operator*= (double s)
 Scale all entries by s: A_scaled = s*A. More...
 
void Threshold (double threshold=0.0)
 Remove values smaller in absolute value than some threshold. More...
 
void DropSmallEntries (double tol)
 Wrapper for hypre_ParCSRMatrixDropSmallEntries in different versions of hypre. Drop off-diagonal entries that are smaller than tol * l2 norm of its row. More...
 
void EliminateZeroRows ()
 If a row contains only zeros, set its diagonal to 1. More...
 
void EliminateRowsCols (const Array< int > &rows_cols, const HypreParVector &X, HypreParVector &B)
 
HypreParMatrixEliminateRowsCols (const Array< int > &rows_cols)
 
HypreParMatrixEliminateCols (const Array< int > &cols)
 
void EliminateRows (const Array< int > &rows)
 Eliminate rows from the diagonal and off-diagonal blocks of the matrix. More...
 
void EliminateBC (const HypreParMatrix &Ae, const Array< int > &ess_dof_list, const Vector &X, Vector &B) const
 Eliminate essential BC specified by ess_dof_list from the solution X to the r.h.s. B. More...
 
void EliminateBC (const Array< int > &ess_dofs, DiagonalPolicy diag_policy)
 Eliminate essential (Dirichlet) boundary conditions. More...
 
void HostRead () const
 Update the internal hypre_ParCSRMatrix object, A, to be on host. More...
 
void HostReadWrite ()
 Update the internal hypre_ParCSRMatrix object, A, to be on host. More...
 
void HostWrite ()
 Update the internal hypre_ParCSRMatrix object, A, to be on host. More...
 
void HypreRead () const
 Update the internal hypre_ParCSRMatrix object, A, to be in hypre memory space. More...
 
void HypreReadWrite ()
 Update the internal hypre_ParCSRMatrix object, A, to be in hypre memory space. More...
 
void HypreWrite ()
 Update the internal hypre_ParCSRMatrix object, A, to be in hypre memory space. More...
 
Memory< HYPRE_Int > & GetDiagMemoryI ()
 
Memory< HYPRE_Int > & GetDiagMemoryJ ()
 
Memory< double > & GetDiagMemoryData ()
 
const Memory< HYPRE_Int > & GetDiagMemoryI () const
 
const Memory< HYPRE_Int > & GetDiagMemoryJ () const
 
const Memory< double > & GetDiagMemoryData () const
 
void Print (const char *fname, HYPRE_Int offi=0, HYPRE_Int offj=0) const
 Prints the locally owned rows in parallel. More...
 
void Read (MPI_Comm comm, const char *fname)
 Reads the matrix from a file. More...
 
void Read_IJMatrix (MPI_Comm comm, const char *fname)
 Read a matrix saved as a HYPRE_IJMatrix. More...
 
void PrintCommPkg (std::ostream &out=mfem::out) const
 Print information about the hypre_ParCSRCommPkg of the HypreParMatrix. More...
 
void PrintHash (std::ostream &out) const
 Print sizes and hashes for all data arrays of the HypreParMatrix from the local MPI rank. More...
 
virtual ~HypreParMatrix ()
 Calls hypre's destroy function. More...
 
Type GetType () const
 
virtual void Mult (const Vector &x, Vector &y) const=0
 Operator application: y=A(x). More...
 
virtual void MultTranspose (const Vector &x, Vector &y) const
 Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an error. More...
 
- Public Member Functions inherited from mfem::Operator
void InitTVectors (const Operator *Po, const Operator *Ri, const Operator *Pi, Vector &x, Vector &b, Vector &X, Vector &B) const
 Initializes memory for true vectors of linear system. More...
 
 Operator (int s=0)
 Construct a square Operator with given size s (default 0). More...
 
 Operator (int h, int w)
 Construct an Operator with the given height (output size) and width (input size). More...
 
int Height () const
 Get the height (size of output) of the Operator. Synonym with NumRows(). More...
 
int NumRows () const
 Get the number of rows (size of output) of the Operator. Synonym with Height(). More...
 
int Width () const
 Get the width (size of input) of the Operator. Synonym with NumCols(). More...
 
int NumCols () const
 Get the number of columns (size of input) of the Operator. Synonym with Width(). More...
 
virtual void ArrayMult (const Array< const Vector *> &X, Array< Vector *> &Y) const
 Operator application on a matrix: Y=A(X). More...
 
virtual void ArrayMultTranspose (const Array< const Vector *> &X, Array< Vector *> &Y) const
 Action of the transpose operator on a matrix: Y=A^t(X). More...
 
virtual void ArrayAddMult (const Array< const Vector *> &X, Array< Vector *> &Y, const double a=1.0) const
 Operator application on a matrix: Y+=A(X) (default) or Y+=a*A(X). More...
 
virtual void ArrayAddMultTranspose (const Array< const Vector *> &X, Array< Vector *> &Y, const double a=1.0) const
 Operator transpose application on a matrix: Y+=A^t(X) (default) or Y+=a*A^t(X). More...
 
virtual OperatorGetGradient (const Vector &x) const
 Evaluate the gradient operator at the point x. The default behavior in class Operator is to generate an error. More...
 
virtual const OperatorGetProlongation () const
 Prolongation operator from linear algebra (linear system) vectors, to input vectors for the operator. NULL means identity. More...
 
virtual const OperatorGetRestriction () const
 Restriction operator from input vectors for the operator to linear algebra (linear system) vectors. NULL means identity. More...
 
virtual const OperatorGetOutputProlongation () const
 Prolongation operator from linear algebra (linear system) vectors, to output vectors for the operator. NULL means identity. More...
 
virtual const OperatorGetOutputRestrictionTranspose () const
 Transpose of GetOutputRestriction, directly available in this form to facilitate matrix-free RAP-type operators. More...
 
virtual const OperatorGetOutputRestriction () const
 Restriction operator from output vectors for the operator to linear algebra (linear system) vectors. NULL means identity. More...
 
void FormLinearSystem (const Array< int > &ess_tdof_list, Vector &x, Vector &b, Operator *&A, Vector &X, Vector &B, int copy_interior=0)
 Form a constrained linear system using a matrix-free approach. More...
 
void FormRectangularLinearSystem (const Array< int > &trial_tdof_list, const Array< int > &test_tdof_list, Vector &x, Vector &b, Operator *&A, Vector &X, Vector &B)
 Form a column-constrained linear system using a matrix-free approach. More...
 
virtual void RecoverFEMSolution (const Vector &X, const Vector &b, Vector &x)
 Reconstruct a solution vector x (e.g. a GridFunction) from the solution X of a constrained linear system obtained from Operator::FormLinearSystem() or Operator::FormRectangularLinearSystem(). More...
 
void FormSystemOperator (const Array< int > &ess_tdof_list, Operator *&A)
 Return in A a parallel (on truedofs) version of this square operator. More...
 
void FormRectangularSystemOperator (const Array< int > &trial_tdof_list, const Array< int > &test_tdof_list, Operator *&A)
 Return in A a parallel (on truedofs) version of this rectangular operator (including constraints). More...
 
void FormDiscreteOperator (Operator *&A)
 Return in A a parallel (on truedofs) version of this rectangular operator. More...
 
void PrintMatlab (std::ostream &out, int n, int m=0) const
 Prints operator with input size n and output size m in Matlab format. More...
 
virtual void PrintMatlab (std::ostream &out) const
 Prints operator in Matlab format. More...
 
virtual ~Operator ()
 Virtual destructor. More...
 
Type GetType () const
 Return the type ID of the Operator class. More...
 

Additional Inherited Members

- Public Types inherited from mfem::Operator
enum  DiagonalPolicy { DIAG_ZERO, DIAG_ONE, DIAG_KEEP }
 Defines operator diagonal policy upon elimination of rows and/or columns. More...
 
enum  Type {
  ANY_TYPE, MFEM_SPARSEMAT, Hypre_ParCSR, PETSC_MATAIJ,
  PETSC_MATIS, PETSC_MATSHELL, PETSC_MATNEST, PETSC_MATHYPRE,
  PETSC_MATGENERIC, Complex_Operator, MFEM_ComplexSparseMat, Complex_Hypre_ParCSR,
  Complex_DenseMat, MFEM_Block_Matrix, MFEM_Block_Operator
}
 Enumeration defining IDs for some classes derived from Operator. More...
 
- Protected Member Functions inherited from mfem::Operator
void FormConstrainedSystemOperator (const Array< int > &ess_tdof_list, ConstrainedOperator *&Aout)
 see FormSystemOperator() More...
 
void FormRectangularConstrainedSystemOperator (const Array< int > &trial_tdof_list, const Array< int > &test_tdof_list, RectangularConstrainedOperator *&Aout)
 see FormRectangularSystemOperator() More...
 
OperatorSetupRAP (const Operator *Pi, const Operator *Po)
 Returns RAP Operator of this, using input/output Prolongation matrices Pi corresponds to "P", Po corresponds to "Rt". More...
 
- Protected Attributes inherited from mfem::Operator
int height
 Dimension of the output / number of rows in the matrix. More...
 
int width
 Dimension of the input / number of columns in the matrix. More...
 

Detailed Description

Wrapper for hypre's ParCSR matrix class.

Definition at line 343 of file hypre.hpp.

Constructor & Destructor Documentation

◆ HypreParMatrix() [1/11]

mfem::HypreParMatrix::HypreParMatrix ( )

An empty matrix to be used as a reference to an existing matrix.

Definition at line 602 of file hypre.cpp.

◆ HypreParMatrix() [2/11]

mfem::HypreParMatrix::HypreParMatrix ( hypre_ParCSRMatrix *  a,
bool  owner = true 
)
inlineexplicit

Converts hypre's format to HypreParMatrix.

If owner is false, ownership of a is not transferred

Definition at line 426 of file hypre.hpp.

◆ HypreParMatrix() [3/11]

mfem::HypreParMatrix::HypreParMatrix ( MPI_Comm  comm,
HYPRE_BigInt  glob_size,
HYPRE_BigInt row_starts,
SparseMatrix diag 
)

Creates block-diagonal square parallel matrix.

Diagonal is given by diag which must be in CSR format (finalized). The new HypreParMatrix does not take ownership of any of the input arrays. See here for a description of the row partitioning array row_starts.

Warning
The ordering of the columns in each row in *diag may be changed by this constructor to ensure that the first entry in each row is the diagonal one. This is expected by most hypre functions.

Definition at line 770 of file hypre.cpp.

◆ HypreParMatrix() [4/11]

mfem::HypreParMatrix::HypreParMatrix ( MPI_Comm  comm,
HYPRE_BigInt  global_num_rows,
HYPRE_BigInt  global_num_cols,
HYPRE_BigInt row_starts,
HYPRE_BigInt col_starts,
SparseMatrix diag 
)

Creates block-diagonal rectangular parallel matrix.

Diagonal is given by diag which must be in CSR format (finalized). The new HypreParMatrix does not take ownership of any of the input arrays. See here for a description of the partitioning arrays row_starts and col_starts.

Definition at line 812 of file hypre.cpp.

◆ HypreParMatrix() [5/11]

mfem::HypreParMatrix::HypreParMatrix ( MPI_Comm  comm,
HYPRE_BigInt  global_num_rows,
HYPRE_BigInt  global_num_cols,
HYPRE_BigInt row_starts,
HYPRE_BigInt col_starts,
SparseMatrix diag,
SparseMatrix offd,
HYPRE_BigInt cmap,
bool  own_diag_offd = false 
)

Creates general (rectangular) parallel matrix.

The new HypreParMatrix does not take ownership of any of the input arrays, if own_diag_offd is false (default). If own_diag_offd is true, ownership of diag and offd is transferred to the HypreParMatrix.

See here for a description of the partitioning arrays row_starts and col_starts.

Definition at line 858 of file hypre.cpp.

◆ HypreParMatrix() [6/11]

mfem::HypreParMatrix::HypreParMatrix ( MPI_Comm  comm,
HYPRE_BigInt  global_num_rows,
HYPRE_BigInt  global_num_cols,
HYPRE_BigInt row_starts,
HYPRE_BigInt col_starts,
HYPRE_Int *  diag_i,
HYPRE_Int *  diag_j,
double *  diag_data,
HYPRE_Int *  offd_i,
HYPRE_Int *  offd_j,
double *  offd_data,
HYPRE_Int  offd_num_cols,
HYPRE_BigInt offd_col_map,
bool  hypre_arrays = false 
)

Creates general (rectangular) parallel matrix.

The new HypreParMatrix takes ownership of all input arrays, except col_starts and row_starts. See here for a description of the partitioning arrays row_starts and col_starts.

If hypre_arrays is false, all arrays (except row_starts and col_starts) are assumed to be allocated according to the MemoryType returned by Device::GetHostMemoryType(). If hypre_arrays is true, then the same arrays are assumed to be allocated by hypre as host arrays.

Definition at line 913 of file hypre.cpp.

◆ HypreParMatrix() [7/11]

mfem::HypreParMatrix::HypreParMatrix ( MPI_Comm  comm,
HYPRE_BigInt row_starts,
HYPRE_BigInt col_starts,
SparseMatrix a 
)

Creates a parallel matrix from SparseMatrix on processor 0.

See here for a description of the partitioning arrays row_starts and col_starts.

Definition at line 987 of file hypre.cpp.

◆ HypreParMatrix() [8/11]

mfem::HypreParMatrix::HypreParMatrix ( MPI_Comm  comm,
HYPRE_BigInt  global_num_rows,
HYPRE_BigInt  global_num_cols,
HYPRE_BigInt row_starts,
HYPRE_BigInt col_starts,
Table diag 
)

Creates boolean block-diagonal rectangular parallel matrix.

The new HypreParMatrix does not take ownership of any of the input arrays. See here for a description of the partitioning arrays row_starts and col_starts.

Definition at line 1031 of file hypre.cpp.

◆ HypreParMatrix() [9/11]

mfem::HypreParMatrix::HypreParMatrix ( MPI_Comm  comm,
int  id,
int  np,
HYPRE_BigInt row,
HYPRE_BigInt col,
HYPRE_Int *  i_diag,
HYPRE_Int *  j_diag,
HYPRE_Int *  i_offd,
HYPRE_Int *  j_offd,
HYPRE_BigInt cmap,
HYPRE_Int  cmap_size 
)

Creates boolean rectangular parallel matrix.

The new HypreParMatrix takes ownership of the arrays i_diag, j_diag, i_offd, j_offd, and cmap; does not take ownership of the arrays row and col. See here for a description of the partitioning arrays row and col.

Definition at line 1077 of file hypre.cpp.

◆ HypreParMatrix() [10/11]

mfem::HypreParMatrix::HypreParMatrix ( MPI_Comm  comm,
int  nrows,
HYPRE_BigInt  glob_nrows,
HYPRE_BigInt  glob_ncols,
int *  I,
HYPRE_BigInt J,
double *  data,
HYPRE_BigInt rows,
HYPRE_BigInt cols 
)

Creates a general parallel matrix from a local CSR matrix on each processor described by the I, J and data arrays.

The local matrix should be of size (local) nrows by (global) glob_ncols. The new parallel matrix contains copies of all input arrays (so they can be deleted). See here for a description of the partitioning arrays rows and cols.

Definition at line 1166 of file hypre.cpp.

◆ HypreParMatrix() [11/11]

mfem::HypreParMatrix::HypreParMatrix ( const HypreParMatrix P)

Copy constructor for a ParCSR matrix which creates a deep copy of structure and data from P.

Definition at line 1310 of file hypre.cpp.

◆ ~HypreParMatrix()

virtual mfem::HypreParMatrix::~HypreParMatrix ( )
inlinevirtual

Calls hypre's destroy function.

Definition at line 896 of file hypre.hpp.

Member Function Documentation

◆ AbsMult()

void mfem::HypreParMatrix::AbsMult ( double  a,
const Vector x,
double  b,
Vector y 
) const

Computes y = a * |A| * x + b * y, using entry-wise absolute values of the matrix A.

Definition at line 1881 of file hypre.cpp.

◆ AbsMultTranspose()

void mfem::HypreParMatrix::AbsMultTranspose ( double  a,
const Vector x,
double  b,
Vector y 
) const

Computes y = a * |At| * x + b * y, using entry-wise absolute values of the transpose of the matrix A.

Definition at line 1898 of file hypre.cpp.

◆ Add()

HypreParMatrix& mfem::HypreParMatrix::Add ( const double  beta,
const HypreParMatrix B 
)
inline

Perform the operation *this += beta*B, assuming that both matrices use the same row and column partitions and the same col_map_offd arrays, or B has an empty off-diagonal block. We also assume that the sparsity pattern of *this contains that of B. For a more general case consider the stand-alone function ParAdd described below.

Definition at line 763 of file hypre.hpp.

◆ AddMult()

void mfem::HypreParMatrix::AddMult ( const Vector x,
Vector y,
const double  a = 1.0 
) const
inlineoverridevirtual

Operator application: y+=A(x) (default) or y+=a*A(x).

Reimplemented from mfem::Operator.

Definition at line 704 of file hypre.hpp.

◆ AddMultTranspose()

void mfem::HypreParMatrix::AddMultTranspose ( const Vector x,
Vector y,
const double  a = 1.0 
) const
inlineoverridevirtual

Operator transpose application: y+=A^t(x) (default) or y+=a*A^t(x).

Reimplemented from mfem::Operator.

Definition at line 706 of file hypre.hpp.

◆ AssembleDiagonal()

void mfem::HypreParMatrix::AssembleDiagonal ( Vector diag) const
inlineoverridevirtual

Return the diagonal of the matrix (Operator interface).

Reimplemented from mfem::Operator.

Definition at line 602 of file hypre.hpp.

◆ BooleanMult()

void mfem::HypreParMatrix::BooleanMult ( int  alpha,
const int *  x,
int  beta,
int *  y 
)
inline

The "Boolean" analog of y = alpha * A * x + beta * y, where elements in the sparsity pattern of the matrix are treated as "true".

Definition at line 723 of file hypre.hpp.

◆ BooleanMultTranspose()

void mfem::HypreParMatrix::BooleanMultTranspose ( int  alpha,
const int *  x,
int  beta,
int *  y 
)
inline

The "Boolean" analog of y = alpha * A^T * x + beta * y, where elements in the sparsity pattern of the matrix are treated as "true".

Definition at line 733 of file hypre.hpp.

◆ ColPart() [1/2]

HYPRE_BigInt* mfem::HypreParMatrix::ColPart ( )
inline

Returns the column partitioning.

See here for a description of the partitioning array.

Definition at line 573 of file hypre.hpp.

◆ ColPart() [2/2]

const HYPRE_BigInt* mfem::HypreParMatrix::ColPart ( ) const
inline

Returns the column partitioning (const version)

See here for a description of the partitioning array.

Definition at line 581 of file hypre.hpp.

◆ CopyColStarts()

void mfem::HypreParMatrix::CopyColStarts ( )

If the HypreParMatrix does not own the col-starts array, make a copy of it that the HypreParMatrix will own. If the row-starts array is the same as the col-starts array, row-starts is also replaced.

Definition at line 1437 of file hypre.cpp.

◆ CopyRowStarts()

void mfem::HypreParMatrix::CopyRowStarts ( )

If the HypreParMatrix does not own the row-starts array, make a copy of it that the HypreParMatrix will own. If the col-starts array is the same as the row-starts array, col-starts is also replaced.

Definition at line 1397 of file hypre.cpp.

◆ DropSmallEntries()

void mfem::HypreParMatrix::DropSmallEntries ( double  tol)

Wrapper for hypre_ParCSRMatrixDropSmallEntries in different versions of hypre. Drop off-diagonal entries that are smaller than tol * l2 norm of its row.

For HYPRE versions < 2.14, this method just calls Threshold() with threshold = tol * max(l2 row norm).

Definition at line 2230 of file hypre.cpp.

◆ EliminateBC() [1/2]

void mfem::HypreParMatrix::EliminateBC ( const HypreParMatrix Ae,
const Array< int > &  ess_dof_list,
const Vector X,
Vector B 
) const

Eliminate essential BC specified by ess_dof_list from the solution X to the r.h.s. B.

This matrix is the matrix with eliminated BC, while Ae is such that (A+Ae) is the original (Neumann) matrix before elimination.

Definition at line 2328 of file hypre.cpp.

◆ EliminateBC() [2/2]

void mfem::HypreParMatrix::EliminateBC ( const Array< int > &  ess_dofs,
DiagonalPolicy  diag_policy 
)

Eliminate essential (Dirichlet) boundary conditions.

Parameters
[in]ess_dofsindices of the degrees of freedom belonging to the essential boundary conditions.
[in]diag_policypolicy for diagonal entries.

Definition at line 2386 of file hypre.cpp.

◆ EliminateCols()

HypreParMatrix * mfem::HypreParMatrix::EliminateCols ( const Array< int > &  cols)

Eliminate columns from the matrix and store the eliminated elements in a new matrix Ae (returned) so that the modified matrix and Ae sum to the original matrix.

Definition at line 2301 of file hypre.cpp.

◆ EliminateRows()

void mfem::HypreParMatrix::EliminateRows ( const Array< int > &  rows)

Eliminate rows from the diagonal and off-diagonal blocks of the matrix.

Definition at line 2315 of file hypre.cpp.

◆ EliminateRowsCols() [1/2]

void mfem::HypreParMatrix::EliminateRowsCols ( const Array< int > &  rows_cols,
const HypreParVector X,
HypreParVector B 
)

Eliminate rows and columns from the matrix, and rows from the vector B. Modify B with the BC values in X.

Definition at line 2276 of file hypre.cpp.

◆ EliminateRowsCols() [2/2]

HypreParMatrix * mfem::HypreParMatrix::EliminateRowsCols ( const Array< int > &  rows_cols)

Eliminate rows and columns from the matrix and store the eliminated elements in a new matrix Ae (returned), so that the modified matrix and Ae sum to the original matrix.

Definition at line 2287 of file hypre.cpp.

◆ EliminateZeroRows()

void mfem::HypreParMatrix::EliminateZeroRows ( )
inline

If a row contains only zeros, set its diagonal to 1.

Definition at line 799 of file hypre.hpp.

◆ EnsureMultTranspose()

void mfem::HypreParMatrix::EnsureMultTranspose ( ) const

Ensure the action of the transpose is performed fast.

When HYPRE is built for GPUs, this method will construct and store the transposes of the 'diag' and 'offd' CSR matrices. When HYPRE is not built for GPUs, this method is a no-op.

This method is automatically called by MultTranspose().

If the matrix is modified the old transpose blocks can be deleted by calling ResetTranspose().

Definition at line 1707 of file hypre.cpp.

◆ ExtractSubmatrix()

HypreParMatrix * mfem::HypreParMatrix::ExtractSubmatrix ( const Array< int > &  indices,
double  threshold = 0.0 
) const

Returns principle submatrix given by array of indices of connections with relative size > threshold in *this.

Definition at line 1630 of file hypre.cpp.

◆ GetBlocks()

void mfem::HypreParMatrix::GetBlocks ( Array2D< HypreParMatrix *> &  blocks,
bool  interleaved_rows = false,
bool  interleaved_cols = false 
) const

Split the matrix into M x N equally sized blocks of parallel matrices. The size of 'blocks' must already be set to M x N.

Definition at line 1587 of file hypre.cpp.

◆ GetColStarts()

HYPRE_BigInt* mfem::HypreParMatrix::GetColStarts ( ) const
inline

Return the parallel column partitioning array.

See here for a description of the partitioning array.

Definition at line 650 of file hypre.hpp.

◆ GetComm()

MPI_Comm mfem::HypreParMatrix::GetComm ( ) const
inline

MPI communicator.

Definition at line 534 of file hypre.hpp.

◆ GetDiag() [1/2]

void mfem::HypreParMatrix::GetDiag ( Vector diag) const

Get the local diagonal of the matrix.

Definition at line 1481 of file hypre.cpp.

◆ GetDiag() [2/2]

void mfem::HypreParMatrix::GetDiag ( SparseMatrix diag) const

Get the local diagonal block. NOTE: 'diag' will not own any data.

Definition at line 1552 of file hypre.cpp.

◆ GetDiagMemoryData() [1/2]

Memory<double>& mfem::HypreParMatrix::GetDiagMemoryData ( )
inline

Definition at line 872 of file hypre.hpp.

◆ GetDiagMemoryData() [2/2]

const Memory<double>& mfem::HypreParMatrix::GetDiagMemoryData ( ) const
inline

Definition at line 876 of file hypre.hpp.

◆ GetDiagMemoryI() [1/2]

Memory<HYPRE_Int>& mfem::HypreParMatrix::GetDiagMemoryI ( )
inline

Definition at line 870 of file hypre.hpp.

◆ GetDiagMemoryI() [2/2]

const Memory<HYPRE_Int>& mfem::HypreParMatrix::GetDiagMemoryI ( ) const
inline

Definition at line 874 of file hypre.hpp.

◆ GetDiagMemoryJ() [1/2]

Memory<HYPRE_Int>& mfem::HypreParMatrix::GetDiagMemoryJ ( )
inline

Definition at line 871 of file hypre.hpp.

◆ GetDiagMemoryJ() [2/2]

const Memory<HYPRE_Int>& mfem::HypreParMatrix::GetDiagMemoryJ ( ) const
inline

Definition at line 875 of file hypre.hpp.

◆ GetGlobalNumCols()

HYPRE_BigInt mfem::HypreParMatrix::GetGlobalNumCols ( ) const
inline

Return the global number of columns.

Definition at line 639 of file hypre.hpp.

◆ GetGlobalNumRows()

HYPRE_BigInt mfem::HypreParMatrix::GetGlobalNumRows ( ) const
inline

Return the global number of rows.

Definition at line 635 of file hypre.hpp.

◆ GetMemoryClass()

MemoryClass mfem::HypreParMatrix::GetMemoryClass ( ) const
inlineoverridevirtual

Return the MemoryClass preferred by the Operator.

This is the MemoryClass that will be used to access the input and output vectors in the Mult() and MultTranspose() methods.

For example, classes using the mfem::forall macro for implementation can return the value returned by Device::GetMemoryClass().

The default implementation of this method in class Operator returns MemoryClass::HOST.

Reimplemented from mfem::Operator.

Definition at line 652 of file hypre.hpp.

◆ GetNumCols()

int mfem::HypreParMatrix::GetNumCols ( ) const
inline

Returns the number of columns in the diagonal block of the ParCSRMatrix.

Definition at line 628 of file hypre.hpp.

◆ GetNumRows()

int mfem::HypreParMatrix::GetNumRows ( ) const
inline

Returns the number of rows in the diagonal block of the ParCSRMatrix.

Definition at line 621 of file hypre.hpp.

◆ GetOffd()

void mfem::HypreParMatrix::GetOffd ( SparseMatrix offd,
HYPRE_BigInt *&  cmap 
) const

Get the local off-diagonal block. NOTE: 'offd' will not own any data.

Definition at line 1557 of file hypre.cpp.

◆ GetRowStarts()

HYPRE_BigInt* mfem::HypreParMatrix::GetRowStarts ( ) const
inline

Return the parallel row partitioning array.

See here for a description of the partitioning array.

Definition at line 645 of file hypre.hpp.

◆ GetType()

Type mfem::HypreParMatrix::GetType ( ) const
inline

Definition at line 898 of file hypre.hpp.

◆ HostRead()

void mfem::HypreParMatrix::HostRead ( ) const
inline

Update the internal hypre_ParCSRMatrix object, A, to be on host.

After this call A's diagonal and off-diagonal should not be modified until after a suitable call to {Host,Hypre}{Write,ReadWrite}.

Definition at line 837 of file hypre.hpp.

◆ HostReadWrite()

void mfem::HypreParMatrix::HostReadWrite ( )
inline

Update the internal hypre_ParCSRMatrix object, A, to be on host.

After this call A's diagonal and off-diagonal can be modified on host and subsequent calls to Hypre{Read,Write,ReadWrite} will require a deep copy of the data if hypre is built with device support.

Definition at line 843 of file hypre.hpp.

◆ HostWrite()

void mfem::HypreParMatrix::HostWrite ( )
inline

Update the internal hypre_ParCSRMatrix object, A, to be on host.

Similar to HostReadWrite(), except that the data will never be copied from device to host to ensure host contains the correct current data.

Definition at line 848 of file hypre.hpp.

◆ HypreRead()

void mfem::HypreParMatrix::HypreRead ( ) const
inline

Update the internal hypre_ParCSRMatrix object, A, to be in hypre memory space.

After this call A's diagonal and off-diagonal should not be modified until after a suitable call to {Host,Hypre}{Write,ReadWrite}.

Definition at line 854 of file hypre.hpp.

◆ HypreReadWrite()

void mfem::HypreParMatrix::HypreReadWrite ( )
inline

Update the internal hypre_ParCSRMatrix object, A, to be in hypre memory space.

After this call A's diagonal and off-diagonal can be modified in hypre memory space and subsequent calls to Host{Read,Write,ReadWrite} will require a deep copy of the data if hypre is built with device support.

Definition at line 861 of file hypre.hpp.

◆ HypreWrite()

void mfem::HypreParMatrix::HypreWrite ( )
inline

Update the internal hypre_ParCSRMatrix object, A, to be in hypre memory space.

Similar to HostReadWrite(), except that the data will never be copied from host to hypre memory space to ensure the latter contains the correct current data.

Definition at line 868 of file hypre.hpp.

◆ InvScaleRows()

void mfem::HypreParMatrix::InvScaleRows ( const Vector s)

Scale the local row i by 1./s(i)

Definition at line 2056 of file hypre.cpp.

◆ LeftDiagMult()

HypreParMatrix * mfem::HypreParMatrix::LeftDiagMult ( const SparseMatrix D,
HYPRE_BigInt row_starts = NULL 
) const

Multiply the HypreParMatrix on the left by a block-diagonal parallel matrix D and return the result as a new HypreParMatrix.

If D has a different number of rows than A (this matrix), D's row starts array needs to be given (as returned by the methods GetDofOffsets/GetTrueDofOffsets of ParFiniteElementSpace). The new matrix D*A uses copies of the row- and column-starts arrays, so "this" matrix and row_starts can be deleted.

Note
This operation is local and does not require communication.

Definition at line 1915 of file hypre.cpp.

◆ M()

HYPRE_BigInt mfem::HypreParMatrix::M ( ) const
inline

Returns the global number of rows.

Definition at line 583 of file hypre.hpp.

◆ MakeRef()

void mfem::HypreParMatrix::MakeRef ( const HypreParMatrix master)

Make this HypreParMatrix a reference to 'master'.

Definition at line 1335 of file hypre.cpp.

◆ MergeDiagAndOffd()

void mfem::HypreParMatrix::MergeDiagAndOffd ( SparseMatrix merged)

Get a single SparseMatrix containing all rows from this processor, merged from the diagonal and off-diagonal blocks stored by the HypreParMatrix.

Note
The number of columns in the SparseMatrix will be the global number of columns in the parallel matrix, so using this method may result in an integer overflow in the column indices.

Definition at line 1563 of file hypre.cpp.

◆ Mult() [1/5]

HYPRE_Int mfem::HypreParMatrix::Mult ( HypreParVector x,
HypreParVector y,
double  alpha = 1.0,
double  beta = 0.0 
) const

Computes y = alpha * A * x + beta * y.

Definition at line 1736 of file hypre.cpp.

◆ Mult() [2/5]

HYPRE_Int mfem::HypreParMatrix::Mult ( HYPRE_ParVector  x,
HYPRE_ParVector  y,
double  alpha = 1.0,
double  beta = 0.0 
) const

Computes y = alpha * A * x + beta * y.

Definition at line 1865 of file hypre.cpp.

◆ Mult() [3/5]

void mfem::HypreParMatrix::Mult ( double  a,
const Vector x,
double  b,
Vector y 
) const

Definition at line 1744 of file hypre.cpp.

◆ Mult() [4/5]

void mfem::HypreParMatrix::Mult ( const Vector x,
Vector y 
) const
inlineoverridevirtual

Operator application: y=A(x).

Implements mfem::Operator.

Definition at line 694 of file hypre.hpp.

◆ Mult() [5/5]

virtual void mfem::Operator::Mult

Operator application: y=A(x).

◆ MultTranspose() [1/4]

HYPRE_Int mfem::HypreParMatrix::MultTranspose ( HypreParVector x,
HypreParVector y,
double  alpha = 1.0,
double  beta = 0.0 
) const

Computes y = alpha * A^t * x + beta * y.

If the matrix is modified, call ResetTranspose() and optionally EnsureMultTranspose() to make sure this method uses the correct updated transpose.

Definition at line 1872 of file hypre.cpp.

◆ MultTranspose() [2/4]

void mfem::HypreParMatrix::MultTranspose ( double  a,
const Vector x,
double  b,
Vector y 
) const

Computes y = alpha * A^t * x + beta * y.

If the matrix is modified, call ResetTranspose() and optionally EnsureMultTranspose() to make sure this method uses the correct updated transpose.

Definition at line 1801 of file hypre.cpp.

◆ MultTranspose() [3/4]

void mfem::HypreParMatrix::MultTranspose ( const Vector x,
Vector y 
) const
inlineoverridevirtual

Computes y = A^t * x.

If the matrix is modified, call ResetTranspose() and optionally EnsureMultTranspose() to make sure this method uses the correct updated transpose.

Reimplemented from mfem::Operator.

Definition at line 701 of file hypre.hpp.

◆ MultTranspose() [4/4]

virtual void mfem::Operator::MultTranspose
inline

Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an error.

Definition at line 93 of file operator.hpp.

◆ N()

HYPRE_BigInt mfem::HypreParMatrix::N ( ) const
inline

Returns the global number of columns.

Definition at line 585 of file hypre.hpp.

◆ NNZ()

HYPRE_BigInt mfem::HypreParMatrix::NNZ ( ) const
inline

Returns the global number of nonzeros.

Definition at line 565 of file hypre.hpp.

◆ operator HYPRE_ParCSRMatrix()

mfem::HypreParMatrix::operator HYPRE_ParCSRMatrix ( )
inline

Typecasting to hypre's HYPRE_ParCSRMatrix, a.k.a. void *.

Definition at line 540 of file hypre.hpp.

◆ operator hypre_ParCSRMatrix *()

mfem::HypreParMatrix::operator hypre_ParCSRMatrix * ( ) const
inline

Typecasting to hypre's hypre_ParCSRMatrix*.

Definition at line 537 of file hypre.hpp.

◆ operator*=()

void mfem::HypreParMatrix::operator*= ( double  s)

Scale all entries by s: A_scaled = s*A.

Definition at line 2102 of file hypre.cpp.

◆ operator+=()

HypreParMatrix& mfem::HypreParMatrix::operator+= ( const HypreParMatrix B)
inline

Perform the operation *this += B, assuming that both matrices use the same row and column partitions and the same col_map_offd arrays, or B has an empty off-diagonal block. We also assume that the sparsity pattern of *this contains that of B.

Definition at line 756 of file hypre.hpp.

◆ operator=()

HypreParMatrix& mfem::HypreParMatrix::operator= ( double  value)
inline

Initialize all entries with value.

Definition at line 742 of file hypre.hpp.

◆ OwnsColMap()

signed char mfem::HypreParMatrix::OwnsColMap ( ) const
inline

Get colmap ownership flag.

Definition at line 553 of file hypre.hpp.

◆ OwnsDiag()

signed char mfem::HypreParMatrix::OwnsDiag ( ) const
inline

Get diag ownership flag.

Definition at line 549 of file hypre.hpp.

◆ OwnsOffd()

signed char mfem::HypreParMatrix::OwnsOffd ( ) const
inline

Get offd ownership flag.

Definition at line 551 of file hypre.hpp.

◆ Print()

void mfem::HypreParMatrix::Print ( const char *  fname,
HYPRE_Int  offi = 0,
HYPRE_Int  offj = 0 
) const

Prints the locally owned rows in parallel.

Definition at line 2530 of file hypre.cpp.

◆ PrintCommPkg()

void mfem::HypreParMatrix::PrintCommPkg ( std::ostream &  out = mfem::out) const

Print information about the hypre_ParCSRCommPkg of the HypreParMatrix.

Definition at line 2574 of file hypre.cpp.

◆ PrintHash()

void mfem::HypreParMatrix::PrintHash ( std::ostream &  out) const

Print sizes and hashes for all data arrays of the HypreParMatrix from the local MPI rank.

This is a compact text representation of the local data of the HypreParMatrix that can be used to compare matrices from different runs without the need to save the whole matrix.

Definition at line 2611 of file hypre.cpp.

◆ Read()

void mfem::HypreParMatrix::Read ( MPI_Comm  comm,
const char *  fname 
)

Reads the matrix from a file.

Definition at line 2538 of file hypre.cpp.

◆ Read_IJMatrix()

void mfem::HypreParMatrix::Read_IJMatrix ( MPI_Comm  comm,
const char *  fname 
)

Read a matrix saved as a HYPRE_IJMatrix.

Definition at line 2553 of file hypre.cpp.

◆ ResetTranspose()

void mfem::HypreParMatrix::ResetTranspose ( ) const

Reset (destroy) the internal transpose matrix that is created by EnsureMultTranspose() and MultTranspose().

If the matrix is modified, this method should be called to delete the out-of-date transpose that is stored internally.

Definition at line 1717 of file hypre.cpp.

◆ RowPart() [1/2]

HYPRE_BigInt* mfem::HypreParMatrix::RowPart ( )
inline

Returns the row partitioning.

See here for a description of the partitioning array.

Definition at line 569 of file hypre.hpp.

◆ RowPart() [2/2]

const HYPRE_BigInt* mfem::HypreParMatrix::RowPart ( ) const
inline

Returns the row partitioning (const version)

See here for a description of the partitioning array.

Definition at line 577 of file hypre.hpp.

◆ ScaleRows()

void mfem::HypreParMatrix::ScaleRows ( const Vector s)

Scale the local row i by s(i).

Definition at line 2017 of file hypre.cpp.

◆ SetOwnerFlags()

void mfem::HypreParMatrix::SetOwnerFlags ( signed char  diag,
signed char  offd,
signed char  colmap 
)

Explicitly set the three ownership flags, see docs for diagOwner etc.

Definition at line 1372 of file hypre.cpp.

◆ StealData()

hypre_ParCSRMatrix * mfem::HypreParMatrix::StealData ( )

Changes the ownership of the matrix.

Definition at line 1353 of file hypre.cpp.

◆ Threshold()

void mfem::HypreParMatrix::Threshold ( double  threshold = 0.0)

Remove values smaller in absolute value than some threshold.

Definition at line 2145 of file hypre.cpp.

◆ Transpose()

HypreParMatrix * mfem::HypreParMatrix::Transpose ( ) const

Returns the transpose of *this.

Definition at line 1611 of file hypre.cpp.

◆ WrapHypreParCSRMatrix()

void mfem::HypreParMatrix::WrapHypreParCSRMatrix ( hypre_ParCSRMatrix *  a,
bool  owner = true 
)

Converts hypre's format to HypreParMatrix.

If owner is false, ownership of a is not transferred

Definition at line 608 of file hypre.cpp.


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