libMesh
Public Member Functions | Static Public Member Functions | Protected Types | Protected Member Functions | Protected Attributes | Static Protected Attributes | Private Attributes | List of all members
libMesh::PetscMatrix< T > Class Template Referencefinal

This class provides a nice interface to the PETSc C-based data structures for parallel, sparse matrices. More...

#include <petsc_matrix.h>

Inheritance diagram for libMesh::PetscMatrix< T >:
[legend]

Public Member Functions

 PetscMatrix (const Parallel::Communicator &comm_in)
 Constructor; initializes the matrix to be empty, without any structure, i.e. More...
 
 PetscMatrix (Mat m, const Parallel::Communicator &comm_in)
 Constructor. More...
 
 PetscMatrix (PetscMatrix &&)=delete
 This class manages a C-style struct (Mat) manually, so we don't want to allow any automatic copy/move functions to be generated, and we can't default the destructor. More...
 
 PetscMatrix (const PetscMatrix &)=delete
 
PetscMatrixoperator= (const PetscMatrix &)=delete
 
PetscMatrixoperator= (PetscMatrix &&)=delete
 
virtual ~PetscMatrix ()
 
void set_matrix_type (PetscMatrixType mat_type)
 
virtual void init (const numeric_index_type m, const numeric_index_type n, const numeric_index_type m_l, const numeric_index_type n_l, const numeric_index_type nnz=30, const numeric_index_type noz=10, const numeric_index_type blocksize=1) override
 Initialize SparseMatrix with the specified sizes. More...
 
void init (const numeric_index_type m, const numeric_index_type n, const numeric_index_type m_l, const numeric_index_type n_l, const std::vector< numeric_index_type > &n_nz, const std::vector< numeric_index_type > &n_oz, const numeric_index_type blocksize=1)
 Initialize a PETSc matrix. More...
 
virtual void init (ParallelType=PARALLEL) override
 Initialize this matrix using the sparsity structure computed by dof_map. More...
 
void update_preallocation_and_zero ()
 Update the sparsity pattern based on dof_map, and set the matrix to zero. More...
 
void reset_preallocation ()
 Reset matrix to use the original nonzero pattern provided by users. More...
 
virtual void clear () noexcept override
 clear() is called from the destructor, so it should not throw. More...
 
virtual void zero () override
 Set all entries to 0. More...
 
virtual std::unique_ptr< SparseMatrix< T > > zero_clone () const override
 
virtual std::unique_ptr< SparseMatrix< T > > clone () const override
 
virtual void zero_rows (std::vector< numeric_index_type > &rows, T diag_value=0.0) override
 Sets all row entries to 0 then puts diag_value in the diagonal entry. More...
 
virtual void close () override
 Calls the SparseMatrix's internal assembly routines, ensuring that the values are consistent across processors. More...
 
virtual void flush () override
 For PETSc matrix , this function is similar to close but without shrinking memory. More...
 
virtual numeric_index_type m () const override
 
virtual numeric_index_type local_m () const final
 Get the number of rows owned by this process. More...
 
virtual numeric_index_type n () const override
 
virtual numeric_index_type local_n () const final
 Get the number of columns owned by this process. More...
 
void get_local_size (numeric_index_type &m, numeric_index_type &n) const
 Get the number of rows and columns owned by this process. More...
 
virtual numeric_index_type row_start () const override
 
virtual numeric_index_type row_stop () const override
 
virtual numeric_index_type col_start () const override
 
virtual numeric_index_type col_stop () const override
 
virtual void set (const numeric_index_type i, const numeric_index_type j, const T value) override
 Set the element (i,j) to value. More...
 
virtual void add (const numeric_index_type i, const numeric_index_type j, const T value) override
 Add value to the element (i,j). More...
 
virtual void add_matrix (const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols) override
 Add the full matrix dm to the SparseMatrix. More...
 
virtual void add_matrix (const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &dof_indices) override
 Same as add_matrix, but assumes the row and column maps are the same. More...
 
virtual void add_block_matrix (const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &brows, const std::vector< numeric_index_type > &bcols) override
 Add the full matrix dm to the SparseMatrix. More...
 
virtual void add_block_matrix (const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &dof_indices) override
 Same as add_block_matrix(), but assumes the row and column maps are the same. More...
 
virtual void add (const T a, const SparseMatrix< T > &X) override
 Compute A += a*X for scalar a, matrix X. More...
 
virtual void matrix_matrix_mult (SparseMatrix< T > &X, SparseMatrix< T > &Y, bool reuse=false) override
 Compute Y = A*X for matrix X. More...
 
virtual void add_sparse_matrix (const SparseMatrix< T > &spm, const std::map< numeric_index_type, numeric_index_type > &row_ltog, const std::map< numeric_index_type, numeric_index_type > &col_ltog, const T scalar) override
 Add scalar* spm to the rows and cols of this matrix (A): A(rows[i], cols[j]) += scalar * spm(i,j) More...
 
virtual T operator() (const numeric_index_type i, const numeric_index_type j) const override
 
virtual Real l1_norm () const override
 
virtual Real linfty_norm () const override
 
virtual bool closed () const override
 
virtual void set_destroy_mat_on_exit (bool destroy=true)
 If set to false, we don't delete the Mat on destruction and allow instead for PETSc to manage it. More...
 
virtual void print_personal (std::ostream &os=libMesh::out) const override
 Print the contents of the matrix to the screen with the PETSc viewer. More...
 
virtual void print_matlab (const std::string &name="") const override
 Print the contents of the matrix in Matlab's sparse matrix format. More...
 
virtual void print_petsc_binary (const std::string &filename) override
 Write the contents of the matrix to a file in PETSc's binary sparse matrix format. More...
 
virtual void print_petsc_hdf5 (const std::string &filename) override
 Write the contents of the matrix to a file in PETSc's HDF5 sparse matrix format. More...
 
virtual void read_petsc_binary (const std::string &filename) override
 Read the contents of the matrix from a file in PETSc's binary sparse matrix format. More...
 
virtual void read_petsc_hdf5 (const std::string &filename) override
 Read the contents of the matrix from a file in PETSc's HDF5 sparse matrix format. More...
 
virtual void get_diagonal (NumericVector< T > &dest) const override
 Copies the diagonal part of the matrix into dest. More...
 
virtual void get_transpose (SparseMatrix< T > &dest) const override
 Copies the transpose of the matrix into dest, which may be *this. More...
 
void swap (PetscMatrix< T > &)
 Swaps the internal data pointers of two PetscMatrices, no actual values are swapped. More...
 
Mat mat ()
 
virtual void get_row (numeric_index_type i, std::vector< numeric_index_type > &indices, std::vector< T > &values) const override
 Get a row from the matrix. More...
 
virtual void create_submatrix_nosort (SparseMatrix< T > &submatrix, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols) const override
 Similar to the create_submatrix function, this function creates a submatrix which is defined by the indices given in the rows and cols vectors. More...
 
virtual bool initialized () const
 
void attach_dof_map (const DofMap &dof_map)
 Set a pointer to the DofMap to use. More...
 
void attach_sparsity_pattern (const SparsityPattern::Build &sp)
 Set a pointer to a sparsity pattern to use. More...
 
virtual bool need_full_sparsity_pattern () const
 
virtual void update_sparsity_pattern (const SparsityPattern::Graph &)
 Updates the matrix sparsity pattern. More...
 
virtual std::size_t n_nonzeros () const
 
void print (std::ostream &os=libMesh::out, const bool sparse=false) const
 Print the contents of the matrix to the screen in a uniform style, regardless of matrix/solver package being used. More...
 
template<>
void print (std::ostream &os, const bool sparse) const
 
virtual void read_matlab (const std::string &filename)
 Read the contents of the matrix from Matlab's sparse matrix format. More...
 
virtual void create_submatrix (SparseMatrix< T > &submatrix, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols) const
 This function creates a matrix called "submatrix" which is defined by the row and column indices given in the "rows" and "cols" entries. More...
 
virtual void reinit_submatrix (SparseMatrix< T > &submatrix, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols) const
 This function is similar to the one above, but it allows you to reuse the existing sparsity pattern of "submatrix" instead of reallocating it again. More...
 
void vector_mult (NumericVector< T > &dest, const NumericVector< T > &arg) const
 Multiplies the matrix by the NumericVector arg and stores the result in NumericVector dest. More...
 
void vector_mult_add (NumericVector< T > &dest, const NumericVector< T > &arg) const
 Multiplies the matrix by the NumericVector arg and adds the result to the NumericVector dest. More...
 
const Parallel::Communicatorcomm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Static Public Member Functions

static std::unique_ptr< SparseMatrix< T > > build (const Parallel::Communicator &comm, const SolverPackage solver_package=libMesh::default_solver_package(), const MatrixBuildType matrix_build_type=MatrixBuildType::AUTOMATIC)
 Builds a SparseMatrix<T> using the linear solver package specified by solver_package. More...
 
static std::string get_info ()
 Gets a string containing the reference information. More...
 
static void print_info (std::ostream &out_stream=libMesh::out)
 Prints the reference information, by default to libMesh::out. More...
 
static unsigned int n_objects ()
 Prints the number of outstanding (created, but not yet destroyed) objects. More...
 
static void enable_print_counter_info ()
 Methods to enable/disable the reference counter output from print_info() More...
 
static void disable_print_counter_info ()
 

Protected Types

typedef std::map< std::string, std::pair< unsigned int, unsigned int > > Counts
 Data structure to log the information. More...
 

Protected Member Functions

virtual void _get_submatrix (SparseMatrix< T > &submatrix, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols, const bool reuse_submatrix) const override
 This function either creates or re-initializes a matrix called submatrix which is defined by the indices given in the rows and cols vectors. More...
 
void _petsc_viewer (const std::string &filename, PetscViewerType viewertype, PetscFileMode filemode)
 
void increment_constructor_count (const std::string &name) noexcept
 Increments the construction counter. More...
 
void increment_destructor_count (const std::string &name) noexcept
 Increments the destruction counter. More...
 

Protected Attributes

DofMap const * _dof_map
 The DofMap object associated with this object. More...
 
SparsityPattern::Build const * _sp
 The sparsity pattern associated with this object. More...
 
bool _is_initialized
 Flag indicating whether or not the matrix has been initialized. More...
 
const Parallel::Communicator_communicator
 

Static Protected Attributes

static Counts _counts
 Actually holds the data. More...
 
static Threads::atomic< unsigned int_n_objects
 The number of objects. More...
 
static Threads::spin_mutex _mutex
 Mutual exclusion object to enable thread-safe reference counting. More...
 
static bool _enable_print_counter = true
 Flag to control whether reference count information is printed when print_info is called. More...
 

Private Attributes

Mat _mat
 PETSc matrix datatype to store values. More...
 
bool _destroy_mat_on_exit
 This boolean value should only be set to false for the constructor which takes a PETSc Mat object. More...
 
PetscMatrixType _mat_type
 
std::mutex _petsc_matrix_mutex
 
Threads::spin_mutex _petsc_matrix_mutex
 

Detailed Description

template<typename T>
class libMesh::PetscMatrix< T >

This class provides a nice interface to the PETSc C-based data structures for parallel, sparse matrices.

All overridden virtual functions are documented in sparse_matrix.h.

Author
Benjamin S. Kirk
Date
2002 SparseMatrix interface to PETSc Mat.

Definition at line 93 of file petsc_matrix.h.

Member Typedef Documentation

◆ Counts

typedef std::map<std::string, std::pair<unsigned int, unsigned int> > libMesh::ReferenceCounter::Counts
protectedinherited

Data structure to log the information.

The log is identified by the class name.

Definition at line 119 of file reference_counter.h.

Constructor & Destructor Documentation

◆ PetscMatrix() [1/4]

template<typename T >
libMesh::PetscMatrix< T >::PetscMatrix ( const Parallel::Communicator comm_in)
explicit

Constructor; initializes the matrix to be empty, without any structure, i.e.

the matrix is not usable at all. This constructor is therefore only useful for matrices which are members of a class. All other matrices should be created at a point in the data flow where all necessary information is available.

You have to initialize the matrix before usage with init(...).

Definition at line 83 of file petsc_matrix.C.

83  :
84  SparseMatrix<T>(comm_in),
85  _mat(nullptr),
87  _mat_type(AIJ)
88 {}
bool _destroy_mat_on_exit
This boolean value should only be set to false for the constructor which takes a PETSc Mat object...
Definition: petsc_matrix.h:373
Generic sparse matrix.
Definition: vector_fe_ex5.C:45
PetscMatrixType _mat_type
Definition: petsc_matrix.h:375
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ PetscMatrix() [2/4]

template<typename T >
libMesh::PetscMatrix< T >::PetscMatrix ( Mat  m,
const Parallel::Communicator comm_in 
)
explicit

Constructor.

Creates a PetscMatrix assuming you already have a valid Mat object. In this case, m is NOT destroyed by the PetscMatrix destructor when this object goes out of scope. This allows ownership of m to remain with the original creator, and to simply provide additional functionality with the PetscMatrix.

Definition at line 95 of file petsc_matrix.C.

96  :
97  SparseMatrix<T>(comm_in),
98  _destroy_mat_on_exit(false),
99  _mat_type(AIJ)
100 {
101  this->_mat = mat_in;
102  this->_is_initialized = true;
103 }
bool _destroy_mat_on_exit
This boolean value should only be set to false for the constructor which takes a PETSc Mat object...
Definition: petsc_matrix.h:373
Generic sparse matrix.
Definition: vector_fe_ex5.C:45
bool _is_initialized
Flag indicating whether or not the matrix has been initialized.
PetscMatrixType _mat_type
Definition: petsc_matrix.h:375
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ PetscMatrix() [3/4]

template<typename T>
libMesh::PetscMatrix< T >::PetscMatrix ( PetscMatrix< T > &&  )
delete

This class manages a C-style struct (Mat) manually, so we don't want to allow any automatic copy/move functions to be generated, and we can't default the destructor.

◆ PetscMatrix() [4/4]

template<typename T>
libMesh::PetscMatrix< T >::PetscMatrix ( const PetscMatrix< T > &  )
delete

◆ ~PetscMatrix()

template<typename T >
libMesh::PetscMatrix< T >::~PetscMatrix ( )
virtual

Definition at line 109 of file petsc_matrix.C.

110 {
111  this->clear();
112 }
virtual void clear() noexcept override
clear() is called from the destructor, so it should not throw.
Definition: petsc_matrix.C:503

Member Function Documentation

◆ _get_submatrix()

template<typename T>
void libMesh::PetscMatrix< T >::_get_submatrix ( SparseMatrix< T > &  submatrix,
const std::vector< numeric_index_type > &  rows,
const std::vector< numeric_index_type > &  cols,
const bool  reuse_submatrix 
) const
overrideprotectedvirtual

This function either creates or re-initializes a matrix called submatrix which is defined by the indices given in the rows and cols vectors.

This function is implemented in terms of MatGetSubMatrix(). The reuse_submatrix parameter determines whether or not PETSc will treat submatrix as one which has already been used (had memory allocated) or as a new matrix.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 885 of file petsc_matrix.C.

889 {
890  if (!this->closed())
891  {
892  libmesh_deprecated();
893  libmesh_warning("The matrix must be assembled before calling PetscMatrix::create_submatrix().\n"
894  "Please update your code, as this warning will become an error in a future release.");
895  const_cast<PetscMatrix<T> *>(this)->close();
896  }
897 
898  // Make sure the SparseMatrix passed in is really a PetscMatrix
899  PetscMatrix<T> * petsc_submatrix = cast_ptr<PetscMatrix<T> *>(&submatrix);
900 
901  // If we're not reusing submatrix and submatrix is already initialized
902  // then we need to clear it, otherwise we get a memory leak.
903  if (!reuse_submatrix && submatrix.initialized())
904  submatrix.clear();
905 
906  // Construct row and column index sets.
907  PetscErrorCode ierr=0;
908 
909  WrappedPetsc<IS> isrow;
910  ierr = ISCreateGeneral(this->comm().get(),
911  cast_int<PetscInt>(rows.size()),
912  numeric_petsc_cast(rows.data()),
913  PETSC_USE_POINTER,
914  isrow.get()); LIBMESH_CHKERR(ierr);
915 
916  WrappedPetsc<IS> iscol;
917  ierr = ISCreateGeneral(this->comm().get(),
918  cast_int<PetscInt>(cols.size()),
919  numeric_petsc_cast(cols.data()),
920  PETSC_USE_POINTER,
921  iscol.get()); LIBMESH_CHKERR(ierr);
922 
923  // Extract submatrix
924  ierr = LibMeshCreateSubMatrix(_mat,
925  isrow,
926  iscol,
927  (reuse_submatrix ? MAT_REUSE_MATRIX : MAT_INITIAL_MATRIX),
928  &(petsc_submatrix->_mat)); LIBMESH_CHKERR(ierr);
929 
930  // Specify that the new submatrix is initialized and close it.
931  petsc_submatrix->_is_initialized = true;
932  petsc_submatrix->close();
933 }
virtual bool initialized() const
const Parallel::Communicator & comm() const
virtual void clear()=0
Restores the SparseMatrix<T> to a pristine state.
PetscInt * numeric_petsc_cast(const numeric_index_type *p)
bool _is_initialized
Flag indicating whether or not the matrix has been initialized.
virtual bool closed() const override
This class provides a nice interface to the PETSc C-based data structures for parallel, sparse matrices.
Definition: petsc_matrix.h:93
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367
virtual void close() override
Calls the SparseMatrix&#39;s internal assembly routines, ensuring that the values are consistent across p...

◆ _petsc_viewer()

template<typename T >
void libMesh::PetscMatrix< T >::_petsc_viewer ( const std::string &  filename,
PetscViewerType  viewertype,
PetscFileMode  filemode 
)
protected

Definition at line 738 of file petsc_matrix.C.

741 {
742  parallel_object_only();
743 
744  PetscErrorCode ierr=0;
745 
746  // We'll get matrix sizes from the file, but we need to at least
747  // have a Mat object
748  if (!this->initialized())
749  {
750  ierr = MatCreate(this->comm().get(), &_mat);
751  LIBMESH_CHKERR(ierr);
752  this->_is_initialized = true;
753  }
754 
755  PetscViewer viewer;
756  ierr = PetscViewerCreate(this->comm().get(), &viewer);
757  LIBMESH_CHKERR(ierr);
758  ierr = PetscViewerSetType(viewer, viewertype);
759  LIBMESH_CHKERR(ierr);
760  ierr = PetscViewerSetFromOptions(viewer);
761  LIBMESH_CHKERR(ierr);
762  ierr = PetscViewerFileSetMode(viewer, filemode);
763  LIBMESH_CHKERR(ierr);
764  ierr = PetscViewerFileSetName(viewer, filename.c_str());
765  LIBMESH_CHKERR(ierr);
766  if (filemode == FILE_MODE_READ)
767  ierr = MatLoad(_mat, viewer);
768  else
769  ierr = MatView(_mat, viewer);
770  LIBMESH_CHKERR(ierr);
771  ierr = PetscViewerDestroy(&viewer);
772  LIBMESH_CHKERR(ierr);
773 }
virtual bool initialized() const
const Parallel::Communicator & comm() const
bool _is_initialized
Flag indicating whether or not the matrix has been initialized.
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ add() [1/2]

