libMesh::PetscMatrix< T > Class Template Reference

SparseMatrix interface to PETSc Mat. More...

#include <petsc_linear_solver.h>

Inheritance diagram for libMesh::PetscMatrix< T >:

Public Member Functions

 PetscMatrix (const Parallel::Communicator &comm_in LIBMESH_CAN_DEFAULT_TO_COMMWORLD)
 
 PetscMatrix (Mat m, const Parallel::Communicator &comm_in LIBMESH_CAN_DEFAULT_TO_COMMWORLD)
 
 ~PetscMatrix ()
 
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) libmesh_override
 
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)
 
virtual void init () libmesh_override
 
void update_preallocation_and_zero ()
 
virtual void clear () libmesh_override
 
virtual void zero () libmesh_override
 
virtual void zero_rows (std::vector< numeric_index_type > &rows, T diag_value=0.0) libmesh_override
 
virtual void close () libmesh_override
 
virtual void flush () libmesh_override
 
virtual numeric_index_type m () const libmesh_override
 
virtual numeric_index_type n () const libmesh_override
 
virtual numeric_index_type row_start () const libmesh_override
 
virtual numeric_index_type row_stop () const libmesh_override
 
virtual void set (const numeric_index_type i, const numeric_index_type j, const T value) libmesh_override
 
virtual void add (const numeric_index_type i, const numeric_index_type j, const T value) libmesh_override
 
virtual void add_matrix (const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols) libmesh_override
 
virtual void add_matrix (const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &dof_indices) libmesh_override
 
virtual void add_block_matrix (const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &brows, const std::vector< numeric_index_type > &bcols) libmesh_override
 
virtual void add_block_matrix (const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &dof_indices) libmesh_override
 
virtual void add (const T a, SparseMatrix< T > &X) libmesh_override
 
virtual T operator() (const numeric_index_type i, const numeric_index_type j) const libmesh_override
 
virtual Real l1_norm () const libmesh_override
 
virtual Real linfty_norm () const libmesh_override
 
virtual bool closed () const libmesh_override
 
virtual void print_personal (std::ostream &os=libMesh::out) const libmesh_override
 
virtual void print_matlab (const std::string &name="") const libmesh_override
 
virtual void get_diagonal (NumericVector< T > &dest) const libmesh_override
 
virtual void get_transpose (SparseMatrix< T > &dest) const libmesh_override
 
void swap (PetscMatrix< T > &)
 
Mat mat ()
 
virtual bool initialized () const
 
void attach_dof_map (const DofMap &dof_map)
 
virtual bool need_full_sparsity_pattern () const
 
virtual void update_sparsity_pattern (const SparsityPattern::Graph &)
 
void print (std::ostream &os=libMesh::out, const bool sparse=false) const
 
template<>
void print (std::ostream &os, const bool sparse) const
 
virtual void create_submatrix (SparseMatrix< T > &submatrix, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols) const
 
virtual void reinit_submatrix (SparseMatrix< T > &submatrix, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols) const
 
void vector_mult (NumericVector< T > &dest, const NumericVector< T > &arg) const
 
void vector_mult_add (NumericVector< T > &dest, const NumericVector< T > &arg) const
 
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())
 
static std::string get_info ()
 
static void print_info (std::ostream &out=libMesh::out)
 
static unsigned int n_objects ()
 
static void enable_print_counter_info ()
 
static void disable_print_counter_info ()
 

Protected Types

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

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 libmesh_override
 
void increment_constructor_count (const std::string &name)
 
void increment_destructor_count (const std::string &name)
 

Protected Attributes

DofMap const * _dof_map
 
bool _is_initialized
 
const Parallel::Communicator_communicator
 

Static Protected Attributes

static Counts _counts
 
static Threads::atomic< unsigned int > _n_objects
 
static Threads::spin_mutex _mutex
 
static bool _enable_print_counter = true
 

Private Attributes

Mat _mat
 
bool _destroy_mat_on_exit
 

Detailed Description

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

SparseMatrix interface to PETSc Mat.

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

Definition at line 71 of file petsc_linear_solver.h.

Member Typedef Documentation

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

template<typename T>
libMesh::PetscMatrix< T >::PetscMatrix ( const Parallel::Communicator &comm_in  LIBMESH_CAN_DEFAULT_TO_COMMWORLD)
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(...).

template<typename T>
libMesh::PetscMatrix< T >::PetscMatrix ( Mat  m,
const Parallel::Communicator &comm_in  LIBMESH_CAN_DEFAULT_TO_COMMWORLD 
)
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.

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

Destructor. Free all memory, but do not release the memory of the sparsity structure.

Definition at line 113 of file petsc_matrix.C.

114 {
115  this->clear();
116 }
virtual void clear() libmesh_override
Definition: petsc_matrix.C:466

Member Function Documentation

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
protectedvirtual

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 789 of file petsc_matrix.C.

References libMesh::SparseMatrix< T >::_is_initialized, libMesh::PetscMatrix< T >::_mat, libMesh::SparseMatrix< T >::clear(), libMesh::PetscMatrix< T >::close(), libMesh::closed(), ierr, libMesh::SparseMatrix< T >::initialized(), libMesh::numeric_petsc_cast(), and PETSC_USE_POINTER.

Referenced by libMesh::PetscMatrix< T >::mat().

793 {
794  if (!this->closed())
795  {
796  libmesh_deprecated();
797  libmesh_warning("The matrix must be assembled before calling PetscMatrix::create_submatrix().\n"
798  "Please update your code, as this warning will become an error in a future release.");
799  const_cast<PetscMatrix<T> *>(this)->close();
800  }
801 
802  // Make sure the SparseMatrix passed in is really a PetscMatrix
803  PetscMatrix<T> * petsc_submatrix = cast_ptr<PetscMatrix<T> *>(&submatrix);
804 
805  // If we're not reusing submatrix and submatrix is already initialized
806  // then we need to clear it, otherwise we get a memory leak.
807  if (!reuse_submatrix && submatrix.initialized())
808  submatrix.clear();
809 
810  // Construct row and column index sets.
811  PetscErrorCode ierr=0;
812  IS isrow, iscol;
813 
814  ierr = ISCreateLibMesh(this->comm().get(),
815  rows.size(),
816  numeric_petsc_cast(&rows[0]),
818  &isrow); LIBMESH_CHKERR(ierr);
819 
820  ierr = ISCreateLibMesh(this->comm().get(),
821  cols.size(),
822  numeric_petsc_cast(&cols[0]),
824  &iscol); LIBMESH_CHKERR(ierr);
825 
826  // Extract submatrix
827  ierr = LibMeshCreateSubMatrix(_mat,
828  isrow,
829  iscol,
830 #if PETSC_RELEASE_LESS_THAN(3,0,1)
831  PETSC_DECIDE,
832 #endif
833  (reuse_submatrix ? MAT_REUSE_MATRIX : MAT_INITIAL_MATRIX),
834  &(petsc_submatrix->_mat)); LIBMESH_CHKERR(ierr);
835 
836  // Specify that the new submatrix is initialized and close it.
837  petsc_submatrix->_is_initialized = true;
838  petsc_submatrix->close();
839 
840  // Clean up PETSc data structures
841  ierr = LibMeshISDestroy(&isrow); LIBMESH_CHKERR(ierr);
842  ierr = LibMeshISDestroy(&iscol); LIBMESH_CHKERR(ierr);
843 }
unsigned int size() const
Definition: parallel.h:725
virtual bool closed() const libmesh_override
virtual void clear()=0
PetscInt * numeric_petsc_cast(const numeric_index_type *p)
virtual bool initialized() const
virtual void close() libmesh_override
Definition: petsc_matrix.C:897
PetscErrorCode ierr
const Parallel::Communicator & comm() const
SparseMatrix interface to PETSc Mat.
template<typename T >
void libMesh::PetscMatrix< T >::add ( const numeric_index_type  i,
const numeric_index_type  j,
const T  value 
)
virtual

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 1014 of file petsc_matrix.C.

References ierr, libMesh::initialized(), and value.

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

1017 {
1018  libmesh_assert (this->initialized());
1019 
1020  PetscErrorCode ierr=0;
1021  PetscInt i_val=i, j_val=j;
1022 
1023  PetscScalar petsc_value = static_cast<PetscScalar>(value);
1024  ierr = MatSetValues(_mat, 1, &i_val, 1, &j_val,
1025  &petsc_value, ADD_VALUES);
1026  LIBMESH_CHKERR(ierr);
1027 }
virtual bool initialized() const
PetscErrorCode ierr
static const bool value
Definition: xdr_io.C:108
template<typename T >
void libMesh::PetscMatrix< T >::add ( const T  a,
SparseMatrix< T > &  X 
)
virtual

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 1045 of file petsc_matrix.C.

References libMesh::PetscMatrix< T >::_mat, libMesh::PetscMatrix< T >::closed(), ierr, libMesh::initialized(), libMesh::SparseMatrix< T >::m(), and libMesh::SparseMatrix< T >::n().

