Ifpack2 Templated Preconditioning Package
Version 1.0
|
#include <Ifpack2_BlockTriDiContainer_decl.hpp>
Classes | |
struct | ApplyParameters |
Input arguments to applyInverseJacobi More... | |
Public Member Functions | |
Constructor and destructor | |
BlockTriDiContainer (const Teuchos::RCP< const row_matrix_type > &matrix, const Teuchos::Array< Teuchos::Array< local_ordinal_type > > &partitions, const Teuchos::RCP< const import_type > &importer, int OverlapLevel, scalar_type DampingFactor) | |
Constructor. More... | |
BlockTriDiContainer (const Teuchos::RCP< const row_matrix_type > &matrix, const Teuchos::Array< Teuchos::Array< local_ordinal_type > > &partitions, bool overlapCommAndComp=false, bool useSequentialMethod=false) | |
Constructor for bare construction. More... | |
~BlockTriDiContainer () override | |
Destructor (declared virtual for memory safety of derived classes). More... | |
Get and set methods | |
bool | isInitialized () const override |
Whether the container has been successfully initialized. More... | |
bool | isComputed () const override |
Whether the container has been successfully computed. More... | |
void | setParameters (const Teuchos::ParameterList &List) override |
Set all necessary parameters. More... | |
void | clearBlocks () override |
Mathematical functions | |
void | initialize () override |
Do all set-up operations that only require matrix structure. More... | |
void | compute () override |
Extract the local tridiagonal block and prepare the solver. More... | |
void | applyInverseJacobi (const mv_type &X, mv_type &Y, bool zeroStartingSolution=false, int numSweeps=1) const override |
Compute Y := (1 - a) Y + a D^{-1} (X - R*Y) . More... | |
ComputeParameters | createDefaultComputeParameters () const |
Create a ComputeParameters struct with default values. More... | |
void | compute (const ComputeParameters &input) |
Extract the local tridiagonal block and prepare the solver. More... | |
ApplyParameters | createDefaultApplyParameters () const |
Create an ApplyParameters struct with default values. More... | |
int | applyInverseJacobi (const mv_type &X, mv_type &Y, const ApplyParameters &input) const |
Compute Y := D^{-1} (X - R*Y) . More... | |
const Teuchos::ArrayRCP< const magnitude_type > | getNorms0 () const |
If a norm-based method was used, return a buffer containing the norms, ordered by degrees of freedom, then RHS; otherwise, return a null ArrayRCp. More... | |
const Teuchos::ArrayRCP< const magnitude_type > | getNormsFinal () const |
If a norm-based method was used, return a buffer containing the norms, ordered by degrees of freedom, then RHS; otherwise, return a null ArrayRCp. More... | |
void | apply (host_view_type &X, host_view_type &Y, int blockIndex, int stride, Teuchos::ETransp mode=Teuchos::NO_TRANS, scalar_type alpha=Teuchos::ScalarTraits< scalar_type >::one(), scalar_type beta=Teuchos::ScalarTraits< scalar_type >::zero()) const override |
void | weightedApply (host_view_type &X, host_view_type &Y, host_view_type &W, int blockIndex, int stride, Teuchos::ETransp mode=Teuchos::NO_TRANS, scalar_type alpha=Teuchos::ScalarTraits< scalar_type >::one(), scalar_type beta=Teuchos::ScalarTraits< scalar_type >::zero()) const override |
Miscellaneous methods | |
std::ostream & | print (std::ostream &os) const override |
Print information about this object to the given output stream. More... | |
Implementation of Teuchos::Describable | |
std::string | description () const override |
A one-line description of this object. More... | |
void | describe (Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const override |
Print the object with some verbosity level to the given FancyOStream. More... | |
Public Member Functions inherited from Ifpack2::Container< MatrixType > | |
Container (const Teuchos::RCP< const row_matrix_type > &matrix, const Teuchos::Array< Teuchos::Array< local_ordinal_type > > &partitions, const Teuchos::RCP< const import_type > &importer, int OverlapLevel, scalar_type DampingFactor) | |
Constructor. More... | |
Container (const Teuchos::RCP< const row_matrix_type > &matrix, const Teuchos::Array< local_ordinal_type > &localRows) | |
Constructor for single block. More... | |
virtual | ~Container () |
Destructor. More... | |
Teuchos::ArrayView< const local_ordinal_type > | getLocalRows (int blockIndex) const |
Local indices of the rows of the input matrix that belong to this block. More... | |
void | setBlockSizes (const Teuchos::Array< Teuchos::Array< local_ordinal_type > > &partitions) |
Initialize arrays with information about block sizes. More... | |
void | applyMV (mv_type &X, mv_type &Y) const |
Wrapper for apply with MVs, used in unit tests (never called by BlockRelaxation) More... | |
void | weightedApplyMV (mv_type &X, mv_type &Y, vector_type &W) |
Wrapper for weightedApply with MVs, used in unit tests (never called by BlockRelaxation) More... | |
Static Public Member Functions | |
static std::string | getName () |
Get the name of this container type for Details::constructContainer() More... | |
Static Public Member Functions inherited from Ifpack2::Container< MatrixType > | |
static std::string | getName () |
Additional Inherited Members | |
Protected Types inherited from Ifpack2::Container< MatrixType > | |
typedef MatrixType::scalar_type | scalar_type |
typedef MatrixType::local_ordinal_type | local_ordinal_type |
typedef MatrixType::global_ordinal_type | global_ordinal_type |
typedef MatrixType::node_type | node_type |
typedef Tpetra::MultiVector < scalar_type, local_ordinal_type, global_ordinal_type, node_type > | mv_type |
typedef Tpetra::Vector < scalar_type, local_ordinal_type, global_ordinal_type, node_type > | vector_type |
typedef Tpetra::Map < local_ordinal_type, global_ordinal_type, node_type > | map_type |
typedef Teuchos::ScalarTraits < scalar_type > | STS |
typedef Tpetra::Import < local_ordinal_type, global_ordinal_type, node_type > | import_type |
typedef Partitioner < Tpetra::RowGraph < local_ordinal_type, global_ordinal_type, node_type > > | partitioner_type |
typedef Tpetra::Experimental::BlockCrsMatrix < scalar_type, local_ordinal_type, global_ordinal_type, node_type > | block_crs_matrix_type |
typedef Tpetra::RowMatrix < scalar_type, local_ordinal_type, global_ordinal_type, node_type > | row_matrix_type |
typedef Kokkos::Details::ArithTraits < scalar_type >::val_type | impl_scalar_type |
Internal representation of Scalar in Kokkos::View. More... | |
Protected Attributes inherited from Ifpack2::Container< MatrixType > | |
Teuchos::RCP< const row_matrix_type > | inputMatrix_ |
The input matrix to the constructor. More... | |
int | numBlocks_ |
The number of blocks (partitions) in the container. More... | |
Teuchos::Array < local_ordinal_type > | partitions_ |
Local indices of the rows of the input matrix that belong to this block. More... | |
Teuchos::Array < local_ordinal_type > | blockRows_ |
Number of rows in each block. More... | |
Teuchos::Array < local_ordinal_type > | partitionIndices_ |
Starting index in partitions_ of local row indices for each block. More... | |
Teuchos::RCP< vector_type > | Diag_ |
Diagonal elements. More... | |
bool | IsParallel_ |
Whether the problem is distributed across multiple MPI processes. More... | |
int | OverlapLevel_ |
Number of rows of overlap for adjacent blocks. More... | |
scalar_type | DampingFactor_ |
Damping factor, passed to apply() as alpha. More... | |
Teuchos::RCP< const Tpetra::Import < local_ordinal_type, global_ordinal_type, node_type > > | Importer_ |
Importer for importing off-process elements of MultiVectors. More... | |
local_ordinal_type | NumLocalRows_ |
Number of local rows in input matrix. More... | |
global_ordinal_type | NumGlobalRows_ |
Number of global rows in input matrix. More... | |
global_ordinal_type | NumGlobalNonzeros_ |
Number of nonzeros in input matrix. More... | |
bool | hasBlockCrs_ |
Whether the input matrix is a BlockCRS matrix. More... | |
int | bcrsBlockSize_ |
If hasBlockCrs_, the number of DOFs per vertex. Otherwise 1. More... | |
Partial specialization with SIMD<ScalarType> internally used. The code does not support Sacado UQ types. As UQ type also uses SIMD like structure, it needs an additional specialization.
Ifpack2::BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::BlockTriDiContainer | ( | const Teuchos::RCP< const row_matrix_type > & | matrix, |
const Teuchos::Array< Teuchos::Array< local_ordinal_type > > & | partitions, | ||
const Teuchos::RCP< const import_type > & | importer, | ||
int | OverlapLevel, | ||
scalar_type | DampingFactor | ||
) |
Constructor.
matrix [in] The original input matrix. This Container will construct a local diagonal block from the rows given by localRows
.
partitions | [in] The set of (local) rows assigned to this container. partitions[i] == j , where i (from 0 to getNumRows() - 1 ) indicates the SparseContainer's row, and j indicates the local row in the calling process. partitions.size() gives the number of rows in the local matrix on each process. This may be different on different processes. If partitions is empty, then a local block-Jacobi preconditioner is formed; this is equivalent to every part of partitions having size 1. |
Ifpack2::BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::BlockTriDiContainer | ( | const Teuchos::RCP< const row_matrix_type > & | matrix, |
const Teuchos::Array< Teuchos::Array< local_ordinal_type > > & | partitions, | ||
bool | overlapCommAndComp = false , |
||
bool | useSequentialMethod = false |
||
) |
Constructor for bare construction.
Constructor for users who want to use BlockTriDiContainer directly. See main constructor for documentation. This constructor removes the Container-general arguments that are not used in BlockTriDiContainer.
overlapCommAndComp | [in] Overlap communication and computation. This is not always better; it depends on (at least) the MPI implementation and the machine architecture. Defaults to false. It has to be specified at construction because the core data structures depend on whether overlapCommAndComp is true or false. |
useSequentialMethod | [in] This is for development and testing purposes only. |
|
override |
Destructor (declared virtual for memory safety of derived classes).
|
overridevirtual |
Whether the container has been successfully initialized.
Implements Ifpack2::Container< MatrixType >.
|
overridevirtual |
Whether the container has been successfully computed.
Implements Ifpack2::Container< MatrixType >.
|
overridevirtual |
Set all necessary parameters.
Implements Ifpack2::Container< MatrixType >.
|
overridevirtual |
Do all set-up operations that only require matrix structure.
Implements Ifpack2::Container< MatrixType >.
|
overridevirtual |
Extract the local tridiagonal block and prepare the solver.
Implements Ifpack2::Container< MatrixType >.
|
overridevirtual |
Compute Y := (1 - a) Y + a D^{-1} (X - R*Y)
.
Reimplemented from Ifpack2::Container< MatrixType >.
BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::ComputeParameters Ifpack2::BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::createDefaultComputeParameters | ( | ) | const |
Create a ComputeParameters struct with default values.
void Ifpack2::BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::compute | ( | const ComputeParameters & | input | ) |
Extract the local tridiagonal block and prepare the solver.
This version of applyInverseJacobi
is meant to be called by direct users of this class, rather than by BlockRelaxation
.
addRadiallyToDiagonal | [in] Add a constant to each diagonal entry of the matrix. It is added according to the (possibly complex) sign of the digonal entry at the time that the entry is the pivot in the LU factorization, such that the entry is moved radially outward in the complex plane. N.B. that this constant modifies the matrix in the linear equation, not simply the diagonal preconditioner. |
BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::ApplyParameters Ifpack2::BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::createDefaultApplyParameters | ( | ) | const |
Create an ApplyParameters struct with default values.
int Ifpack2::BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::applyInverseJacobi | ( | const mv_type & | X, |
mv_type & | Y, | ||
const ApplyParameters & | input | ||
) | const |
Compute Y := D^{-1} (X - R*Y)
.
This version of applyInverseJacobi
is meant to be called by direct users of this class, rather than by BlockRelaxation
. It supports 2-norm-based termination. It returns the number of sweeps performed.
const Teuchos::ArrayRCP< const typename BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::magnitude_type > Ifpack2::BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::getNorms0 | ( | ) | const |
If a norm-based method was used, return a buffer containing the norms, ordered by degrees of freedom, then RHS; otherwise, return a null ArrayRCp.
const Teuchos::ArrayRCP< const typename BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::magnitude_type > Ifpack2::BlockTriDiContainer< MatrixType, BlockTriDiContainerDetails::ImplSimdTag >::getNormsFinal | ( | ) | const |
If a norm-based method was used, return a buffer containing the norms, ordered by degrees of freedom, then RHS; otherwise, return a null ArrayRCp.
|
overridevirtual |
Compute Y := (1 - a) Y + a D^{-1} (X - R*Y)
. Not supported. Call applyInverseJacobi
instead.
Implements Ifpack2::Container< MatrixType >.
|
overridevirtual |
Compute Y := alpha * diag(D) * M^{-1} (diag(D) * X) + beta*Y
. Not supported.
Implements Ifpack2::Container< MatrixType >.
|
overridevirtual |
Print information about this object to the given output stream.
operator<< uses this method.
Implements Ifpack2::Container< MatrixType >.
|
overridevirtual |
A one-line description of this object.
Reimplemented from Teuchos::Describable.
|
overridevirtual |
Print the object with some verbosity level to the given FancyOStream.
Reimplemented from Teuchos::Describable.
|
static |
Get the name of this container type for Details::constructContainer()