Reference documentation for deal.II version GIT relicensing-214-g6e74dec06b 2024-03-27 18:10:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Classes | Public Types | Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | Private Types | Private Member Functions | Private Attributes | Static Private Attributes | Friends | List of all members
PETScWrappers::MPI::SparseMatrix Class Reference

#include <deal.II/lac/petsc_sparse_matrix.h>

Inheritance diagram for PETScWrappers::MPI::SparseMatrix:
Inheritance graph
[legend]

Classes

struct  Traits
 

Public Types

using size_type = types::global_dof_index
 
using const_iterator = MatrixIterators::const_iterator
 
using value_type = PetscScalar
 

Public Member Functions

 SparseMatrix ()
 
 SparseMatrix (const Mat &)
 
 ~SparseMatrix () override
 
template<typename SparsityPatternType >
 SparseMatrix (const MPI_Comm communicator, const SparsityPatternType &sparsity_pattern, const std::vector< size_type > &local_rows_per_process, const std::vector< size_type > &local_columns_per_process, const unsigned int this_process, const bool preset_nonzero_locations=true)
 
SparseMatrixoperator= (const value_type d)
 
void copy_from (const SparseMatrix &other)
 
template<typename SparsityPatternType >
void reinit (const MPI_Comm communicator, const SparsityPatternType &sparsity_pattern, const std::vector< size_type > &local_rows_per_process, const std::vector< size_type > &local_columns_per_process, const unsigned int this_process, const bool preset_nonzero_locations=true)
 
template<typename SparsityPatternType >
void reinit (const IndexSet &local_partitioning, const SparsityPatternType &sparsity_pattern, const MPI_Comm communicator)
 
template<typename SparsityPatternType >
void reinit (const IndexSet &local_rows, const IndexSet &local_columns, const SparsityPatternType &sparsity_pattern, const MPI_Comm communicator)
 
void reinit (const SparseMatrix &other)
 
template<typename SparsityPatternType >
void reinit (const IndexSet &local_rows, const IndexSet &local_active_rows, const IndexSet &local_columns, const IndexSet &local_active_columns, const SparsityPatternType &sparsity_pattern, const MPI_Comm communicator)
 
PetscScalar matrix_norm_square (const Vector &v) const
 
PetscScalar matrix_scalar_product (const Vector &u, const Vector &v) const
 
IndexSet locally_owned_domain_indices () const
 
IndexSet locally_owned_range_indices () const
 
void mmult (SparseMatrix &C, const SparseMatrix &B, const MPI::Vector &V=MPI::Vector()) const
 
void Tmmult (SparseMatrix &C, const SparseMatrix &B, const MPI::Vector &V=MPI::Vector()) const
 
void reinit (Mat A)
 
void clear ()
 
void set (const size_type i, const size_type j, const PetscScalar value)
 
void set (const std::vector< size_type > &indices, const FullMatrix< PetscScalar > &full_matrix, const bool elide_zero_values=false)
 
void set (const std::vector< size_type > &row_indices, const std::vector< size_type > &col_indices, const FullMatrix< PetscScalar > &full_matrix, const bool elide_zero_values=false)
 
void set (const size_type row, const std::vector< size_type > &col_indices, const std::vector< PetscScalar > &values, const bool elide_zero_values=false)
 
void set (const size_type row, const size_type n_cols, const size_type *col_indices, const PetscScalar *values, const bool elide_zero_values=false)
 
void add (const size_type i, const size_type j, const PetscScalar value)
 
void add (const std::vector< size_type > &indices, const FullMatrix< PetscScalar > &full_matrix, const bool elide_zero_values=true)
 
void add (const std::vector< size_type > &row_indices, const std::vector< size_type > &col_indices, const FullMatrix< PetscScalar > &full_matrix, const bool elide_zero_values=true)
 
void add (const size_type row, const std::vector< size_type > &col_indices, const std::vector< PetscScalar > &values, const bool elide_zero_values=true)
 
void add (const size_type row, const size_type n_cols, const size_type *col_indices, const PetscScalar *values, const bool elide_zero_values=true, const bool col_indices_are_sorted=false)
 
MatrixBaseadd (const PetscScalar factor, const MatrixBase &other)
 
void clear_row (const size_type row, const PetscScalar new_diag_value=0)
 
void clear_rows (const ArrayView< const size_type > &rows, const PetscScalar new_diag_value=0)
 
void clear_rows_columns (const std::vector< size_type > &row_and_column_indices, const PetscScalar new_diag_value=0)
 
void compress (const VectorOperation::values operation)
 
PetscScalar operator() (const size_type i, const size_type j) const
 
PetscScalar el (const size_type i, const size_type j) const
 
PetscScalar diag_element (const size_type i) const
 
size_type m () const
 
size_type n () const
 
size_type local_size () const
 
std::pair< size_type, size_typelocal_range () const
 
bool in_local_range (const size_type index) const
 
size_type local_domain_size () const
 
std::pair< size_type, size_typelocal_domain () const
 
MPI_Comm get_mpi_communicator () const
 
std::uint64_t n_nonzero_elements () const
 
size_type row_length (const size_type row) const
 
PetscReal l1_norm () const
 
PetscReal linfty_norm () const
 
PetscReal frobenius_norm () const
 
PetscScalar matrix_norm_square (const VectorBase &v) const
 
PetscScalar matrix_scalar_product (const VectorBase &u, const VectorBase &v) const
 
PetscScalar trace () const
 
MatrixBaseoperator*= (const PetscScalar factor)
 
MatrixBaseoperator/= (const PetscScalar factor)
 
void vmult (VectorBase &dst, const VectorBase &src) const
 
void Tvmult (VectorBase &dst, const VectorBase &src) const
 
void vmult_add (VectorBase &dst, const VectorBase &src) const
 
void Tvmult_add (VectorBase &dst, const VectorBase &src) const
 
PetscScalar residual (VectorBase &dst, const VectorBase &x, const VectorBase &b) const
 
const_iterator begin () const
 
const_iterator begin (const size_type r) const
 
const_iterator end () const
 
const_iterator end (const size_type r) const
 
 operator Mat () const
 
Mat & petsc_matrix ()
 
void transpose ()
 
PetscBool is_symmetric (const double tolerance=1.e-12)
 
PetscBool is_hermitian (const double tolerance=1.e-12)
 
void write_ascii (const PetscViewerFormat format=PETSC_VIEWER_DEFAULT)
 
void print (std::ostream &out, const bool alternative_output=false) const
 
std::size_t memory_consumption () const
 
template<class Archive >
void serialize (Archive &ar, const unsigned int version)
 
Subscriptor functionality

Classes derived from Subscriptor provide a facility to subscribe to this object. This is mostly used by the SmartPointer class.

void subscribe (std::atomic< bool > *const validity, const std::string &identifier="") const
 
void unsubscribe (std::atomic< bool > *const validity, const std::string &identifier="") const
 
unsigned int n_subscriptions () const
 