1046 {
1047  libmesh_assert (this->initialized());
1048 
1049  // sanity check. but this cannot avoid
1050  // crash due to incompatible sparsity structure...
1051  libmesh_assert_equal_to (this->m(), X_in.m());
1052  libmesh_assert_equal_to (this->n(), X_in.n());
1053 
1054  PetscScalar a = static_cast<PetscScalar> (a_in);
1055  PetscMatrix<T> * X = cast_ptr<PetscMatrix<T> *> (&X_in);
1056 
1057  libmesh_assert (X);
1058 
1059  PetscErrorCode ierr=0;
1060 
1061  // the matrix from which we copy the values has to be assembled/closed
1062  libmesh_assert(X->closed());
1063 
1064  semiparallel_only();
1065 
1066  ierr = MatAXPY(_mat, a, X->_mat, DIFFERENT_NONZERO_PATTERN);
1067  LIBMESH_CHKERR(ierr);
1068 }
virtual numeric_index_type m() const libmesh_override
Definition: petsc_matrix.C:932
virtual bool closed() const libmesh_override
virtual bool initialized() const
virtual numeric_index_type n() const libmesh_override
Definition: petsc_matrix.C:948
PetscErrorCode ierr
SparseMatrix interface to PETSc Mat.
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 
)
virtual

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 745 of file petsc_matrix.C.

References libMesh::DenseMatrix< T >::get_values(), ierr, libMesh::initialized(), libMesh::DenseMatrixBase< T >::m(), libMesh::DenseMatrixBase< T >::n(), and libMesh::numeric_petsc_cast().

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

748 {
749  libmesh_assert (this->initialized());
750 
751  const numeric_index_type n_brows =
752  cast_int<numeric_index_type>(brows.size());
753  const numeric_index_type n_bcols =
754  cast_int<numeric_index_type>(bcols.size());
755 
756  PetscErrorCode ierr=0;
757 
758 #ifndef NDEBUG
759  const numeric_index_type n_rows =
760  cast_int<numeric_index_type>(dm.m());
761  const numeric_index_type n_cols =
762  cast_int<numeric_index_type>(dm.n());
763  const numeric_index_type blocksize = n_rows / n_brows;
764 
765  libmesh_assert_equal_to (n_cols / n_bcols, blocksize);
766  libmesh_assert_equal_to (blocksize*n_brows, n_rows);
767  libmesh_assert_equal_to (blocksize*n_bcols, n_cols);
768 
769  PetscInt petsc_blocksize;
770  ierr = MatGetBlockSize(_mat, &petsc_blocksize);
771  LIBMESH_CHKERR(ierr);
772  libmesh_assert_equal_to (blocksize, static_cast<numeric_index_type>(petsc_blocksize));
773 #endif
774 
775  // These casts are required for PETSc <= 2.1.5
776  ierr = MatSetValuesBlocked(_mat,
777  n_brows, numeric_petsc_cast(&brows[0]),
778  n_bcols, numeric_petsc_cast(&bcols[0]),
779  const_cast<PetscScalar *>(&dm.get_values()[0]),
780  ADD_VALUES);
781  LIBMESH_CHKERR(ierr);
782 }
unsigned int n() const
PetscInt * numeric_petsc_cast(const numeric_index_type *p)
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:92
unsigned int m() const
std::vector< T > & get_values()
Definition: dense_matrix.h:325
PetscErrorCode ierr
template<typename T>
virtual void libMesh::PetscMatrix< T >::add_block_matrix ( const DenseMatrix< T > &  dm,
const std::vector< numeric_index_type > &  dof_indices 
)
inlinevirtual

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 182 of file petsc_matrix.h.

References libMesh::PetscMatrix< T >::add(), libMesh::PetscMatrix< T >::add_block_matrix(), libMesh::PetscMatrix< T >::closed(), libMesh::PetscMatrix< T >::get_diagonal(), libMesh::PetscMatrix< T >::get_transpose(), libMesh::PetscMatrix< T >::l1_norm(), libMesh::PetscMatrix< T >::linfty_norm(), libMesh::Quality::name(), libMesh::PetscMatrix< T >::operator()(), libMesh::out, libMesh::PetscMatrix< T >::print_matlab(), libMesh::PetscMatrix< T >::print_personal(), libMesh::Real, and libMesh::PetscMatrix< T >::swap().

184  { 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) libmesh_override
Definition: petsc_matrix.C:745
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 
)
virtual

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 716 of file petsc_matrix.C.

References libMesh::DenseMatrix< T >::get_values(), ierr, libMesh::initialized(), libMesh::DenseMatrixBase< T >::m(), libMesh::DenseMatrixBase< T >::n(), and libMesh::numeric_petsc_cast().

719 {
720  libmesh_assert (this->initialized());
721 
722  const numeric_index_type n_rows = dm.m();
723  const numeric_index_type n_cols = dm.n();
724 
725  libmesh_assert_equal_to (rows.size(), n_rows);
726  libmesh_assert_equal_to (cols.size(), n_cols);
727 
728  PetscErrorCode ierr=0;
729 
730  // These casts are required for PETSc <= 2.1.5
731  ierr = MatSetValues(_mat,
732  n_rows, numeric_petsc_cast(&rows[0]),
733  n_cols, numeric_petsc_cast(&cols[0]),
734  const_cast<PetscScalar *>(&dm.get_values()[0]),
735  ADD_VALUES);
736  LIBMESH_CHKERR(ierr);
737 }
unsigned int n() const
PetscInt * numeric_petsc_cast(const numeric_index_type *p)
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:92
unsigned int m() const
std::vector< T > & get_values()
Definition: dense_matrix.h:325
PetscErrorCode ierr
template<typename T >
void libMesh::PetscMatrix< T >::add_matrix ( const DenseMatrix< T > &  dm,
const std::vector< numeric_index_type > &  dof_indices 
)
virtual

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 1032 of file petsc_matrix.C.

1034 {
1035  this->add_matrix (dm, dof_indices, dof_indices);
1036 }
virtual void add_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols) libmesh_override
Definition: petsc_matrix.C:716
template<typename T>
void libMesh::SparseMatrix< T >::attach_dof_map ( const DofMap dof_map)
inlineinherited

Get a pointer to the DofMap to use.

Definition at line 106 of file sparse_matrix.h.

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

107  { _dof_map = &dof_map; }
DofMap const * _dof_map
template<typename T >
std::unique_ptr< SparseMatrix< T > > libMesh::SparseMatrix< T >::build ( const Parallel::Communicator comm,
const SolverPackage  solver_package = libMesh::default_solver_package() 
)
staticinherited

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

Definition at line 136 of file sparse_matrix.C.

References libMesh::ParallelObject::comm(), libMesh::EIGEN_SOLVERS, libMesh::LASPACK_SOLVERS, libMesh::libmesh_ignore(), libMesh::PETSC_SOLVERS, and libMesh::TRILINOS_SOLVERS.

Referenced by libMesh::ImplicitSystem::add_matrix(), libMesh::EigenSystem::init_data(), and libMesh::EigenSystem::init_matrices().

138 {
139  // Avoid unused parameter warnings when no solver packages are enabled.
141 
142  // Build the appropriate vector
143  switch (solver_package)
144  {
145 
146 #ifdef LIBMESH_HAVE_LASPACK
147  case LASPACK_SOLVERS:
148  return libmesh_make_unique<LaspackMatrix<T>>(comm);
149 #endif
150 
151 
152 #ifdef LIBMESH_HAVE_PETSC
153  case PETSC_SOLVERS:
154  return libmesh_make_unique<PetscMatrix<T>>(comm);
155 #endif
156 
157 
158 #ifdef LIBMESH_TRILINOS_HAVE_EPETRA
159  case TRILINOS_SOLVERS:
160  return libmesh_make_unique<EpetraMatrix<T>>(comm);
161 #endif
162 
163 
164 #ifdef LIBMESH_HAVE_EIGEN
165  case EIGEN_SOLVERS:
166  return libmesh_make_unique<EigenSparseMatrix<T>>(comm);
167 #endif
168 
169  default:
170  libmesh_error_msg("ERROR: Unrecognized solver package: " << solver_package);
171  }
172 }
EIGEN_SOLVERS
Definition: libmesh.C:261
TRILINOS_SOLVERS
Definition: libmesh.C:259
LASPACK_SOLVERS
Definition: libmesh.C:263
void libmesh_ignore(const T &)
const Parallel::Communicator & comm() const
template<typename T >
void libMesh::PetscMatrix< T >::clear ( )
virtual

Restores the SparseMatrix<T> to a pristine state.

Implements libMesh::SparseMatrix< T >.

Definition at line 466 of file petsc_matrix.C.

References libMesh::libMeshPrivateData::_is_initialized, ierr, and libMesh::initialized().

467 {
468  PetscErrorCode ierr=0;
469 
470  if ((this->initialized()) && (this->_destroy_mat_on_exit))
471  {
472  semiparallel_only();
473 
474  ierr = LibMeshMatDestroy (&_mat);
475  LIBMESH_CHKERR(ierr);
476 
477  this->_is_initialized = false;
478  }
479 }
virtual bool initialized() const
PetscErrorCode ierr
template<typename T >
void libMesh::PetscMatrix< T >::close ( )
virtual

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

Implements libMesh::SparseMatrix< T >.

Definition at line 897 of file petsc_matrix.C.

References ierr.

Referenced by libMesh::PetscMatrix< T >::_get_submatrix(), libMesh::PetscLinearSolver< T >::adjoint_solve(), libMesh::PetscMatrix< T >::get_transpose(), libMesh::PetscLinearSolver< T >::solve(), libMesh::SlepcEigenSolver< T >::solve_generalized(), and libMesh::SlepcEigenSolver< T >::solve_standard().

