MFEM  v4.3.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
Public Member Functions | Protected Member Functions | Protected Attributes | Friends | List of all members
mfem::PetscParMatrix Class Reference

Wrapper for PETSc's matrix class. More...

#include <petsc.hpp>

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

Public Member Functions

 PetscParMatrix ()
 Create an empty matrix to be used as a reference to an existing matrix. More...
 
 PetscParMatrix (petsc::Mat a, bool ref=false)
 Creates PetscParMatrix out of PETSc's Mat. More...
 
 PetscParMatrix (const PetscParMatrix *pa, Operator::Type tid)
 Convert a PetscParMatrix pa with a new PETSc format tid. Note that if pa is already a PetscParMatrix of the same type as tid, the resulting PetscParMatrix will share the same Mat object. More...
 
 PetscParMatrix (const PetscParMatrix &A, const Array< PetscInt > &rows, const Array< PetscInt > &cols)
 Creates a PetscParMatrix extracting the submatrix of A with rows row indices and cols column indices. More...
 
 PetscParMatrix (const HypreParMatrix *ha, Operator::Type tid=Operator::PETSC_MATAIJ)
 Convert a HypreParMatrix ha to a PetscParMatrix in the given PETSc format tid. More...
 
 PetscParMatrix (const SparseMatrix *sa, Operator::Type tid=Operator::PETSC_MATAIJ)
 Convert a SparseMatrix ha to a PetscParMatrix in the given PETSc format tid. More...
 
 PetscParMatrix (MPI_Comm comm, const Operator *op, Operator::Type tid=Operator::PETSC_MATSHELL)
 Convert an mfem::Operator into a PetscParMatrix in the given PETSc format tid. More...
 
 PetscParMatrix (MPI_Comm comm, PetscInt glob_size, PetscInt *row_starts, SparseMatrix *diag, Operator::Type tid)
 Creates block-diagonal square parallel matrix. More...
 
 PetscParMatrix (MPI_Comm comm, PetscInt global_num_rows, PetscInt global_num_cols, PetscInt *row_starts, PetscInt *col_starts, SparseMatrix *diag, Operator::Type tid)
 Creates block-diagonal rectangular parallel matrix. More...
 
virtual ~PetscParMatrix ()
 Calls PETSc's destroy function. More...
 
void SetMat (petsc::Mat newA)
 Replace the inner Mat Object. The reference count of newA is increased. More...
 
void Mult (double a, const Vector &x, double b, Vector &y) const
 Matvec: y = a A x + b y. More...
 
void MultTranspose (double a, const Vector &x, double b, Vector &y) const
 Matvec transpose: y = a A^T x + b y. More...
 
virtual void Mult (const Vector &x, Vector &y) const
 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...
 
MPI_Comm GetComm () const
 Get the associated MPI communicator. More...
 
 operator petsc::Mat () const
 Typecasting to PETSc's Mat type. More...
 
 operator PetscObject () const
 Typecasting to PETSc object. More...
 
PetscInt GetRowStart () const
 Returns the global index of the first local row. More...
 
PetscInt GetColStart () const
 Returns the global index of the first local column. More...
 
PetscInt GetNumRows () const
 Returns the local number of rows. More...
 
PetscInt GetNumCols () const
 Returns the local number of columns. More...
 
PetscInt M () const
 Returns the global number of rows. More...
 
PetscInt N () const
 Returns the global number of columns. More...
 
PetscInt GetGlobalNumRows () const
 Returns the global number of rows. More...
 
PetscInt GetGlobalNumCols () const
 Returns the global number of columns. More...
 
PetscInt NNZ () const
 Returns the number of nonzeros. More...
 
PetscParVectorGetX () const
 Returns the inner vector in the domain of A (it creates it if needed) More...
 
PetscParVectorGetY () const
 Returns the inner vector in the range of A (it creates it if needed) More...
 
PetscParMatrixTranspose (bool action=false)
 Returns the transpose of the PetscParMatrix. More...
 