template<typename T>
void libMesh::PetscMatrix< T >::add ( const numeric_index_type  i,
const numeric_index_type  j,
const T  value 
)
overridevirtual

Add value to the element (i,j).

Throws an error if the entry does not exist. Zero values can be "added" to non-existent entries.

Implements libMesh::SparseMatrix< T >.

Definition at line 1234 of file petsc_matrix.C.

1237 {
1238  libmesh_assert (this->initialized());
1239 
1240  PetscErrorCode ierr=0;
1241  PetscInt i_val=i, j_val=j;
1242 
1243  PetscScalar petsc_value = static_cast<PetscScalar>(value);
1244  ierr = MatSetValues(_mat, 1, &i_val, 1, &j_val,
1245  &petsc_value, ADD_VALUES);
1246  LIBMESH_CHKERR(ierr);
1247 }
virtual bool initialized() const
libmesh_assert(ctx)
static const bool value
Definition: xdr_io.C:54
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ add() [2/2]

template<typename T>
void libMesh::PetscMatrix< T >::add ( const T  a,
const SparseMatrix< T > &  X 
)
overridevirtual

Compute A += a*X for scalar a, matrix X.

Note
The matrices A and X need to have the same nonzero pattern, otherwise PETSc will crash!
It is advisable to not only allocate appropriate memory with init(), but also explicitly zero the terms of this whenever you add a non-zero value to X.
X will be closed, if not already done, before performing any work.

Implements libMesh::SparseMatrix< T >.

Definition at line 1265 of file petsc_matrix.C.

1266 {
1267  libmesh_assert (this->initialized());
1268 
1269  // sanity check. but this cannot avoid
1270  // crash due to incompatible sparsity structure...
1271  libmesh_assert_equal_to (this->m(), X_in.m());
1272  libmesh_assert_equal_to (this->n(), X_in.n());
1273 
1274  PetscScalar a = static_cast<PetscScalar> (a_in);
1275  const PetscMatrix<T> * X = cast_ptr<const PetscMatrix<T> *> (&X_in);
1276 
1277  libmesh_assert (X);
1278 
1279  PetscErrorCode ierr=0;
1280 
1281  // the matrix from which we copy the values has to be assembled/closed
1282  libmesh_assert(X->closed());
1283 
1284  semiparallel_only();
1285 
1286  ierr = MatAXPY(_mat, a, X->_mat, DIFFERENT_NONZERO_PATTERN);
1287  LIBMESH_CHKERR(ierr);
1288 }
virtual bool initialized() const
virtual numeric_index_type n() const override
libmesh_assert(ctx)
virtual numeric_index_type m() const override
virtual bool closed() const override
This class provides a nice interface to the PETSc C-based data structures for parallel, sparse matrices.
Definition: petsc_matrix.h:93
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ add_block_matrix() [1/2]

template<typename T>
void libMesh::PetscMatrix< T >::add_block_matrix ( const DenseMatrix< T > &  dm,
const std::vector< numeric_index_type > &  brows,
const std::vector< numeric_index_type > &  bcols 
)
overridevirtual

Add the full matrix dm to the SparseMatrix.

This is useful for adding an element matrix at assembly time. The matrix is assumed blocked, and brow, bcol correspond to the block row and column indices.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 841 of file petsc_matrix.C.

Referenced by libMesh::PetscMatrix< libMesh::Number >::add_block_matrix().

844 {
845  libmesh_assert (this->initialized());
846 
847  const numeric_index_type n_brows =
848  cast_int<numeric_index_type>(brows.size());
849  const numeric_index_type n_bcols =
850  cast_int<numeric_index_type>(bcols.size());
851 
852  PetscErrorCode ierr=0;
853 
854 #ifndef NDEBUG
855  const numeric_index_type n_rows =
856  cast_int<numeric_index_type>(dm.m());
857  const numeric_index_type n_cols =
858  cast_int<numeric_index_type>(dm.n());
859  const numeric_index_type blocksize = n_rows / n_brows;
860 
861  libmesh_assert_equal_to (n_cols / n_bcols, blocksize);
862  libmesh_assert_equal_to (blocksize*n_brows, n_rows);
863  libmesh_assert_equal_to (blocksize*n_bcols, n_cols);
864 
865  PetscInt petsc_blocksize;
866  ierr = MatGetBlockSize(_mat, &petsc_blocksize);
867  LIBMESH_CHKERR(ierr);
868  libmesh_assert_equal_to (blocksize, static_cast<numeric_index_type>(petsc_blocksize));
869 #endif
870 
871  // These casts are required for PETSc <= 2.1.5
872  ierr = MatSetValuesBlocked(_mat,
873  n_brows, numeric_petsc_cast(brows.data()),
874  n_bcols, numeric_petsc_cast(bcols.data()),
875  pPS(const_cast<T*>(dm.get_values().data())),
876  ADD_VALUES);
877  LIBMESH_CHKERR(ierr);
878 }
virtual bool initialized() const
unsigned int m() const
PetscScalar * pPS(T *ptr)
Definition: petsc_macro.h:172
PetscInt * numeric_petsc_cast(const numeric_index_type *p)
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
std::vector< T > & get_values()
Definition: dense_matrix.h:382
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367
unsigned int n() const

◆ add_block_matrix() [2/2]

template<typename T>
virtual void libMesh::PetscMatrix< T >::add_block_matrix ( const DenseMatrix< T > &  dm,
const std::vector< numeric_index_type > &  dof_indices 
)
inlineoverridevirtual

Same as add_block_matrix(), but assumes the row and column maps are the same.

Thus the matrix dm must be square.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 236 of file petsc_matrix.h.

238  { this->add_block_matrix (dm, dof_indices, dof_indices); }
virtual void add_block_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &brows, const std::vector< numeric_index_type > &bcols) override
Add the full matrix dm to the SparseMatrix.
Definition: petsc_matrix.C:841

◆ add_matrix() [1/2]

template<typename T>
void libMesh::PetscMatrix< T >::add_matrix ( const DenseMatrix< T > &  dm,
const std::vector< numeric_index_type > &  rows,
const std::vector< numeric_index_type > &  cols 
)
overridevirtual

Add the full matrix dm to the SparseMatrix.

This is useful for adding an element matrix at assembly time.

Implements libMesh::SparseMatrix< T >.

Definition at line 814 of file petsc_matrix.C.

817 {
818  libmesh_assert (this->initialized());
819 
820  const numeric_index_type n_rows = dm.m();
821  const numeric_index_type n_cols = dm.n();
822 
823  libmesh_assert_equal_to (rows.size(), n_rows);
824  libmesh_assert_equal_to (cols.size(), n_cols);
825 
826  PetscErrorCode ierr=0;
827  ierr = MatSetValues(_mat,
828  n_rows, numeric_petsc_cast(rows.data()),
829  n_cols, numeric_petsc_cast(cols.data()),
830  pPS(const_cast<T*>(dm.get_values().data())),
831  ADD_VALUES);
832  LIBMESH_CHKERR(ierr);
833 }
virtual bool initialized() const
unsigned int m() const
PetscScalar * pPS(T *ptr)
Definition: petsc_macro.h:172
PetscInt * numeric_petsc_cast(const numeric_index_type *p)
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
std::vector< T > & get_values()
Definition: dense_matrix.h:382
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367
unsigned int n() const

◆ add_matrix() [2/2]

template<typename T>
void libMesh::PetscMatrix< T >::add_matrix ( const DenseMatrix< T > &  dm,
const std::vector< numeric_index_type > &  dof_indices 
)
overridevirtual

Same as add_matrix, but assumes the row and column maps are the same.

Thus the matrix dm must be square.

Implements libMesh::SparseMatrix< T >.

Definition at line 1252 of file petsc_matrix.C.

1254 {
1255  this->add_matrix (dm, dof_indices, dof_indices);
1256 }
virtual void add_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols) override
Add the full matrix dm to the SparseMatrix.
Definition: petsc_matrix.C:814

◆ add_sparse_matrix()

template<typename T>
void libMesh::PetscMatrix< T >::add_sparse_matrix ( const SparseMatrix< T > &  spm,
const std::map< numeric_index_type, numeric_index_type > &  row_ltog,
const std::map< numeric_index_type, numeric_index_type > &  col_ltog,
const T  scalar 
)
overridevirtual

Add scalar* spm to the rows and cols of this matrix (A): A(rows[i], cols[j]) += scalar * spm(i,j)

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 1328 of file petsc_matrix.C.

1332 {
1333  // size of spm is usually greater than row_ltog and col_ltog in parallel as the indices are owned by the processor
1334  // also, we should allow adding certain parts of spm to _mat
1335  libmesh_assert_greater_equal(spm.m(), row_ltog.size());
1336  libmesh_assert_greater_equal(spm.n(), col_ltog.size());
1337 
1338  // make sure matrix has larger size than spm
1339  libmesh_assert_greater_equal(this->m(), spm.m());
1340  libmesh_assert_greater_equal(this->n(), spm.n());
1341 
1342  if (!this->closed())
1343  this->close();
1344 
1345  PetscErrorCode ierr = 0;
1346 
1347  auto pscm = cast_ptr<const PetscMatrix<T> *>(&spm);
1348 
1349  PetscInt ncols = 0;
1350 
1351  const PetscInt * lcols;
1352  const PetscScalar * vals;
1353 
1354  std::vector<PetscInt> gcols;
1355  std::vector<PetscScalar> values;
1356 
1357  for (auto ltog : row_ltog)
1358  {
1359  PetscInt grow[] = {static_cast<PetscInt>(ltog.second)}; // global row index
1360 
1361  ierr = MatGetRow(pscm->_mat, static_cast<PetscInt>(ltog.first), &ncols, &lcols, &vals);
1362  LIBMESH_CHKERR(ierr);
1363 
1364  // get global indices (gcols) from lcols, increment values = vals*scalar
1365  gcols.resize(ncols);
1366  values.resize(ncols);
1367  for (auto i : index_range(gcols))
1368  {
1369  gcols[i] = libmesh_map_find(col_ltog, lcols[i]);
1370  values[i] = PS(scalar) * vals[i];
1371  }
1372 
1373  ierr = MatSetValues(_mat, 1, grow, ncols, gcols.data(), values.data(), ADD_VALUES);
1374  LIBMESH_CHKERR(ierr);
1375  ierr = MatRestoreRow(pscm->_mat, static_cast<PetscInt>(ltog.first), &ncols, &lcols, &vals);
1376  LIBMESH_CHKERR(ierr);
1377  }
1378  // Note: We are not closing the matrix because it is expensive to do so when adding multiple sparse matrices.
1379  // Remember to manually close the matrix once all changes to the matrix have been made.
1380 }
virtual numeric_index_type n() const override
PetscScalar PS(T val)
Definition: petsc_macro.h:166
virtual numeric_index_type m() const =0
virtual numeric_index_type m() const override
virtual bool closed() const override
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367
virtual void close() override
Calls the SparseMatrix&#39;s internal assembly routines, ensuring that the values are consistent across p...
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111
virtual numeric_index_type n() const =0

◆ attach_dof_map()

template<typename T >
void libMesh::SparseMatrix< T >::attach_dof_map ( const DofMap dof_map)
inherited

Set a pointer to the DofMap to use.

If a separate sparsity pattern is not being used, use the one from the DofMap.

The lifetime of dof_map must exceed the lifetime of this.

Definition at line 63 of file sparse_matrix.C.

Referenced by libMesh::__libmesh_tao_hessian(), and libMesh::DofMap::update_sparsity_pattern().

64 {
65  _dof_map = &dof_map;
66  if (!_sp)
67  _sp = dof_map.get_sparsity_pattern();
68 }
SparsityPattern::Build const * _sp
The sparsity pattern associated with this object.
const SparsityPattern::Graph & get_sparsity_pattern() const
Rows of sparse matrix indices, indexed by the offset from the first DoF on this processor.
DofMap const * _dof_map
The DofMap object associated with this object.

◆ attach_sparsity_pattern()

template<typename T >
void libMesh::SparseMatrix< T >::attach_sparsity_pattern ( const SparsityPattern::Build sp)
inherited

Set a pointer to a sparsity pattern to use.

Useful in cases where a matrix requires a wider (or for efficiency narrower) pattern than most matrices in the system, or in cases where no system sparsity pattern is being calculated by the DofMap.

The lifetime of sp must exceed the lifetime of this.

Definition at line 73 of file sparse_matrix.C.

Referenced by libMesh::DofMap::update_sparsity_pattern().

74 {
75  _sp = &sp;
76 }
SparsityPattern::Build const * _sp
The sparsity pattern associated with this object.

◆ build()

template<typename T >
std::unique_ptr< SparseMatrix< T > > libMesh::SparseMatrix< T >::build ( const Parallel::Communicator comm,
const SolverPackage  solver_package = libMesh::default_solver_package(),
const MatrixBuildType  matrix_build_type = MatrixBuildType::AUTOMATIC 
)
staticinherited

Builds a SparseMatrix<T> using the linear solver package specified by solver_package.

Definition at line 156 of file sparse_matrix.C.

Referenced by libMesh::TransientRBConstruction::allocate_data_structures(), libMesh::RBConstruction::allocate_data_structures(), libMesh::DofMap::process_mesh_constraint_rows(), ConstraintOperatorTest::test1DCoarseningNewNodes(), ConstraintOperatorTest::test1DCoarseningOperator(), ConstraintOperatorTest::testCoreform(), SystemsTest::testProjectMatrix1D(), SystemsTest::testProjectMatrix2D(), and SystemsTest::testProjectMatrix3D().

159 {
160  // Avoid unused parameter warnings when no solver packages are enabled.
162 
163  if (matrix_build_type == MatrixBuildType::DIAGONAL)
164  return std::make_unique<DiagonalMatrix<T>>(comm);
165 
166  // Build the appropriate vector
167  switch (solver_package)
168  {
169 
170 #ifdef LIBMESH_HAVE_LASPACK
171  case LASPACK_SOLVERS:
172  return std::make_unique<LaspackMatrix<T>>(comm);
173 #endif
174 
175 
176 #ifdef LIBMESH_HAVE_PETSC
177  case PETSC_SOLVERS:
178  return std::make_unique<PetscMatrix<T>>(comm);
179 #endif
180 
181 
182 #ifdef LIBMESH_TRILINOS_HAVE_EPETRA
183  case TRILINOS_SOLVERS:
184  return std::make_unique<EpetraMatrix<T>>(comm);
185 #endif
186 
187 
188 #ifdef LIBMESH_HAVE_EIGEN
189  case EIGEN_SOLVERS:
190  return std::make_unique<EigenSparseMatrix<T>>(comm);
191 #endif
192 
193  default:
194  libmesh_error_msg("ERROR: Unrecognized solver package: " << solver_package);
195  }
196 }
const Parallel::Communicator & comm() const
void libmesh_ignore(const Args &...)

◆ clear()

template<typename T >
void libMesh::PetscMatrix< T >::clear ( )
overridevirtualnoexcept

clear() is called from the destructor, so it should not throw.

Implements libMesh::SparseMatrix< T >.

Definition at line 503 of file petsc_matrix.C.

504 {
505  if ((this->initialized()) && (this->_destroy_mat_on_exit))
506  {
507  exceptionless_semiparallel_only();
508 
509  // If we encounter an error here, print a warning but otherwise
510  // keep going since we may be recovering from an exception.
511  PetscErrorCode ierr = MatDestroy (&_mat);
512  if (ierr)
513  libmesh_warning("Warning: MatDestroy returned a non-zero error code which we ignored.");
514 
515  this->_is_initialized = false;
516  }
517 }
virtual bool initialized() const
bool _destroy_mat_on_exit
This boolean value should only be set to false for the constructor which takes a PETSc Mat object...
Definition: petsc_matrix.h:373
bool _is_initialized
Flag indicating whether or not the matrix has been initialized.
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ clone()

template<typename T >
std::unique_ptr< SparseMatrix< T > > libMesh::PetscMatrix< T >::clone ( ) const
overridevirtual
Returns
A smart pointer to a copy of this matrix.
Note
This must be overridden in the derived classes.

Implements libMesh::SparseMatrix< T >.

Definition at line 483 of file petsc_matrix.C.

484 {
485  libmesh_error_msg_if(!this->closed(), "Matrix must be closed before it can be cloned!");
486 
487  // Copy the nonzero pattern and numerical values
488  Mat copy;
489  PetscErrorCode ierr = MatDuplicate(_mat, MAT_COPY_VALUES, &copy);
490  LIBMESH_CHKERR(ierr);
491 
492  // Call wrapping PetscMatrix constructor, have it take over
493  // ownership.
494  auto ret = std::make_unique<PetscMatrix<T>>(copy, this->comm());
495  ret->set_destroy_mat_on_exit(true);
496 
497  // Work around an issue on older compilers. We are able to simply
498  // "return ret;" on newer compilers
499  return std::unique_ptr<SparseMatrix<T>>(ret.release());
500 }
const Parallel::Communicator & comm() const
virtual bool closed() const override
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ close()

template<typename T >
void libMesh::PetscMatrix< T >::close ( )
overridevirtual

Calls the SparseMatrix's internal assembly routines, ensuring that the values are consistent across processors.

Implements libMesh::SparseMatrix< T >.

Definition at line 1056 of file petsc_matrix.C.

Referenced by libMesh::PetscMatrix< libMesh::Number >::_get_submatrix(), libMesh::PetscMatrix< libMesh::Number >::create_submatrix_nosort(), libMesh::PetscMatrix< libMesh::Number >::get_transpose(), libMesh::libmesh_petsc_snes_jacobian(), libMesh::PetscLinearSolver< Number >::solve_common(), libMesh::SlepcEigenSolver< libMesh::Number >::solve_generalized(), and libMesh::SlepcEigenSolver< libMesh::Number >::solve_standard().

1057 {
1058  semiparallel_only();
1059 
1060  // BSK - 1/19/2004
1061  // strictly this check should be OK, but it seems to
1062  // fail on matrix-free matrices. Do they falsely
1063  // state they are assembled? Check with the developers...
1064  // if (this->closed())
1065  // return;
1066 
1067  MatAssemblyBeginEnd(this->comm(), _mat, MAT_FINAL_ASSEMBLY);
1068 }
const Parallel::Communicator & comm() const
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ closed()

template<typename T >
bool libMesh::PetscMatrix< T >::closed ( ) const
overridevirtual
Returns
true if the matrix has been assembled.

Implements libMesh::SparseMatrix< T >.

Definition at line 1497 of file petsc_matrix.C.

Referenced by libMesh::PetscMatrix< libMesh::Number >::add().

1498 {
1499  libmesh_assert (this->initialized());
1500 
1501  PetscErrorCode ierr=0;
1502  PetscBool assembled;
1503 
1504  ierr = MatAssembled(_mat, &assembled);
1505  LIBMESH_CHKERR(ierr);
1506 
1507  return (assembled == PETSC_TRUE);
1508 }
virtual bool initialized() const
libmesh_assert(ctx)
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ col_start()