template<typename StreamType >
void list_subscribers (StreamType &stream) const
 
void list_subscribers () const
 

Static Public Member Functions

static ::ExceptionBaseExcLocalRowsTooLarge (int arg1, int arg2)
 
static ::ExceptionBaseExcSourceEqualsDestination ()
 
static ::ExceptionBaseExcWrongMode (int arg1, int arg2)
 
static ::ExceptionBaseExcInUse (int arg1, std::string arg2, std::string arg3)
 
static ::ExceptionBaseExcNoSubscriber (std::string arg1, std::string arg2)
 

Protected Member Functions

void prepare_action (const VectorOperation::values new_action)
 
void assert_is_compressed ()
 
void prepare_add ()
 
void prepare_set ()
 
void mmult (MatrixBase &C, const MatrixBase &B, const VectorBase &V) const
 
void Tmmult (MatrixBase &C, const MatrixBase &B, const VectorBase &V) const
 

Protected Attributes

Mat matrix
 
VectorOperation::values last_action
 

Private Types

using map_value_type = decltype(counter_map)::value_type
 
using map_iterator = decltype(counter_map)::iterator
 

Private Member Functions

template<typename SparsityPatternType >
void do_reinit (const MPI_Comm comm, const SparsityPatternType &sparsity_pattern, const std::vector< size_type > &local_rows_per_process, const std::vector< size_type > &local_columns_per_process, const unsigned int this_process, const bool preset_nonzero_locations)
 
template<typename SparsityPatternType >
void do_reinit (const MPI_Comm comm, const IndexSet &local_rows, const IndexSet &local_columns, const SparsityPatternType &sparsity_pattern)
 
template<typename SparsityPatternType >
void do_reinit (const MPI_Comm comm, const IndexSet &local_rows, const IndexSet &local_active_rows, const IndexSet &local_columns, const IndexSet &local_active_columns, const SparsityPatternType &sparsity_pattern)
 
void check_no_subscribers () const noexcept
 

Private Attributes

std::vector< PetscInt > column_indices
 
std::vector< PetscScalar > column_values
 
std::atomic< unsigned intcounter
 
std::map< std::string, unsigned intcounter_map
 
std::vector< std::atomic< bool > * > validity_pointers
 
const std::type_info * object_info
 

Static Private Attributes

static std::mutex mutex
 

Friends

class BlockMatrixBase< SparseMatrix >
 

Detailed Description

Implementation of a parallel sparse matrix class based on PETSc, with rows of the matrix distributed across an MPI network. All the functionality is actually in the base class, except for the calls to generate a parallel sparse matrix. This is possible since PETSc only works on an abstract matrix type and internally distributes to functions that do the actual work depending on the actual matrix type (much like using virtual functions). Only the functions creating a matrix of specific type differ, and are implemented in this particular class.

There are a number of comments on the communication model as well as access to individual elements in the documentation to the parallel vector class. These comments apply here as well.

Partitioning of matrices

PETSc partitions parallel matrices so that each MPI process "owns" a certain number of rows (i.e. only this process stores the respective entries in these rows). The number of rows each process owns has to be passed to the constructors and reinit() functions via the argument local_rows. The individual values passed as local_rows on all the MPI processes of course have to add up to the global number of rows of the matrix.

In addition to this, PETSc also partitions the rectangular chunk of the matrix it owns (i.e. the local_rows times n() elements in the matrix), so that matrix vector multiplications can be performed efficiently. This column-partitioning therefore has to match the partitioning of the vectors with which the matrix is multiplied, just as the row-partitioning has to match the partitioning of destination vectors. This partitioning is passed to the constructors and reinit() functions through the local_columns variable, which again has to add up to the global number of columns in the matrix. The name local_columns may be named inappropriately since it does not reflect that only these columns are stored locally, but it reflects the fact that these are the columns for which the elements of incoming vectors are stored locally.

To make things even more complicated, PETSc needs a very good estimate of the number of elements to be stored in each row to be efficient. Otherwise it spends most of the time with allocating small chunks of memory, a process that can slow down programs to a crawl if it happens to often. As if a good estimate of the number of entries per row isn't even, it even needs to split this as follows: for each row it owns, it needs an estimate for the number of elements in this row that fall into the columns that are set apart for this process (see above), and the number of elements that are in the rest of the columns.

Since in general this information is not readily available, most of the initializing functions of this class assume that all of the number of elements you give as an argument to n_nonzero_per_row or by row_lengths fall into the columns "owned" by this process, and none into the other ones. This is a fair guess for most of the rows, since in a good domain partitioning, nodes only interact with nodes that are within the same subdomain. It does not hold for nodes on the interfaces of subdomain, however, and for the rows corresponding to these nodes, PETSc will have to allocate additional memory, a costly process.

The only way to avoid this is to tell PETSc where the actual entries of the matrix will be. For this, there are constructors and reinit() functions of this class that take a DynamicSparsityPattern object containing all this information. While in the general case it is sufficient if the constructors and reinit() functions know the number of local rows and columns, the functions getting a sparsity pattern also need to know the number of local rows (local_rows_per_process) and columns (local_columns_per_process) for all other processes, in order to compute which parts of the matrix are which. Thus, it is not sufficient to just count the number of degrees of freedom that belong to a particular process, but you have to have the numbers for all processes available at all processes.

Definition at line 366 of file petsc_sparse_matrix.h.

Member Typedef Documentation

◆ size_type

Declare type for container size.

Definition at line 372 of file petsc_sparse_matrix.h.

◆ const_iterator

Declare an alias for the iterator class.

Definition at line 289 of file petsc_matrix_base.h.

◆ value_type

using PETScWrappers::MatrixBase::value_type = PetscScalar
inherited

Declare an alias in analogy to all the other container classes.

Definition at line 299 of file petsc_matrix_base.h.

◆ map_value_type

using Subscriptor::map_value_type = decltype(counter_map)::value_type
privateinherited

The data type used in counter_map.

Definition at line 229 of file subscriptor.h.

◆ map_iterator

using Subscriptor::map_iterator = decltype(counter_map)::iterator
privateinherited

The iterator type used in counter_map.

Definition at line 234 of file subscriptor.h.

Constructor & Destructor Documentation

◆ SparseMatrix() [1/3]

Default constructor. Create an empty matrix.

Definition at line 33 of file petsc_parallel_sparse_matrix.cc.

◆ SparseMatrix() [2/3]

SparseMatrix< number >::SparseMatrix ( const Mat &  A)
explicit

Initialize a SparseMatrix from a PETSc Mat object. Note that we do not copy the matrix. The Mat object is referenced by the newly created instance of the class using PetscObjectReference. This is in line with the PETSc approach to object ownership, which mimics std::shared_ptr.

Definition at line 44 of file petsc_parallel_sparse_matrix.cc.

◆ ~SparseMatrix()

SparseMatrix< number >::~SparseMatrix ( )
override

Destructor to free the PETSc object.