898 {
899  semiparallel_only();
900 
901  // BSK - 1/19/2004
902  // strictly this check should be OK, but it seems to
903  // fail on matrix-free matrices. Do they falsely
904  // state they are assembled? Check with the developers...
905  // if (this->closed())
906  // return;
907 
908  PetscErrorCode ierr=0;
909 
910  ierr = MatAssemblyBegin (_mat, MAT_FINAL_ASSEMBLY);
911  LIBMESH_CHKERR(ierr);
912  ierr = MatAssemblyEnd (_mat, MAT_FINAL_ASSEMBLY);
913  LIBMESH_CHKERR(ierr);
914 }
PetscErrorCode ierr
template<typename T >
bool libMesh::PetscMatrix< T >::closed ( ) const
virtual
Returns
true if the matrix has been assembled.

Implements libMesh::SparseMatrix< T >.

Definition at line 1134 of file petsc_matrix.C.

References ierr, and libMesh::initialized().

Referenced by libMesh::PetscMatrix< T >::add(), and libMesh::PetscMatrix< T >::add_block_matrix().

1135 {
1136  libmesh_assert (this->initialized());
1137 
1138  PetscErrorCode ierr=0;
1139  PetscBool assembled;
1140 
1141  ierr = MatAssembled(_mat, &assembled);
1142  LIBMESH_CHKERR(ierr);
1143 
1144  return (assembled == PETSC_TRUE);
1145 }
virtual bool initialized() const
PetscTruth PetscBool
Definition: petsc_macro.h:64
PetscErrorCode ierr
const Parallel::Communicator& libMesh::ParallelObject::comm ( ) const
inlineinherited
Returns
A reference to the Parallel::Communicator object used by this mesh.