template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::col_start ( ) const
overridevirtual
Returns
The index of the first matrix column owned by this processor.

Implements libMesh::SparseMatrix< T >.

Definition at line 1184 of file petsc_matrix.C.

1185 {
1186  libmesh_assert (this->initialized());
1187 
1188  PetscInt start=0, stop=0;
1189  PetscErrorCode ierr=0;
1190 
1191  ierr = MatGetOwnershipRangeColumn(_mat, &start, &stop);
1192  LIBMESH_CHKERR(ierr);
1193 
1194  return static_cast<numeric_index_type>(start);
1195 }
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
void stop(const char *file, int line, const char *date, const char *time)
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ col_stop()

template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::col_stop ( ) const
overridevirtual
Returns
The index of the last matrix column (+1) owned by this processor.

Implements libMesh::SparseMatrix< T >.

Definition at line 1200 of file petsc_matrix.C.

1201 {
1202  libmesh_assert (this->initialized());
1203 
1204  PetscInt start=0, stop=0;
1205  PetscErrorCode ierr=0;
1206 
1207  ierr = MatGetOwnershipRangeColumn(_mat, &start, &stop);
1208  LIBMESH_CHKERR(ierr);
1209 
1210  return static_cast<numeric_index_type>(stop);
1211 }
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
void stop(const char *file, int line, const char *date, const char *time)
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ comm()

const Parallel::Communicator& libMesh::ParallelObject::comm ( ) const
inlineinherited
Returns
A reference to the Parallel::Communicator object used by this mesh.

Definition at line 97 of file parallel_object.h.

References libMesh::ParallelObject::_communicator.