Definition at line 48 of file petsc_parallel_sparse_matrix.cc.

◆ SparseMatrix() [3/3]

template<typename SparsityPatternType >
SparseMatrix< SparsityPatternType >::SparseMatrix ( const MPI_Comm  communicator,
const SparsityPatternType &  sparsity_pattern,
const std::vector< size_type > &  local_rows_per_process,
const std::vector< size_type > &  local_columns_per_process,
const unsigned int  this_process,
const bool  preset_nonzero_locations = true 
)

Initialize using the given sparsity pattern with communication happening over the provided communicator.

For the meaning of the local_rows_per_process and local_columns_per_process parameters, see the class documentation.

Note that PETSc can be very slow if you do not provide it with a good estimate of the lengths of rows. Using the present function is a very efficient way to do this, as it uses the exact number of nonzero entries for each row of the matrix by using the given sparsity pattern argument. If the preset_nonzero_locations flag is true, this function in addition not only sets the correct row sizes up front, but also pre-allocated the correct nonzero entries in the matrix.

PETsc allows to later add additional nonzero entries to a matrix, by simply writing to these elements. However, this will then lead to additional memory allocations which are very inefficient and will greatly slow down your program. It is therefore significantly more efficient to get memory allocation right from the start.

Definition at line 58 of file petsc_parallel_sparse_matrix.cc.

Member Function Documentation

◆ operator=()

SparseMatrix & SparseMatrix< number >::operator= ( const value_type  d)

This operator assigns a scalar to a matrix. Since this does usually not make much sense (should we set all matrix entries to this value? Only the nonzero entries of the sparsity pattern?), this operation is only allowed if the actual value to be assigned is zero. This operator only exists to allow for the obvious notation matrix=0, which sets all elements of the matrix to zero, but keep the sparsity pattern previously used.

Definition at line 112 of file petsc_parallel_sparse_matrix.cc.

◆ copy_from()

void SparseMatrix< number >::copy_from ( const SparseMatrix other)

Make a copy of the PETSc matrix other. It is assumed that both matrices have the same SparsityPattern.

Definition at line 119 of file petsc_parallel_sparse_matrix.cc.

◆ reinit() [1/6]

template<typename SparsityPatternType >
void SparseMatrix< SparsityPatternType >::reinit ( const MPI_Comm  communicator,
const SparsityPatternType &  sparsity_pattern,
const std::vector< size_type > &  local_rows_per_process,
const std::vector< size_type > &  local_columns_per_process,
const unsigned int  this_process,
const bool  preset_nonzero_locations = true 
)

Initialize using the given sparsity pattern with communication happening over the provided communicator.

Note that PETSc can be very slow if you do not provide it with a good estimate of the lengths of rows. Using the present function is a very efficient way to do this, as it uses the exact number of nonzero entries for each row of the matrix by using the given sparsity pattern argument. If the preset_nonzero_locations flag is true, this function in addition not only sets the correct row sizes up front, but also pre-allocated the correct nonzero entries in the matrix.

PETsc allows to later add additional nonzero entries to a matrix, by simply writing to these elements. However, this will then lead to additional memory allocations which are very inefficient and will greatly slow down your program. It is therefore significantly more efficient to get memory allocation right from the start.

Definition at line 133 of file petsc_parallel_sparse_matrix.cc.

◆ reinit() [2/6]

template<typename SparsityPatternType >
void SparseMatrix< SparsityPatternType >::reinit ( const IndexSet local_partitioning,
const SparsityPatternType &  sparsity_pattern,
const MPI_Comm  communicator 
)

Create a square matrix where the size() of the IndexSet determines the global number of rows and columns and the entries of the IndexSet give the rows and columns for the calling processor. Note that only ascending, 1:1 IndexSets are supported.

Definition at line 158 of file petsc_parallel_sparse_matrix.cc.

◆ reinit() [3/6]

template<typename SparsityPatternType >
void SparseMatrix< SparsityPatternType >::reinit ( const IndexSet local_rows,
const IndexSet local_columns,
const SparsityPatternType &  sparsity_pattern,
const MPI_Comm  communicator 
)

Create a matrix where the size() of the IndexSets determine the global number of rows and columns and the entries of the IndexSet give the rows and columns for the calling processor. Note that only ascending, 1:1 IndexSets are supported.

Definition at line 167 of file petsc_parallel_sparse_matrix.cc.

◆ reinit() [4/6]

void SparseMatrix< number >::reinit ( const SparseMatrix other)

Initialize this matrix to have the same structure as other. This will not copy the values of the other matrix, but you can use copy_from() for this.

Definition at line 77 of file petsc_parallel_sparse_matrix.cc.

◆ reinit() [5/6]

template<typename SparsityPatternType >
void SparseMatrix< SparsityPatternType >::reinit ( const IndexSet local_rows,
const IndexSet local_active_rows,
const IndexSet local_columns,
const IndexSet local_active_columns,
const SparsityPatternType &  sparsity_pattern,
const MPI_Comm  communicator 
)

Create a matrix where the size of the IndexSets determine the global number of rows and columns and the entries of the IndexSet give the rows and columns for the calling processor. Note that only ascending, 1:1 IndexSets are supported. The additional call to the local to global mappings is required to create the matrix of type IS (see DoFTools::extract_locally_active_dofs). This is required by the BDDC preconditioner.

Definition at line 91 of file petsc_parallel_sparse_matrix.cc.

◆ matrix_norm_square() [1/2]

PetscScalar SparseMatrix< number >::matrix_norm_square ( const Vector v) const

Return the square of the norm of the vector \(v\) with respect to the norm induced by this matrix, i.e. \(\left(v^\ast,Mv\right)\). This is useful, e.g. in the finite element context, where the \(L_2\) norm of a function equals the matrix norm with respect to the mass matrix of the vector representing the nodal values of the finite element function.

Obviously, the matrix needs to be quadratic for this operation.

The implementation of this function is not as efficient as the one in the MatrixBase class used in deal.II (i.e. the original one, not the PETSc wrapper class) since PETSc doesn't support this operation and needs a temporary vector.

Definition at line 805 of file petsc_parallel_sparse_matrix.cc.

◆ matrix_scalar_product() [1/2]

PetscScalar SparseMatrix< number >::matrix_scalar_product ( const Vector u,
const Vector v 
) const

Compute the matrix scalar product \(\left(u^\ast,Mv\right)\).

The implementation of this function is not as efficient as the one in the MatrixBase class used in deal.II (i.e. the original one, not the PETSc wrapper class) since PETSc doesn't support this operation and needs a temporary vector.

Definition at line 814 of file petsc_parallel_sparse_matrix.cc.

◆ locally_owned_domain_indices()

IndexSet SparseMatrix< number >::locally_owned_domain_indices ( ) const

Return the partitioning of the domain space of this matrix, i.e., the partitioning of the vectors this matrix has to be multiplied with.

Definition at line 823 of file petsc_parallel_sparse_matrix.cc.

◆ locally_owned_range_indices()