Definition at line 87 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_petsc_snes_fd_residual(), libMesh::__libmesh_petsc_snes_jacobian(), libMesh::__libmesh_petsc_snes_mffd_residual(), libMesh::__libmesh_petsc_snes_postcheck(), libMesh::__libmesh_petsc_snes_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::ParmetisPartitioner::_do_repartition(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult_add(), libMesh::EquationSystems::_read_impl(), libMesh::MeshRefinement::_refine_elements(), libMesh::MeshRefinement::_smooth_flags(), libMesh::ImplicitSystem::add_matrix(), libMesh::System::add_vector(), libMesh::EigenSparseLinearSolver< T >::adjoint_solve(), libMesh::UnstructuredMesh::all_second_order(), libMesh::MeshTools::Modification::all_tri(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assemble_qoi(), libMesh::MeshCommunication::assign_global_indices(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::DofMap::attach_matrix(), libMesh::Parallel::BinSorter< KeyType, IdxType >::binsort(), libMesh::Parallel::Sort< KeyType, IdxType >::binsort(), libMesh::MeshCommunication::broadcast(), libMesh::SparseMatrix< T >::build(), libMesh::MeshTools::Generation::build_extrusion(), libMesh::Parallel::Histogram< KeyType, IdxType >::build_histogram(), libMesh::PetscNonlinearSolver< T >::build_mat_null_space(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_solution_vector(), libMesh::MeshBase::cache_elem_dims(), libMesh::System::calculate_norm(), libMesh::DofMap::check_dirichlet_bcid_consistency(), libMesh::DistributedVector< T >::clone(), libMesh::EigenSparseVector< T >::clone(), libMesh::LaspackVector< T >::clone(), libMesh::EpetraVector< T >::clone(), libMesh::PetscVector< T >::clone(), libMesh::EpetraVector< T >::close(), libMesh::Parallel::Sort< KeyType, IdxType >::communicate_bins(), 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::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::ExodusII_IO::copy_elemental_solution(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::MeshTools::create_bounding_box(), 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::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), DMlibMeshFunction(), DMlibMeshJacobian(), DMlibMeshSetSystem_libMesh(), DMVariableBounds_libMesh(), libMesh::MeshRefinement::eliminate_unrefined_patches(), libMesh::EpetraVector< T >::EpetraVector(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::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_nelem_target(), libMesh::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::DofMap::get_info(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::EquationSystems::get_solution(), libMesh::LocationMap< T >::init(), libMesh::TopologyMap::init(), libMesh::TimeSolver::init(), libMesh::TaoOptimizationSolver< T >::init(), libMesh::PetscNonlinearSolver< T >::init(), libMesh::DistributedVector< T >::init(), libMesh::EpetraVector< T >::init(), libMesh::PetscVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::EigenSystem::init_data(), libMesh::EigenSystem::init_matrices(), libMesh::ParmetisPartitioner::initialize(), libMesh::OptimizationSystem::initialize_equality_constraints_storage(), libMesh::OptimizationSystem::initialize_inequality_constraints_storage(), 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_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::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_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::DistributedVector< T >::max(), libMesh::FEMSystem::mesh_position_set(), libMesh::MeshSerializer::MeshSerializer(), libMesh::DistributedVector< T >::min(), libMesh::DistributedMesh::n_active_elem(), libMesh::MeshTools::n_active_levels(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::CondensedEigenSystem::n_global_non_condensed_dofs(), libMesh::MeshTools::n_levels(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::MeshTools::n_p_levels(), libMesh::BoundaryInfo::n_shellface_conds(), 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::MeshTools::paranoid_n_levels(), libMesh::Partitioner::partition(), libMesh::LinearPartitioner::partition_range(), libMesh::MetisPartitioner::partition_range(), 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::SparseMatrix< T >::print(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::MeshBase::recalculate_n_partitions(), libMesh::MeshCommunication::redistribute(), libMesh::MeshRefinement::refine_and_coarsen_elements(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::CheckpointIO::select_split_config(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::Partitioner::set_parent_processor_ids(), libMesh::PetscDiffSolver::setup_petsc_data(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::Parallel::Sort< KeyType, IdxType >::sort(), libMesh::split_mesh(), libMesh::MeshBase::subdomain_ids(), libMesh::BoundaryInfo::sync(), libMesh::Parallel::sync_element_data_by_parent_id(), libMesh::Parallel::sync_node_data_by_element_id(), libMesh::MeshRefinement::test_level_one(), libMesh::MeshRefinement::test_unflagged(), libMesh::MeshTools::total_weight(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::DistributedVector< T >::zero_clone(), libMesh::EigenSparseVector< T >::zero_clone(), libMesh::LaspackVector< T >::zero_clone(), libMesh::EpetraVector< T >::zero_clone(), and libMesh::PetscVector< T >::zero_clone().

88  { return _communicator; }
const Parallel::Communicator & _communicator
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.

Definition at line 351 of file sparse_matrix.h.

354  {
355  this->_get_submatrix(submatrix,
356  rows,
357  cols,
358  false); // false means DO NOT REUSE submatrix
359  }
virtual void _get_submatrix(SparseMatrix< T > &, const std::vector< numeric_index_type > &, const std::vector< numeric_index_type > &, const bool) const
void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited
void libMesh::ReferenceCounter::enable_print_counter_info ( )
staticinherited

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

Definition at line 101 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

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

102 {
103  _enable_print_counter = true;
104  return;
105 }
template<typename T >
void libMesh::PetscMatrix< T >::flush ( )
virtual

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 917 of file petsc_matrix.C.

References ierr.

918 {
919  semiparallel_only();
920 
921  PetscErrorCode ierr=0;
922 
923  ierr = MatAssemblyBegin (_mat, MAT_FLUSH_ASSEMBLY);
924  LIBMESH_CHKERR(ierr);
925  ierr = MatAssemblyEnd (_mat, MAT_FLUSH_ASSEMBLY);
926  LIBMESH_CHKERR(ierr);
927 }
PetscErrorCode ierr
template<typename T >
void libMesh::PetscMatrix< T >::get_diagonal ( NumericVector< T > &  dest) const
virtual

Copies the diagonal part of the matrix into dest.

Implements libMesh::SparseMatrix< T >.

Definition at line 848 of file petsc_matrix.C.

References ierr, and libMesh::PetscVector< T >::vec().

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

849 {
850  // Make sure the NumericVector passed in is really a PetscVector
851  PetscVector<T> & petsc_dest = cast_ref<PetscVector<T> &>(dest);
852 
853  // Needs a const_cast since PETSc does not work with const.
854  PetscErrorCode ierr =
855  MatGetDiagonal(const_cast<PetscMatrix<T> *>(this)->mat(),petsc_dest.vec()); LIBMESH_CHKERR(ierr);
856 }
NumericVector interface to PETSc Vec.
Definition: petsc_vector.h:64
PetscErrorCode ierr
SparseMatrix interface to PETSc Mat.
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 (Counts::iterator it = _counts.begin();
59  it != _counts.end(); ++it)
60  {
61  const std::string name(it->first);
62  const unsigned int creations = it->second.first;
63  const unsigned int destructions = it->second.second;
64 
65  oss << "| " << name << " reference count information:\n"
66  << "| Creations: " << creations << '\n'
67  << "| Destructions: " << destructions << '\n';
68  }
69 
70  oss << " ---------------------------------------------------------------------------- \n";
71 
72  return oss.str();
73 
74 #else
75 
76  return "";
77 
78 #endif
79 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
template<typename T >
void libMesh::PetscMatrix< T >::get_transpose ( SparseMatrix< T > &  dest) const
virtual

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

Implements libMesh::SparseMatrix< T >.

Definition at line 861 of file petsc_matrix.C.

References libMesh::SparseMatrix< T >::_is_initialized, libMesh::PetscMatrix< T >::_mat, libMesh::SparseMatrix< T >::clear(), libMesh::PetscMatrix< T >::close(), and ierr.

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

862 {
863  // Make sure the SparseMatrix passed in is really a PetscMatrix
864  PetscMatrix<T> & petsc_dest = cast_ref<PetscMatrix<T> &>(dest);
865 
866  // If we aren't reusing the matrix then need to clear dest,
867  // otherwise we get a memory leak
868  if (&petsc_dest != this)
869  dest.clear();
870 
871  PetscErrorCode ierr;
872 #if PETSC_VERSION_LESS_THAN(3,0,0)
873  if (&petsc_dest == this)
874  ierr = MatTranspose(_mat,PETSC_NULL);
875  else
876  ierr = MatTranspose(_mat,&petsc_dest._mat);
877  LIBMESH_CHKERR(ierr);
878 #else
879  // FIXME - we can probably use MAT_REUSE_MATRIX in more situations
880  if (&petsc_dest == this)
881  ierr = MatTranspose(_mat,MAT_REUSE_MATRIX,&petsc_dest._mat);
882  else
883  ierr = MatTranspose(_mat,MAT_INITIAL_MATRIX,&petsc_dest._mat);
884  LIBMESH_CHKERR(ierr);
885 #endif
886 
887  // Specify that the transposed matrix is initialized and close it.
888  petsc_dest._is_initialized = true;
889  petsc_dest.close();
890 }
virtual void clear()=0
virtual void close() libmesh_override
Definition: petsc_matrix.C:897
PetscErrorCode ierr
SparseMatrix interface to PETSc Mat.
void libMesh::ReferenceCounter::increment_constructor_count ( const std::string &  name)
inlineprotectedinherited

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

Definition at line 185 of file reference_counter.h.

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

Referenced by libMesh::ReferenceCounter::n_objects(), and libMesh::ReferenceCountedObject< RBParametrized >::ReferenceCountedObject().

186 {
187  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
188  std::pair<unsigned int, unsigned int> & p = _counts[name];
189 
190  p.first++;
191 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
spin_mutex spin_mtx
Definition: threads.C:29
void libMesh::ReferenceCounter::increment_destructor_count ( const std::string &  name)
inlineprotectedinherited

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

Definition at line 198 of file reference_counter.h.

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

Referenced by libMesh::ReferenceCounter::n_objects(), and libMesh::ReferenceCountedObject< RBParametrized >::~ReferenceCountedObject().

199 {
200  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
201  std::pair<unsigned int, unsigned int> & p = _counts[name];
202 
203  p.second++;
204 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
spin_mutex spin_mtx
Definition: threads.C:29
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 
)
virtual

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 120 of file petsc_matrix.C.

References libMesh::libMeshPrivateData::_is_initialized, ierr, libMesh::initialized(), libMesh::libmesh_ignore(), and libMesh::zero.

127 {
128  // So compilers don't warn when !LIBMESH_ENABLE_BLOCKED_STORAGE
129  libmesh_ignore(blocksize_in);
130 
131  // Clear initialized matrices
132  if (this->initialized())
133  this->clear();
134 
135  this->_is_initialized = true;
136 
137 
138  PetscErrorCode ierr = 0;
139  PetscInt m_global = static_cast<PetscInt>(m_in);
140  PetscInt n_global = static_cast<PetscInt>(n_in);
141  PetscInt m_local = static_cast<PetscInt>(m_l);
142  PetscInt n_local = static_cast<PetscInt>(n_l);
143  PetscInt n_nz = static_cast<PetscInt>(nnz);
144  PetscInt n_oz = static_cast<PetscInt>(noz);
145 
146  ierr = MatCreate(this->comm().get(), &_mat);
147  LIBMESH_CHKERR(ierr);
148  ierr = MatSetSizes(_mat, m_local, n_local, m_global, n_global);
149  LIBMESH_CHKERR(ierr);
150 
151 #ifdef LIBMESH_ENABLE_BLOCKED_STORAGE
152  PetscInt blocksize = static_cast<PetscInt>(blocksize_in);
153  if (blocksize > 1)
154  {
155  // specified blocksize, bs>1.
156  // double check sizes.
157  libmesh_assert_equal_to (m_local % blocksize, 0);
158  libmesh_assert_equal_to (n_local % blocksize, 0);
159  libmesh_assert_equal_to (m_global % blocksize, 0);
160  libmesh_assert_equal_to (n_global % blocksize, 0);
161  libmesh_assert_equal_to (n_nz % blocksize, 0);
162  libmesh_assert_equal_to (n_oz % blocksize, 0);
163 
164  ierr = MatSetType(_mat, MATBAIJ); // Automatically chooses seqbaij or mpibaij
165  LIBMESH_CHKERR(ierr);
166  ierr = MatSetBlockSize(_mat, blocksize);
167  LIBMESH_CHKERR(ierr);
168  ierr = MatSeqBAIJSetPreallocation(_mat, blocksize, n_nz/blocksize, PETSC_NULL);
169  LIBMESH_CHKERR(ierr);
170  ierr = MatMPIBAIJSetPreallocation(_mat, blocksize,
171  n_nz/blocksize, PETSC_NULL,
172  n_oz/blocksize, PETSC_NULL);
173  LIBMESH_CHKERR(ierr);
174  }
175  else
176 #endif
177  {
178  ierr = MatSetType(_mat, MATAIJ); // Automatically chooses seqaij or mpiaij
179  LIBMESH_CHKERR(ierr);
180  ierr = MatSeqAIJSetPreallocation(_mat, n_nz, PETSC_NULL);
181  LIBMESH_CHKERR(ierr);
182  ierr = MatMPIAIJSetPreallocation(_mat, n_nz, PETSC_NULL, n_oz, PETSC_NULL);
183  LIBMESH_CHKERR(ierr);
184  }
185 
186  // Make it an error for PETSc to allocate new nonzero entries during assembly
187 #if PETSC_VERSION_LESS_THAN(3,0,0)
188  ierr = MatSetOption(_mat, MAT_NEW_NONZERO_ALLOCATION_ERR);
189 #else
190  ierr = MatSetOption(_mat, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
191 #endif
192  LIBMESH_CHKERR(ierr);
193 
194  // Is prefix information available somewhere? Perhaps pass in the system name?
195  ierr = MatSetOptionsPrefix(_mat, "");
196  LIBMESH_CHKERR(ierr);
197  ierr = MatSetFromOptions(_mat);
198  LIBMESH_CHKERR(ierr);
199 
200  this->zero ();
201 }
virtual bool initialized() const
virtual void clear() libmesh_override
Definition: petsc_matrix.C:466
void libmesh_ignore(const T &)
PetscErrorCode ierr
const Parallel::Communicator & comm() const
virtual void zero() libmesh_override
Definition: petsc_matrix.C:413
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 205 of file petsc_matrix.C.

References libMesh::libMeshPrivateData::_is_initialized, ierr, libMesh::initialized(), libMesh::libmesh_ignore(), libmesh_nullptr, libMesh::numeric_petsc_cast(), and libMesh::zero.

212 {
213  // So compilers don't warn when !LIBMESH_ENABLE_BLOCKED_STORAGE
214  libmesh_ignore(blocksize_in);
215 
216  // Clear initialized matrices
217  if (this->initialized())
218  this->clear();
219 
220  this->_is_initialized = true;
221 
222  // Make sure the sparsity pattern isn't empty unless the matrix is 0x0
223  libmesh_assert_equal_to (n_nz.size(), m_l);
224  libmesh_assert_equal_to (n_oz.size(), m_l);
225 
226  PetscErrorCode ierr = 0;
227  PetscInt m_global = static_cast<PetscInt>(m_in);
228  PetscInt n_global = static_cast<PetscInt>(n_in);
229  PetscInt m_local = static_cast<PetscInt>(m_l);
230  PetscInt n_local = static_cast<PetscInt>(n_l);
231 
232  ierr = MatCreate(this->comm().get(), &_mat);
233  LIBMESH_CHKERR(ierr);
234  ierr = MatSetSizes(_mat, m_local, n_local, m_global, n_global);
235  LIBMESH_CHKERR(ierr);
236 
237 #ifdef LIBMESH_ENABLE_BLOCKED_STORAGE
238  PetscInt blocksize = static_cast<PetscInt>(blocksize_in);
239  if (blocksize > 1)
240  {
241  // specified blocksize, bs>1.
242  // double check sizes.
243  libmesh_assert_equal_to (m_local % blocksize, 0);
244  libmesh_assert_equal_to (n_local % blocksize, 0);
245  libmesh_assert_equal_to (m_global % blocksize, 0);
246  libmesh_assert_equal_to (n_global % blocksize, 0);
247 
248  ierr = MatSetType(_mat, MATBAIJ); // Automatically chooses seqbaij or mpibaij
249  LIBMESH_CHKERR(ierr);
250  ierr = MatSetBlockSize(_mat, blocksize);
251  LIBMESH_CHKERR(ierr);
252 
253  // transform the per-entry n_nz and n_oz arrays into their block counterparts.
254  std::vector<numeric_index_type> b_n_nz, b_n_oz;
255 
256  transform_preallocation_arrays (blocksize,
257  n_nz, n_oz,
258  b_n_nz, b_n_oz);
259 
260  ierr = MatSeqBAIJSetPreallocation (_mat,
261  blocksize,
262  0,
263  numeric_petsc_cast(b_n_nz.empty() ? libmesh_nullptr : &b_n_nz[0]));
264  LIBMESH_CHKERR(ierr);
265 
266  ierr = MatMPIBAIJSetPreallocation (_mat,
267  blocksize,
268  0,
269  numeric_petsc_cast(b_n_nz.empty() ? libmesh_nullptr : &b_n_nz[0]),
270  0,
271  numeric_petsc_cast(b_n_oz.empty() ? libmesh_nullptr : &b_n_oz[0]));
272  LIBMESH_CHKERR(ierr);
273  }
274  else
275 #endif
276  {
277 
278  ierr = MatSetType(_mat, MATAIJ); // Automatically chooses seqaij or mpiaij
279  LIBMESH_CHKERR(ierr);
280  ierr = MatSeqAIJSetPreallocation (_mat,
281  0,
282  numeric_petsc_cast(n_nz.empty() ? libmesh_nullptr : &n_nz[0]));
283  LIBMESH_CHKERR(ierr);
284  ierr = MatMPIAIJSetPreallocation (_mat,
285  0,
286  numeric_petsc_cast(n_nz.empty() ? libmesh_nullptr : &n_nz[0]),
287  0,
288  numeric_petsc_cast(n_oz.empty() ? libmesh_nullptr : &n_oz[0]));
289  LIBMESH_CHKERR(ierr);
290  }
291 
292  // Is prefix information available somewhere? Perhaps pass in the system name?
293  ierr = MatSetOptionsPrefix(_mat, "");
294  LIBMESH_CHKERR(ierr);
295  ierr = MatSetFromOptions(_mat);
296  LIBMESH_CHKERR(ierr);
297 
298 
299  this->zero();
300 }
const class libmesh_nullptr_t libmesh_nullptr
PetscInt * numeric_petsc_cast(const numeric_index_type *p)
virtual bool initialized() const
virtual void clear() libmesh_override
Definition: petsc_matrix.C:466
void libmesh_ignore(const T &)
PetscErrorCode ierr
const Parallel::Communicator & comm() const
virtual void zero() libmesh_override
Definition: petsc_matrix.C:413
template<typename T >
void libMesh::PetscMatrix< T >::init ( )
virtual

Initialize this matrix using the sparsity structure computed by dof_map.

Implements libMesh::SparseMatrix< T >.

Definition at line 304 of file petsc_matrix.C.

References libMesh::libMeshPrivateData::_is_initialized, ierr, libMesh::initialized(), libmesh_nullptr, libMesh::numeric_petsc_cast(), libMesh::processor_id(), and libMesh::zero.

305 {
306  libmesh_assert(this->_dof_map);
307 
308  // Clear initialized matrices
309  if (this->initialized())
310  this->clear();
311 
312  this->_is_initialized = true;
313 
314 
315  const numeric_index_type my_m = this->_dof_map->n_dofs();
316  const numeric_index_type my_n = my_m;
317  const numeric_index_type n_l = this->_dof_map->n_dofs_on_processor(this->processor_id());
318  const numeric_index_type m_l = n_l;
319 
320 
321  const std::vector<numeric_index_type> & n_nz = this->_dof_map->get_n_nz();
322  const std::vector<numeric_index_type> & n_oz = this->_dof_map->get_n_oz();
323 
324  // Make sure the sparsity pattern isn't empty unless the matrix is 0x0
325  libmesh_assert_equal_to (n_nz.size(), m_l);
326  libmesh_assert_equal_to (n_oz.size(), m_l);
327 
328  PetscErrorCode ierr = 0;
329  PetscInt m_global = static_cast<PetscInt>(my_m);
330  PetscInt n_global = static_cast<PetscInt>(my_n);
331  PetscInt m_local = static_cast<PetscInt>(m_l);
332  PetscInt n_local = static_cast<PetscInt>(n_l);
333 
334  ierr = MatCreate(this->comm().get(), &_mat);
335  LIBMESH_CHKERR(ierr);
336  ierr = MatSetSizes(_mat, m_local, n_local, m_global, n_global);
337  LIBMESH_CHKERR(ierr);
338 
339 #ifdef LIBMESH_ENABLE_BLOCKED_STORAGE
340  PetscInt blocksize = static_cast<PetscInt>(this->_dof_map->block_size());
341  if (blocksize > 1)
342  {
343  // specified blocksize, bs>1.
344  // double check sizes.
345  libmesh_assert_equal_to (m_local % blocksize, 0);
346  libmesh_assert_equal_to (n_local % blocksize, 0);
347  libmesh_assert_equal_to (m_global % blocksize, 0);
348  libmesh_assert_equal_to (n_global % blocksize, 0);
349 
350  ierr = MatSetType(_mat, MATBAIJ); // Automatically chooses seqbaij or mpibaij
351  LIBMESH_CHKERR(ierr);
352  ierr = MatSetBlockSize(_mat, blocksize);
353  LIBMESH_CHKERR(ierr);
354 
355  // transform the per-entry n_nz and n_oz arrays into their block counterparts.
356  std::vector<numeric_index_type> b_n_nz, b_n_oz;
357 
358  transform_preallocation_arrays (blocksize,
359  n_nz, n_oz,
360  b_n_nz, b_n_oz);
361 
362  ierr = MatSeqBAIJSetPreallocation (_mat,
363  blocksize,
364  0,
365  numeric_petsc_cast(b_n_nz.empty() ? libmesh_nullptr : &b_n_nz[0]));
366  LIBMESH_CHKERR(ierr);
367 
368  ierr = MatMPIBAIJSetPreallocation (_mat,
369  blocksize,
370  0,
371  numeric_petsc_cast(b_n_nz.empty() ? libmesh_nullptr : &b_n_nz[0]),
372  0,
373  numeric_petsc_cast(b_n_oz.empty() ? libmesh_nullptr : &b_n_oz[0]));
374  LIBMESH_CHKERR(ierr);
375  }
376  else
377 #endif
378  {
379  // no block storage case
380  ierr = MatSetType(_mat, MATAIJ); // Automatically chooses seqaij or mpiaij
381  LIBMESH_CHKERR(ierr);
382 
383  ierr = MatSeqAIJSetPreallocation (_mat,
384  0,
385  numeric_petsc_cast(n_nz.empty() ? libmesh_nullptr : &n_nz[0]));
386  LIBMESH_CHKERR(ierr);
387  ierr = MatMPIAIJSetPreallocation (_mat,
388  0,
389  numeric_petsc_cast(n_nz.empty() ? libmesh_nullptr : &n_nz[0]),
390  0,
391  numeric_petsc_cast(n_oz.empty() ? libmesh_nullptr : &n_oz[0]));
392  LIBMESH_CHKERR(ierr);
393  }
394 
395  // Is prefix information available somewhere? Perhaps pass in the system name?
396  ierr = MatSetOptionsPrefix(_mat, "");
397  LIBMESH_CHKERR(ierr);
398  ierr = MatSetFromOptions(_mat);
399  LIBMESH_CHKERR(ierr);
400 
401  this->zero();
402 }
const class libmesh_nullptr_t libmesh_nullptr
unsigned int block_size() const
Definition: dof_map.h:507
dof_id_type n_dofs() const
Definition: dof_map.h:519
dof_id_type n_dofs_on_processor(const processor_id_type proc) const
Definition: dof_map.h:535
PetscInt * numeric_petsc_cast(const numeric_index_type *p)
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:92
virtual void clear() libmesh_override
Definition: petsc_matrix.C:466
const std::vector< dof_id_type > & get_n_oz() const
Definition: dof_map.h:417
DofMap const * _dof_map
PetscErrorCode ierr
const Parallel::Communicator & comm() const
virtual void zero() libmesh_override
Definition: petsc_matrix.C:413
const std::vector< dof_id_type > & get_n_nz() const
Definition: dof_map.h:404
processor_id_type processor_id() const
template<typename T >
Real libMesh::PetscMatrix< T >::l1_norm ( ) const
virtual
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 484 of file petsc_matrix.C.

References libMesh::closed(), ierr, libMesh::initialized(), libMesh::Real, and value.

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

485 {
486  libmesh_assert (this->initialized());
487 
488  semiparallel_only();
489 
490  PetscErrorCode ierr=0;
491  PetscReal petsc_value;
492  Real value;
493 
494  libmesh_assert (this->closed());
495 
496  ierr = MatNorm(_mat, NORM_1, &petsc_value);
497  LIBMESH_CHKERR(ierr);
498 
499  value = static_cast<Real>(petsc_value);
500 
501  return value;
502 }
virtual bool closed() const libmesh_override
virtual bool initialized() const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
PetscErrorCode ierr
static const bool value
Definition: xdr_io.C:108
template<typename T >
Real libMesh::PetscMatrix< T >::linfty_norm ( ) const
virtual
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 507 of file petsc_matrix.C.

References libMesh::closed(), ierr, libMesh::initialized(), libMesh::Real, and value.

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

508 {
509  libmesh_assert (this->initialized());
510 
511  semiparallel_only();
512 
513  PetscErrorCode ierr=0;
514  PetscReal petsc_value;
515  Real value;
516 
517  libmesh_assert (this->closed());
518 
519  ierr = MatNorm(_mat, NORM_INFINITY, &petsc_value);
520  LIBMESH_CHKERR(ierr);
521 
522  value = static_cast<Real>(petsc_value);
523 
524  return value;
525 }
virtual bool closed() const libmesh_override
virtual bool initialized() const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
PetscErrorCode ierr
static const bool value
Definition: xdr_io.C:108
template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::m ( ) const
virtual
Returns
The row-dimension of the matrix.

Implements libMesh::SparseMatrix< T >.

Definition at line 932 of file petsc_matrix.C.

References ierr, and libMesh::initialized().

933 {
934  libmesh_assert (this->initialized());
935 
936  PetscInt petsc_m=0, petsc_n=0;
937  PetscErrorCode ierr=0;
938 
939  ierr = MatGetSize (_mat, &petsc_m, &petsc_n);
940  LIBMESH_CHKERR(ierr);
941 
942  return static_cast<numeric_index_type>(petsc_m);
943 }
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:92
PetscErrorCode ierr
template<typename T>
Mat libMesh::PetscMatrix< T >::mat ( )
inline
template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::n ( ) const
virtual
Returns
The column-dimension of the matrix.

Implements libMesh::SparseMatrix< T >.

Definition at line 948 of file petsc_matrix.C.

References ierr, and libMesh::initialized().

949 {
950  libmesh_assert (this->initialized());
951 
952  PetscInt petsc_m=0, petsc_n=0;
953  PetscErrorCode ierr=0;
954 
955  ierr = MatGetSize (_mat, &petsc_m, &petsc_n);
956  LIBMESH_CHKERR(ierr);
957 
958  return static_cast<numeric_index_type>(petsc_n);
959 }
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:92
PetscErrorCode ierr
static unsigned int libMesh::ReferenceCounter::n_objects ( )
inlinestaticinherited
processor_id_type libMesh::ParallelObject::n_processors ( ) const
inlineinherited
Returns
The number of processors in the group.

Definition at line 93 of file parallel_object.h.

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

Referenced by libMesh::ParmetisPartitioner::_do_repartition(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::DistributedMesh::add_elem(), libMesh::DistributedMesh::add_node(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assembly(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::AztecLinearSolver< T >::AztecLinearSolver(), libMesh::MeshCommunication::broadcast(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::DistributedMesh::clear(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::UnstructuredMesh::create_pid_mesh(), libMesh::MeshTools::create_processor_bounding_box(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DistributedMesh::DistributedMesh(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), libMesh::MeshBase::get_info(), libMesh::EquationSystems::get_solution(), libMesh::DistributedVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::DistributedMesh::insert_elem(), 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::SparsityPattern::Build::parallel_sync(), libMesh::Partitioner::partition(), libMesh::MeshBase::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::PetscLinearSolver< T >::PetscLinearSolver(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::SparseMatrix< T >::print(), libMesh::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::CheckpointIO::read(), libMesh::CheckpointIO::read_connectivity(), libMesh::XdrIO::read_header(), libMesh::CheckpointIO::read_nodes(), libMesh::MeshCommunication::redistribute(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::Partitioner::repartition(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::Parallel::Sort< KeyType, IdxType >::sort(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), and libMesh::XdrIO::write_serialized_nodesets().

94  { return cast_int<processor_id_type>(_communicator.size()); }
unsigned int size() const
Definition: parallel.h:725
const Parallel::Communicator & _communicator
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 118 of file sparse_matrix.h.

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

119  { return false; }
template<typename T >
T libMesh::PetscMatrix< T >::operator() ( const numeric_index_type  i,
const numeric_index_type  j 
) const
virtual
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 1074 of file petsc_matrix.C.

References libMesh::closed(), ierr, libMesh::initialized(), and value.

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

1076 {
1077  libmesh_assert (this->initialized());
1078 
1079  // PETSc 2.2.1 & newer
1080  const PetscScalar * petsc_row;
1081  const PetscInt * petsc_cols;
1082 
1083  // If the entry is not in the sparse matrix, it is 0.
1084  T value=0.;
1085 
1086  PetscErrorCode
1087  ierr=0;
1088  PetscInt
1089  ncols=0,
1090  i_val=static_cast<PetscInt>(i_in),
1091  j_val=static_cast<PetscInt>(j_in);
1092 
1093 
1094  // the matrix needs to be closed for this to work
1095  // this->close();
1096  // but closing it is a semiparallel operation; we want operator()
1097  // to run on one processor.
1098  libmesh_assert(this->closed());
1099 
1100  ierr = MatGetRow(_mat, i_val, &ncols, &petsc_cols, &petsc_row);
1101  LIBMESH_CHKERR(ierr);
1102 
1103  // Perform a binary search to find the contiguous index in
1104  // petsc_cols (resp. petsc_row) corresponding to global index j_val
1105  std::pair<const PetscInt *, const PetscInt *> p =
1106  std::equal_range (&petsc_cols[0], &petsc_cols[0] + ncols, j_val);
1107 
1108  // Found an entry for j_val
1109  if (p.first != p.second)
1110  {
1111  // The entry in the contiguous row corresponding
1112  // to the j_val column of interest
1113  const std::size_t j =
1114  std::distance (const_cast<PetscInt *>(&petsc_cols[0]),
1115  const_cast<PetscInt *>(p.first));
1116 
1117  libmesh_assert_less (static_cast<PetscInt>(j), ncols);
1118  libmesh_assert_equal_to (petsc_cols[j], j_val);
1119 
1120  value = static_cast<T> (petsc_row[j]);
1121  }
1122 
1123  ierr = MatRestoreRow(_mat, i_val,
1124  &ncols, &petsc_cols, &petsc_row);
1125  LIBMESH_CHKERR(ierr);
1126 
1127  return value;
1128 }
virtual bool closed() const libmesh_override
virtual bool initialized() const
PetscErrorCode ierr
static const bool value
Definition: xdr_io.C:108
template<>
void libMesh::SparseMatrix< Complex >::print ( std::ostream &  os,
const bool  sparse 
) const
inherited

Definition at line 102 of file sparse_matrix.C.

References libMesh::SparseMatrix< T >::m(), and libMesh::SparseMatrix< T >::n().

103 {
104  // std::complex<>::operator<<() is defined, but use this form
105 
106  if (sparse)
107  {
108  libmesh_not_implemented();
109  }
110 
111  os << "Real part:" << std::endl;
112  for (numeric_index_type i=0; i<this->m(); i++)
113  {
114  for (numeric_index_type j=0; j<this->n(); j++)
115  os << std::setw(8) << (*this)(i,j).real() << " ";
116  os << std::endl;
117  }
118 
119  os << std::endl << "Imaginary part:" << std::endl;
120  for (numeric_index_type i=0; i<this->m(); i++)
121  {
122  for (numeric_index_type j=0; j<this->n(); j++)
123  os << std::setw(8) << (*this)(i,j).imag() << " ";
124  os << std::endl;
125  }
126 }
virtual numeric_index_type n() const =0
virtual numeric_index_type m() const =0
dof_id_type numeric_index_type
Definition: id_types.h:92
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 206 of file sparse_matrix.C.

References libMesh::SparseMatrix< T >::_dof_map, libMesh::ParallelObject::comm(), libMesh::DofMap::end_dof(), libMesh::DofMap::first_dof(), libMesh::SparseMatrix< T >::initialized(), libMesh::SparseMatrix< T >::m(), libMesh::SparseMatrix< T >::n(), libMesh::ParallelObject::n_processors(), libMesh::ParallelObject::processor_id(), libMesh::Parallel::Communicator::receive(), and libMesh::Parallel::Communicator::send().

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

207 {
208  parallel_object_only();
209 
210  libmesh_assert (this->initialized());
211 
212  if (!this->_dof_map)
213  libmesh_error_msg("Error! Trying to print a matrix with no dof_map set!");
214 
215  // We'll print the matrix from processor 0 to make sure
216  // it's serialized properly
217  if (this->processor_id() == 0)
218  {
219  libmesh_assert_equal_to (this->_dof_map->first_dof(), 0);
220  for (numeric_index_type i=this->_dof_map->first_dof();
221  i!=this->_dof_map->end_dof(); ++i)
222  {
223  if (sparse)
224  {
225  for (numeric_index_type j=0; j<this->n(); j++)
226  {
227  T c = (*this)(i,j);
228  if (c != static_cast<T>(0.0))
229  {
230  os << i << " " << j << " " << c << std::endl;
231  }
232  }
233  }
234  else
235  {
236  for (numeric_index_type j=0; j<this->n(); j++)
237  os << (*this)(i,j) << " ";
238  os << std::endl;
239  }
240  }
241 
242  std::vector<numeric_index_type> ibuf, jbuf;
243  std::vector<T> cbuf;
244  numeric_index_type currenti = this->_dof_map->end_dof();
245  for (processor_id_type p=1; p < this->n_processors(); ++p)
246  {
247  this->comm().receive(p, ibuf);
248  this->comm().receive(p, jbuf);
249  this->comm().receive(p, cbuf);
250  libmesh_assert_equal_to (ibuf.size(), jbuf.size());
251  libmesh_assert_equal_to (ibuf.size(), cbuf.size());
252 
253  if (ibuf.empty())
254  continue;
255  libmesh_assert_greater_equal (ibuf.front(), currenti);
256  libmesh_assert_greater_equal (ibuf.back(), ibuf.front());
257 
258  std::size_t currentb = 0;
259  for (;currenti <= ibuf.back(); ++currenti)
260  {
261  if (sparse)
262  {
263  for (numeric_index_type j=0; j<this->n(); j++)
264  {
265  if (currentb < ibuf.size() &&
266  ibuf[currentb] == currenti &&
267  jbuf[currentb] == j)
268  {
269  os << currenti << " " << j << " " << cbuf[currentb] << std::endl;
270  currentb++;
271  }
272  }
273  }
274  else
275  {
276  for (numeric_index_type j=0; j<this->n(); j++)
277  {
278  if (currentb < ibuf.size() &&
279  ibuf[currentb] == currenti &&
280  jbuf[currentb] == j)
281  {
282  os << cbuf[currentb] << " ";
283  currentb++;
284  }
285  else
286  os << static_cast<T>(0.0) << " ";
287  }
288  os << std::endl;
289  }
290  }
291  }
292  if (!sparse)
293  {
294  for (; currenti != this->m(); ++currenti)
295  {
296  for (numeric_index_type j=0; j<this->n(); j++)
297  os << static_cast<T>(0.0) << " ";
298  os << std::endl;
299  }
300  }
301  }
302  else
303  {
304  std::vector<numeric_index_type> ibuf, jbuf;
305  std::vector<T> cbuf;
306 
307  // We'll assume each processor has access to entire
308  // matrix rows, so (*this)(i,j) is valid if i is a local index.
309  for (numeric_index_type i=this->_dof_map->first_dof();
310  i!=this->_dof_map->end_dof(); ++i)
311  {
312  for (numeric_index_type j=0; j<this->n(); j++)
313  {
314  T c = (*this)(i,j);
315  if (c != static_cast<T>(0.0))
316  {
317  ibuf.push_back(i);
318  jbuf.push_back(j);
319  cbuf.push_back(c);
320  }
321  }
322  }
323  this->comm().send(0,ibuf);
324  this->comm().send(0,jbuf);
325  this->comm().send(0,cbuf);
326  }
327 }
virtual numeric_index_type n() const =0
virtual numeric_index_type m() const =0
processor_id_type n_processors() const
uint8_t processor_id_type
Definition: id_types.h:99
virtual bool initialized() const
dof_id_type first_dof(const processor_id_type proc) const
Definition: dof_map.h:544
dof_id_type end_dof(const processor_id_type proc) const
Definition: dof_map.h:586
dof_id_type numeric_index_type
Definition: id_types.h:92
void send(const unsigned int dest_processor_id, const T &buf, const MessageTag &tag=no_tag) const
DofMap const * _dof_map
const Parallel::Communicator & comm() const
Status receive(const unsigned int dest_processor_id, T &buf, const MessageTag &tag=any_tag) const
processor_id_type processor_id() const
void libMesh::ReferenceCounter::print_info ( std::ostream &  out = libMesh::out)
staticinherited

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

Definition at line 88 of file reference_counter.C.

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

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

89 {
91  out_stream << ReferenceCounter::get_info();
92 }
static std::string get_info()
template<typename T >
void libMesh::PetscMatrix< T >::print_matlab ( const std::string &  = "") const
virtual

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.

Create an ASCII file containing the matrix if a filename was provided.

Otherwise the matrix will be dumped to the screen.

Destroy the viewer.

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 530 of file petsc_matrix.C.

References libMesh::closed(), ierr, and libMesh::initialized().

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

531 {
532  libmesh_assert (this->initialized());
533 
534  semiparallel_only();
535 
536  if (!this->closed())
537  {
538  libmesh_deprecated();
539  libmesh_warning("The matrix must be assembled before calling PetscMatrix::print_matlab().\n"
540  "Please update your code, as this warning will become an error in a future release.");
541  const_cast<PetscMatrix<T> *>(this)->close();
542  }
543 
544  PetscErrorCode ierr=0;
545  PetscViewer petsc_viewer;
546 
547 
548  ierr = PetscViewerCreate (this->comm().get(),
549  &petsc_viewer);
550  LIBMESH_CHKERR(ierr);
551 
556  if (name != "")
557  {
558  ierr = PetscViewerASCIIOpen( this->comm().get(),
559  name.c_str(),
560  &petsc_viewer);
561  LIBMESH_CHKERR(ierr);
562 
563 #if PETSC_VERSION_LESS_THAN(3,7,0)
564  ierr = PetscViewerSetFormat (petsc_viewer,
565  PETSC_VIEWER_ASCII_MATLAB);
566 #else
567  ierr = PetscViewerPushFormat (petsc_viewer,
568  PETSC_VIEWER_ASCII_MATLAB);
569 #endif
570 
571  LIBMESH_CHKERR(ierr);
572 
573  ierr = MatView (_mat, petsc_viewer);
574  LIBMESH_CHKERR(ierr);
575  }
576 
580  else
581  {
582 #if PETSC_VERSION_LESS_THAN(3,7,0)
583  ierr = PetscViewerSetFormat (PETSC_VIEWER_STDOUT_WORLD,
584  PETSC_VIEWER_ASCII_MATLAB);
585 #else
586  ierr = PetscViewerPushFormat (PETSC_VIEWER_STDOUT_WORLD,
587  PETSC_VIEWER_ASCII_MATLAB);
588 #endif
589 
590  LIBMESH_CHKERR(ierr);
591 
592  ierr = MatView (_mat, PETSC_VIEWER_STDOUT_WORLD);
593  LIBMESH_CHKERR(ierr);
594  }
595 
596 
600  ierr = LibMeshPetscViewerDestroy (&petsc_viewer);
601  LIBMESH_CHKERR(ierr);
602 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
virtual bool closed() const libmesh_override
virtual bool initialized() const
virtual void close() libmesh_override
Definition: petsc_matrix.C:897
PetscErrorCode ierr
const Parallel::Communicator & comm() const
SparseMatrix interface to PETSc Mat.
template<typename T >
void libMesh::PetscMatrix< T >::print_personal ( std::ostream &  os = libMesh::out) const
virtual

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 609 of file petsc_matrix.C.

References libMesh::closed(), ierr, libMesh::initialized(), and libMesh::processor_id().

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

610 {
611  libmesh_assert (this->initialized());
612 
613  // Routine must be called in parallel on parallel matrices
614  // and serial on serial matrices.
615  semiparallel_only();
616 
617  // #ifndef NDEBUG
618  // if (os != std::cout)
619  // libMesh::err << "Warning! PETSc can only print to std::cout!" << std::endl;
620  // #endif
621 
622  // Matrix must be in an assembled state to be printed
623  if (!this->closed())
624  {
625  libmesh_deprecated();
626  libmesh_warning("The matrix must be assembled before calling PetscMatrix::print_personal().\n"
627  "Please update your code, as this warning will become an error in a future release.");
628  const_cast<PetscMatrix<T> *>(this)->close();
629  }
630 
631  PetscErrorCode ierr=0;
632 
633  // Print to screen if ostream is stdout
634  if (os.rdbuf() == std::cout.rdbuf())
635  {
636  ierr = MatView(_mat, PETSC_VIEWER_STDOUT_SELF);
637  LIBMESH_CHKERR(ierr);
638  }
639 
640  // Otherwise, print to the requested file, in a roundabout way...
641  else
642  {
643  // We will create a temporary filename, and file, for PETSc to
644  // write to.
645  std::string temp_filename;
646 
647  {
648  // Template for temporary filename
649  char c[] = "temp_petsc_matrix.XXXXXX";
650 
651  // Generate temporary, unique filename only on processor 0. We will
652  // use this filename for PetscViewerASCIIOpen, before copying it into
653  // the user's stream
654  if (this->processor_id() == 0)
655  {
656  int fd = mkstemp(c);
657 
658  // Check to see that mkstemp did not fail.
659  if (fd == -1)
660  libmesh_error_msg("mkstemp failed in PetscMatrix::print_personal()");
661 
662  // mkstemp returns a file descriptor for an open file,
663  // so let's close it before we hand it to PETSc!
664  ::close (fd);
665  }
666 
667  // Store temporary filename as string, makes it easier to broadcast
668  temp_filename = c;
669  }
670 
671  // Now broadcast the filename from processor 0 to all processors.
672  this->comm().broadcast(temp_filename);
673 
674  // PetscViewer object for passing to MatView
675  PetscViewer petsc_viewer;
676 
677  // This PETSc function only takes a string and handles the opening/closing
678  // of the file internally. Since print_personal() takes a reference to
679  // an ostream, we have to do an extra step... print_personal() should probably
680  // have a version that takes a string to get rid of this problem.
681  ierr = PetscViewerASCIIOpen( this->comm().get(),
682  temp_filename.c_str(),
683  &petsc_viewer);
684  LIBMESH_CHKERR(ierr);
685 
686  // Probably don't need to set the format if it's default...
687  // ierr = PetscViewerSetFormat (petsc_viewer,
688  // PETSC_VIEWER_DEFAULT);
689  // LIBMESH_CHKERR(ierr);
690 
691  // Finally print the matrix using the viewer
692  ierr = MatView (_mat, petsc_viewer);
693  LIBMESH_CHKERR(ierr);
694 
695  if (this->processor_id() == 0)
696  {
697  // Now the inefficient bit: open temp_filename as an ostream and copy the contents
698  // into the user's desired ostream. We can't just do a direct file copy, we don't even have the filename!
699  std::ifstream input_stream(temp_filename.c_str());
700  os << input_stream.rdbuf(); // The "most elegant" way to copy one stream into another.
701  // os.close(); // close not defined in ostream
702 
703  // Now remove the temporary file
704  input_stream.close();
705  std::remove(temp_filename.c_str());
706  }
707  }
708 }
virtual bool closed() const libmesh_override
virtual bool initialized() const
virtual void close() libmesh_override
Definition: petsc_matrix.C:897
void broadcast(T &data, const unsigned int root_id=0) const
PetscErrorCode ierr
const Parallel::Communicator & comm() const
SparseMatrix interface to PETSc Mat.
processor_id_type processor_id() const
processor_id_type libMesh::ParallelObject::processor_id ( ) const
inlineinherited
Returns
The rank of this processor in the group.

Definition at line 99 of file parallel_object.h.

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

Referenced by libMesh::BoundaryInfo::_find_id_maps(), libMesh::EquationSystems::_read_impl(), libMesh::DistributedMesh::add_elem(), libMesh::BoundaryInfo::add_elements(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::DistributedMesh::add_node(), libMesh::MeshRefinement::add_node(), libMesh::MeshTools::Modification::all_tri(), libMesh::FEMSystem::assembly(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::MeshCommunication::broadcast(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::Nemesis_IO_Helper::build_element_and_node_maps(), libMesh::ParmetisPartitioner::build_graph(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::DofMap::build_sparsity(), libMesh::DistributedMesh::clear(), 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::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::MeshTools::correct_node_proc_ids(), libMesh::ExodusII_IO_Helper::create(), libMesh::DistributedMesh::delete_elem(), libMesh::DistributedMesh::delete_node(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DistributedMesh::DistributedMesh(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshFunction::find_element(), libMesh::MeshFunction::find_elements(), libMesh::UnstructuredMesh::find_neighbors(), libMesh::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), 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::Nemesis_IO_Helper::get_node_cmap(), libMesh::Nemesis_IO_Helper::get_node_map(), libMesh::Nemesis_IO_Helper::get_ns_param_global(), libMesh::EquationSystems::get_solution(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), libMesh::SparsityPattern::Build::handle_vi_vj(), libMesh::DistributedVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), 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::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_object_ids(), libMesh::DofMap::local_variable_indices(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshBase::n_active_local_elem(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), 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::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::SparsityPattern::Build::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::SparsityPattern::Build::parallel_sync(), libMesh::MetisPartitioner::partition_range(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::SparseMatrix< T >::print(), libMesh::NumericVector< T >::print_global(), 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::ExodusII_IO_Helper::read_elem_num_map(), libMesh::ExodusII_IO_Helper::read_global_values(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::ExodusII_IO_Helper::read_node_num_map(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::MeshCommunication::redistribute(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::CheckpointIO::select_split_config(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::MeshTools::total_weight(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::Parallel::Packing< Node * >::unpack(), libMesh::Parallel::Packing< Elem * >::unpack(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::MeshTools::weight(), 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_elements(), libMesh::ExodusII_IO::write_global_data(), libMesh::ExodusII_IO_Helper::write_global_values(), 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::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::ExodusII_IO_Helper::write_nodal_values(), libMesh::Nemesis_IO_Helper::write_nodesets(), libMesh::ExodusII_IO_Helper::write_nodesets(), libMesh::XdrIO::write_serialized_bc_names(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::XdrIO::write_serialized_subdomain_names(), 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().

100  { return cast_int<processor_id_type>(_communicator.rank()); }
const Parallel::Communicator & _communicator
unsigned int rank() const
Definition: parallel.h:723
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 367 of file sparse_matrix.h.

370  {
371  this->_get_submatrix(submatrix,
372  rows,
373  cols,
374  true); // true means REUSE submatrix
375  }
virtual void _get_submatrix(SparseMatrix< T > &, const std::vector< numeric_index_type > &, const std::vector< numeric_index_type > &, const bool) const
template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::row_start ( ) const
virtual
Returns
The index of the first matrix row stored on this processor.

Implements libMesh::SparseMatrix< T >.

Definition at line 964 of file petsc_matrix.C.

References ierr, libMesh::initialized(), and libMesh::MacroFunctions::stop().

965 {
966  libmesh_assert (this->initialized());
967 
968  PetscInt start=0, stop=0;
969  PetscErrorCode ierr=0;
970 
971  ierr = MatGetOwnershipRange(_mat, &start, &stop);
972  LIBMESH_CHKERR(ierr);
973 
974  return static_cast<numeric_index_type>(start);
975 }
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:92
void stop(const char *file, int line, const char *date, const char *time)
PetscErrorCode ierr
template<typename T >
numeric_index_type libMesh::PetscMatrix< T >::row_stop ( ) const
virtual
Returns
The index of the last matrix row (+1) stored on this processor.

Implements libMesh::SparseMatrix< T >.

Definition at line 980 of file petsc_matrix.C.

References ierr, libMesh::initialized(), and libMesh::MacroFunctions::stop().

981 {
982  libmesh_assert (this->initialized());
983 
984  PetscInt start=0, stop=0;
985  PetscErrorCode ierr=0;
986 
987  ierr = MatGetOwnershipRange(_mat, &start, &stop);
988  LIBMESH_CHKERR(ierr);
989 
990  return static_cast<numeric_index_type>(stop);
991 }
virtual bool initialized() const
dof_id_type numeric_index_type
Definition: id_types.h:92
void stop(const char *file, int line, const char *date, const char *time)
PetscErrorCode ierr
template<typename T >
void libMesh::PetscMatrix< T >::set ( const numeric_index_type  i,
const numeric_index_type  j,
const T  value 
)
virtual

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 996 of file petsc_matrix.C.

References ierr, libMesh::initialized(), and value.

999 {
1000  libmesh_assert (this->initialized());
1001 
1002  PetscErrorCode ierr=0;
1003  PetscInt i_val=i, j_val=j;
1004 
1005  PetscScalar petsc_value = static_cast<PetscScalar>(value);
1006  ierr = MatSetValues(_mat, 1, &i_val, 1, &j_val,
1007  &petsc_value, INSERT_VALUES);
1008  LIBMESH_CHKERR(ierr);
1009 }
virtual bool initialized() const
PetscErrorCode ierr
static const bool value
Definition: xdr_io.C:108
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 1150 of file petsc_matrix.C.

References libMesh::PetscMatrix< T >::_destroy_mat_on_exit, libMesh::PetscMatrix< T >::_mat, and swap().

Referenced by libMesh::__libmesh_petsc_diff_solver_jacobian(), and libMesh::PetscMatrix< T >::add_block_matrix().

1151 {
1152  std::swap(_mat, m_in._mat);
1154 }
void swap(Iterator &lhs, Iterator &rhs)
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 406 of file petsc_matrix.C.

407 {
408  libmesh_not_implemented();
409 }
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 126 of file sparse_matrix.h.

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

126 {}
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 176 of file sparse_matrix.C.

References libMesh::SparseMatrix< T >::vector_mult_add(), and libMesh::NumericVector< T >::zero().

Referenced by libMesh::ImplicitSystem::qoi_parameter_hessian(), and libMesh::ImplicitSystem::qoi_parameter_hessian_vector_product().

178 {
179  dest.zero();
180  this->vector_mult_add(dest,arg);
181 }
void vector_mult_add(NumericVector< T > &dest, const NumericVector< T > &arg) const
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 186 of file sparse_matrix.C.

References libMesh::NumericVector< T >::add_vector().

Referenced by libMesh::SparseMatrix< T >::vector_mult(), and libMesh::ImplicitSystem::weighted_sensitivity_adjoint_solve().

188 {
189  /* This functionality is actually implemented in the \p
190  NumericVector class. */
191  dest.add_vector(arg,*this);
192 }
template<typename T >
void libMesh::PetscMatrix< T >::zero ( )
virtual

Set all entries to 0.

Implements libMesh::SparseMatrix< T >.

Definition at line 413 of file petsc_matrix.C.

References ierr, and libMesh::initialized().

414 {
415  libmesh_assert (this->initialized());
416 
417  semiparallel_only();
418 
419  PetscErrorCode ierr=0;
420 
421  PetscInt m_l, n_l;
422 
423  ierr = MatGetLocalSize(_mat,&m_l,&n_l);
424  LIBMESH_CHKERR(ierr);
425 
426  if (n_l)
427  {
428  ierr = MatZeroEntries(_mat);
429  LIBMESH_CHKERR(ierr);
430  }
431 }
virtual bool initialized() const
PetscErrorCode ierr
template<typename T >
void libMesh::PetscMatrix< T >::zero_rows ( std::vector< numeric_index_type > &  rows,
diag_value = 0.0 
)
virtual

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

Reimplemented from libMesh::SparseMatrix< T >.

Definition at line 434 of file petsc_matrix.C.

References ierr, libMesh::initialized(), and libMesh::numeric_petsc_cast().

435 {
436  libmesh_assert (this->initialized());
437 
438  semiparallel_only();
439 
440  PetscErrorCode ierr=0;
441 
442 #if PETSC_RELEASE_LESS_THAN(3,1,1)
443  if (!rows.empty())
444  ierr = MatZeroRows(_mat, rows.size(),
445  numeric_petsc_cast(&rows[0]), diag_value);
446  else
447  ierr = MatZeroRows(_mat, 0, PETSC_NULL, diag_value);
448 #else
449  // As of petsc-dev at the time of 3.1.0, MatZeroRows now takes two additional
450  // optional arguments. The optional arguments (x,b) can be used to specify the
451  // solutions for the zeroed rows (x) and right hand side (b) to update.
452  // Could be useful for setting boundary conditions...
453  if (!rows.empty())
454  ierr = MatZeroRows(_mat, cast_int<PetscInt>(rows.size()),
455  numeric_petsc_cast(&rows[0]), diag_value,
456  PETSC_NULL, PETSC_NULL);
457  else
458  ierr = MatZeroRows(_mat, 0, PETSC_NULL, diag_value, PETSC_NULL,
459  PETSC_NULL);
460 #endif
461 
462  LIBMESH_CHKERR(ierr);
463 }
PetscInt * numeric_petsc_cast(const numeric_index_type *p)
virtual bool initialized() const
PetscErrorCode ierr

Member Data Documentation

ReferenceCounter::Counts libMesh::ReferenceCounter::_counts
staticprotectedinherited
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 268 of file petsc_matrix.h.

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

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

The DofMap object associated with this object.

Definition at line 422 of file sparse_matrix.h.

Referenced by libMesh::SparseMatrix< T >::print().

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().

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

Flag indicating whether or not the matrix has been initialized.

Definition at line 427 of file sparse_matrix.h.

Referenced by libMesh::PetscMatrix< T >::_get_submatrix(), libMesh::EigenSparseMatrix< T >::clear(), libMesh::LaspackMatrix< T >::clear(), and libMesh::PetscMatrix< T >::get_transpose().

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.

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().


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