void Print (const char *fname=NULL, bool binary=false) const
 Prints the matrix (to stdout if fname is NULL) More...
 
void operator*= (double s)
 Scale all entries by s: A_scaled = s*A. More...
 
void EliminateRowsCols (const Array< int > &rows_cols, const PetscParVector &X, PetscParVector &B, double diag=1.)
 Eliminate rows and columns from the matrix, and rows from the vector B. Modify B with the BC values in X. Put diag on the diagonal corresponding to eliminated entries. More...
 
void EliminateRowsCols (const Array< int > &rows_cols, const HypreParVector &X, HypreParVector &B, double diag=1.)
 
PetscParMatrixEliminateRowsCols (const Array< int > &rows_cols)
 Eliminate rows and columns from the matrix and store the eliminated elements in a new matrix Ae (returned). More...
 
void ScaleRows (const Vector &s)
 Scale the local row i by s(i). More...
 
void ScaleCols (const Vector &s)
 Scale the local col i by s(i). More...
 
void Shift (double s)
 Shift diagonal by a constant. More...
 
void Shift (const Vector &s)
 Shift diagonal by a vector. More...
 
void EliminateRows (const Array< int > &rows)
 Eliminate only the rows from the matrix. More...
 
void MakeRef (const PetscParMatrix &master)
 Makes this object a reference to another PetscParMatrix. More...
 
petsc::Mat ReleaseMat (bool dereference)
 Release the PETSc Mat object. If dereference is true, decrement the refcount of the Mat object. More...
 
Type GetType () const
 
Assignment operators
PetscParMatrixoperator= (const PetscParMatrix &B)
 
PetscParMatrixoperator= (const HypreParMatrix &B)
 
PetscParMatrixoperator+= (const PetscParMatrix &B)
 
PetscParMatrixoperator-= (const PetscParMatrix &B)
 
- 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 MemoryClass GetMemoryClass () const
 Return the MemoryClass preferred by the Operator. 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 void AssembleDiagonal (Vector &diag) const
 Computes the diagonal entries into diag. Typically, this operation only makes sense for linear Operators. In some cases, only an approximation of the diagonal is computed. 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=0, int m=0) const
 Prints operator with input size n and output size m in Matlab format. More...
 
virtual ~Operator ()
 Virtual destructor. More...
 
Type GetType () const
 Return the type ID of the Operator class. More...
 

Protected Member Functions

void Init ()
 Initialize with defaults. Does not initialize inherited members. More...
 
void Destroy ()
 Delete all owned data. Does not perform re-initialization with defaults. More...
 
void MakeWrapper (MPI_Comm comm, const Operator *op, petsc::Mat *B)
 Creates a wrapper around a mfem::Operator op using PETSc's MATSHELL object and returns the Mat in B. More...
 
void ConvertOperator (MPI_Comm comm, const Operator &op, petsc::Mat *B, Operator::Type tid)
 
- 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

petsc::Mat A
 The actual PETSc object. More...
 
PetscParVectorX
 Auxiliary vectors for typecasting. More...
 
PetscParVectorY
 
- 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...
 

Friends

class PetscLinearSolver
 
class PetscPreconditioner
 

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
}
 Enumeration defining IDs for some classes derived from Operator. More...
 

Detailed Description

Wrapper for PETSc's matrix class.

Definition at line 307 of file petsc.hpp.

Constructor & Destructor Documentation

mfem::PetscParMatrix::PetscParMatrix ( )

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

Definition at line 956 of file petsc.cpp.

mfem::PetscParMatrix::PetscParMatrix ( petsc::Mat  a,
bool  ref = false 
)

Creates PetscParMatrix out of PETSc's Mat.

Parameters
[in]aThe PETSc Mat object.
[in]refIf true, we increase the reference count of a.

Definition at line 1703 of file petsc.cpp.

mfem::PetscParMatrix::PetscParMatrix ( const PetscParMatrix pa,
Operator::Type  tid 
)
explicit