IndexSet SparseMatrix< number >::locally_owned_range_indices ( ) const

Return the partitioning of the range space of this matrix, i.e., the partitioning of the vectors that result from matrix-vector products.

Definition at line 849 of file petsc_parallel_sparse_matrix.cc.

◆ mmult() [1/2]

void SparseMatrix< number >::mmult ( SparseMatrix C,
const SparseMatrix B,
const MPI::Vector V = MPI::Vector() 
) const

Perform the matrix-matrix multiplication \(C = AB\), or, \(C = A \text{diag}(V) B\) given a compatible vector \(V\).

This function calls MatrixBase::mmult() to do the actual work.

Definition at line 875 of file petsc_parallel_sparse_matrix.cc.

◆ Tmmult() [1/2]

void SparseMatrix< number >::Tmmult ( SparseMatrix C,
const SparseMatrix B,
const MPI::Vector V = MPI::Vector() 
) const

Perform the matrix-matrix multiplication with the transpose of this, i.e., \(C = A^T B\), or, \(C = A^T \text{diag}(V) B\) given a compatible vector \(V\).

This function calls MatrixBase::Tmmult() to do the actual work.

Definition at line 886 of file petsc_parallel_sparse_matrix.cc.

◆ do_reinit() [1/3]

template<typename SparsityPatternType >
void SparseMatrix< SparsityPatternType >::do_reinit ( const MPI_Comm  comm,
const SparsityPatternType &  sparsity_pattern,
const std::vector< size_type > &  local_rows_per_process,
const std::vector< size_type > &  local_columns_per_process,
const unsigned int  this_process,
const bool  preset_nonzero_locations 
)
private

Same as previous functions.

Definition at line 341 of file petsc_parallel_sparse_matrix.cc.

◆ do_reinit() [2/3]

template<typename SparsityPatternType >
void SparseMatrix< SparsityPatternType >::do_reinit ( const MPI_Comm  comm,
const IndexSet local_rows,
const IndexSet local_columns,
const SparsityPatternType &  sparsity_pattern 
)
private

Same as previous functions.

Definition at line 183 of file petsc_parallel_sparse_matrix.cc.

◆ do_reinit() [3/3]

template<typename SparsityPatternType >
void SparseMatrix< SparsityPatternType >::do_reinit ( const MPI_Comm  comm,
const IndexSet local_rows,
const IndexSet local_active_rows,
const IndexSet local_columns,
const IndexSet local_active_columns,
const SparsityPatternType &  sparsity_pattern 
)
private

Same as previous functions, but here we consider active dofs for matrices of IS type.

Definition at line 467 of file petsc_parallel_sparse_matrix.cc.

◆ reinit() [6/6]

void PETScWrappers::MatrixBase::reinit ( Mat  A)
inherited

This method associates the PETSc Mat to the instance of the class. This is particularly useful when performing PETSc to Deal.II operations since it allows to reuse the Deal.II MatrixBase and the PETSc Mat without incurring in memory copies.

Definition at line 92 of file petsc_matrix_base.cc.

◆ clear()

void PETScWrappers::MatrixBase::clear ( )
inherited

Release all memory and return to a state just like after having called the default constructor.

Definition at line 112 of file petsc_matrix_base.cc.

◆ set() [1/5]

void PETScWrappers::MatrixBase::set ( const size_type  i,
const size_type  j,
const PetscScalar  value 
)
inherited

Set the element (i,j) to value.

If the present object (from a derived class of this one) happens to be a sparse matrix, then this function adds a new entry to the matrix if it didn't exist before, very much in contrast to the SparseMatrix class which throws an error if the entry does not exist. If value is not a finite number an exception is thrown.

◆ set() [2/5]

void PETScWrappers::MatrixBase::set ( const std::vector< size_type > &  indices,
const FullMatrix< PetscScalar > &  full_matrix,
const bool  elide_zero_values = false 
)
inherited

Set all elements given in a FullMatrix<double> into the sparse matrix locations given by indices. In other words, this function writes the elements in full_matrix into the calling matrix, using the local-to-global indexing specified by indices for both the rows and the columns of the matrix. This function assumes a quadratic sparse matrix and a quadratic full_matrix, the usual situation in FE calculations.

If the present object (from a derived class of this one) happens to be a sparse matrix, then this function adds some new entries to the matrix if they didn't exist before, very much in contrast to the SparseMatrix class which throws an error if the entry does not exist.

The optional parameter elide_zero_values can be used to specify whether zero values should be inserted anyway or they should be filtered away. The default value is false, i.e., even zero values are inserted/replaced.

◆ set() [3/5]

void PETScWrappers::MatrixBase::set ( const std::vector< size_type > &  row_indices,
const std::vector< size_type > &  col_indices,
const FullMatrix< PetscScalar > &  full_matrix,
const bool  elide_zero_values = false 
)
inherited

Same function as before, but now including the possibility to use rectangular full_matrices and different local-to-global indexing on rows and columns, respectively.

◆ set() [4/5]

void PETScWrappers::MatrixBase::set ( const size_type  row,
const std::vector< size_type > &  col_indices,
const std::vector< PetscScalar > &  values,
const bool  elide_zero_values = false 
)
inherited

Set several elements in the specified row of the matrix with column indices as given by col_indices to the respective value.

If the present object (from a derived class of this one) happens to be a sparse matrix, then this function adds some new entries to the matrix if they didn't exist before, very much in contrast to the SparseMatrix class which throws an error if the entry does not exist.

The optional parameter elide_zero_values can be used to specify whether zero values should be inserted anyway or they should be filtered away. The default value is false, i.e., even zero values are inserted/replaced.

◆ set() [5/5]

void PETScWrappers::MatrixBase::set ( const size_type  row,
const size_type  n_cols,
const size_type col_indices,
const PetscScalar *  values,
const bool  elide_zero_values = false 
)
inherited

Set several elements to values given by values in a given row in columns given by col_indices into the sparse matrix.

If the present object (from a derived class of this one) happens to be a sparse matrix, then this function adds some new entries to the matrix if they didn't exist before, very much in contrast to the SparseMatrix class which throws an error if the entry does not exist.

The optional parameter elide_zero_values can be used to specify whether zero values should be inserted anyway or they should be filtered away. The default value is false, i.e., even zero values are inserted/replaced.

◆ add() [1/6]

void PETScWrappers::MatrixBase::add ( const size_type  i,
const size_type  j,
const PetscScalar  value 
)
inherited

Add value to the element (i,j).

If the present object (from a derived class of this one) happens to be a sparse matrix, then this function adds a new entry to the matrix if it didn't exist before, very much in contrast to the SparseMatrix class which throws an error if the entry does not exist. If value is not a finite number an exception is thrown.

◆ add() [2/6]

void PETScWrappers::MatrixBase::add ( const std::vector< size_type > &  indices,
const FullMatrix< PetscScalar > &  full_matrix,
const bool  elide_zero_values = true 
)
inherited