Referenced by libMesh::__libmesh_petsc_diff_solver_jacobian(), libMesh::__libmesh_petsc_diff_solver_monitor(), libMesh::__libmesh_petsc_diff_solver_residual(), libMesh::__libmesh_tao_equality_constraints(), libMesh::__libmesh_tao_equality_constraints_jacobian(), libMesh::__libmesh_tao_gradient(), libMesh::__libmesh_tao_hessian(), libMesh::__libmesh_tao_inequality_constraints(), libMesh::__libmesh_tao_inequality_constraints_jacobian(), libMesh::__libmesh_tao_objective(), libMesh::MeshRefinement::_coarsen_elements(), libMesh::ExactSolution::_compute_error(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::Partitioner::_find_global_index_by_pid_map(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::SlepcEigenSolver< libMesh::Number >::_petsc_shell_matrix_get_diagonal(), libMesh::PetscLinearSolver< Number >::_petsc_shell_matrix_get_diagonal(), libMesh::SlepcEigenSolver< libMesh::Number >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< Number >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< Number >::_petsc_shell_matrix_mult_add(), libMesh::MeshRefinement::_refine_elements(), libMesh::MeshRefinement::_smooth_flags(), libMesh::DofMap::add_constraints_to_send_list(), add_cube_convex_hull_to_mesh(), libMesh::PetscDMWrapper::add_dofs_helper(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::TransientRBConstruction::add_IC_to_RB_space(), libMesh::EigenSystem::add_matrices(), libMesh::System::add_matrix(), libMesh::RBConstruction::add_scaled_matrix_and_vector(), libMesh::System::add_variable(), libMesh::System::add_variables(), libMesh::System::add_vector(), libMesh::MeshTools::Modification::all_tri(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::DofMap::allgather_recursive_constraints(), libMesh::TransientRBConstruction::allocate_data_structures(), libMesh::RBConstruction::allocate_data_structures(), libMesh::TransientRBConstruction::assemble_affine_expansion(), libMesh::FEMSystem::assemble_qoi(), libMesh::Nemesis_IO::assert_symmetric_cmaps(), libMesh::MeshCommunication::assign_global_indices(), libMesh::Partitioner::assign_partitioning(), libMesh::MeshTools::Generation::build_extrusion(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::EquationSystems::build_parallel_solution_vector(), libMesh::PetscDMWrapper::build_section(), libMesh::PetscDMWrapper::build_sf(), libMesh::MeshBase::cache_elem_data(), libMesh::System::calculate_norm(), libMesh::DofMap::check_dirichlet_bcid_consistency(), libMesh::RBConstruction::compute_Fq_representor_innerprods(), libMesh::RBConstruction::compute_max_error_bound(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::RBConstruction::compute_output_dual_innerprods(), libMesh::RBConstruction::compute_residual_dual_norm_slow(), libMesh::RBSCMConstruction::compute_SCM_bounds_on_training_set(), libMesh::DofMap::computed_sparsity_already(), libMesh::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::ContinuationSystem::ContinuationSystem(), libMesh::MeshBase::copy_constraint_rows(), libMesh::ExodusII_IO::copy_elemental_solution(), libMesh::ExodusII_IO::copy_nodal_solution(), libMesh::ExodusII_IO::copy_scalar_solution(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::MeshTools::create_bounding_box(), libMesh::DofMap::create_dof_constraints(), libMesh::MeshTools::create_nodal_bounding_box(), libMesh::MeshRefinement::create_parent_error_vector(), libMesh::MeshTools::create_processor_bounding_box(), libMesh::MeshTools::create_subdomain_bounding_box(), libMesh::PetscMatrix< libMesh::Number >::create_submatrix_nosort(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::RBEIMEvaluation::distribute_bfs(), libMesh::DofMap::distribute_dofs(), DMlibMeshFunction(), DMlibMeshJacobian(), DMlibMeshSetSystem_libMesh(), DMVariableBounds_libMesh(), libMesh::DTKSolutionTransfer::DTKSolutionTransfer(), libMesh::MeshRefinement::eliminate_unrefined_patches(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_interiors(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_nodes(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_sides(), libMesh::TransientRBConstruction::enrich_RB_space(), libMesh::EpetraVector< T >::EpetraVector(), AssembleOptimization::equality_constraints(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::MeshRefinement::flag_elements_by_elem_fraction(), libMesh::MeshRefinement::flag_elements_by_error_fraction(), libMesh::MeshRefinement::flag_elements_by_error_tolerance(), libMesh::MeshRefinement::flag_elements_by_mean_stddev(), libMesh::MeshRefinement::flag_elements_by_nelem_target(), libMesh::RBEIMEvaluation::gather_bfs(), libMesh::DofMap::gather_constraints(), libMesh::MeshfreeInterpolation::gather_remote_data(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::RBEIMEvaluation::get_eim_basis_function_node_value(), libMesh::RBEIMEvaluation::get_eim_basis_function_side_value(), libMesh::RBEIMEvaluation::get_eim_basis_function_value(), libMesh::MeshBase::get_info(), libMesh::System::get_info(), libMesh::DofMap::get_info(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::RBEIMConstruction::get_max_abs_value(), libMesh::RBEIMConstruction::get_node_max_abs_value(), libMesh::RBEIMEvaluation::get_parametrized_function_node_value(), libMesh::RBEIMEvaluation::get_parametrized_function_side_value(), libMesh::RBEIMEvaluation::get_parametrized_function_value(), libMesh::RBEIMConstruction::get_random_point(), AssembleOptimization::inequality_constraints(), AssembleOptimization::inequality_constraints_jacobian(), libMesh::LocationMap< T >::init(), libMesh::TimeSolver::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_and_attach_petscdm(), libMesh::ExodusII_IO_Helper::initialize(), libMesh::OptimizationSystem::initialize_equality_constraints_storage(), libMesh::OptimizationSystem::initialize_inequality_constraints_storage(), libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set(), libMesh::RBEIMConstruction::inner_product(), integrate_function(), libMesh::MeshTools::libmesh_assert_consistent_distributed(), libMesh::MeshTools::libmesh_assert_consistent_distributed_nodes(), libMesh::MeshTools::libmesh_assert_contiguous_dof_ids(), libMesh::MeshTools::libmesh_assert_equal_connectivity(), libMesh::MeshTools::libmesh_assert_equal_points(), libMesh::MeshTools::libmesh_assert_parallel_consistent_new_node_procids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_topology_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_boundary_ids(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_flags(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_object_ids(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_p_levels(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::MeshTools::libmesh_assert_valid_unique_ids(), libMesh::libmesh_petsc_linesearch_shellfunc(), libMesh::libmesh_petsc_preconditioner_apply(), libMesh::libmesh_petsc_recalculate_monitor(), libMesh::libmesh_petsc_snes_fd_residual(), libMesh::libmesh_petsc_snes_jacobian(), libMesh::libmesh_petsc_snes_mffd_interface(), libMesh::libmesh_petsc_snes_mffd_residual(), libMesh::libmesh_petsc_snes_postcheck(), libMesh::libmesh_petsc_snes_precheck(), libMesh::libmesh_petsc_snes_residual(), libMesh::libmesh_petsc_snes_residual_helper(), libMesh::MeshRefinement::limit_level_mismatch_at_edge(), libMesh::MeshRefinement::limit_level_mismatch_at_node(), libMesh::MeshRefinement::limit_overrefined_boundary(), libMesh::MeshRefinement::limit_underrefined_boundary(), libMesh::LinearImplicitSystem::LinearImplicitSystem(), main(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshCommunication::make_elems_parallel_consistent(), libMesh::MeshRefinement::make_flags_parallel_consistent(), libMesh::MeshCommunication::make_new_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_new_nodes_parallel_consistent(), libMesh::MeshCommunication::make_node_bcids_parallel_consistent(), libMesh::MeshCommunication::make_node_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_unique_ids_parallel_consistent(), libMesh::MeshCommunication::make_nodes_parallel_consistent(), libMesh::MeshCommunication::make_p_levels_parallel_consistent(), libMesh::MeshRefinement::make_refinement_compatible(), libMesh::TransientRBConstruction::mass_matrix_scaled_matvec(), libMesh::FEMSystem::mesh_position_set(), libMesh::TriangulatorInterface::MeshedHole::MeshedHole(), LinearElasticityWithContact::move_mesh(), libMesh::DistributedMesh::n_active_elem(), libMesh::MeshTools::n_active_levels(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::DofMap::n_constrained_dofs(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::CondensedEigenSystem::n_global_non_condensed_dofs(), libMesh::MeshTools::n_levels(), MixedOrderTest::n_neighbor_links(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::SparsityPattern::Build::n_nonzeros(), libMesh::MeshTools::n_p_levels(), libMesh::BoundaryInfo::n_shellface_conds(), libMesh::RBEIMEvaluation::node_distribute_bfs(), libMesh::RBEIMEvaluation::node_gather_bfs(), libMesh::RBEIMConstruction::node_inner_product(), libMesh::MeshBase::operator==(), libMesh::DistributedMesh::parallel_max_elem_id(), libMesh::DistributedMesh::parallel_max_node_id(), libMesh::ReplicatedMesh::parallel_max_unique_id(), libMesh::DistributedMesh::parallel_max_unique_id(), libMesh::DistributedMesh::parallel_n_elem(), libMesh::DistributedMesh::parallel_n_nodes(), libMesh::SparsityPattern::Build::parallel_sync(), libMesh::BoundaryInfo::parallel_sync_node_ids(), libMesh::BoundaryInfo::parallel_sync_side_ids(), libMesh::MeshTools::paranoid_n_levels(), libMesh::Partitioner::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::petsc_auto_fieldsplit(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshBase::prepare_for_use(), libMesh::DofMap::print_dof_constraints(), libMesh::DofMap::process_mesh_constraint_rows(), libMesh::Partitioner::processor_pairs_to_interface_nodes(), libMesh::InterMeshProjection::project_system_vectors(), FEMParameters::read(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::EquationSystems::read(), libMesh::ExodusII_IO::read_header(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::System::read_header(), libMesh::RBEIMEvaluation::read_in_interior_basis_functions(), libMesh::RBEIMEvaluation::read_in_node_basis_functions(), libMesh::RBEIMEvaluation::read_in_side_basis_functions(), libMesh::RBEvaluation::read_in_vectors_from_multiple_files(), libMesh::System::read_legacy_data(), libMesh::TransientRBConstruction::read_riesz_representors_from_files(), libMesh::RBConstruction::read_riesz_representors_from_files(), libMesh::System::read_SCALAR_dofs(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::System::read_serialized_vector(), libMesh::Nemesis_IO_Helper::read_var_names_impl(), libMesh::MeshBase::recalculate_n_partitions(), libMesh::MeshRefinement::refine_and_coarsen_elements(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::DistributedMesh::renumber_nodes_and_elements(), LinearElasticityWithContact::residual_and_jacobian(), OverlappingAlgebraicGhostingTest::run_ghosting_test(), OverlappingCouplingGhostingTest::run_sparsity_pattern_test(), scale_mesh_and_plot(), libMesh::DofMap::scatter_constraints(), libMesh::CheckpointIO::select_split_config(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::send_and_insert_dof_values(), libMesh::TransientRBConstruction::set_error_temporal_data(), libMesh::Partitioner::set_interface_node_processor_ids_BFS(), libMesh::Partitioner::set_interface_node_processor_ids_linear(), libMesh::Partitioner::set_interface_node_processor_ids_petscpartitioner(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::Partitioner::set_parent_processor_ids(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::PetscDiffSolver::setup_petsc_data(), libMesh::RBEIMEvaluation::side_distribute_bfs(), libMesh::RBEIMEvaluation::side_gather_bfs(), libMesh::RBEIMConstruction::side_inner_product(), libMesh::Partitioner::single_partition(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::split_mesh(), libMesh::RBEIMConstruction::store_eim_solutions_for_training_set(), libMesh::MeshBase::subdomain_ids(), libMesh::BoundaryInfo::sync(), ConstraintOperatorTest::test1DCoarseningNewNodes(), ConstraintOperatorTest::test1DCoarseningOperator(), libMesh::MeshRefinement::test_level_one(), MeshfunctionDFEM::test_mesh_function_dfem(), MeshfunctionDFEM::test_mesh_function_dfem_grad(), MeshFunctionTest::test_p_level(), libMesh::MeshRefinement::test_unflagged(), DofMapTest::testBadElemFECombo(), SystemsTest::testBlockRestrictedVarNDofs(), BoundaryInfoTest::testBoundaryOnChildrenErrors(), ConstraintOperatorTest::testCoreform(), MeshInputTest::testExodusIGASidesets(), MeshTriangulationTest::testFoundCenters(), PointLocatorTest::testLocator(), BoundaryInfoTest::testMesh(), PointLocatorTest::testPlanar(), MeshTriangulationTest::testPoly2TriRefinementBase(), SystemsTest::testProjectCubeWithMeshFunction(), BoundaryInfoTest::testRenumber(), CheckpointIOTest::testSplitter(), MeshInputTest::testTetgenIO(), MeshTriangulationTest::testTriangulatorInterp(), MeshTriangulationTest::testTriangulatorMeshedHoles(), MeshTriangulationTest::testTriangulatorRoundHole(), libMesh::MeshTools::total_weight(), libMesh::RBConstruction::train_reduced_basis_with_POD(), libMesh::MeshFunctionSolutionTransfer::transfer(), libMesh::MeshfreeSolutionTransfer::transfer(), libMesh::Poly2TriTriangulator::triangulate(), libMesh::TransientRBConstruction::truth_assembly(), libMesh::RBConstruction::truth_assembly(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::TransientRBConstruction::update_RB_initial_condition_all_N(), libMesh::TransientRBConstruction::update_RB_system_matrices(), libMesh::RBConstruction::update_RB_system_matrices(), libMesh::TransientRBConstruction::update_residual_terms(), libMesh::RBConstruction::update_residual_terms(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::VTKIO::write_nodal_data(), libMesh::RBEIMEvaluation::write_out_interior_basis_functions(), libMesh::RBEIMEvaluation::write_out_node_basis_functions(), libMesh::RBEIMEvaluation::write_out_side_basis_functions(), libMesh::RBEvaluation::write_out_vectors(), libMesh::TransientRBConstruction::write_riesz_representors_to_files(), libMesh::RBConstruction::write_riesz_representors_to_files(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::RBDataSerialization::RBEvaluationSerialization::write_to_file(), libMesh::RBDataSerialization::TransientRBEvaluationSerialization::write_to_file(), libMesh::RBDataSerialization::RBEIMEvaluationSerialization::write_to_file(), and libMesh::RBDataSerialization::RBSCMEvaluationSerialization::write_to_file().

98  { return _communicator; }
const Parallel::Communicator & _communicator

◆ create_submatrix()

template<typename T>
virtual void libMesh::SparseMatrix< T >::create_submatrix ( SparseMatrix< T > &  submatrix,
const std::vector< numeric_index_type > &  rows,
const std::vector< numeric_index_type > &  cols 
) const
inlinevirtualinherited

This function creates a matrix called "submatrix" which is defined by the row and column indices given in the "rows" and "cols" entries.

Currently this operation is only defined for the PetscMatrix type. Note: The rows and cols vectors need to be sorted; Use the nosort version below if rows and cols vectors are not sorted; The rows and cols only contain indices that are owned by this processor.

Definition at line 469 of file sparse_matrix.h.

Referenced by libMesh::libmesh_petsc_DMCreateInterpolation(), and libMesh::CondensedEigenSystem::solve().

472  {
473  this->_get_submatrix(submatrix,
474  rows,
475  cols,
476  false); // false means DO NOT REUSE submatrix
477  }
virtual void _get_submatrix(SparseMatrix< T > &, const std::vector< numeric_index_type > &, const std::vector< numeric_index_type > &, const bool) const
Protected implementation of the create_submatrix and reinit_submatrix routines.

◆ create_submatrix_nosort()

template<typename T>
void libMesh::PetscMatrix< T >::create_submatrix_nosort ( SparseMatrix< T > &  submatrix,
const std::vector< numeric_index_type > &  rows,
const std::vector< numeric_index_type > &  cols 
) const
overridevirtual

Similar to the create_submatrix function, this function creates a submatrix which is defined by the indices given in the rows and cols vectors.

Note: Both rows and cols can be unsorted; Use the above function for better efficiency if your indices are sorted; rows and cols can contain indices that are owned by other processors.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 936 of file petsc_matrix.C.

939 {
940  if (!this->closed())
941  {
942  libmesh_deprecated();
943  libmesh_warning("The matrix must be assembled before calling PetscMatrix::create_submatrix_nosort().\n"
944  "Please update your code, as this warning will become an error in a future release.");
945  const_cast<PetscMatrix<T> *>(this)->close();
946  }
947 
948  // Make sure the SparseMatrix passed in is really a PetscMatrix
949  PetscMatrix<T> * petsc_submatrix = cast_ptr<PetscMatrix<T> *>(&submatrix);
950 
951  PetscErrorCode ierr=0;
952 
953  ierr=MatZeroEntries(petsc_submatrix->_mat);
954  LIBMESH_CHKERR(ierr);
955 
956  PetscInt pc_ncols = 0;
957  const PetscInt * pc_cols;
958  const PetscScalar * pc_vals;
959 
960  // // data for creating the submatrix
961  std::vector<PetscInt> sub_cols;
962  std::vector<PetscScalar> sub_vals;
963 
964  for (auto i : index_range(rows))
965  {
966  PetscInt sub_rid[] = {static_cast<PetscInt>(i)};
967  PetscInt rid = static_cast<PetscInt>(rows[i]);
968  // only get value from local rows, and set values to the corresponding columns in the submatrix
969  if (rows[i]>= this->row_start() && rows[i]< this->row_stop())
970  {
971  // get one row of data from the original matrix
972  ierr = MatGetRow(_mat, rid, &pc_ncols, &pc_cols, &pc_vals);
973  LIBMESH_CHKERR(ierr);
974  // extract data from certain cols, save the indices and entries sub_cols and sub_vals
975  for (auto j : index_range(cols))
976  {
977  for (unsigned int idx = 0; idx< static_cast<unsigned int>(pc_ncols); idx++)
978  {
979  if (pc_cols[idx] == static_cast<PetscInt>(cols[j]))
980  {
981  sub_cols.push_back(static_cast<PetscInt>(j));
982  sub_vals.push_back(pc_vals[idx]);
983  }
984  }
985  }
986  // set values
987  ierr = MatSetValues(petsc_submatrix->_mat,
988  1,
989  sub_rid,
990  static_cast<PetscInt>(sub_vals.size()),
991  sub_cols.data(),
992  sub_vals.data(),
993  INSERT_VALUES);
994  LIBMESH_CHKERR(ierr);
995  ierr = MatRestoreRow(_mat, rid, &pc_ncols, &pc_cols, &pc_vals);
996  LIBMESH_CHKERR(ierr);
997  // clear data for this row
998  sub_cols.clear();
999  sub_vals.clear();
1000  }
1001  }
1002  MatAssemblyBeginEnd(petsc_submatrix->comm(), petsc_submatrix->_mat, MAT_FINAL_ASSEMBLY);
1003  // Specify that the new submatrix is initialized and close it.
1004  petsc_submatrix->_is_initialized = true;
1005  petsc_submatrix->close();
1006 }
const Parallel::Communicator & comm() const
bool _is_initialized
Flag indicating whether or not the matrix has been initialized.
virtual numeric_index_type row_stop() const override
virtual bool closed() const override
virtual numeric_index_type row_start() const override
This class provides a nice interface to the PETSc C-based data structures for parallel, sparse matrices.
Definition: petsc_matrix.h:93
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367
virtual void close() override
Calls the SparseMatrix&#39;s internal assembly routines, ensuring that the values are consistent across p...
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111
unsigned int idx(const ElemType type, const unsigned int nx, const unsigned int i, const unsigned int j)
A useful inline function which replaces the macros used previously.

◆ disable_print_counter_info()

void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited

Definition at line 100 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

101 {
102  _enable_print_counter = false;
103  return;
104 }
static bool _enable_print_counter
Flag to control whether reference count information is printed when print_info is called...

◆ enable_print_counter_info()

void libMesh::ReferenceCounter::enable_print_counter_info ( )
staticinherited

Methods to enable/disable the reference counter output from print_info()

Definition at line 94 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

95 {
96  _enable_print_counter = true;
97  return;
98 }
static bool _enable_print_counter
Flag to control whether reference count information is printed when print_info is called...

◆ flush()

template<typename T >
void libMesh::PetscMatrix< T >::flush ( )
overridevirtual

For PETSc matrix , this function is similar to close but without shrinking memory.

This is useful when we want to switch between ADD_VALUES and INSERT_VALUES. close should be called before using the matrix.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 1071 of file petsc_matrix.C.

1072 {
1073  semiparallel_only();
1074 
1075  MatAssemblyBeginEnd(this->comm(), _mat, MAT_FLUSH_ASSEMBLY);
1076 }
const Parallel::Communicator & comm() const
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ get_diagonal()

template<typename T>
void libMesh::PetscMatrix< T >::get_diagonal ( NumericVector< T > &  dest) const
overridevirtual

Copies the diagonal part of the matrix into dest.

Implements libMesh::SparseMatrix< T >.

Definition at line 1010 of file petsc_matrix.C.

1011 {
1012  // Make sure the NumericVector passed in is really a PetscVector
1013  PetscVector<T> & petsc_dest = cast_ref<PetscVector<T> &>(dest);
1014 
1015  // Needs a const_cast since PETSc does not work with const.
1016  PetscErrorCode ierr =
1017  MatGetDiagonal(const_cast<PetscMatrix<T> *>(this)->mat(),petsc_dest.vec()); LIBMESH_CHKERR(ierr);
1018 }
This class provides a nice interface to PETSc&#39;s Vec object.
Definition: petsc_vector.h:73
This class provides a nice interface to the PETSc C-based data structures for parallel, sparse matrices.
Definition: petsc_matrix.h:93

◆ get_info()

std::string libMesh::ReferenceCounter::get_info ( )
staticinherited

Gets a string containing the reference information.

Definition at line 47 of file reference_counter.C.

References libMesh::ReferenceCounter::_counts, and libMesh::Quality::name().

Referenced by libMesh::ReferenceCounter::print_info().

48 {
49 #if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
50 
51  std::ostringstream oss;
52 
53  oss << '\n'
54  << " ---------------------------------------------------------------------------- \n"
55  << "| Reference count information |\n"
56  << " ---------------------------------------------------------------------------- \n";
57 
58  for (const auto & [name, cd] : _counts)
59  oss << "| " << name << " reference count information:\n"
60  << "| Creations: " << cd.first << '\n'
61  << "| Destructions: " << cd.second << '\n';
62 
63  oss << " ---------------------------------------------------------------------------- \n";
64 
65  return oss.str();
66 
67 #else
68 
69  return "";
70 
71 #endif
72 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
static Counts _counts
Actually holds the data.

◆ get_local_size()

template<typename T >
void libMesh::PetscMatrix< T >::get_local_size ( numeric_index_type m,
numeric_index_type n 
) const

Get the number of rows and columns owned by this process.

Parameters
iRow size
jColumn size

Definition at line 1135 of file petsc_matrix.C.

1137 {
1138  libmesh_assert (this->initialized());
1139 
1140  PetscInt petsc_m = 0, petsc_n = 0;
1141 
1142  auto ierr = MatGetLocalSize (_mat, &petsc_m, &petsc_n);
1143  LIBMESH_CHKERR(ierr);
1144 
1145  m = static_cast<numeric_index_type>(petsc_m);
1146  n = static_cast<numeric_index_type>(petsc_n);
1147 }
virtual bool initialized() const
virtual numeric_index_type n() const override
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
virtual numeric_index_type m() const override
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ get_row()

template<typename T>
void libMesh::PetscMatrix< T >::get_row ( numeric_index_type  i,
std::vector< numeric_index_type > &  indices,
std::vector< T > &  values 
) const
overridevirtual

Get a row from the matrix.

Parameters
iThe matrix row to get
indicesA container that will be filled with the column indices corresponding to (possibly) non-zero values
valuesA container holding the column values

Implements libMesh::SparseMatrix< T >.

Definition at line 1440 of file petsc_matrix.C.

1443 {
1444  libmesh_assert (this->initialized());
1445 
1446  const PetscScalar * petsc_row;
1447  const PetscInt * petsc_cols;
1448 
1449  PetscErrorCode ierr=0;
1450  PetscInt
1451  ncols=0,
1452  i_val = static_cast<PetscInt>(i_in);
1453 
1454  // the matrix needs to be closed for this to work
1455  // this->close();
1456  // but closing it is a semiparallel operation; we want operator()
1457  // to run on one processor.
1458  libmesh_assert(this->closed());
1459 
1460  // PETSc makes no effort at being thread safe. Helgrind complains about
1461  // possible data races even just in PetscFunctionBegin (due to things
1462  // like stack counter incrementing). Perhaps we could ignore
1463  // this, but there are legitimate data races for Mat data members like
1464  // mat->getrowactive between MatGetRow and MatRestoreRow. Moreover,
1465  // there could be a write into mat->rowvalues during MatGetRow from
1466  // one thread while we are attempting to read from mat->rowvalues
1467  // (through petsc_cols) during data copy in another thread. So
1468  // the safe thing to do is to lock the whole method
1469 
1470 #ifdef LIBMESH_HAVE_CXX11_THREAD
1471  std::lock_guard<std::mutex>
1472 #else
1473  Threads::spin_mutex::scoped_lock
1474 #endif
1475  lock(_petsc_matrix_mutex);
1476 
1477  ierr = MatGetRow(_mat, i_val, &ncols, &petsc_cols, &petsc_row);
1478  LIBMESH_CHKERR(ierr);
1479 
1480  // Copy the data
1481  indices.resize(static_cast<std::size_t>(ncols));
1482  values.resize(static_cast<std::size_t>(ncols));
1483 
1484  for (auto i : index_range(indices))
1485  {
1486  indices[i] = static_cast<numeric_index_type>(petsc_cols[i]);
1487  values[i] = static_cast<T>(petsc_row[i]);
1488  }
1489 
1490  ierr = MatRestoreRow(_mat, i_val,
1491  &ncols, &petsc_cols, &petsc_row);
1492  LIBMESH_CHKERR(ierr);
1493 }
virtual bool initialized() const
std::mutex _petsc_matrix_mutex
Definition: petsc_matrix.h:378
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
virtual bool closed() const override
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111

◆ get_transpose()

template<typename T>
void libMesh::PetscMatrix< T >::get_transpose ( SparseMatrix< T > &  dest) const
overridevirtual

Copies the transpose of the matrix into dest, which may be *this.

Implements libMesh::SparseMatrix< T >.

Definition at line 1023 of file petsc_matrix.C.

1024 {
1025  // Make sure the SparseMatrix passed in is really a PetscMatrix
1026  PetscMatrix<T> & petsc_dest = cast_ref<PetscMatrix<T> &>(dest);
1027 
1028  // If we aren't reusing the matrix then need to clear dest,
1029  // otherwise we get a memory leak
1030  if (&petsc_dest != this)
1031  dest.clear();
1032 
1033  PetscErrorCode ierr;
1034  if (&petsc_dest == this)
1035  // The MAT_REUSE_MATRIX flag was replaced by MAT_INPLACE_MATRIX
1036  // in PETSc 3.7.0
1037 #if PETSC_VERSION_LESS_THAN(3,7,0)
1038  ierr = MatTranspose(_mat,MAT_REUSE_MATRIX,&petsc_dest._mat);
1039 #else
1040  ierr = MatTranspose(_mat, MAT_INPLACE_MATRIX, &petsc_dest._mat);
1041 #endif
1042  else
1043  ierr = MatTranspose(_mat,MAT_INITIAL_MATRIX,&petsc_dest._mat);
1044  LIBMESH_CHKERR(ierr);
1045 
1046  // Specify that the transposed matrix is initialized and close it.
1047  petsc_dest._is_initialized = true;
1048  petsc_dest.close();
1049 }
virtual void clear()=0
Restores the SparseMatrix<T> to a pristine state.
bool _is_initialized
Flag indicating whether or not the matrix has been initialized.
This class provides a nice interface to the PETSc C-based data structures for parallel, sparse matrices.
Definition: petsc_matrix.h:93
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367
virtual void close() override
Calls the SparseMatrix&#39;s internal assembly routines, ensuring that the values are consistent across p...

◆ increment_constructor_count()

void libMesh::ReferenceCounter::increment_constructor_count ( const std::string &  name)
inlineprotectednoexceptinherited

Increments the construction counter.

Should be called in the constructor of any derived class that will be reference counted.

Definition at line 183 of file reference_counter.h.

References libMesh::err, libMesh::BasicOStreamProxy< charT, traits >::get(), libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCountedObject< RBParametrized >::ReferenceCountedObject().

184 {
185  libmesh_try
186  {
187  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
188  std::pair<unsigned int, unsigned int> & p = _counts[name];
189  p.first++;
190  }
191  libmesh_catch (...)
192  {
193  auto stream = libMesh::err.get();
194  stream->exceptions(stream->goodbit); // stream must not throw
195  libMesh::err << "Encountered unrecoverable error while calling "
196  << "ReferenceCounter::increment_constructor_count() "
197  << "for a(n) " << name << " object." << std::endl;
198  std::terminate();
199  }
200 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
OStreamProxy err
static Counts _counts
Actually holds the data.
streamT * get()
Rather than implement every ostream/ios/ios_base function, we&#39;ll be lazy and make esoteric uses go th...
spin_mutex spin_mtx
A convenient spin mutex object which can be used for obtaining locks.
Definition: threads.C:30

◆ increment_destructor_count()

void libMesh::ReferenceCounter::increment_destructor_count ( const std::string &  name)
inlineprotectednoexceptinherited

Increments the destruction counter.

Should be called in the destructor of any derived class that will be reference counted.

Definition at line 207 of file reference_counter.h.

References libMesh::err, libMesh::BasicOStreamProxy< charT, traits >::get(), libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCountedObject< RBParametrized >::~ReferenceCountedObject().

208 {
209  libmesh_try
210  {
211  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
212  std::pair<unsigned int, unsigned int> & p = _counts[name];
213  p.second++;
214  }
215  libmesh_catch (...)
216  {
217  auto stream = libMesh::err.get();
218  stream->exceptions(stream->goodbit); // stream must not throw
219  libMesh::err << "Encountered unrecoverable error while calling "
220  << "ReferenceCounter::increment_destructor_count() "
221  << "for a(n) " << name << " object." << std::endl;
222  std::terminate();
223  }
224 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
OStreamProxy err
static Counts _counts
Actually holds the data.
streamT * get()
Rather than implement every ostream/ios/ios_base function, we&#39;ll be lazy and make esoteric uses go th...
spin_mutex spin_mtx
A convenient spin mutex object which can be used for obtaining locks.
Definition: threads.C:30

◆ init() [1/3]

template<typename T >
void libMesh::PetscMatrix< T >::init ( const numeric_index_type  m,
const numeric_index_type  n,
const numeric_index_type  m_l,
const numeric_index_type  n_l,
const numeric_index_type  nnz = 30,
const numeric_index_type  noz = 10,
const numeric_index_type  blocksize = 1 
)
overridevirtual

Initialize SparseMatrix with the specified sizes.

Parameters
mThe global number of rows.
nThe global number of columns.
m_lThe local number of rows.
n_lThe local number of columns.
nnzThe number of on-diagonal nonzeros per row (defaults to 30).
nozThe number of off-diagonal nonzeros per row (defaults to 10).
blocksizeOptional value indicating dense coupled blocks for systems with multiple variables all of the same type.

Implements libMesh::SparseMatrix< T >.

Definition at line 121 of file petsc_matrix.C.

128 {
129  // So compilers don't warn when !LIBMESH_ENABLE_BLOCKED_STORAGE
130  libmesh_ignore(blocksize_in);
131 
132  // Clear initialized matrices
133  if (this->initialized())
134  this->clear();
135 
136  this->_is_initialized = true;
137 
138 
139  PetscErrorCode ierr = 0;
140  PetscInt m_global = static_cast<PetscInt>(m_in);
141  PetscInt n_global = static_cast<PetscInt>(n_in);
142  PetscInt m_local = static_cast<PetscInt>(m_l);
143  PetscInt n_local = static_cast<PetscInt>(n_l);
144  PetscInt n_nz = static_cast<PetscInt>(nnz);
145  PetscInt n_oz = static_cast<PetscInt>(noz);
146 
147  ierr = MatCreate(this->comm().get(), &_mat);
148  LIBMESH_CHKERR(ierr);
149  ierr = MatSetSizes(_mat, m_local, n_local, m_global, n_global);
150  LIBMESH_CHKERR(ierr);
151  PetscInt blocksize = static_cast<PetscInt>(blocksize_in);
152  ierr = MatSetBlockSize(_mat,blocksize);
153  LIBMESH_CHKERR(ierr);
154 
155 #ifdef LIBMESH_ENABLE_BLOCKED_STORAGE
156  if (blocksize > 1)
157  {
158  // specified blocksize, bs>1.
159  // double check sizes.
160  libmesh_assert_equal_to (m_local % blocksize, 0);
161  libmesh_assert_equal_to (n_local % blocksize, 0);
162  libmesh_assert_equal_to (m_global % blocksize, 0);
163  libmesh_assert_equal_to (n_global % blocksize, 0);
164  libmesh_assert_equal_to (n_nz % blocksize, 0);
165  libmesh_assert_equal_to (n_oz % blocksize, 0);
166 
167  ierr = MatSetType(_mat, MATBAIJ); // Automatically chooses seqbaij or mpibaij
168  LIBMESH_CHKERR(ierr);
169 
170  // MatSetFromOptions needs to happen before Preallocation routines
171  // since MatSetFromOptions can change matrix type and remove incompatible
172  // preallocation
173  ierr = MatSetOptionsPrefix(_mat, "");
174  LIBMESH_CHKERR(ierr);
175  ierr = MatSetFromOptions(_mat);
176  LIBMESH_CHKERR(ierr);
177  ierr = MatSeqBAIJSetPreallocation(_mat, blocksize, n_nz/blocksize, NULL);
178  LIBMESH_CHKERR(ierr);
179  ierr = MatMPIBAIJSetPreallocation(_mat, blocksize,
180  n_nz/blocksize, NULL,
181  n_oz/blocksize, NULL);
182  LIBMESH_CHKERR(ierr);
183  }
184  else
185 #endif
186  {
187  switch (_mat_type) {
188  case AIJ:
189  ierr = MatSetType(_mat, MATAIJ); // Automatically chooses seqaij or mpiaij
190  LIBMESH_CHKERR(ierr);
191 
192  // MatSetFromOptions needs to happen before Preallocation routines
193  // since MatSetFromOptions can change matrix type and remove incompatible
194  // preallocation
195  ierr = MatSetOptionsPrefix(_mat, "");
196  LIBMESH_CHKERR(ierr);
197  ierr = MatSetFromOptions(_mat);
198  LIBMESH_CHKERR(ierr);
199  ierr = MatSeqAIJSetPreallocation(_mat, n_nz, NULL);
200  LIBMESH_CHKERR(ierr);
201  ierr = MatMPIAIJSetPreallocation(_mat, n_nz, NULL, n_oz, NULL);
202  LIBMESH_CHKERR(ierr);
203  break;
204 
205  case HYPRE:
206 #if !PETSC_VERSION_LESS_THAN(3,9,4) && LIBMESH_HAVE_PETSC_HYPRE
207  ierr = MatSetType(_mat, MATHYPRE);
208 
209  // MatSetFromOptions needs to happen before Preallocation routines
210  // since MatSetFromOptions can change matrix type and remove incompatible
211  // preallocation
212  LIBMESH_CHKERR(ierr);
213  ierr = MatSetOptionsPrefix(_mat, "");
214  LIBMESH_CHKERR(ierr);
215  ierr = MatSetFromOptions(_mat);
216  LIBMESH_CHKERR(ierr);
217  ierr = MatHYPRESetPreallocation(_mat, n_nz, NULL, n_oz, NULL);
218  LIBMESH_CHKERR(ierr);
219 #else
220  libmesh_error_msg("PETSc 3.9.4 or higher with hypre is required for MatHypre");
221 #endif
222  break;
223 
224  default: libmesh_error_msg("Unsupported petsc matrix type");
225  }
226  }
227 
228  // Make it an error for PETSc to allocate new nonzero entries during assembly
229  ierr = MatSetOption(_mat, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
230  LIBMESH_CHKERR(ierr);
231 
232  this->zero ();
233 }
virtual bool initialized() const
virtual void zero() override
Set all entries to 0.
Definition: petsc_matrix.C:417
const Parallel::Communicator & comm() const
virtual void clear() noexcept override
clear() is called from the destructor, so it should not throw.
Definition: petsc_matrix.C:503
void libmesh_ignore(const Args &...)
bool _is_initialized
Flag indicating whether or not the matrix has been initialized.
PetscMatrixType _mat_type
Definition: petsc_matrix.h:375
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ init() [2/3]

template<typename T >
void libMesh::PetscMatrix< T >::init ( const numeric_index_type  m,
const numeric_index_type  n,
const numeric_index_type  m_l,
const numeric_index_type  n_l,
const std::vector< numeric_index_type > &  n_nz,
const std::vector< numeric_index_type > &  n_oz,
const numeric_index_type  blocksize = 1 
)

Initialize a PETSc matrix.

Parameters
mThe global number of rows.
nThe global number of columns.
m_lThe local number of rows.
n_lThe local number of columns.
n_nzarray containing the number of nonzeros in each row of the DIAGONAL portion of the local submatrix.
n_ozArray containing the number of nonzeros in each row of the OFF-DIAGONAL portion of the local submatrix.
blocksizeOptional value indicating dense coupled blocks for systems with multiple variables all of the same type.

Definition at line 237 of file petsc_matrix.C.

244 {
245  // So compilers don't warn when !LIBMESH_ENABLE_BLOCKED_STORAGE
246  libmesh_ignore(blocksize_in);
247 
248  PetscInt blocksize = static_cast<PetscInt>(blocksize_in);
249 
250  // Clear initialized matrices
251  if (this->initialized())
252  this->clear();
253 
254  this->_is_initialized = true;
255 
256  // Make sure the sparsity pattern isn't empty unless the matrix is 0x0
257  libmesh_assert_equal_to (n_nz.size(), m_l);
258  libmesh_assert_equal_to (n_oz.size(), m_l);
259 
260  PetscErrorCode ierr = 0;
261  PetscInt m_global = static_cast<PetscInt>(m_in);
262  PetscInt n_global = static_cast<PetscInt>(n_in);
263  PetscInt m_local = static_cast<PetscInt>(m_l);
264  PetscInt n_local = static_cast<PetscInt>(n_l);
265 
266  ierr = MatCreate(this->comm().get(), &_mat);
267  LIBMESH_CHKERR(ierr);
268  ierr = MatSetSizes(_mat, m_local, n_local, m_global, n_global);
269  LIBMESH_CHKERR(ierr);
270  ierr = MatSetBlockSize(_mat,blocksize);
271  LIBMESH_CHKERR(ierr);
272 
273 #ifdef LIBMESH_ENABLE_BLOCKED_STORAGE
274  if (blocksize > 1)
275  {
276  // specified blocksize, bs>1.
277  // double check sizes.
278  libmesh_assert_equal_to (m_local % blocksize, 0);
279  libmesh_assert_equal_to (n_local % blocksize, 0);
280  libmesh_assert_equal_to (m_global % blocksize, 0);
281  libmesh_assert_equal_to (n_global % blocksize, 0);
282 
283  ierr = MatSetType(_mat, MATBAIJ); // Automatically chooses seqbaij or mpibaij
284  LIBMESH_CHKERR(ierr);
285 
286  // MatSetFromOptions needs to happen before Preallocation routines
287  // since MatSetFromOptions can change matrix type and remove incompatible
288  // preallocation
289  LIBMESH_CHKERR(ierr);
290  ierr = MatSetOptionsPrefix(_mat, "");
291  LIBMESH_CHKERR(ierr);
292  ierr = MatSetFromOptions(_mat);
293  LIBMESH_CHKERR(ierr);
294  // transform the per-entry n_nz and n_oz arrays into their block counterparts.
295  std::vector<numeric_index_type> b_n_nz, b_n_oz;
296 
297  transform_preallocation_arrays (blocksize,
298  n_nz, n_oz,
299  b_n_nz, b_n_oz);
300 
301  ierr = MatSeqBAIJSetPreallocation (_mat,
302  blocksize,
303  0,
304  numeric_petsc_cast(b_n_nz.empty() ? nullptr : b_n_nz.data()));
305  LIBMESH_CHKERR(ierr);
306 
307  ierr = MatMPIBAIJSetPreallocation (_mat,
308  blocksize,
309  0,
310  numeric_petsc_cast(b_n_nz.empty() ? nullptr : b_n_nz.data()),
311  0,
312  numeric_petsc_cast(b_n_oz.empty() ? nullptr : b_n_oz.data()));
313  LIBMESH_CHKERR(ierr);
314  }
315  else
316 #endif
317  {
318  switch (_mat_type) {
319  case AIJ:
320  ierr = MatSetType(_mat, MATAIJ); // Automatically chooses seqaij or mpiaij
321  LIBMESH_CHKERR(ierr);
322 
323  // MatSetFromOptions needs to happen before Preallocation routines
324  // since MatSetFromOptions can change matrix type and remove incompatible
325  // preallocation
326  LIBMESH_CHKERR(ierr);
327  ierr = MatSetOptionsPrefix(_mat, "");
328  LIBMESH_CHKERR(ierr);
329  ierr = MatSetFromOptions(_mat);
330  LIBMESH_CHKERR(ierr);
331  ierr = MatSeqAIJSetPreallocation (_mat,
332  0,
333  numeric_petsc_cast(n_nz.empty() ? nullptr : n_nz.data()));
334  LIBMESH_CHKERR(ierr);
335  ierr = MatMPIAIJSetPreallocation (_mat,
336  0,
337  numeric_petsc_cast(n_nz.empty() ? nullptr : n_nz.data()),
338  0,
339  numeric_petsc_cast(n_oz.empty() ? nullptr : n_oz.data()));
340  break;
341 
342  case HYPRE:
343 #if !PETSC_VERSION_LESS_THAN(3,9,4) && LIBMESH_HAVE_PETSC_HYPRE
344  ierr = MatSetType(_mat, MATHYPRE);
345  LIBMESH_CHKERR(ierr);
346 
347  // MatSetFromOptions needs to happen before Preallocation routines
348  // since MatSetFromOptions can change matrix type and remove incompatible
349  // preallocation
350  LIBMESH_CHKERR(ierr);
351  ierr = MatSetOptionsPrefix(_mat, "");
352  LIBMESH_CHKERR(ierr);
353  ierr = MatSetFromOptions(_mat);
354  LIBMESH_CHKERR(ierr);
355  ierr = MatHYPRESetPreallocation (_mat,
356  0,
357  numeric_petsc_cast(n_nz.empty() ? nullptr : n_nz.data()),
358  0,
359  numeric_petsc_cast(n_oz.empty() ? nullptr : n_oz.data()));
360  LIBMESH_CHKERR(ierr);
361 #else
362  libmesh_error_msg("PETSc 3.9.4 or higher with hypre is required for MatHypre");
363 #endif
364  break;
365 
366  default: libmesh_error_msg("Unsupported petsc matrix type");
367  }
368 
369  }
370 
371  // Make it an error for PETSc to allocate new nonzero entries during assembly
372  ierr = MatSetOption(_mat, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
373  LIBMESH_CHKERR(ierr);
374  this->zero();
375 }
virtual bool initialized() const
virtual void zero() override
Set all entries to 0.
Definition: petsc_matrix.C:417
const Parallel::Communicator & comm() const
PetscInt * numeric_petsc_cast(const numeric_index_type *p)
virtual void clear() noexcept override
clear() is called from the destructor, so it should not throw.
Definition: petsc_matrix.C:503
void libmesh_ignore(const Args &...)
bool _is_initialized
Flag indicating whether or not the matrix has been initialized.
PetscMatrixType _mat_type
Definition: petsc_matrix.h:375
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ init() [3/3]

template<typename T >
void libMesh::PetscMatrix< T >::init ( ParallelType  type = PARALLEL)
overridevirtual

Initialize this matrix using the sparsity structure computed by dof_map.

Parameters
typeThe serial/parallel/ghosted type of the matrix

Implements libMesh::SparseMatrix< T >.

Definition at line 379 of file petsc_matrix.C.

380 {
381  libmesh_assert(this->_dof_map);
382 
383  const numeric_index_type m_in = this->_dof_map->n_dofs();
384  const numeric_index_type m_l = this->_dof_map->n_dofs_on_processor(this->processor_id());
385 
386  const std::vector<numeric_index_type> & n_nz = this->_sp->get_n_nz();
387  const std::vector<numeric_index_type> & n_oz = this->_sp->get_n_oz();
388 
389  PetscInt blocksize = static_cast<PetscInt>(this->_dof_map->block_size());
390 
391  this->init(m_in, m_in, m_l, m_l, n_nz, n_oz, blocksize);
392 }
dof_id_type n_dofs_on_processor(const processor_id_type proc) const
Definition: dof_map.h:675
virtual void init(const numeric_index_type m, const numeric_index_type n, const numeric_index_type m_l, const numeric_index_type n_l, const numeric_index_type nnz=30, const numeric_index_type noz=10, const numeric_index_type blocksize=1) override
Initialize SparseMatrix with the specified sizes.
Definition: petsc_matrix.C:121
const std::vector< dof_id_type > & get_n_oz() const
The number of off-processor nonzeros in my portion of the global matrix.
unsigned int block_size() const
Definition: dof_map.h:651
SparsityPattern::Build const * _sp
The sparsity pattern associated with this object.
dof_id_type n_dofs() const
Definition: dof_map.h:659
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
const std::vector< dof_id_type > & get_n_nz() const
The number of on-processor nonzeros in my portion of the global matrix.
DofMap const * _dof_map
The DofMap object associated with this object.
processor_id_type processor_id() const

◆ initialized()

template<typename T>
virtual bool libMesh::SparseMatrix< T >::initialized ( ) const
inlinevirtualinherited
Returns
true if the matrix has been initialized, false otherwise.

Definition at line 109 of file sparse_matrix.h.

Referenced by libMesh::PetscMatrix< libMesh::Number >::_get_submatrix(), libMesh::ImplicitSystem::assemble(), and libMesh::System::init_matrices().

109 { return _is_initialized; }
bool _is_initialized
Flag indicating whether or not the matrix has been initialized.

◆ l1_norm()

template<typename T >
Real libMesh::PetscMatrix< T >::l1_norm ( ) const
overridevirtual
Returns
The \( \ell_1 \)-norm of the matrix, that is the max column sum: \( |M|_1 = \max_{j} \sum_{i} |M_{ij}| \)

This is the natural matrix norm that is compatible with the \( \ell_1 \)-norm for vectors, i.e. \( |Mv|_1 \leq |M|_1 |v|_1 \). (cf. Haemmerlin-Hoffmann : Numerische Mathematik)

Implements libMesh::SparseMatrix< T >.

Definition at line 522 of file petsc_matrix.C.

523 {
524  libmesh_assert (this->initialized());
525 
526  semiparallel_only();
527 
528  PetscErrorCode ierr=0;
529  PetscReal petsc_value;
530  Real value;
531 
532  libmesh_assert (this->closed());
533 
534  ierr = MatNorm(_mat, NORM_1, &petsc_value);
535  LIBMESH_CHKERR(ierr);
536 
537  value = static_cast<Real>(petsc_value);
538 
539  return value;
540 }
virtual bool initialized() const
libmesh_assert(ctx)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual bool closed() const override
static const bool value
Definition: xdr_io.C:54
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ linfty_norm()

template<typename T >
Real libMesh::PetscMatrix< T >::linfty_norm ( ) const
overridevirtual
Returns
The \( \ell_{\infty} \)-norm of the matrix, that is the max row sum:

\( |M|_{\infty} = \max_{i} \sum_{j} |M_{ij}| \)

This is the natural matrix norm that is compatible to the \( \ell_{\infty} \)-norm of vectors, i.e. \( |Mv|_{\infty} \leq |M|_{\infty} |v|_{\infty} \). (cf. Haemmerlin-Hoffmann : Numerische Mathematik)

Implements libMesh::SparseMatrix< T >.

Definition at line 545 of file petsc_matrix.C.

546 {
547  libmesh_assert (this->initialized());
548 
549  semiparallel_only();
550 
551  PetscErrorCode ierr=0;
552  PetscReal petsc_value;
553  Real value;
554 
555  libmesh_assert (this->closed());
556 
557  ierr = MatNorm(_mat, NORM_INFINITY, &petsc_value);
558  LIBMESH_CHKERR(ierr);
559 
560  value = static_cast<Real>(petsc_value);
561 
562  return value;
563 }
virtual bool initialized() const
libmesh_assert(ctx)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual bool closed() const override
static const bool value
Definition: xdr_io.C:54
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ local_m()

template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::local_m ( ) const
finalvirtual

Get the number of rows owned by this process.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 1095 of file petsc_matrix.C.

1096 {
1097  libmesh_assert (this->initialized());
1098 
1099  PetscInt m = 0;
1100 
1101  auto ierr = MatGetLocalSize (_mat, &m, NULL);
1102  LIBMESH_CHKERR(ierr);
1103 
1104  return static_cast<numeric_index_type>(m);
1105 }
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
virtual numeric_index_type m() const override
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ local_n()

template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::local_n ( ) const
finalvirtual

Get the number of columns owned by this process.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 1122 of file petsc_matrix.C.

1123 {
1124  libmesh_assert (this->initialized());
1125 
1126  PetscInt n = 0;
1127 
1128  auto ierr = MatGetLocalSize (_mat, NULL, &n);
1129  LIBMESH_CHKERR(ierr);
1130 
1131  return static_cast<numeric_index_type>(n);
1132 }
virtual bool initialized() const
virtual numeric_index_type n() const override
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ m()

template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::m ( ) const
overridevirtual
Returns
The row-dimension of the matrix.

Implements libMesh::SparseMatrix< T >.

Definition at line 1081 of file petsc_matrix.C.

1082 {
1083  libmesh_assert (this->initialized());
1084 
1085  PetscInt petsc_m=0, petsc_n=0;
1086  PetscErrorCode ierr=0;
1087 
1088  ierr = MatGetSize (_mat, &petsc_m, &petsc_n);
1089  LIBMESH_CHKERR(ierr);
1090 
1091  return static_cast<numeric_index_type>(petsc_m);
1092 }
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ mat()

template<typename T>
Mat libMesh::PetscMatrix< T >::mat ( )
inline

◆ matrix_matrix_mult()

template<typename T>
void libMesh::PetscMatrix< T >::matrix_matrix_mult ( SparseMatrix< T > &  X,
SparseMatrix< T > &  Y,
bool  reuse = false 
)
overridevirtual

Compute Y = A*X for matrix X.

Set reuse = true if this->_mat and X have the same nonzero pattern as before, and Y is obtained from a previous call to this function with reuse = false

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 1292 of file petsc_matrix.C.

1293 {
1294  libmesh_assert (this->initialized());
1295 
1296  // sanity check
1297  // we do not check the Y_out size here as we will initialize & close it at the end
1298  libmesh_assert_equal_to (this->n(), X_in.m());
1299 
1300  const PetscMatrix<T> * X = cast_ptr<const PetscMatrix<T> *> (&X_in);
1301  PetscMatrix<T> * Y = cast_ptr<PetscMatrix<T> *> (&Y_out);
1302 
1303  PetscErrorCode ierr = 0;
1304 
1305  // the matrix from which we copy the values has to be assembled/closed
1306  libmesh_assert(X->closed());
1307 
1308  semiparallel_only();
1309 
1310  if (reuse)
1311  {
1312  ierr = MatMatMult(_mat, X->_mat, MAT_REUSE_MATRIX, PETSC_DEFAULT, &Y->_mat);
1313  LIBMESH_CHKERR(ierr);
1314  }
1315  else
1316  {
1317  Y->clear();
1318  ierr = MatMatMult(_mat, X->_mat, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &Y->_mat);
1319  LIBMESH_CHKERR(ierr);
1320  }
1321  // Specify that the new matrix is initialized
1322  // We do not close it here as `MatMatMult` ensures Y being closed
1323  Y->_is_initialized = true;
1324 }
virtual bool initialized() const
virtual numeric_index_type n() const override
libmesh_assert(ctx)
This class provides a nice interface to the PETSc C-based data structures for parallel, sparse matrices.
Definition: petsc_matrix.h:93
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ n()

template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::n ( ) const
overridevirtual
Returns
The column-dimension of the matrix.

Implements libMesh::SparseMatrix< T >.

Definition at line 1108 of file petsc_matrix.C.

1109 {
1110  libmesh_assert (this->initialized());
1111 
1112  PetscInt petsc_m=0, petsc_n=0;
1113  PetscErrorCode ierr=0;
1114 
1115  ierr = MatGetSize (_mat, &petsc_m, &petsc_n);
1116  LIBMESH_CHKERR(ierr);
1117 
1118  return static_cast<numeric_index_type>(petsc_n);
1119 }
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ n_nonzeros()

template<typename T >
std::size_t libMesh::SparseMatrix< T >::n_nonzeros ( ) const
virtualinherited
Returns
the global number of non-zero entries in the matrix sparsity pattern

Definition at line 229 of file sparse_matrix.C.

Referenced by libMesh::SparseMatrix< ValOut >::n_nonzeros().

230 {
231  if (!_sp)
232  return 0;
233  return _sp->n_nonzeros();
234 }
std::size_t n_nonzeros() const
The total number of nonzeros in the global matrix.
SparsityPattern::Build const * _sp
The sparsity pattern associated with this object.

◆ n_objects()

static unsigned int libMesh::ReferenceCounter::n_objects ( )
inlinestaticinherited

Prints the number of outstanding (created, but not yet destroyed) objects.

Definition at line 85 of file reference_counter.h.

References libMesh::ReferenceCounter::_n_objects.

Referenced by libMesh::LibMeshInit::~LibMeshInit().

86  { return _n_objects; }
static Threads::atomic< unsigned int > _n_objects
The number of objects.

◆ n_processors()

processor_id_type libMesh::ParallelObject::n_processors ( ) const
inlineinherited
Returns
The number of processors in the group.

Definition at line 103 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, libMesh::libmesh_assert(), and TIMPI::Communicator::size().

Referenced by libMesh::Partitioner::_find_global_index_by_pid_map(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::DofMap::add_constraints_to_send_list(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::DistributedMesh::add_elem(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::DistributedMesh::add_node(), libMesh::System::add_vector(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::DofMap::allgather_recursive_constraints(), libMesh::FEMSystem::assembly(), libMesh::Nemesis_IO::assert_symmetric_cmaps(), libMesh::Partitioner::assign_partitioning(), libMesh::AztecLinearSolver< T >::AztecLinearSolver(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::DistributedMesh::clear(), libMesh::DistributedMesh::clear_elems(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::ExodusII_IO::copy_scalar_solution(), libMesh::Nemesis_IO::copy_scalar_solution(), libMesh::UnstructuredMesh::create_pid_mesh(), libMesh::MeshTools::create_processor_bounding_box(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_scalar_dofs(), libMesh::DistributedMesh::DistributedMesh(), libMesh::EnsightIO::EnsightIO(), libMesh::RBEIMEvaluation::gather_bfs(), libMesh::MeshBase::get_info(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_and_attach_petscdm(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::ExodusII_IO_Helper::initialize(), libMesh::DistributedMesh::insert_elem(), libMesh::MeshTools::libmesh_assert_contiguous_dof_ids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_new_node_procids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_topology_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_boundary_ids(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::DofMap::local_variable_indices(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshBase::n_active_elem_on_proc(), libMesh::MeshBase::n_elem_on_proc(), libMesh::MeshBase::n_nodes_on_proc(), libMesh::RBEIMEvaluation::node_gather_bfs(), libMesh::Partitioner::partition(), libMesh::MeshBase::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::DofMap::prepare_send_list(), libMesh::DofMap::print_dof_constraints(), libMesh::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::CheckpointIO::read(), libMesh::CheckpointIO::read_connectivity(), libMesh::XdrIO::read_header(), libMesh::CheckpointIO::read_nodes(), libMesh::System::read_parallel_data(), libMesh::System::read_SCALAR_dofs(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::System::read_serialized_vector(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::Partitioner::repartition(), OverlappingFunctorTest::run_partitioner_test(), libMesh::DofMap::scatter_constraints(), libMesh::DistributedMesh::set_next_unique_id(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), WriteVecAndScalar::setupTests(), libMesh::RBEIMEvaluation::side_gather_bfs(), DistributedMeshTest::testRemoteElemError(), CheckpointIOTest::testSplitter(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO_Helper::write_nodal_coordinates(), libMesh::VTKIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::System::write_parallel_data(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), and libMesh::XdrIO::write_serialized_nodesets().

104  {
105  processor_id_type returnval =
106  cast_int<processor_id_type>(_communicator.size());
107  libmesh_assert(returnval); // We never have an empty comm
108  return returnval;
109  }
const Parallel::Communicator & _communicator
processor_id_type size() const
uint8_t processor_id_type
libmesh_assert(ctx)

◆ need_full_sparsity_pattern()

template<typename T>
virtual bool libMesh::SparseMatrix< T >::need_full_sparsity_pattern ( ) const
inlinevirtualinherited
Returns
true if this sparse matrix format needs to be fed the graph of the sparse matrix.

This is true for LaspackMatrix, but not PetscMatrix. In the case where the full graph is not required, we can efficiently approximate it to provide a good estimate of the required size of the sparse matrix.

Reimplemented in libMesh::EpetraMatrix< T >, and libMesh::LaspackMatrix< T >.

Definition at line 138 of file sparse_matrix.h.

Referenced by libMesh::DofMap::attach_matrix(), and libMesh::DofMap::update_sparsity_pattern().

139  { return false; }

◆ operator()()

template<typename T >
T libMesh::PetscMatrix< T >::operator() ( const numeric_index_type  i,
const numeric_index_type  j 
) const
overridevirtual
Returns
A copy of matrix entry (i,j).
Note
This may be an expensive operation, and you should always be careful where you call this function.

Implements libMesh::SparseMatrix< T >.

Definition at line 1383 of file petsc_matrix.C.

1385 {
1386  libmesh_assert (this->initialized());
1387 
1388  // PETSc 2.2.1 & newer
1389  const PetscScalar * petsc_row;
1390  const PetscInt * petsc_cols;
1391 
1392  // If the entry is not in the sparse matrix, it is 0.
1393  T value=0.;
1394 
1395  PetscErrorCode
1396  ierr=0;
1397  PetscInt
1398  ncols=0,
1399  i_val=static_cast<PetscInt>(i_in),
1400  j_val=static_cast<PetscInt>(j_in);
1401 
1402 
1403  // the matrix needs to be closed for this to work
1404  // this->close();
1405  // but closing it is a semiparallel operation; we want operator()
1406  // to run on one processor.
1407  libmesh_assert(this->closed());
1408 
1409  ierr = MatGetRow(_mat, i_val, &ncols, &petsc_cols, &petsc_row);
1410  LIBMESH_CHKERR(ierr);
1411 
1412  // Perform a binary search to find the contiguous index in
1413  // petsc_cols (resp. petsc_row) corresponding to global index j_val
1414  std::pair<const PetscInt *, const PetscInt *> p =
1415  std::equal_range (petsc_cols, petsc_cols + ncols, j_val);
1416 
1417  // Found an entry for j_val
1418  if (p.first != p.second)
1419  {
1420  // The entry in the contiguous row corresponding
1421  // to the j_val column of interest
1422  const std::size_t j =
1423  std::distance (const_cast<PetscInt *>(petsc_cols),
1424  const_cast<PetscInt *>(p.first));
1425 
1426  libmesh_assert_less (j, ncols);
1427  libmesh_assert_equal_to (petsc_cols[j], j_val);
1428 
1429  value = static_cast<T> (petsc_row[j]);
1430  }
1431 
1432  ierr = MatRestoreRow(_mat, i_val,
1433  &ncols, &petsc_cols, &petsc_row);
1434  LIBMESH_CHKERR(ierr);
1435 
1436  return value;
1437 }
virtual bool initialized() const
Real distance(const Point &p)
libmesh_assert(ctx)
virtual bool closed() const override
static const bool value
Definition: xdr_io.C:54
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ operator=() [1/2]

template<typename T>
PetscMatrix& libMesh::PetscMatrix< T >::operator= ( const PetscMatrix< T > &  )
delete

◆ operator=() [2/2]

template<typename T>
PetscMatrix& libMesh::PetscMatrix< T >::operator= ( PetscMatrix< T > &&  )
delete

◆ print() [1/2]

template<>
void libMesh::SparseMatrix< Complex >::print ( std::ostream &  os,
const bool  sparse 
) const
inherited

Definition at line 122 of file sparse_matrix.C.

123 {
124  // std::complex<>::operator<<() is defined, but use this form
125 
126  if (sparse)
127  {
128  libmesh_not_implemented();
129  }
130 
131  os << "Real part:" << std::endl;
132  for (auto i : make_range(this->m()))
133  {
134  for (auto j : make_range(this->n()))
135  os << std::setw(8) << (*this)(i,j).real() << " ";
136  os << std::endl;
137  }
138 
139  os << std::endl << "Imaginary part:" << std::endl;
140  for (auto i : make_range(this->m()))
141  {
142  for (auto j : make_range(this->n()))
143  os << std::setw(8) << (*this)(i,j).imag() << " ";
144  os << std::endl;
145  }
146 }
boost::multiprecision::float128 real(const boost::multiprecision::float128 in)
virtual numeric_index_type m() const =0
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
boost::multiprecision::float128 imag(const boost::multiprecision::float128)
virtual numeric_index_type n() const =0

◆ print() [2/2]

template<typename T >
void libMesh::SparseMatrix< T >::print ( std::ostream &  os = libMesh::out,
const bool  sparse = false 
) const
inherited

Print the contents of the matrix to the screen in a uniform style, regardless of matrix/solver package being used.

Definition at line 238 of file sparse_matrix.C.

Referenced by libMesh::EigenSparseMatrix< T >::print_personal(), and libMesh::LaspackMatrix< T >::print_personal().

239 {
240  parallel_object_only();
241 
242  libmesh_assert (this->initialized());
243 
244  const numeric_index_type first_dof = this->row_start(),
245  end_dof = this->row_stop();
246 
247  // We'll print the matrix from processor 0 to make sure
248  // it's serialized properly
249  if (this->processor_id() == 0)
250  {
251  libmesh_assert_equal_to (first_dof, 0);
252  for (numeric_index_type i : make_range(end_dof))
253  {
254  if (sparse)
255  {
256  for (auto j : make_range(this->n()))
257  {
258  T c = (*this)(i,j);
259  if (c != static_cast<T>(0.0))
260  {
261  os << i << " " << j << " " << c << std::endl;
262  }
263  }
264  }
265  else
266  {
267  for (auto j : make_range(this->n()))
268  os << (*this)(i,j) << " ";
269  os << std::endl;
270  }
271  }
272 
273  std::vector<numeric_index_type> ibuf, jbuf;
274  std::vector<T> cbuf;
275  numeric_index_type currenti = end_dof;
276  for (auto p : IntRange<processor_id_type>(1, this->n_processors()))
277  {
278  this->comm().receive(p, ibuf);
279  this->comm().receive(p, jbuf);
280  this->comm().receive(p, cbuf);
281  libmesh_assert_equal_to (ibuf.size(), jbuf.size());
282  libmesh_assert_equal_to (ibuf.size(), cbuf.size());
283 
284  if (ibuf.empty())
285  continue;
286  libmesh_assert_greater_equal (ibuf.front(), currenti);
287  libmesh_assert_greater_equal (ibuf.back(), ibuf.front());
288 
289  std::size_t currentb = 0;
290  for (;currenti <= ibuf.back(); ++currenti)
291  {
292  if (sparse)
293  {
294  for (numeric_index_type j=0; j<this->n(); j++)
295  {
296  if (currentb < ibuf.size() &&
297  ibuf[currentb] == currenti &&
298  jbuf[currentb] == j)
299  {
300  os << currenti << " " << j << " " << cbuf[currentb] << std::endl;
301  currentb++;
302  }
303  }
304  }
305  else
306  {
307  for (auto j : make_range(this->n()))
308  {
309  if (currentb < ibuf.size() &&
310  ibuf[currentb] == currenti &&
311  jbuf[currentb] == j)
312  {
313  os << cbuf[currentb] << " ";
314  currentb++;
315  }
316  else
317  os << static_cast<T>(0.0) << " ";
318  }
319  os << std::endl;
320  }
321  }
322  }
323  if (!sparse)
324  {
325  for (; currenti != this->m(); ++currenti)
326  {
327  for (numeric_index_type j=0; j<this->n(); j++)
328  os << static_cast<T>(0.0) << " ";
329  os << std::endl;
330  }
331  }
332  }
333  else
334  {
335  std::vector<numeric_index_type> ibuf, jbuf;
336  std::vector<T> cbuf;
337 
338  // We'll assume each processor has access to entire
339  // matrix rows, so (*this)(i,j) is valid if i is a local index.
340  for (numeric_index_type i : make_range(first_dof, end_dof))
341  {
342  for (auto j : make_range(this->n()))
343  {
344  T c = (*this)(i,j);
345  if (c != static_cast<T>(0.0))
346  {
347  ibuf.push_back(i);
348  jbuf.push_back(j);
349  cbuf.push_back(c);
350  }
351  }
352  }
353  this->comm().send(0,ibuf);
354  this->comm().send(0,jbuf);
355  this->comm().send(0,cbuf);
356  }
357 }
virtual bool initialized() const
const Parallel::Communicator & comm() const
virtual numeric_index_type row_stop() const =0
processor_id_type n_processors() const
Status receive(const unsigned int dest_processor_id, T &buf, const MessageTag &tag=any_tag) const
dof_id_type numeric_index_type
Definition: id_types.h:99
virtual numeric_index_type m() const =0
libmesh_assert(ctx)
virtual numeric_index_type row_start() const =0
void send(const unsigned int dest_processor_id, const T &buf, const MessageTag &tag=no_tag) const
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
processor_id_type processor_id() const
virtual numeric_index_type n() const =0

◆ print_info()

void libMesh::ReferenceCounter::print_info ( std::ostream &  out_stream = libMesh::out)
staticinherited

Prints the reference information, by default to libMesh::out.

Definition at line 81 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter, and libMesh::ReferenceCounter::get_info().

Referenced by libMesh::LibMeshInit::~LibMeshInit().

82 {
84  out_stream << ReferenceCounter::get_info();
85 }
static std::string get_info()
Gets a string containing the reference information.
static bool _enable_print_counter
Flag to control whether reference count information is printed when print_info is called...

◆ print_matlab()

template<typename T >
void libMesh::PetscMatrix< T >::print_matlab ( const std::string &  name = "") const
overridevirtual

Print the contents of the matrix in Matlab's sparse matrix format.

Optionally prints the matrix to the file named name. If name is not specified it is dumped to the screen.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 568 of file petsc_matrix.C.

569 {
570  libmesh_assert (this->initialized());
571 
572  semiparallel_only();
573 
574  if (!this->closed())
575  {
576  libmesh_deprecated();
577  libmesh_warning("The matrix must be assembled before calling PetscMatrix::print_matlab().\n"
578  "Please update your code, as this warning will become an error in a future release.");
579  const_cast<PetscMatrix<T> *>(this)->close();
580  }
581 
582  PetscErrorCode ierr = 0;
583 
584  WrappedPetsc<PetscViewer> petsc_viewer;
585  ierr = PetscViewerCreate (this->comm().get(),
586  petsc_viewer.get());
587  LIBMESH_CHKERR(ierr);
588 
589  // Create an ASCII file containing the matrix
590  // if a filename was provided.
591  if (name != "")
592  {
593  ierr = PetscViewerASCIIOpen( this->comm().get(),
594  name.c_str(),
595  petsc_viewer.get());
596  LIBMESH_CHKERR(ierr);
597 
598 #if PETSC_VERSION_LESS_THAN(3,7,0)
599  ierr = PetscViewerSetFormat (petsc_viewer,
600  PETSC_VIEWER_ASCII_MATLAB);
601 #else
602  ierr = PetscViewerPushFormat (petsc_viewer,
603  PETSC_VIEWER_ASCII_MATLAB);
604 #endif
605 
606  LIBMESH_CHKERR(ierr);
607 
608  ierr = MatView (_mat, petsc_viewer);
609  LIBMESH_CHKERR(ierr);
610  }
611 
612  // Otherwise the matrix will be dumped to the screen.
613  else
614  {
615 #if PETSC_VERSION_LESS_THAN(3,7,0)
616  ierr = PetscViewerSetFormat (PETSC_VIEWER_STDOUT_WORLD,
617  PETSC_VIEWER_ASCII_MATLAB);
618 #else
619  ierr = PetscViewerPushFormat (PETSC_VIEWER_STDOUT_WORLD,
620  PETSC_VIEWER_ASCII_MATLAB);
621 #endif
622 
623  LIBMESH_CHKERR(ierr);
624 
625  ierr = MatView (_mat, PETSC_VIEWER_STDOUT_WORLD);
626  LIBMESH_CHKERR(ierr);
627  }
628 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
virtual bool initialized() const
const Parallel::Communicator & comm() const
libmesh_assert(ctx)
virtual bool closed() const override
This class provides a nice interface to the PETSc C-based data structures for parallel, sparse matrices.
Definition: petsc_matrix.h:93
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367
virtual void close() override
Calls the SparseMatrix&#39;s internal assembly routines, ensuring that the values are consistent across p...

◆ print_personal()

template<typename T >
void libMesh::PetscMatrix< T >::print_personal ( std::ostream &  os = libMesh::out) const
overridevirtual

Print the contents of the matrix to the screen with the PETSc viewer.

This function only allows printing to standard out since we have limited ourselves to one PETSc implementation for writing.

Implements libMesh::SparseMatrix< T >.

Definition at line 635 of file petsc_matrix.C.

636 {
637  libmesh_assert (this->initialized());
638 
639  // Routine must be called in parallel on parallel matrices
640  // and serial on serial matrices.
641  semiparallel_only();
642 
643  // #ifndef NDEBUG
644  // if (os != std::cout)
645  // libMesh::err << "Warning! PETSc can only print to std::cout!" << std::endl;
646  // #endif
647 
648  // Matrix must be in an assembled state to be printed
649  if (!this->closed())
650  {
651  libmesh_deprecated();
652  libmesh_warning("The matrix must be assembled before calling PetscMatrix::print_personal().\n"
653  "Please update your code, as this warning will become an error in a future release.");
654  const_cast<PetscMatrix<T> *>(this)->close();
655  }
656 
657  PetscErrorCode ierr=0;
658 
659  // Print to screen if ostream is stdout
660  if (os.rdbuf() == std::cout.rdbuf())
661  {
662  ierr = MatView(_mat, NULL);
663  LIBMESH_CHKERR(ierr);
664  }
665 
666  // Otherwise, print to the requested file, in a roundabout way...
667  else
668  {
669  // We will create a temporary filename, and file, for PETSc to
670  // write to.
671  std::string temp_filename;
672 
673  {
674  // Template for temporary filename
675  char c[] = "temp_petsc_matrix.XXXXXX";
676 
677  // Generate temporary, unique filename only on processor 0. We will
678  // use this filename for PetscViewerASCIIOpen, before copying it into
679  // the user's stream
680  if (this->processor_id() == 0)
681  {
682  int fd = mkstemp(c);
683 
684  // Check to see that mkstemp did not fail.
685  libmesh_error_msg_if(fd == -1, "mkstemp failed in PetscMatrix::print_personal()");
686 
687  // mkstemp returns a file descriptor for an open file,
688  // so let's close it before we hand it to PETSc!
689  ::close (fd);
690  }
691 
692  // Store temporary filename as string, makes it easier to broadcast
693  temp_filename = c;
694  }
695 
696  // Now broadcast the filename from processor 0 to all processors.
697  this->comm().broadcast(temp_filename);
698 
699  // PetscViewer object for passing to MatView
700  PetscViewer petsc_viewer;
701 
702  // This PETSc function only takes a string and handles the opening/closing
703  // of the file internally. Since print_personal() takes a reference to
704  // an ostream, we have to do an extra step... print_personal() should probably
705  // have a version that takes a string to get rid of this problem.
706  ierr = PetscViewerASCIIOpen( this->comm().get(),
707  temp_filename.c_str(),
708  &petsc_viewer);
709  LIBMESH_CHKERR(ierr);
710 
711  // Probably don't need to set the format if it's default...
712  // ierr = PetscViewerSetFormat (petsc_viewer,
713  // PETSC_VIEWER_DEFAULT);
714  // LIBMESH_CHKERR(ierr);
715 
716  // Finally print the matrix using the viewer
717  ierr = MatView (_mat, petsc_viewer);
718  LIBMESH_CHKERR(ierr);
719 
720  if (this->processor_id() == 0)
721  {
722  // Now the inefficient bit: open temp_filename as an ostream and copy the contents
723  // into the user's desired ostream. We can't just do a direct file copy, we don't even have the filename!
724  std::ifstream input_stream(temp_filename.c_str());
725  os << input_stream.rdbuf(); // The "most elegant" way to copy one stream into another.
726  // os.close(); // close not defined in ostream
727 
728  // Now remove the temporary file
729  input_stream.close();
730  std::remove(temp_filename.c_str());
731  }
732  }
733 }
virtual bool initialized() const
const Parallel::Communicator & comm() const
int mkstemp(char *tmpl)
Definition: win_mkstemp.h:13
libmesh_assert(ctx)
void broadcast(T &data, const unsigned int root_id=0, const bool identical_sizes=false) const
virtual bool closed() const override
This class provides a nice interface to the PETSc C-based data structures for parallel, sparse matrices.
Definition: petsc_matrix.h:93
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367
processor_id_type processor_id() const
virtual void close() override
Calls the SparseMatrix&#39;s internal assembly routines, ensuring that the values are consistent across p...

◆ print_petsc_binary()

template<typename T >
void libMesh::PetscMatrix< T >::print_petsc_binary ( const std::string &  filename)
overridevirtual

Write the contents of the matrix to a file in PETSc's binary sparse matrix format.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 778 of file petsc_matrix.C.

779 {
780  libmesh_assert (this->initialized());
781 
782  this->_petsc_viewer(filename, PETSCVIEWERBINARY, FILE_MODE_WRITE);
783 }
virtual bool initialized() const
libmesh_assert(ctx)
void _petsc_viewer(const std::string &filename, PetscViewerType viewertype, PetscFileMode filemode)
Definition: petsc_matrix.C:738

◆ print_petsc_hdf5()

template<typename T >
void libMesh::PetscMatrix< T >::print_petsc_hdf5 ( const std::string &  filename)
overridevirtual

Write the contents of the matrix to a file in PETSc's HDF5 sparse matrix format.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 788 of file petsc_matrix.C.

789 {
790  libmesh_assert (this->initialized());
791 
792  this->_petsc_viewer(filename, PETSCVIEWERHDF5, FILE_MODE_WRITE);
793 }
virtual bool initialized() const
libmesh_assert(ctx)
void _petsc_viewer(const std::string &filename, PetscViewerType viewertype, PetscFileMode filemode)
Definition: petsc_matrix.C:738

◆ processor_id()

processor_id_type libMesh::ParallelObject::processor_id ( ) const
inlineinherited
Returns
The rank of this processor in the group.

Definition at line 114 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, and TIMPI::Communicator::rank().

Referenced by libMesh::BoundaryInfo::_find_id_maps(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::DistributedMesh::add_elem(), libMesh::BoundaryInfo::add_elements(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::DistributedMesh::add_node(), libMesh::MeshTools::Modification::all_tri(), libMesh::DofMap::allgather_recursive_constraints(), libMesh::FEMSystem::assembly(), libMesh::Nemesis_IO::assert_symmetric_cmaps(), libMesh::Partitioner::assign_partitioning(), libMesh::Nemesis_IO_Helper::build_element_and_node_maps(), libMesh::Partitioner::build_graph(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::EquationSystems::build_parallel_solution_vector(), libMesh::DistributedMesh::clear(), libMesh::DistributedMesh::clear_elems(), libMesh::ExodusII_IO_Helper::close(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::compute_communication_map_parameters(), libMesh::Nemesis_IO_Helper::compute_internal_and_border_elems_and_internal_nodes(), libMesh::RBConstruction::compute_max_error_bound(), libMesh::Nemesis_IO_Helper::compute_node_communication_maps(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::ExodusII_IO::copy_elemental_solution(), libMesh::ExodusII_IO::copy_nodal_solution(), libMesh::ExodusII_IO::copy_scalar_solution(), libMesh::Nemesis_IO::copy_scalar_solution(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::ExodusII_IO_Helper::create(), libMesh::DistributedMesh::delete_elem(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DofMap::distribute_scalar_dofs(), libMesh::DistributedMesh::DistributedMesh(), libMesh::DofMap::end_dof(), libMesh::DofMap::end_old_dof(), libMesh::EnsightIO::EnsightIO(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::find_dofs_to_send(), libMesh::MeshFunction::find_element(), libMesh::MeshFunction::find_elements(), libMesh::UnstructuredMesh::find_neighbors(), libMesh::DofMap::first_dof(), libMesh::DofMap::first_old_dof(), libMesh::RBEIMEvaluation::gather_bfs(), libMesh::Nemesis_IO_Helper::get_cmap_params(), libMesh::Nemesis_IO_Helper::get_eb_info_global(), libMesh::Nemesis_IO_Helper::get_elem_cmap(), libMesh::Nemesis_IO_Helper::get_elem_map(), libMesh::MeshBase::get_info(), libMesh::DofMap::get_info(), libMesh::Nemesis_IO_Helper::get_init_global(), libMesh::Nemesis_IO_Helper::get_init_info(), libMesh::Nemesis_IO_Helper::get_loadbal_param(), libMesh::DofMap::get_local_constraints(), libMesh::Nemesis_IO_Helper::get_node_cmap(), libMesh::Nemesis_IO_Helper::get_node_map(), libMesh::Nemesis_IO_Helper::get_ns_param_global(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), libMesh::SparsityPattern::Build::handle_vi_vj(), libMesh::LaplaceMeshSmoother::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_and_attach_petscdm(), HeatSystem::init_data(), libMesh::ExodusII_IO_Helper::initialize(), libMesh::ExodusII_IO_Helper::initialize_element_variables(), libMesh::ExodusII_IO_Helper::initialize_global_variables(), libMesh::ExodusII_IO_Helper::initialize_nodal_variables(), libMesh::DistributedMesh::insert_elem(), libMesh::DofMap::is_evaluable(), libMesh::SparsityPattern::Build::join(), libMesh::TransientRBEvaluation::legacy_write_offline_data_to_files(), libMesh::RBSCMEvaluation::legacy_write_offline_data_to_files(), libMesh::RBEvaluation::legacy_write_offline_data_to_files(), libMesh::MeshTools::libmesh_assert_consistent_distributed(), libMesh::MeshTools::libmesh_assert_consistent_distributed_nodes(), libMesh::MeshTools::libmesh_assert_contiguous_dof_ids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_object_ids(), libMesh::DofMap::local_variable_indices(), main(), libMesh::MeshRefinement::make_coarsening_compatible(), AugmentSparsityOnInterface::mesh_reinit(), libMesh::TriangulatorInterface::MeshedHole::MeshedHole(), libMesh::MeshBase::n_active_local_elem(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::DofMap::n_local_dofs(), libMesh::System::n_local_dofs(), libMesh::MeshBase::n_local_elem(), libMesh::MeshBase::n_local_nodes(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::BoundaryInfo::n_shellface_conds(), libMesh::RBEIMEvaluation::node_gather_bfs(), libMesh::SparsityPattern::Build::operator()(), libMesh::DistributedMesh::own_node(), libMesh::BoundaryInfo::parallel_sync_node_ids(), libMesh::BoundaryInfo::parallel_sync_side_ids(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::DofMap::print_dof_constraints(), libMesh::DofMap::process_mesh_constraint_rows(), libMesh::Nemesis_IO_Helper::put_cmap_params(), libMesh::Nemesis_IO_Helper::put_elem_cmap(), libMesh::Nemesis_IO_Helper::put_elem_map(), libMesh::Nemesis_IO_Helper::put_loadbal_param(), libMesh::Nemesis_IO_Helper::put_node_cmap(), libMesh::Nemesis_IO_Helper::put_node_map(), libMesh::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::CheckpointIO::read(), libMesh::EquationSystems::read(), libMesh::ExodusII_IO_Helper::read_elem_num_map(), libMesh::ExodusII_IO_Helper::read_global_values(), libMesh::ExodusII_IO::read_header(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), libMesh::DynaIO::read_mesh(), libMesh::ExodusII_IO_Helper::read_node_num_map(), libMesh::System::read_parallel_data(), libMesh::TransientRBConstruction::read_riesz_representors_from_files(), libMesh::RBConstruction::read_riesz_representors_from_files(), libMesh::System::read_SCALAR_dofs(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::System::read_serialized_data(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::System::read_serialized_vector(), libMesh::System::read_serialized_vectors(), libMesh::Nemesis_IO_Helper::read_var_names_impl(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::DistributedMesh::renumber_nodes_and_elements(), libMesh::DofMap::scatter_constraints(), libMesh::CheckpointIO::select_split_config(), libMesh::DistributedMesh::set_next_unique_id(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::RBEIMEvaluation::side_gather_bfs(), ExodusTest< elem_type >::test_read_gold(), ExodusTest< elem_type >::test_write(), MeshInputTest::testAbaqusRead(), MeshInputTest::testCopyElementSolutionImpl(), MeshInputTest::testCopyElementVectorImpl(), MeshInputTest::testCopyNodalSolutionImpl(), DefaultCouplingTest::testCoupling(), PointNeighborCouplingTest::testCoupling(), MeshInputTest::testDynaFileMappings(), MeshInputTest::testDynaNoSplines(), MeshInputTest::testDynaReadElem(), MeshInputTest::testDynaReadPatch(), MeshInputTest::testExodusFileMappings(), MeshInputTest::testExodusIGASidesets(), MeshInputTest::testExodusWriteElementDataFromDiscontinuousNodalData(), MeshInputTest::testLowOrderEdgeBlocks(), SystemsTest::testProjectMatrix1D(), SystemsTest::testProjectMatrix2D(), SystemsTest::testProjectMatrix3D(), BoundaryInfoTest::testShellFaceConstraints(), MeshInputTest::testSingleElementImpl(), WriteVecAndScalar::testSolution(), CheckpointIOTest::testSplitter(), MeshInputTest::testTetgenIO(), libMesh::MeshTools::total_weight(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::DTKAdapter::update_variable_values(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::CheckpointIO::write(), libMesh::EquationSystems::write(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO::write_element_data(), libMesh::ExodusII_IO_Helper::write_element_values(), libMesh::ExodusII_IO_Helper::write_element_values_element_major(), libMesh::ExodusII_IO_Helper::write_elements(), libMesh::ExodusII_IO_Helper::write_elemset_data(), libMesh::ExodusII_IO_Helper::write_elemsets(), libMesh::ExodusII_IO::write_global_data(), libMesh::ExodusII_IO_Helper::write_global_values(), libMesh::System::write_header(), libMesh::ExodusII_IO::write_information_records(), libMesh::ExodusII_IO_Helper::write_information_records(), libMesh::ExodusII_IO_Helper::write_nodal_coordinates(), libMesh::UCDIO::write_nodal_data(), libMesh::VTKIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_common(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::ExodusII_IO_Helper::write_nodal_values(), libMesh::ExodusII_IO_Helper::write_nodeset_data(), libMesh::Nemesis_IO_Helper::write_nodesets(), libMesh::ExodusII_IO_Helper::write_nodesets(), libMesh::RBEIMEvaluation::write_out_interior_basis_functions(), libMesh::RBEIMEvaluation::write_out_node_basis_functions(), libMesh::RBEIMEvaluation::write_out_side_basis_functions(), write_output_solvedata(), libMesh::System::write_parallel_data(), libMesh::RBConstruction::write_riesz_representors_to_files(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bc_names(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::System::write_serialized_data(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::XdrIO::write_serialized_subdomain_names(), libMesh::System::write_serialized_vector(), libMesh::System::write_serialized_vectors(), libMesh::ExodusII_IO_Helper::write_sideset_data(), libMesh::Nemesis_IO_Helper::write_sidesets(), libMesh::ExodusII_IO_Helper::write_sidesets(), libMesh::ExodusII_IO::write_timestep(), libMesh::ExodusII_IO_Helper::write_timestep(), and libMesh::ExodusII_IO::write_timestep_discontinuous().

115  { return cast_int<processor_id_type>(_communicator.rank()); }
processor_id_type rank() const
const Parallel::Communicator & _communicator

◆ read_matlab()

template<typename T >
void libMesh::SparseMatrix< T >::read_matlab ( const std::string &  filename)
virtualinherited

Read the contents of the matrix from Matlab's sparse matrix format.

If the size and sparsity of the matrix in filename appears consistent with the existing sparsity of this then the existing parallel decomposition and sparsity will be retained. If not, then this will be initialized with the sparsity from the file, linearly partitioned onto the number of processors available.

Definition at line 502 of file sparse_matrix.C.

Referenced by ConstraintOperatorTest::test1DCoarseningNewNodes(), and ConstraintOperatorTest::testCoreform().

503 {
504 #ifndef LIBMESH_HAVE_CXX11_REGEX
505  libmesh_not_implemented(); // What is your compiler?!? Email us!
506  libmesh_ignore(filename);
507 #else
508  parallel_object_only();
509 
510  // The sizes we get from the file
511  std::size_t m = 0,
512  n = 0;
513 
514  // We'll read through the file three times: once to get a reliable
515  // value for the matrix size (so we can divvy it up among
516  // processors), then again to get the sparsity to send to each
517  // processor, then a final time to get the entries to send to each
518  // processor.
519  std::unique_ptr<std::ifstream> file;
520 
521  // We'll be using regular expressions to make ourselves slightly
522  // more robust to formatting.
523  const std::regex start_regex // assignment like "zzz = ["
524  ("\\s*\\w+\\s*=\\s*\\[");
525  const std::regex entry_regex // row/col/val like "1 1 -2.0e-4"
526  ("(\\d+)\\s+(\\d+)\\s+([+-]?(\\d+([.]\\d*)?([eE][+-]?\\d+)?|[.]\\d+([eE][+-]?\\d+)?))");
527  const std::regex end_regex // end of assignment
528  ("^[^%]*\\]");
529 
530  // We'll read the matrix on processor 0 rather than try to juggle
531  // parallel I/O. Start by reading to deduce the size.
532  if (this->processor_id() == 0)
533  {
534  file = std::make_unique<std::ifstream>(filename.c_str());
535 
536  // If we have a matrix with all-zero trailing rows, the only
537  // way to get the size is if it ended up in a comment
538  const std::regex size_regex // comment like "% size = 8 8"
539  ("%\\s*[Ss][Ii][Zz][Ee]\\s*=\\s*(\\d+)\\s+(\\d+)");
540 
541  // Get the size
542  bool have_started = false;
543  bool have_ended = false;
544  std::size_t largest_i_seen = 0, largest_j_seen = 0;
545  for (std::string line; std::getline(*file, line);)
546  {
547  std::smatch sm;
548 
549  if (std::regex_search(line, sm, size_regex))
550  {
551  const std::string msize = sm[1];
552  const std::string nsize = sm[2];
553  m = std::stoull(msize);
554  n = std::stoull(nsize);
555  }
556 
557  if (std::regex_search(line, start_regex))
558  have_started = true;
559 
560  if (std::regex_search(line, sm, entry_regex))
561  {
562  libmesh_error_msg_if
563  (!have_started, "Confused by premature entries in matrix file " << filename);
564 
565  const std::string istr = sm[1];
566  const std::string jstr = sm[2];
567 
568  std::size_t i = std::stoull(istr);
569  largest_i_seen = std::max(i, largest_i_seen);
570  std::size_t j = std::stoull(jstr);
571  largest_j_seen = std::max(j, largest_j_seen);
572  }
573 
574  if (std::regex_search(line, end_regex))
575  {
576  have_ended = true;
577  break;
578  }
579  }
580 
581  libmesh_error_msg_if
582  (!have_started, "Confused by missing assignment beginning in matrix file " << filename);
583 
584  libmesh_error_msg_if
585  (!have_ended, "Confused by missing assignment ending in matrix file " << filename);
586 
587  libmesh_error_msg_if
588  (m > largest_i_seen, "Confused by missing final row(s) in matrix file " << filename);
589 
590  libmesh_error_msg_if
591  (m > 0 && m < largest_i_seen, "Confused by extra final row(s) in matrix file " << filename);
592 
593  if (!m)
594  m = largest_i_seen;
595 
596  libmesh_error_msg_if
597  (n > largest_j_seen, "Confused by missing final column(s) in matrix file " << filename);
598 
599  libmesh_error_msg_if
600  (n > 0 && n < largest_j_seen, "Confused by extra final column(s) in matrix file " << filename);
601 
602  if (!n)
603  n = largest_j_seen;
604 
605  this->comm().broadcast(m);
606  this->comm().broadcast(n);
607  }
608  else
609  {
610  this->comm().broadcast(m);
611  this->comm().broadcast(n);
612  }
613 
614  // If we don't already have this size, we'll need to reinit later,
615  // and we'll need to redetermine which rows each processor is in
616  // charge of.
617 
619  new_row_start = this->processor_id() * m / this->n_processors(),
620  new_row_stop = (this->processor_id()+1) * m / this->n_processors();
622  new_col_start = this->processor_id() * n / this->n_processors(),
623  new_col_stop = (this->processor_id()+1) * n / this->n_processors();
624 
625  if (this->initialized() &&
626  m == this->m() &&
627  n == this->n())
628  {
629  new_row_start = this->row_start(),
630  new_row_stop = this->row_stop();
631 
632  new_col_start = this->col_start(),
633  new_col_stop = this->col_stop();
634  }
635 
636  std::vector<numeric_index_type> new_row_starts, new_row_stops,
637  new_col_starts, new_col_stops;
638 
639  this->comm().gather(0, new_row_start, new_row_starts);
640  this->comm().gather(0, new_row_stop, new_row_stops);
641  this->comm().gather(0, new_col_start, new_col_starts);
642  this->comm().gather(0, new_col_stop, new_col_stops);
643 
644  // Reread to deduce the sparsity pattern, or at least the maximum
645  // number of on- and off- diagonal non-zeros per row.
646  numeric_index_type on_diagonal_nonzeros =0,
647  off_diagonal_nonzeros =0;
648 
649  if (this->processor_id() == 0)
650  {
651  file->seekg(0);
652 
653  bool have_started = false;
654 
655  // Data for the row we're working on
656 
657  // Use 1-based indexing for current_row, as in the file
658  numeric_index_type current_row = 1;
659  processor_id_type current_proc = 0;
660  numeric_index_type current_on_diagonal_nonzeros = 0;
661  numeric_index_type current_off_diagonal_nonzeros = 0;
662 
663  for (std::string line; std::getline(*file, line);)
664  {
665  std::smatch sm;
666 
667  if (std::regex_search(line, start_regex))
668  have_started = true;
669 
670  if (have_started && std::regex_search(line, sm, entry_regex))
671  {
672  const std::string istr = sm[1];
673  const std::string jstr = sm[2];
674 
675  const numeric_index_type i = std::stoull(istr);
676  const numeric_index_type j = std::stoull(jstr);
677 
678  libmesh_error_msg_if
679  (!i || !j, "Expected 1-based indexing in matrix file "
680  << filename);
681 
682  libmesh_error_msg_if
683  (i < current_row,
684  "Can't handle out-of-order entries in matrix file "
685  << filename);
686  if (i > current_row)
687  {
688  current_row = i;
689  // +1 for 1-based indexing in file
690  while (current_row >= (new_row_stops[current_proc]+1))
691  ++current_proc;
692  current_on_diagonal_nonzeros = 0;
693  current_off_diagonal_nonzeros = 0;
694  }
695 
696  // +1 for 1-based indexing in file
697  if (j >= (new_col_starts[current_proc]+1) &&
698  j < (new_col_stops[current_proc]+1))
699  {
700  ++current_on_diagonal_nonzeros;
701  on_diagonal_nonzeros =
702  std::max(on_diagonal_nonzeros,
703  current_on_diagonal_nonzeros);
704  }
705  else
706  {
707  ++current_off_diagonal_nonzeros;
708  off_diagonal_nonzeros =
709  std::max(off_diagonal_nonzeros,
710  current_off_diagonal_nonzeros);
711  }
712  }
713 
714  if (std::regex_search(line, end_regex))
715  break;
716  }
717  }
718 
719  this->comm().broadcast(on_diagonal_nonzeros);
720  this->comm().broadcast(off_diagonal_nonzeros);
721 
722  this->init(m, n,
723  new_row_stop-new_row_start,
724  new_col_stop-new_col_start,
725  on_diagonal_nonzeros,
726  off_diagonal_nonzeros);
727 
728  // One last reread to set values
729  if (this->processor_id() == 0)
730  {
731  file->seekg(0);
732 
733  bool have_started = false;
734 
735  for (std::string line; std::getline(*file, line);)
736  {
737  std::smatch sm;
738 
739  if (std::regex_search(line, start_regex))
740  have_started = true;
741 
742  if (have_started && std::regex_search(line, sm, entry_regex))
743  {
744  const std::string istr = sm[1];
745  const std::string jstr = sm[2];
746  const std::string cstr = sm[3];
747 
748  const numeric_index_type i = std::stoull(istr);
749  const numeric_index_type j = std::stoull(jstr);
750 
751  // Try to be compatible with
752  // higher-than-double-precision T
753  std::stringstream ss(cstr);
754  T value;
755  ss >> value;
756 
757  // Convert from 1-based to 0-based indexing
758  this->set(i-1, j-1, value);
759  }
760 
761  if (std::regex_search(line, end_regex))
762  break;
763  }
764  }
765  this->close();
766 #endif
767 }
virtual bool initialized() const
void gather(const unsigned int root_id, const T &send_data, std::vector< T, A > &recv) const
const Parallel::Communicator & comm() const
virtual numeric_index_type row_stop() const =0
uint8_t processor_id_type
processor_id_type n_processors() const
void libmesh_ignore(const Args &...)
dof_id_type numeric_index_type
Definition: id_types.h:99
virtual numeric_index_type m() const =0
virtual numeric_index_type col_stop() const =0
virtual numeric_index_type col_start() const =0
void broadcast(T &data, const unsigned int root_id=0, const bool identical_sizes=false) const
virtual void close()=0
Calls the SparseMatrix&#39;s internal assembly routines, ensuring that the values are consistent across p...
virtual numeric_index_type row_start() const =0
static const bool value
Definition: xdr_io.C:54
processor_id_type processor_id() const
virtual numeric_index_type n() const =0
virtual void init(const numeric_index_type m, const numeric_index_type n, const numeric_index_type m_l, const numeric_index_type n_l, const numeric_index_type nnz=30, const numeric_index_type noz=10, const numeric_index_type blocksize=1)=0
Initialize SparseMatrix with the specified sizes.

◆ read_petsc_binary()

template<typename T >
void libMesh::PetscMatrix< T >::read_petsc_binary ( const std::string &  filename)
overridevirtual

Read the contents of the matrix from a file in PETSc's binary sparse matrix format.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 798 of file petsc_matrix.C.

799 {
800  this->_petsc_viewer(filename, PETSCVIEWERBINARY, FILE_MODE_READ);
801 }
void _petsc_viewer(const std::string &filename, PetscViewerType viewertype, PetscFileMode filemode)
Definition: petsc_matrix.C:738

◆ read_petsc_hdf5()

template<typename T >
void libMesh::PetscMatrix< T >::read_petsc_hdf5 ( const std::string &  filename)
overridevirtual

Read the contents of the matrix from a file in PETSc's HDF5 sparse matrix format.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 806 of file petsc_matrix.C.

807 {
808  this->_petsc_viewer(filename, PETSCVIEWERHDF5, FILE_MODE_READ);
809 }
void _petsc_viewer(const std::string &filename, PetscViewerType viewertype, PetscFileMode filemode)
Definition: petsc_matrix.C:738

◆ reinit_submatrix()

template<typename T>
virtual void libMesh::SparseMatrix< T >::reinit_submatrix ( SparseMatrix< T > &  submatrix,
const std::vector< numeric_index_type > &  rows,
const std::vector< numeric_index_type > &  cols 
) const
inlinevirtualinherited

This function is similar to the one above, but it allows you to reuse the existing sparsity pattern of "submatrix" instead of reallocating it again.

This should hopefully be more efficient if you are frequently extracting submatrices of the same size.

Definition at line 500 of file sparse_matrix.h.

503  {
504  this->_get_submatrix(submatrix,
505  rows,
506  cols,
507  true); // true means REUSE submatrix
508  }
virtual void _get_submatrix(SparseMatrix< T > &, const std::vector< numeric_index_type > &, const std::vector< numeric_index_type > &, const bool) const
Protected implementation of the create_submatrix and reinit_submatrix routines.

◆ reset_preallocation()

template<typename T >
void libMesh::PetscMatrix< T >::reset_preallocation ( )

Reset matrix to use the original nonzero pattern provided by users.

Definition at line 402 of file petsc_matrix.C.

403 {
404 #if !PETSC_VERSION_LESS_THAN(3,9,0)
405  libmesh_assert (this->initialized());
406 
407  auto ierr = MatResetPreallocation(_mat);
408  LIBMESH_CHKERR(ierr);
409 #else
410  libmesh_warning("Your version of PETSc doesn't support resetting of "
411  "preallocation, so we will use your most recent sparsity "
412  "pattern. This may result in a degradation of performance\n");
413 #endif
414 }
virtual bool initialized() const
libmesh_assert(ctx)
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ row_start()

template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::row_start ( ) const
overridevirtual
Returns
The index of the first matrix row stored on this processor.

Implements libMesh::SparseMatrix< T >.

Definition at line 1152 of file petsc_matrix.C.

1153 {
1154  libmesh_assert (this->initialized());
1155 
1156  PetscInt start=0, stop=0;
1157  PetscErrorCode ierr=0;
1158 
1159  ierr = MatGetOwnershipRange(_mat, &start, &stop);
1160  LIBMESH_CHKERR(ierr);
1161 
1162  return static_cast<numeric_index_type>(start);
1163 }
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
void stop(const char *file, int line, const char *date, const char *time)
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ row_stop()

template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::row_stop ( ) const
overridevirtual
Returns
The index of the last matrix row (+1) stored on this processor.

Implements libMesh::SparseMatrix< T >.

Definition at line 1168 of file petsc_matrix.C.

1169 {
1170  libmesh_assert (this->initialized());
1171 
1172  PetscInt start=0, stop=0;
1173  PetscErrorCode ierr=0;
1174 
1175  ierr = MatGetOwnershipRange(_mat, &start, &stop);
1176  LIBMESH_CHKERR(ierr);
1177 
1178  return static_cast<numeric_index_type>(stop);
1179 }
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:99
libmesh_assert(ctx)
void stop(const char *file, int line, const char *date, const char *time)
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ set()

template<typename T>
void libMesh::PetscMatrix< T >::set ( const numeric_index_type  i,
const numeric_index_type  j,
const T  value 
)
overridevirtual

Set the element (i,j) to value.

Throws an error if the entry does not exist. Zero values can be "stored" in non-existent fields.

Implements libMesh::SparseMatrix< T >.

Definition at line 1216 of file petsc_matrix.C.

1219 {
1220  libmesh_assert (this->initialized());
1221 
1222  PetscErrorCode ierr=0;
1223  PetscInt i_val=i, j_val=j;
1224 
1225  PetscScalar petsc_value = static_cast<PetscScalar>(value);
1226  ierr = MatSetValues(_mat, 1, &i_val, 1, &j_val,
1227  &petsc_value, INSERT_VALUES);
1228  LIBMESH_CHKERR(ierr);
1229 }
virtual bool initialized() const
libmesh_assert(ctx)
static const bool value
Definition: xdr_io.C:54
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ set_destroy_mat_on_exit()

template<typename T >
void libMesh::PetscMatrix< T >::set_destroy_mat_on_exit ( bool  destroy = true)
virtual

If set to false, we don't delete the Mat on destruction and allow instead for PETSc to manage it.

Definition at line 1511 of file petsc_matrix.C.

1512 {
1513  this->_destroy_mat_on_exit = destroy;
1514 }
bool _destroy_mat_on_exit
This boolean value should only be set to false for the constructor which takes a PETSc Mat object...
Definition: petsc_matrix.h:373
void destroy(triangulateio &t, IO_Type)
Frees any memory which has been dynamically allocated by Triangle.

◆ set_matrix_type()

template<typename T >
void libMesh::PetscMatrix< T >::set_matrix_type ( PetscMatrixType  mat_type)

Definition at line 115 of file petsc_matrix.C.

116 {
117  _mat_type = mat_type;
118 }
PetscMatrixType _mat_type
Definition: petsc_matrix.h:375

◆ swap()

template<typename T>
void libMesh::PetscMatrix< T >::swap ( PetscMatrix< T > &  m_in)

Swaps the internal data pointers of two PetscMatrices, no actual values are swapped.

Definition at line 1518 of file petsc_matrix.C.

Referenced by libMesh::__libmesh_petsc_diff_solver_jacobian().

1519 {
1520  std::swap(_mat, m_in._mat);
1521  std::swap(_destroy_mat_on_exit, m_in._destroy_mat_on_exit);
1522 }
bool _destroy_mat_on_exit
This boolean value should only be set to false for the constructor which takes a PETSc Mat object...
Definition: petsc_matrix.h:373
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ update_preallocation_and_zero()

template<typename T >
void libMesh::PetscMatrix< T >::update_preallocation_and_zero ( )

Update the sparsity pattern based on dof_map, and set the matrix to zero.

This is useful in cases where the sparsity pattern changes during a computation.

Definition at line 396 of file petsc_matrix.C.

397 {
398  libmesh_not_implemented();
399 }

◆ update_sparsity_pattern()

template<typename T>
virtual void libMesh::SparseMatrix< T >::update_sparsity_pattern ( const SparsityPattern::Graph )
inlinevirtualinherited

Updates the matrix sparsity pattern.

When your SparseMatrix<T> implementation does not need this data, simply do not override this method.

Reimplemented in libMesh::EpetraMatrix< T >, and libMesh::LaspackMatrix< T >.

Definition at line 146 of file sparse_matrix.h.

Referenced by libMesh::DofMap::update_sparsity_pattern().

146 {}

◆ vector_mult()

template<typename T>
void libMesh::SparseMatrix< T >::vector_mult ( NumericVector< T > &  dest,
const NumericVector< T > &  arg 
) const
inherited

Multiplies the matrix by the NumericVector arg and stores the result in NumericVector dest.

Definition at line 200 of file sparse_matrix.C.

Referenced by libMesh::TransientRBConstruction::add_IC_to_RB_space(), libMesh::RBSCMConstruction::Aq_inner_product(), libMesh::RBSCMConstruction::B_inner_product(), libMesh::RBConstruction::compute_Fq_representor_innerprods(), libMesh::RBConstruction::compute_output_dual_innerprods(), libMesh::RBConstruction::compute_residual_dual_norm_slow(), libMesh::TransientRBConstruction::enrich_RB_space(), libMesh::RBConstruction::enrich_RB_space(), AssembleOptimization::gradient(), libMesh::TransientRBConstruction::mass_matrix_scaled_matvec(), AssembleOptimization::objective(), libMesh::RBConstruction::print_basis_function_orthogonality(), libMesh::ImplicitSystem::qoi_parameter_hessian(), libMesh::ImplicitSystem::qoi_parameter_hessian_vector_product(), libMesh::TransientRBConstruction::set_error_temporal_data(), libMesh::RBConstruction::train_reduced_basis_with_POD(), libMesh::TransientRBConstruction::truth_assembly(), libMesh::RBConstruction::truth_solve(), libMesh::TransientRBConstruction::update_RB_system_matrices(), libMesh::RBConstruction::update_RB_system_matrices(), libMesh::TransientRBConstruction::update_residual_terms(), and libMesh::RBConstruction::update_residual_terms().

202 {
203  dest.zero();
204  this->vector_mult_add(dest,arg);
205 }
void vector_mult_add(NumericVector< T > &dest, const NumericVector< T > &arg) const
Multiplies the matrix by the NumericVector arg and adds the result to the NumericVector dest...

◆ vector_mult_add()

template<typename T>
void libMesh::SparseMatrix< T >::vector_mult_add ( NumericVector< T > &  dest,
const NumericVector< T > &  arg 
) const
inherited

Multiplies the matrix by the NumericVector arg and adds the result to the NumericVector dest.

Definition at line 210 of file sparse_matrix.C.

Referenced by libMesh::ImplicitSystem::weighted_sensitivity_adjoint_solve().

212 {
213  /* This functionality is actually implemented in the \p
214  NumericVector class. */
215  dest.add_vector(arg,*this);
216 }

◆ zero()

template<typename T >
void libMesh::PetscMatrix< T >::zero ( )
overridevirtual

Set all entries to 0.

Implements libMesh::SparseMatrix< T >.

Definition at line 417 of file petsc_matrix.C.

418 {
419  libmesh_assert (this->initialized());
420 
421  semiparallel_only();
422 
423  PetscErrorCode ierr=0;
424 
425  PetscInt m_l, n_l;
426 
427  ierr = MatGetLocalSize(_mat,&m_l,&n_l);
428  LIBMESH_CHKERR(ierr);
429 
430  if (n_l)
431  {
432  ierr = MatZeroEntries(_mat);
433  LIBMESH_CHKERR(ierr);
434  }
435 }
virtual bool initialized() const
libmesh_assert(ctx)
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ zero_clone()

template<typename T >
std::unique_ptr< SparseMatrix< T > > libMesh::PetscMatrix< T >::zero_clone ( ) const
overridevirtual
Returns
A smart pointer to a copy of this matrix with the same type, size, and partitioning, but with all zero entries.
Note
This must be overridden in the derived classes.

Implements libMesh::SparseMatrix< T >.

Definition at line 461 of file petsc_matrix.C.

462 {
463  libmesh_error_msg_if(!this->closed(), "Matrix must be closed before it can be cloned!");
464 
465  // Copy the nonzero pattern only
466  Mat copy;
467  PetscErrorCode ierr = MatDuplicate(_mat, MAT_DO_NOT_COPY_VALUES, &copy);
468  LIBMESH_CHKERR(ierr);
469 
470  // Call wrapping PetscMatrix constructor, have it take over
471  // ownership.
472  auto ret = std::make_unique<PetscMatrix<T>>(copy, this->comm());
473  ret->set_destroy_mat_on_exit(true);
474 
475  // Work around an issue on older compilers. We are able to simply
476  // "return ret;" on newer compilers
477  return std::unique_ptr<SparseMatrix<T>>(ret.release());
478 }
const Parallel::Communicator & comm() const
virtual bool closed() const override
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

◆ zero_rows()

template<typename T>
void libMesh::PetscMatrix< T >::zero_rows ( std::vector< numeric_index_type > &  rows,
diag_value = 0.0 
)
overridevirtual

Sets all row entries to 0 then puts diag_value in the diagonal entry.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 438 of file petsc_matrix.C.

439 {
440  libmesh_assert (this->initialized());
441 
442  semiparallel_only();
443 
444  PetscErrorCode ierr=0;
445 
446  // As of petsc-dev at the time of 3.1.0, MatZeroRows now takes two additional
447  // optional arguments. The optional arguments (x,b) can be used to specify the
448  // solutions for the zeroed rows (x) and right hand side (b) to update.
449  // Could be useful for setting boundary conditions...
450  if (!rows.empty())
451  ierr = MatZeroRows(_mat, cast_int<PetscInt>(rows.size()),
452  numeric_petsc_cast(rows.data()), PS(diag_value),
453  NULL, NULL);
454  else
455  ierr = MatZeroRows(_mat, 0, NULL, PS(diag_value), NULL, NULL);
456 
457  LIBMESH_CHKERR(ierr);
458 }
virtual bool initialized() const
PetscInt * numeric_petsc_cast(const numeric_index_type *p)
PetscScalar PS(T val)
Definition: petsc_macro.h:166
libmesh_assert(ctx)
Mat _mat
PETSc matrix datatype to store values.
Definition: petsc_matrix.h:367

Member Data Documentation

◆ _communicator

const Parallel::Communicator& libMesh::ParallelObject::_communicator
protectedinherited

◆ _counts

ReferenceCounter::Counts libMesh::ReferenceCounter::_counts
staticprotectedinherited

Actually holds the data.

Definition at line 124 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::get_info().

◆ _destroy_mat_on_exit

template<typename T>
bool libMesh::PetscMatrix< T >::_destroy_mat_on_exit
private

This boolean value should only be set to false for the constructor which takes a PETSc Mat object.

Definition at line 373 of file petsc_matrix.h.

Referenced by libMesh::PetscMatrix< libMesh::Number >::swap().

◆ _dof_map

template<typename T>
DofMap const* libMesh::SparseMatrix< T >::_dof_map
protectedinherited

The DofMap object associated with this object.

May be queried for degree-of-freedom counts on processors.

Definition at line 567 of file sparse_matrix.h.

◆ _enable_print_counter

bool libMesh::ReferenceCounter::_enable_print_counter = true
staticprotectedinherited

Flag to control whether reference count information is printed when print_info is called.

Definition at line 143 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::disable_print_counter_info(), libMesh::ReferenceCounter::enable_print_counter_info(), and libMesh::ReferenceCounter::print_info().

◆ _is_initialized

template<typename T>
bool libMesh::SparseMatrix< T >::_is_initialized
protectedinherited

◆ _mat

template<typename T>
Mat libMesh::PetscMatrix< T >::_mat
private

◆ _mat_type

template<typename T>
PetscMatrixType libMesh::PetscMatrix< T >::_mat_type
private

Definition at line 375 of file petsc_matrix.h.

◆ _mutex

Threads::spin_mutex libMesh::ReferenceCounter::_mutex
staticprotectedinherited

Mutual exclusion object to enable thread-safe reference counting.

Definition at line 137 of file reference_counter.h.

◆ _n_objects

Threads::atomic< unsigned int > libMesh::ReferenceCounter::_n_objects
staticprotectedinherited

The number of objects.

Print the reference count information when the number returns to 0.

Definition at line 132 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::n_objects(), libMesh::ReferenceCounter::ReferenceCounter(), and libMesh::ReferenceCounter::~ReferenceCounter().

◆ _petsc_matrix_mutex [1/2]

template<typename T>
std::mutex libMesh::PetscMatrix< T >::_petsc_matrix_mutex
mutableprivate

Definition at line 378 of file petsc_matrix.h.

◆ _petsc_matrix_mutex [2/2]

template<typename T>
Threads::spin_mutex libMesh::PetscMatrix< T >::_petsc_matrix_mutex
mutableprivate

Definition at line 380 of file petsc_matrix.h.

◆ _sp

template<typename T>
SparsityPattern::Build const* libMesh::SparseMatrix< T >::_sp
protectedinherited

The sparsity pattern associated with this object.

Should be queried for entry counts (or with need_full_sparsity_pattern, patterns) when needed.

Definition at line 574 of file sparse_matrix.h.


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