Convert a PetscParMatrix pa with a new PETSc format tid. Note that if pa is already a PetscParMatrix of the same type as tid, the resulting PetscParMatrix will share the same Mat object.

Definition at line 981 of file petsc.cpp.

mfem::PetscParMatrix::PetscParMatrix ( const PetscParMatrix A,
const Array< PetscInt > &  rows,
const Array< PetscInt > &  cols 
)

Creates a PetscParMatrix extracting the submatrix of A with rows row indices and cols column indices.

Definition at line 961 of file petsc.cpp.

mfem::PetscParMatrix::PetscParMatrix ( const HypreParMatrix ha,
Operator::Type  tid = Operator::PETSC_MATAIJ 
)
explicit

Convert a HypreParMatrix ha to a PetscParMatrix in the given PETSc format tid.

The supported type ids are: Operator::PETSC_MATAIJ, Operator::PETSC_MATIS, Operator::PETSC_MATSHELL and Operator::PETSC_MATHYPRE ha can be destroyed unless tid == PETSC_MATSHELL or tid == PETSC_MATHYPRE

Definition at line 989 of file petsc.cpp.

mfem::PetscParMatrix::PetscParMatrix ( const SparseMatrix sa,
Operator::Type  tid = Operator::PETSC_MATAIJ 
)
explicit

Convert a SparseMatrix ha to a PetscParMatrix in the given PETSc format tid.

Definition at line 997 of file petsc.cpp.

mfem::PetscParMatrix::PetscParMatrix ( MPI_Comm  comm,
const Operator op,
Operator::Type  tid = Operator::PETSC_MATSHELL 
)

Convert an mfem::Operator into a PetscParMatrix in the given PETSc format tid.

If tid is Operator::PETSC_MATSHELL and op is not a PetscParMatrix, it converts any mfem::Operator op implementing Operator::Mult() and Operator::MultTranspose() into a PetscParMatrix. The Operator op should not be deleted while the constructed PetscParMatrix is used.

Otherwise, it tries to convert the operator in PETSc's classes. op cannot be destroyed if tid == PETSC_MATHYPRE.

In particular, if op is a BlockOperator, then a MATNEST Mat object is created using tid as the type for the blocks. Note that if op is already a PetscParMatrix of the same type as tid, the resulting PetscParMatrix will share the same Mat object

Definition at line 1005 of file petsc.cpp.

mfem::PetscParMatrix::PetscParMatrix ( MPI_Comm  comm,
PetscInt  glob_size,
PetscInt row_starts,
SparseMatrix diag,
Operator::Type  tid 
)

Creates block-diagonal square parallel matrix.

The block-diagonal is given by diag which must be in CSR format (finalized). The new PetscParMatrix does not take ownership of any of the input arrays. The type id tid can be either PETSC_MATAIJ (parallel distributed CSR) or PETSC_MATIS.

Definition at line 1014 of file petsc.cpp.

mfem::PetscParMatrix::PetscParMatrix ( MPI_Comm  comm,
PetscInt  global_num_rows,
PetscInt  global_num_cols,
PetscInt row_starts,
PetscInt col_starts,
SparseMatrix diag,
Operator::Type  tid 
)

Creates block-diagonal rectangular parallel matrix.

The block-diagonal is given by diag which must be in CSR format (finalized). The new PetscParMatrix does not take ownership of any of the input arrays. The type id tid can be either PETSC_MATAIJ (parallel distributed CSR) or PETSC_MATIS.

Definition at line 1027 of file petsc.cpp.

virtual mfem::PetscParMatrix::~PetscParMatrix ( )
inlinevirtual

Calls PETSc's destroy function.

Definition at line 417 of file petsc.hpp.

Member Function Documentation

void mfem::PetscParMatrix::ConvertOperator ( MPI_Comm  comm,
const Operator op,
petsc::Mat B,
Operator::Type  tid 
)
protected