Add all elements given in a FullMatrix<double> into sparse matrix locations given by indices. In other words, this function adds the elements in full_matrix to the respective entries in calling matrix, using the local-to-global indexing specified by indices for both the rows and the columns of the matrix. This function assumes a quadratic sparse matrix and a quadratic full_matrix, the usual situation in FE calculations.

If the present object (from a derived class of this one) happens to be a sparse matrix, then this function adds some new entries to the matrix if they didn't exist before, very much in contrast to the SparseMatrix class which throws an error if the entry does not exist.

The optional parameter elide_zero_values can be used to specify whether zero values should be added anyway or these should be filtered away and only non-zero data is added. The default value is true, i.e., zero values won't be added into the matrix.

◆ add() [3/6]

void PETScWrappers::MatrixBase::add ( const std::vector< size_type > &  row_indices,
const std::vector< size_type > &  col_indices,
const FullMatrix< PetscScalar > &  full_matrix,
const bool  elide_zero_values = true 
)
inherited

Same function as before, but now including the possibility to use rectangular full_matrices and different local-to-global indexing on rows and columns, respectively.

◆ add() [4/6]

void PETScWrappers::MatrixBase::add ( const size_type  row,
const std::vector< size_type > &  col_indices,
const std::vector< PetscScalar > &  values,
const bool  elide_zero_values = true 
)
inherited

Set several elements in the specified row of the matrix with column indices as given by col_indices to the respective value.

If the present object (from a derived class of this one) happens to be a sparse matrix, then this function adds some new entries to the matrix if they didn't exist before, very much in contrast to the SparseMatrix class which throws an error if the entry does not exist.

The optional parameter elide_zero_values can be used to specify whether zero values should be added anyway or these should be filtered away and only non-zero data is added. The default value is true, i.e., zero values won't be added into the matrix.

◆ add() [5/6]

void PETScWrappers::MatrixBase::add ( const size_type  row,
const size_type  n_cols,
const size_type col_indices,
const PetscScalar *  values,
const bool  elide_zero_values = true,
const bool  col_indices_are_sorted = false 
)
inherited

Add an array of values given by values in the given global matrix row at columns specified by col_indices in the sparse matrix.

If the present object (from a derived class of this one) happens to be a sparse matrix, then this function adds some new entries to the matrix if they didn't exist before, very much in contrast to the SparseMatrix class which throws an error if the entry does not exist.

The optional parameter elide_zero_values can be used to specify whether zero values should be added anyway or these should be filtered away and only non-zero data is added. The default value is true, i.e., zero values won't be added into the matrix.

◆ add() [6/6]

MatrixBase & PETScWrappers::MatrixBase::add ( const PetscScalar  factor,
const MatrixBase other 
)
inherited

Add the matrix other scaled by the factor factor to the current matrix.

Definition at line 516 of file petsc_matrix_base.cc.

◆ clear_row()

void PETScWrappers::MatrixBase::clear_row ( const size_type  row,
const PetscScalar  new_diag_value = 0 
)
inherited

Remove all elements from this row by setting them to zero. The function does not modify the number of allocated nonzero entries, it only sets some entries to zero. It may drop them from the sparsity pattern, though (but retains the allocated memory in case new entries are again added later).

This operation is used in eliminating constraints (e.g. due to hanging nodes) and makes sure that we can write this modification to the matrix without having to read entries (such as the locations of non-zero elements) from it – without this operation, removing constraints on parallel matrices is a rather complicated procedure.

The second parameter can be used to set the diagonal entry of this row to a value different from zero. The default is to set it to zero.

Definition at line 147 of file petsc_matrix_base.cc.

◆ clear_rows()

void PETScWrappers::MatrixBase::clear_rows ( const ArrayView< const size_type > &  rows,
const PetscScalar  new_diag_value = 0 
)
inherited

Same as clear_row(), except that it works on a number of rows at once.

The second parameter can be used to set the diagonal entries of all cleared rows to something different from zero. Note that all of these diagonal entries get the same value – if you want different values for the diagonal entries, you have to set them by hand.

Definition at line 155 of file petsc_matrix_base.cc.

◆ clear_rows_columns()

void PETScWrappers::MatrixBase::clear_rows_columns ( const std::vector< size_type > &  row_and_column_indices,
const PetscScalar  new_diag_value = 0 
)
inherited

Same as clear_rows(), except that the function also zeros the columns.

Definition at line 184 of file petsc_matrix_base.cc.

◆ compress()

void PETScWrappers::MatrixBase::compress ( const VectorOperation::values  operation)
inherited

PETSc matrices store their own sparsity patterns. So, in analogy to our own SparsityPattern class, this function compresses the sparsity pattern and allows the resulting matrix to be used in all other operations where before only assembly functions were allowed. This function must therefore be called once you have assembled the matrix.

See Compressing distributed objects for more information.

Definition at line 244 of file petsc_matrix_base.cc.

◆ operator()()

PetscScalar PETScWrappers::MatrixBase::operator() ( const size_type  i,
const size_type  j 
) const
inherited

Return the value of the entry (i,j). This may be an expensive operation and you should always take care where to call this function. In contrast to the respective function in the MatrixBase class, we don't throw an exception if the respective entry doesn't exist in the sparsity pattern of this class, since PETSc does not transmit this information.

This function is therefore exactly equivalent to the el() function.

◆ el()

PetscScalar PETScWrappers::MatrixBase::el ( const size_type  i,
const size_type  j 
) const
inherited

Return the value of the matrix entry (i,j). If this entry does not exist in the sparsity pattern, then zero is returned. While this may be convenient in some cases, note that it is simple to write algorithms that are slow compared to an optimal solution, since the sparsity of the matrix is not used.

Definition at line 216 of file petsc_matrix_base.cc.

◆ diag_element()

PetscScalar PETScWrappers::MatrixBase::diag_element ( const size_type  i) const
inherited

Return the main diagonal element in the ith row. This function throws an error if the matrix is not quadratic.

Since we do not have direct access to the underlying data structure, this function is no faster than the elementwise access using the el() function. However, we provide this function for compatibility with the SparseMatrix class.

Definition at line 232 of file petsc_matrix_base.cc.

◆ m()

MatrixBase::size_type PETScWrappers::MatrixBase::m ( ) const
inherited

Return the number of rows in this matrix.

Definition at line 286 of file petsc_matrix_base.cc.

◆ n()

MatrixBase::size_type PETScWrappers::MatrixBase::n ( ) const
inherited

Return the number of columns in this matrix.

Definition at line 299 of file petsc_matrix_base.cc.

◆ local_size()

MatrixBase::size_type PETScWrappers::MatrixBase::local_size ( ) const
inherited

Return the local dimension of the matrix, i.e. the number of rows stored on the present MPI process. For sequential matrices, this number is the same as m(), but for parallel matrices it may be smaller.

To figure out which elements exactly are stored locally, use local_range().

Definition at line 312 of file petsc_matrix_base.cc.

