Reference documentation for deal.II version 9.4.1
|
#include <deal.II/matrix_free/operators.h>
Public Types | |
using | value_type = typename Base< dim, VectorType, VectorizedArrayType >::value_type |
using | size_type = typename Base< dim, VectorType, VectorizedArrayType >::size_type |
Public Member Functions | |
MassOperator () | |
virtual void | compute_diagonal () override |
void | compute_lumped_diagonal () |
const std::shared_ptr< DiagonalMatrix< VectorType > > & | get_matrix_lumped_diagonal () const |
const std::shared_ptr< DiagonalMatrix< VectorType > > & | get_matrix_lumped_diagonal_inverse () const |
virtual void | clear () |
void | initialize (std::shared_ptr< const MatrixFree< dim, value_type, VectorizedArrayType > > data, const std::vector< unsigned int > &selected_row_blocks=std::vector< unsigned int >(), const std::vector< unsigned int > &selected_column_blocks=std::vector< unsigned int >()) |
void | initialize (std::shared_ptr< const MatrixFree< dim, value_type, VectorizedArrayType > > data, const MGConstrainedDoFs &mg_constrained_dofs, const unsigned int level, const std::vector< unsigned int > &selected_row_blocks=std::vector< unsigned int >()) |
void | initialize (std::shared_ptr< const MatrixFree< dim, value_type, VectorizedArrayType > > data_, const std::vector< MGConstrainedDoFs > &mg_constrained_dofs, const unsigned int level, const std::vector< unsigned int > &selected_row_blocks=std::vector< unsigned int >()) |
size_type | m () const |
size_type | n () const |
void | vmult_interface_down (VectorType &dst, const VectorType &src) const |
void | vmult_interface_up (VectorType &dst, const VectorType &src) const |
void | vmult (VectorType &dst, const VectorType &src) const |
void | Tvmult (VectorType &dst, const VectorType &src) const |
void | vmult_add (VectorType &dst, const VectorType &src) const |
void | Tvmult_add (VectorType &dst, const VectorType &src) const |
value_type | el (const unsigned int row, const unsigned int col) const |
virtual std::size_t | memory_consumption () const |
void | initialize_dof_vector (VectorType &vec) const |
std::shared_ptr< const MatrixFree< dim, value_type, VectorizedArrayType > > | get_matrix_free () const |
const std::shared_ptr< DiagonalMatrix< VectorType > > & | get_matrix_diagonal_inverse () const |
const std::shared_ptr< DiagonalMatrix< VectorType > > & | get_matrix_diagonal () const |
void | precondition_Jacobi (VectorType &dst, const VectorType &src, const value_type omega) const |
Protected Member Functions | |
void | preprocess_constraints (VectorType &dst, const VectorType &src) const |
void | postprocess_constraints (VectorType &dst, const VectorType &src) const |
void | set_constrained_entries_to_one (VectorType &dst) const |
virtual void | Tapply_add (VectorType &dst, const VectorType &src) const |
Protected Attributes | |
std::shared_ptr< const MatrixFree< dim, value_type, VectorizedArrayType > > | data |
std::shared_ptr< DiagonalMatrix< VectorType > > | diagonal_entries |
std::shared_ptr< DiagonalMatrix< VectorType > > | inverse_diagonal_entries |
std::vector< unsigned int > | selected_rows |
std::vector< unsigned int > | selected_columns |
Private Member Functions | |
virtual void | apply_add (VectorType &dst, const VectorType &src) const override |
void | local_apply_cell (const MatrixFree< dim, value_type, VectorizedArrayType > &data, VectorType &dst, const VectorType &src, const std::pair< unsigned int, unsigned int > &cell_range) const |
void | mult_add (VectorType &dst, const VectorType &src, const bool transpose) const |
void | adjust_ghost_range_if_necessary (const VectorType &vec, const bool is_row) const |
Private Attributes | |
std::shared_ptr< DiagonalMatrix< VectorType > > | lumped_diagonal_entries |
std::shared_ptr< DiagonalMatrix< VectorType > > | inverse_lumped_diagonal_entries |
std::vector< std::vector< unsigned int > > | edge_constrained_indices |
std::vector< std::vector< std::pair< value_type, value_type > > > | edge_constrained_values |
bool | have_interface_matrices |
Subscriptor functionality | |
Classes derived from Subscriptor provide a facility to subscribe to this object. This is mostly used by the SmartPointer class. | |
std::atomic< unsigned int > | counter |
std::map< std::string, unsigned int > | counter_map |
std::vector< std::atomic< bool > * > | validity_pointers |
const std::type_info * | object_info |
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 |
template<class Archive > | |
void | serialize (Archive &ar, const unsigned int version) |
void | check_no_subscribers () const noexcept |
using | map_value_type = decltype(counter_map)::value_type |
using | map_iterator = decltype(counter_map)::iterator |
static std::mutex | mutex |
static ::ExceptionBase & | ExcInUse (int arg1, std::string arg2, std::string arg3) |
static ::ExceptionBase & | ExcNoSubscriber (std::string arg1, std::string arg2) |
This class implements the operation of the action of a mass matrix.
Note that this class only supports the non-blocked vector variant of the Base operator because only a single FEEvaluation object is used in the apply function.
Definition at line 744 of file operators.h.
using MatrixFreeOperators::MassOperator< dim, fe_degree, n_q_points_1d, n_components, VectorType, VectorizedArrayType >::value_type = typename Base<dim, VectorType, VectorizedArrayType>::value_type |
Number alias.
Definition at line 750 of file operators.h.
using MatrixFreeOperators::MassOperator< dim, fe_degree, n_q_points_1d, n_components, VectorType, VectorizedArrayType >::size_type = typename Base<dim, VectorType, VectorizedArrayType>::size_type |
size_type needed for preconditioner classes.
Definition at line 756 of file operators.h.
MatrixFreeOperators::MassOperator< dim, fe_degree, n_q_points_1d, n_components, VectorType, VectorizedArrayType >::MassOperator |
Constructor.
Definition at line 1832 of file operators.h.
|
overridevirtual |
Same as the base class.
Implements MatrixFreeOperators::Base< dim, VectorType, VectorizedArrayType >.
Definition at line 1850 of file operators.h.
void MatrixFreeOperators::MassOperator< dim, fe_degree, n_q_points_1d, n_components, VectorType, VectorizedArrayType >::compute_lumped_diagonal |
Compute the lumped mass matrix. This is equal to the mass matrix times a vector of all ones and is equivalent to approximating the mass matrix with a nodal quadrature rule.
The lumped mass matrix is an excellent preconditioner for mass matrices corresponding to FE_Q elements on axis-aligned cells. However, some elements (like FE_SimplexP with degrees higher than 1) have basis functions whose integrals are zero or negative (and therefore their lumped mass matrix entries are zero or negative). For such elements a lumped mass matrix is a very poor approximation of the operator - the diagonal should be used instead. If you are interested in using mass lumping with simplices then use FE_SimplexP_Bubbles instead of FE_SimplexP.
Definition at line 1929 of file operators.h.
const std::shared_ptr< DiagonalMatrix< VectorType > > & MatrixFreeOperators::MassOperator< dim, fe_degree, n_q_points_1d, n_components, VectorType, VectorizedArrayType >::get_matrix_lumped_diagonal |
Get read access to the lumped diagonal of this operator.
Definition at line 2008 of file operators.h.
const std::shared_ptr< DiagonalMatrix< VectorType > > & MatrixFreeOperators::MassOperator< dim, fe_degree, n_q_points_1d, n_components, VectorType, VectorizedArrayType >::get_matrix_lumped_diagonal_inverse |
Get read access to the inverse lumped diagonal of this operator.
Definition at line 1986 of file operators.h.
|
overrideprivatevirtual |
Applies the mass matrix operation on an input vector. It is assumed that the passed input and output vector are correctly initialized using initialize_dof_vector().
Implements MatrixFreeOperators::Base< dim, VectorType, VectorizedArrayType >.
Definition at line 2030 of file operators.h.
|
private |
For this operator, there is just a cell contribution.
Definition at line 2051 of file operators.h.
|
virtualinherited |
Release all memory and return to a state just like after having called the default constructor.
Reimplemented in MatrixFreeOperators::LaplaceOperator< dim, fe_degree, n_q_points_1d, n_components, VectorType, VectorizedArrayType >.
Definition at line 1188 of file operators.h.
|
inherited |
Initialize operator on fine scale.
The optional selection vector allows to choose only some components from the underlying MatrixFree object, e.g. just a single one. The entry selected_row_blocks
[i] in the vector chooses the DoFHandler and AffineConstraints object that was given as the selected_row_blocks
[i]-th argument to the MatrixFree::reinit() call. Different arguments for rows and columns also make it possible to select non-diagonal blocks or rectangular blocks. If the row vector is empty, all components are selected, otherwise its size must be smaller or equal to MatrixFree::n_components() and all indices need to be unique and within the range of 0 and MatrixFree::n_components(). If the column selection vector is empty, it is taken the same as the row selection, defining a diagonal block.
Definition at line 1238 of file operators.h.
|
inherited |
Initialize operator on a level level
for a single FiniteElement.
The optional selection vector allows to choose only some components from the underlying MatrixFree object, e.g. just a single one. The entry selected_row_blocks
[i] in the vector chooses the DoFHandler and AffineConstraints object that was given as the selected_row_blocks
[i]-th argument to the MatrixFree::reinit() call. Since a multigrid operator is always associated to inverting a matrix and thus represents a diagonal block, the same vector for rows and columns is used as opposed to the non-level initialization function. If empty, all components are selected.
Definition at line 1291 of file operators.h.
|
inherited |
Initialize operator on a level level
for multiple FiniteElement objects.
The optional selection vector allows to choose only some components from the underlying MatrixFree object, e.g. just a single one. The entry selected_row_blocks
[i] in the vector chooses the DoFHandler and AffineConstraints object that was given as the selected_row_blocks
[i]-th argument to the MatrixFree::reinit() call. Since a multigrid operator is always associated to inverting a matrix and thus represents a diagonal block, the same vector for rows and columns is used as opposed to the non-level initialization function. If empty, all components are selected.
Definition at line 1307 of file operators.h.
|
inherited |
Return the dimension of the codomain (or range) space.
Definition at line 1160 of file operators.h.
|
inherited |
Return the dimension of the domain space.
Definition at line 1174 of file operators.h.
|
inherited |
vmult operator for interface.
Definition at line 1557 of file operators.h.
|
inherited |
vmult operator for interface.
Definition at line 1611 of file operators.h.
|
inherited |
Matrix-vector multiplication.
Definition at line 1396 of file operators.h.
|
inherited |
Transpose matrix-vector multiplication.
Definition at line 1652 of file operators.h.
|
inherited |
Adding Matrix-vector multiplication.
Definition at line 1409 of file operators.h.
|
inherited |
Adding transpose matrix-vector multiplication.
Definition at line 1420 of file operators.h.
|
inherited |
Return the value of the matrix entry (row,col). In matrix-free context this function is valid only for row==col when diagonal is initialized.
Definition at line 1198 of file operators.h.
|
virtualinherited |
Determine an estimate for the memory consumption (in bytes) of this object.
Definition at line 1666 of file operators.h.
|
inherited |
A wrapper for initialize_dof_vector() of MatrixFree object.
Definition at line 1213 of file operators.h.
|
inherited |
Get read access to the MatrixFree object stored with this operator.
Definition at line 1680 of file operators.h.
|
inherited |
Get read access to the inverse diagonal of this operator.
Definition at line 1689 of file operators.h.
|
inherited |
Get read access to the diagonal of this operator.
Definition at line 1702 of file operators.h.
|
inherited |
Apply the Jacobi preconditioner, which multiplies every element of the src
vector by the inverse of the respective diagonal element and multiplies the result with the relaxation factor omega
.
Definition at line 1724 of file operators.h.
|
protectedinherited |
Perform necessary operations related to constraints before calling apply_add() or Tapply_add() inside mult_add().
Definition at line 1473 of file operators.h.
|
protectedinherited |
Perform necessary operations related to constraints after calling apply_add() or Tapply_add() inside mult_add().
Definition at line 1523 of file operators.h.
|
protectedinherited |
Set constrained entries (both from hanging nodes and edge constraints) of dst
to one.
Definition at line 1377 of file operators.h.
|
protectedvirtualinherited |
Apply transpose operator to src
and add result in dst
.
Default implementation is to call apply_add().
Definition at line 1713 of file operators.h.
|
privateinherited |
Function which implements vmult_add (transpose
= false) and Tvmult_add (transpose
= true).
Definition at line 1503 of file operators.h.
|
privateinherited |
Adjust the ghost range of the vectors to the storage requirements of the underlying MatrixFree class. This is used inside the mult_add() as well as vmult_interface_up() and vmult_interface_down() methods in order to ensure that the cell loops will be able to access the ghost indices with the correct local indices.
Definition at line 1431 of file operators.h.
|
private |
A shared pointer to a diagonal matrix that stores the lumped diagonal elements as a vector.
Definition at line 823 of file operators.h.
|
private |
A shared pointer to a diagonal matrix that stores the inverse of lumped diagonal elements as a vector.
Definition at line 829 of file operators.h.
|
protectedinherited |
MatrixFree object to be used with this operator.
Definition at line 438 of file operators.h.
|
protectedinherited |
A shared pointer to a diagonal matrix that stores the diagonal elements as a vector.
Definition at line 444 of file operators.h.
|
protectedinherited |
A shared pointer to a diagonal matrix that stores the inverse of diagonal elements as a vector.
Definition at line 450 of file operators.h.
|
protectedinherited |
A vector which defines the selection of sub-components of MatrixFree for the rows of the matrix representation.
Definition at line 456 of file operators.h.
|
protectedinherited |
A vector which defines the selection of sub-components of MatrixFree for the columns of the matrix representation.
Definition at line 462 of file operators.h.
|
privateinherited |
Indices of DoFs on edge in case the operator is used in GMG context.
Definition at line 468 of file operators.h.
|
mutableprivateinherited |
Auxiliary vector.
Definition at line 474 of file operators.h.
|
privateinherited |
A flag which determines whether or not this operator has interface matrices in GMG context.
Definition at line 480 of file operators.h.