Convert an mfem::Operator into a Mat B; op can be destroyed unless tid == PETSC_MATSHELL or tid == PETSC_MATHYPRE if op is a BlockOperator, the operator type is relevant to the individual blocks

Definition at line 1294 of file petsc.cpp.

void mfem::PetscParMatrix::Destroy ( )
protected

Delete all owned data. Does not perform re-initialization with defaults.

Definition at line 1690 of file petsc.cpp.

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

Eliminate only the rows from the matrix.

Definition at line 2200 of file petsc.cpp.

void mfem::PetscParMatrix::EliminateRowsCols ( const Array< int > &  rows_cols,
const PetscParVector X,
PetscParVector B,
double  diag = 1. 
)

Eliminate rows and columns from the matrix, and rows from the vector B. Modify B with the BC values in X. Put diag on the diagonal corresponding to eliminated entries.

Definition at line 2171 of file petsc.cpp.

void mfem::PetscParMatrix::EliminateRowsCols ( const Array< int > &  rows_cols,
const HypreParVector X,
HypreParVector B,
double  diag = 1. 
)

Definition at line 2163 of file petsc.cpp.

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

Eliminate rows and columns from the matrix and store the eliminated elements in a new matrix Ae (returned).

The sum of the modified matrix and the returned matrix, Ae, is equal to the original matrix.

Definition at line 2152 of file petsc.cpp.

PetscInt mfem::PetscParMatrix::GetColStart ( ) const

Returns the global index of the first local column.

Definition at line 907 of file petsc.cpp.

MPI_Comm mfem::PetscParMatrix::GetComm ( ) const

Get the associated MPI communicator.

Definition at line 1247 of file petsc.cpp.

PetscInt mfem::PetscParMatrix::GetGlobalNumCols ( ) const
inline

Returns the global number of columns.

Definition at line 473 of file petsc.hpp.

PetscInt mfem::PetscParMatrix::GetGlobalNumRows ( ) const
inline

Returns the global number of rows.

Definition at line 470 of file petsc.hpp.

PetscInt mfem::PetscParMatrix::GetNumCols ( ) const

Returns the local number of columns.

Definition at line 921 of file petsc.cpp.

PetscInt mfem::PetscParMatrix::GetNumRows ( ) const

Returns the local number of rows.

Definition at line 914 of file petsc.cpp.

PetscInt mfem::PetscParMatrix::GetRowStart ( ) const

Returns the global index of the first local row.

Definition at line 900 of file petsc.cpp.

Operator::Type mfem::PetscParMatrix::GetType ( ) const

Definition at line 2227 of file petsc.cpp.

PetscParVector * mfem::PetscParMatrix::GetX ( ) const

Returns the inner vector in the domain of A (it creates it if needed)

Definition at line 1871 of file petsc.cpp.

PetscParVector * mfem::PetscParMatrix::GetY ( ) const

Returns the inner vector in the range of A (it creates it if needed)

Definition at line 1881 of file petsc.cpp.

void mfem::PetscParMatrix::Init ( )
protected

Initialize with defaults. Does not initialize inherited members.

Definition at line 949 of file petsc.cpp.

PetscInt mfem::PetscParMatrix::M ( ) const

Returns the global number of rows.

Definition at line 928 of file petsc.cpp.

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

Makes this object a reference to another PetscParMatrix.

Definition at line 1861 of file petsc.cpp.

void mfem::PetscParMatrix::MakeWrapper ( MPI_Comm  comm,
const Operator op,
petsc::Mat B 
)
protected

Creates a wrapper around a mfem::Operator op using PETSc's MATSHELL object and returns the Mat in B.

This does not take any reference to op, that should not be destroyed until B is needed.

Definition at line 1260 of file petsc.cpp.

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

Matvec: y = a A x + b y.

Definition at line 1910 of file petsc.cpp.

virtual void mfem::PetscParMatrix::Mult ( const Vector x,
Vector y 
) const
inlinevirtual

Operator application: y=A(x).