◆ local_range()

std::pair< MatrixBase::size_type, MatrixBase::size_type > PETScWrappers::MatrixBase::local_range ( ) const
inherited

Return a pair of indices indicating which rows of this matrix are stored locally. The first number is the index of the first row stored, the second the index of the one past the last one that is stored locally. If this is a sequential matrix, then the result will be the pair (0,m()), otherwise it will be a pair (i,i+n), where n=local_size().

Definition at line 325 of file petsc_matrix_base.cc.

◆ in_local_range()

bool PETScWrappers::MatrixBase::in_local_range ( const size_type  index) const
inherited

Return whether index is in the local range or not, see also local_range().

◆ local_domain_size()

MatrixBase::size_type PETScWrappers::MatrixBase::local_domain_size ( ) const
inherited

Return the local number of columns stored on the present MPI process.

To figure out which elements exactly are stored locally, use local_domain().

Definition at line 339 of file petsc_matrix_base.cc.

◆ local_domain()

std::pair< MatrixBase::size_type, MatrixBase::size_type > PETScWrappers::MatrixBase::local_domain ( ) const
inherited

Return a pair of indices indicating which columns of this matrix are stored locally. The first number is the index of the first column stored, the second the index of the one past the last one that is stored locally.

Definition at line 352 of file petsc_matrix_base.cc.

◆ get_mpi_communicator()

MPI_Comm PETScWrappers::MatrixBase::get_mpi_communicator ( ) const
inherited

Return the underlying MPI communicator.

◆ n_nonzero_elements()

std::uint64_t PETScWrappers::MatrixBase::n_nonzero_elements ( ) const
inherited

Return the number of nonzero elements of this matrix. Actually, it returns the number of entries in the sparsity pattern; if any of the entries should happen to be zero, it is counted anyway.

Definition at line 368 of file petsc_matrix_base.cc.

◆ row_length()

MatrixBase::size_type PETScWrappers::MatrixBase::row_length ( const size_type  row) const
inherited

Number of entries in a specific row.

Definition at line 382 of file petsc_matrix_base.cc.

◆ l1_norm()

PetscReal PETScWrappers::MatrixBase::l1_norm ( ) const
inherited

Return the l1-norm of the matrix, that is \(|M|_1=max_{all columns j}\sum_{all rows i} |M_ij|\), (max. sum of columns). This is the natural matrix norm that is compatible to the l1-norm for vectors, i.e. \(|Mv|_1\leq |M|_1 |v|_1\). (cf. Haemmerlin-Hoffmann: Numerische Mathematik)

Definition at line 418 of file petsc_matrix_base.cc.

◆ linfty_norm()

PetscReal PETScWrappers::MatrixBase::linfty_norm ( ) const
inherited

Return the linfty-norm of the matrix, that is \(|M|_infty=max_{all rows i}\sum_{all columns j} |M_ij|\), (max. sum of rows). This is the natural matrix norm that is compatible to the linfty-norm of vectors, i.e. \(|Mv|_infty \leq |M|_infty |v|_infty\). (cf. Haemmerlin-Hoffmann: Numerische Mathematik)

Definition at line 431 of file petsc_matrix_base.cc.

◆ frobenius_norm()

PetscReal PETScWrappers::MatrixBase::frobenius_norm ( ) const
inherited

Return the frobenius norm of the matrix, i.e. the square root of the sum of squares of all entries in the matrix.

Definition at line 444 of file petsc_matrix_base.cc.

◆ matrix_norm_square() [2/2]

PetscScalar PETScWrappers::MatrixBase::matrix_norm_square ( const VectorBase v) const
inherited

Return the square of the norm of the vector \(v\) with respect to the norm induced by this matrix, i.e. \(\left(v,Mv\right)\). This is useful, e.g. in the finite element context, where the \(L_2\) norm of a function equals the matrix norm with respect to the mass matrix of the vector representing the nodal values of the finite element function.

Obviously, the matrix needs to be quadratic for this operation.

The implementation of this function is not as efficient as the one in the MatrixBase class used in deal.II (i.e. the original one, not the PETSc wrapper class) since PETSc doesn't support this operation and needs a temporary vector.

Note that if the current object represents a parallel distributed matrix (of type PETScWrappers::MPI::SparseMatrix), then the given vector has to be a distributed vector as well. Conversely, if the matrix is not distributed, then neither may the vector be.

Definition at line 456 of file petsc_matrix_base.cc.

◆ matrix_scalar_product() [2/2]

PetscScalar PETScWrappers::MatrixBase::matrix_scalar_product ( const VectorBase u,
const VectorBase v 
) const
inherited

Compute the matrix scalar product \(\left(u,Mv\right)\).

The implementation of this function is not as efficient as the one in the MatrixBase class used in deal.II (i.e. the original one, not the PETSc wrapper class) since PETSc doesn't support this operation and needs a temporary vector.

Note that if the current object represents a parallel distributed matrix (of type PETScWrappers::MPI::SparseMatrix), then both vectors have to be distributed vectors as well. Conversely, if the matrix is not distributed, then neither of the vectors may be.

Definition at line 467 of file petsc_matrix_base.cc.

◆ trace()

PetscScalar PETScWrappers::MatrixBase::trace ( ) const
inherited

Return the trace of the matrix, i.e. the sum of all diagonal entries in the matrix.

Definition at line 480 of file petsc_matrix_base.cc.

◆ operator*=()

MatrixBase & PETScWrappers::MatrixBase::operator*= ( const PetscScalar  factor)
inherited

Multiply the entire matrix by a fixed factor.

Definition at line 493 of file petsc_matrix_base.cc.

◆ operator/=()

MatrixBase & PETScWrappers::MatrixBase::operator/= ( const PetscScalar  factor)
inherited

Divide the entire matrix by a fixed factor.

Definition at line 504 of file petsc_matrix_base.cc.

◆ vmult()

void PETScWrappers::MatrixBase::vmult ( VectorBase dst,
const VectorBase src 
) const
inherited

Matrix-vector multiplication: let dst = M*src with M being this matrix.

Source and destination must not be the same vector.

Note that if the current object represents a parallel distributed matrix (of type PETScWrappers::MPI::SparseMatrix), then both vectors have to be distributed vectors as well. Conversely, if the matrix is not distributed, then neither of the vectors may be.

Definition at line 527 of file petsc_matrix_base.cc.

◆ Tvmult()

void PETScWrappers::MatrixBase::Tvmult ( VectorBase dst,
const VectorBase src 
) const
inherited

Matrix-vector multiplication: let dst = MT*src with M being this matrix. This function does the same as vmult() but takes the transposed matrix.

Source and destination must not be the same vector.

Note that if the current object represents a parallel distributed matrix (of type PETScWrappers::MPI::SparseMatrix), then both vectors have to be distributed vectors as well. Conversely, if the matrix is not distributed, then neither of the vectors may be.

Definition at line 538 of file petsc_matrix_base.cc.

◆ vmult_add()