Implements mfem::Operator.

Definition at line 436 of file petsc.hpp.

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

Matvec transpose: y = a A^T x + b y.

Definition at line 1927 of file petsc.cpp.

virtual void mfem::PetscParMatrix::MultTranspose ( const Vector x,
Vector y 
) const
inlinevirtual

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

Reimplemented from mfem::Operator.

Definition at line 439 of file petsc.hpp.

PetscInt mfem::PetscParMatrix::N ( ) const

Returns the global number of columns.

Definition at line 935 of file petsc.cpp.

PetscInt mfem::PetscParMatrix::NNZ ( ) const

Returns the number of nonzeros.

Differently from HYPRE, this call is collective on the communicator, as this number is not stored inside PETSc, but needs to be computed.

Definition at line 942 of file petsc.cpp.

mfem::PetscParMatrix::operator petsc::Mat ( ) const
inline

Typecasting to PETSc's Mat type.

Definition at line 446 of file petsc.hpp.

mfem::PetscParMatrix::operator PetscObject ( ) const
inline

Typecasting to PETSc object.

Definition at line 449 of file petsc.hpp.

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

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

Definition at line 1905 of file petsc.cpp.

PetscParMatrix & mfem::PetscParMatrix::operator+= ( const PetscParMatrix B)

Definition at line 1078 of file petsc.cpp.

PetscParMatrix & mfem::PetscParMatrix::operator-= ( const PetscParMatrix B)

Definition at line 1093 of file petsc.cpp.

PetscParMatrix & mfem::PetscParMatrix::operator= ( const PetscParMatrix B)

Definition at line 1062 of file petsc.cpp.

PetscParMatrix & mfem::PetscParMatrix::operator= ( const HypreParMatrix B)

Definition at line 1041 of file petsc.cpp.

void mfem::PetscParMatrix::Print ( const char *  fname = NULL,
bool  binary = false 
) const

Prints the matrix (to stdout if fname is NULL)

Definition at line 1945 of file petsc.cpp.

Mat mfem::PetscParMatrix::ReleaseMat ( bool  dereference)

Release the PETSc Mat object. If dereference is true, decrement the refcount of the Mat object.

Definition at line 2214 of file petsc.cpp.

void mfem::PetscParMatrix::ScaleCols ( const Vector s)

Scale the local col i by s(i).

Definition at line 1981 of file petsc.cpp.

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

Scale the local row i by s(i).

Definition at line 1970 of file petsc.cpp.

void mfem::PetscParMatrix::SetMat ( petsc::Mat  newA)

Replace the inner Mat Object. The reference count of newA is increased.

Definition at line 1715 of file petsc.cpp.

void mfem::PetscParMatrix::Shift ( double  s)

Shift diagonal by a constant.

Definition at line 1992 of file petsc.cpp.

void mfem::PetscParMatrix::Shift ( const Vector s)

Shift diagonal by a vector.

Definition at line 1997 of file petsc.cpp.

PetscParMatrix * mfem::PetscParMatrix::Transpose ( bool  action = false)

Returns the transpose of the PetscParMatrix.

If action is false, the new matrix is constructed with the PETSc function MatTranspose(). If action is true, then the matrix is not actually transposed. Instead, an object that behaves like the transpose is returned.

Definition at line 1891 of file petsc.cpp.

Friends And Related Function Documentation

friend class PetscLinearSolver
friend

Definition at line 336 of file petsc.hpp.

friend class PetscPreconditioner
friend

Definition at line 337 of file petsc.hpp.

Member Data Documentation

petsc::Mat mfem::PetscParMatrix::A
protected

The actual PETSc object.

Definition at line 311 of file petsc.hpp.

PetscParVector* mfem::PetscParMatrix::X
mutableprotected

Auxiliary vectors for typecasting.

Definition at line 314 of file petsc.hpp.

PetscParVector * mfem::PetscParMatrix::Y
mutableprotected

Definition at line 314 of file petsc.hpp.


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