void PETScWrappers::MatrixBase::vmult_add ( VectorBase dst,
const VectorBase src 
) const
inherited

Adding Matrix-vector multiplication. Add M*src on dst with M being this matrix.

Source and destination must not be the same vector.

Note that if the current object represents a parallel distributed matrix (of type PETScWrappers::MPI::SparseMatrix), then both vectors have to be distributed vectors as well. Conversely, if the matrix is not distributed, then neither of the vectors may be.

Definition at line 549 of file petsc_matrix_base.cc.

◆ Tvmult_add()

void PETScWrappers::MatrixBase::Tvmult_add ( VectorBase dst,
const VectorBase src 
) const
inherited

Adding Matrix-vector multiplication. Add MT*src to dst with M being this matrix. This function does the same as vmult_add() but takes the transposed matrix.

Source and destination must not be the same vector.

Note that if the current object represents a parallel distributed matrix (of type PETScWrappers::MPI::SparseMatrix), then both vectors have to be distributed vectors as well. Conversely, if the matrix is not distributed, then neither of the vectors may be.

Definition at line 560 of file petsc_matrix_base.cc.

◆ residual()

PetscScalar PETScWrappers::MatrixBase::residual ( VectorBase dst,
const VectorBase x,
const VectorBase b 
) const
inherited

Compute the residual of an equation Mx=b, where the residual is defined to be r=b-Mx. Write the residual into dst. The l2 norm of the residual vector is returned.

Source x and destination dst must not be the same vector.

Note that if the current object represents a parallel distributed matrix (of type PETScWrappers::MPI::SparseMatrix), then all vectors have to be distributed vectors as well. Conversely, if the matrix is not distributed, then neither of the vectors may be.

Definition at line 660 of file petsc_matrix_base.cc.

◆ begin() [1/2]

const_iterator PETScWrappers::MatrixBase::begin ( ) const
inherited

Iterator starting at the first entry. This can only be called on a processor owning the entire matrix. In all other cases refer to the version of begin() taking a row number as an argument.

◆ begin() [2/2]

const_iterator PETScWrappers::MatrixBase::begin ( const size_type  r) const
inherited

Iterator starting at the first entry of row r.

Note that if the given row is empty, i.e. does not contain any nonzero entries, then the iterator returned by this function equals end(r). Note also that the iterator may not be dereferenceable in that case.

◆ end() [1/2]

const_iterator PETScWrappers::MatrixBase::end ( ) const
inherited

Final iterator. This can only be called on a processor owning the entire matrix. In all other cases refer to the version of end() taking a row number as an argument.

◆ end() [2/2]

const_iterator PETScWrappers::MatrixBase::end ( const size_type  r) const
inherited

Final iterator of row r. It points to the first element past the end of line r, or past the end of the entire sparsity pattern.

Note that the end iterator is not necessarily dereferenceable. This is in particular the case if it is the end iterator for the last row of a matrix.

◆ operator Mat()

PETScWrappers::MatrixBase::operator Mat ( ) const
inherited

Conversion operator to gain access to the underlying PETSc type. If you do this, you cut this class off some information it may need, so this conversion operator should only be used if you know what you do. In particular, it should only be used for read-only operations into the matrix.

Definition at line 676 of file petsc_matrix_base.cc.

◆ petsc_matrix()

Mat & PETScWrappers::MatrixBase::petsc_matrix ( )
inherited

Return a reference to the underlying PETSc type. It can be used to modify the underlying data, so use it only when you know what you are doing.

Definition at line 682 of file petsc_matrix_base.cc.

◆ transpose()

void PETScWrappers::MatrixBase::transpose ( )
inherited

Make an in-place transpose of a matrix.

Definition at line 688 of file petsc_matrix_base.cc.

◆ is_symmetric()

PetscBool PETScWrappers::MatrixBase::is_symmetric ( const double  tolerance = 1.e-12)
inherited

Test whether a matrix is symmetric. Default tolerance is \(1000\times32\)-bit machine precision.

Definition at line 700 of file petsc_matrix_base.cc.

◆ is_hermitian()

PetscBool PETScWrappers::MatrixBase::is_hermitian ( const double  tolerance = 1.e-12)
inherited

Test whether a matrix is Hermitian, i.e. it is the complex conjugate of its transpose. Default tolerance is \(1000\times32\)-bit machine precision.

Definition at line 710 of file petsc_matrix_base.cc.

◆ write_ascii()

void PETScWrappers::MatrixBase::write_ascii ( const PetscViewerFormat  format = PETSC_VIEWER_DEFAULT)
inherited

Print the PETSc matrix object values using PETSc internal matrix viewer function MatView. The default format prints the non- zero matrix elements. For other valid view formats, consult http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatView.html

Definition at line 722 of file petsc_matrix_base.cc.

◆ print()

void PETScWrappers::MatrixBase::print ( std::ostream &  out,
const bool  alternative_output = false 
) const
inherited

Print the elements of a matrix to the given output stream.

Parameters
[in,out]outThe output stream to which to write.
[in]alternative_outputThis argument is ignored. It exists for compatibility with similar functions in other matrix classes.

Definition at line 738 of file petsc_matrix_base.cc.

◆ memory_consumption()

std::size_t PETScWrappers::MatrixBase::memory_consumption ( ) const
inherited

Return the number bytes consumed by this matrix on this CPU.

Definition at line 785 of file petsc_matrix_base.cc.

◆ prepare_action()

void PETScWrappers::MatrixBase::prepare_action ( const VectorOperation::values  new_action)
protectedinherited

Ensure that the add/set mode that is required for actions following this call is compatible with the current mode. Should be called from all internal functions accessing matrix elements.

◆ assert_is_compressed()

void PETScWrappers::MatrixBase::assert_is_compressed ( )
protectedinherited

Internal function that checks that there are no pending insert/add operations. Throws an exception otherwise. Useful before calling any PETSc internal functions modifying the matrix.

◆ prepare_add()

void PETScWrappers::MatrixBase::prepare_add ( )
protectedinherited

For some matrix storage formats, in particular for the PETSc distributed blockmatrices, set and add operations on individual elements can not be freely mixed. Rather, one has to synchronize operations when one wants to switch from setting elements to adding to elements. BlockMatrixBase automatically synchronizes the access by calling this helper function for each block. This function ensures that the matrix is in a state that allows adding elements; if it previously already was in this state, the function does nothing.

◆ prepare_set()

void PETScWrappers::MatrixBase::prepare_set ( )
protectedinherited

Same as prepare_add() but prepare the matrix for setting elements if the representation of elements in this class requires such an operation.

◆ mmult() [2/2]

void PETScWrappers::MatrixBase::mmult ( MatrixBase C,
const MatrixBase B,
const VectorBase V 
) const
protectedinherited

Base function to perform the matrix-matrix multiplication \(C = AB\), or, if a vector \(V\) whose size is compatible with B is given, \(C = A \text{diag}(V) B\), where \(\text{diag}(V)\) defines a diagonal matrix with the vector entries.

This function assumes that the calling matrix \(A\) and \(B\) have compatible sizes. The size of \(C\) will be set within this function.

The content as well as the sparsity pattern of the matrix \(C\) will be reset by this function, so make sure that the sparsity pattern is not used somewhere else in your program. This is an expensive operation, so think twice before you use this function.

Definition at line 644 of file petsc_matrix_base.cc.

◆ Tmmult() [2/2]

void PETScWrappers::MatrixBase::Tmmult ( MatrixBase C,
const MatrixBase B,
const VectorBase V 
) const
protectedinherited

Base function to perform the matrix-matrix multiplication with the transpose of this, i.e., \(C = A^T B\), or, if an optional vector \(V\) whose size is compatible with \(B\) is given, \(C = A^T \text{diag}(V) B\), where \(\text{diag}(V)\) defines a diagonal matrix with the vector entries.

This function assumes that the calling matrix \(A\) and \(B\) have compatible sizes. The size of \(C\) will be set within this function.

The content as well as the sparsity pattern of the matrix \(C\) will be changed by this function, so make sure that the sparsity pattern is not used somewhere else in your program. This is an expensive operation, so think twice before you use this function.

Definition at line 652 of file petsc_matrix_base.cc.

◆ subscribe()

void Subscriptor::subscribe ( std::atomic< bool > *const  validity,
const std::string &  identifier = "" 
) const
inherited

Subscribes a user of the object by storing the pointer validity. The subscriber may be identified by text supplied as identifier.

Definition at line 135 of file subscriptor.cc.

◆ unsubscribe()

void Subscriptor::unsubscribe ( std::atomic< bool > *const  validity,
const std::string &  identifier = "" 
) const
inherited

Unsubscribes a user from the object.

Note
The identifier and the validity pointer must be the same as the one supplied to subscribe().

Definition at line 155 of file subscriptor.cc.

◆ n_subscriptions()

unsigned int Subscriptor::n_subscriptions ( ) const
inlineinherited

Return the present number of subscriptions to this object. This allows to use this class for reference counted lifetime determination where the last one to unsubscribe also deletes the object.

Definition at line 300 of file subscriptor.h.

◆ list_subscribers() [1/2]

template<typename StreamType >
void Subscriptor::list_subscribers ( StreamType &  stream) const
inlineinherited

List the subscribers to the input stream.

Definition at line 317 of file subscriptor.h.

◆ list_subscribers() [2/2]

void Subscriptor::list_subscribers ( ) const
inherited

List the subscribers to deallog.

Definition at line 203 of file subscriptor.cc.

◆ serialize()

template<class Archive >
void Subscriptor::serialize ( Archive &  ar,
const unsigned int  version 
)
inlineinherited

Read or write the data of this object to or from a stream for the purpose of serialization using the BOOST serialization library.

This function does not actually serialize any of the member variables of this class. The reason is that what this class stores is only who subscribes to this object, but who does so at the time of storing the contents of this object does not necessarily have anything to do with who subscribes to the object when it is restored. Consequently, we do not want to overwrite the subscribers at the time of restoring, and then there is no reason to write the subscribers out in the first place.

Definition at line 309 of file subscriptor.h.

◆ check_no_subscribers()

void Subscriptor::check_no_subscribers ( ) const
privatenoexceptinherited

Check that there are no objects subscribing to this object. If this check passes then it is safe to destroy the current object. It this check fails then this function will either abort or print an error message to deallog (by using the AssertNothrow mechanism), but will not throw an exception.

Note
Since this function is just a consistency check it does nothing in release mode.
If this function is called when there is an uncaught exception then, rather than aborting, this function prints an error message to the standard error stream and returns.

Definition at line 52 of file subscriptor.cc.

Friends And Related Symbol Documentation

◆ BlockMatrixBase< SparseMatrix >

friend class BlockMatrixBase< SparseMatrix >
friend

Definition at line 651 of file petsc_sparse_matrix.h.

Member Data Documentation

◆ matrix

Mat PETScWrappers::MatrixBase::matrix
protectedinherited

A generic matrix object in PETSc. The actual type, a sparse matrix, is set in the constructor.

Definition at line 1011 of file petsc_matrix_base.h.

◆ last_action

VectorOperation::values PETScWrappers::MatrixBase::last_action
protectedinherited

Store whether the last action was a write or add operation.

Definition at line 1016 of file petsc_matrix_base.h.

◆ column_indices

std::vector<PetscInt> PETScWrappers::MatrixBase::column_indices
mutableprivateinherited

An internal array of integer values that is used to store the column indices when adding/inserting local data into the (large) sparse matrix.

This variable does not store any "state" of the matrix object. Rather, it is only used as a temporary buffer by some of the member functions of this class. As with all mutable member variables, the use of this variable is not thread-safe unless guarded by a mutex. However, since PETSc matrix operations are not thread-safe anyway, there is no need to attempt to make things thread-safe, and so there is no mutex associated with this variable.

Definition at line 1106 of file petsc_matrix_base.h.

◆ column_values

std::vector<PetscScalar> PETScWrappers::MatrixBase::column_values
mutableprivateinherited

An internal array of double values that is used to store the column indices when adding/inserting local data into the (large) sparse matrix.

The same comment as for the column_indices variable above applies.

Definition at line 1116 of file petsc_matrix_base.h.

◆ counter

std::atomic<unsigned int> Subscriptor::counter
mutableprivateinherited

Store the number of objects which subscribed to this object. Initially, this number is zero, and upon destruction it shall be zero again (i.e. all objects which subscribed should have unsubscribed again).

The creator (and owner) of an object is counted in the map below if HE manages to supply identification.

We use the mutable keyword in order to allow subscription to constant objects also.

This counter may be read from and written to concurrently in multithreaded code: hence we use the std::atomic class template.

Definition at line 218 of file subscriptor.h.

◆ counter_map

std::map<std::string, unsigned int> Subscriptor::counter_map
mutableprivateinherited

In this map, we count subscriptions for each different identification string supplied to subscribe().

Definition at line 224 of file subscriptor.h.

◆ validity_pointers

std::vector<std::atomic<bool> *> Subscriptor::validity_pointers
mutableprivateinherited

In this vector, we store pointers to the validity bool in the SmartPointer objects that subscribe to this class.

Definition at line 240 of file subscriptor.h.

◆ object_info

const std::type_info* Subscriptor::object_info
mutableprivateinherited

Pointer to the typeinfo object of this object, from which we can later deduce the class name. Since this information on the derived class is neither available in the destructor, nor in the constructor, we obtain it in between and store it here.

Definition at line 248 of file subscriptor.h.

◆ mutex

std::mutex Subscriptor::mutex
staticprivateinherited

A mutex used to ensure data consistency when accessing the mutable members of this class. This lock is used in the subscribe() and unsubscribe() functions, as well as in list_subscribers().

Definition at line 271 of file subscriptor.h.


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