Ifpack2 Templated Preconditioning Package
Version 1.0
|
Sparse matrix (Tpetra::RowMatrix subclass) with ghost rows. More...
#include <Ifpack2_OverlappingRowMatrix_decl.hpp>
Public Types | |
Typedefs | |
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 Teuchos::ScalarTraits < scalar_type >::magnitudeType | magnitude_type |
typedef MatrixType::global_inds_host_view_type | global_inds_host_view_type |
typedef MatrixType::local_inds_host_view_type | local_inds_host_view_type |
typedef MatrixType::values_host_view_type | values_host_view_type |
typedef MatrixType::nonconst_global_inds_host_view_type | nonconst_global_inds_host_view_type |
typedef MatrixType::nonconst_local_inds_host_view_type | nonconst_local_inds_host_view_type |
typedef MatrixType::nonconst_values_host_view_type | nonconst_values_host_view_type |
using | row_matrix_type = Tpetra::RowMatrix< scalar_type, local_ordinal_type, global_ordinal_type, node_type > |
using | crs_matrix_type = Tpetra::CrsMatrix< scalar_type, local_ordinal_type, global_ordinal_type, node_type > |
typedef MatrixType::node_type::device_type | device_type |
typedef device_type::execution_space | execution_space |
typedef MatrixType::local_inds_device_view_type | local_inds_device_view_type |
typedef MatrixType::global_inds_device_view_type | global_inds_device_view_type |
typedef MatrixType::values_device_view_type | values_device_view_type |
typedef row_matrix_type::mag_type | mag_type |
Public Types inherited from Ifpack2::Details::RowMatrix< MatrixType > | |
using | scalar_type = typename MatrixType::scalar_type |
using | local_ordinal_type = typename MatrixType::local_ordinal_type |
using | global_ordinal_type = typename MatrixType::global_ordinal_type |
using | node_type = typename MatrixType::node_type |
Public Member Functions | |
Constructors and destructor | |
OverlappingRowMatrix (const Teuchos::RCP< const row_matrix_type > &A, const int overlapLevel) | |
~OverlappingRowMatrix ()=default | |
Destructor. More... | |
Matrix query methods | |
virtual Teuchos::RCP< const Teuchos::Comm< int > > | getComm () const |
The communicator over which the matrix is distributed. More... | |
virtual Teuchos::RCP< const Tpetra::Map < local_ordinal_type, global_ordinal_type, node_type > > | getRowMap () const |
The Map that describes the distribution of rows over processes. More... | |
virtual Teuchos::RCP< const Tpetra::Map < local_ordinal_type, global_ordinal_type, node_type > > | getColMap () const |
The Map that describes the distribution of columns over processes. More... | |
virtual Teuchos::RCP< const Tpetra::Map < local_ordinal_type, global_ordinal_type, node_type > > | getDomainMap () const |
The Map that describes the domain of this matrix. More... | |
virtual Teuchos::RCP< const Tpetra::Map < local_ordinal_type, global_ordinal_type, node_type > > | getRangeMap () const |
The Map that describes the range of this matrix. More... | |
virtual Teuchos::RCP< const Tpetra::RowGraph < local_ordinal_type, global_ordinal_type, node_type > > | getGraph () const |
This matrix's graph. More... | |
virtual global_size_t | getGlobalNumRows () const |
The global number of rows in this matrix. More... | |
virtual global_size_t | getGlobalNumCols () const |
The global number of columns in this matrix. More... | |
virtual size_t | getLocalNumRows () const |
The number of rows owned by the calling process. More... | |
virtual size_t | getLocalNumCols () const |
The number of columns owned by the calling process. More... | |
virtual global_ordinal_type | getIndexBase () const |
The index base for global indices for this matrix. More... | |
virtual global_size_t | getGlobalNumEntries () const |
The global number of entries in this matrix. More... | |
virtual size_t | getLocalNumEntries () const |
The number of entries in this matrix owned by the calling process. More... | |
virtual size_t | getNumEntriesInGlobalRow (global_ordinal_type globalRow) const |
The number of entries in the given global row that are owned by the calling process. More... | |
virtual size_t | getNumEntriesInLocalRow (local_ordinal_type localRow) const |
The number of entries in the given local row that are owned by the calling process. More... | |
virtual size_t | getGlobalMaxNumRowEntries () const |
The maximum number of entries in any row on any process. More... | |
virtual size_t | getLocalMaxNumRowEntries () const |
The maximum number of entries in any row on the calling process. More... | |
virtual local_ordinal_type | getBlockSize () const |
The number of degrees of freedom per mesh point. More... | |
virtual bool | hasColMap () const |
Whether this matrix has a column Map. More... | |
virtual bool | isLocallyIndexed () const |
Whether this matrix is locally indexed. More... | |
virtual bool | isGloballyIndexed () const |
Whether this matrix is globally indexed. More... | |
virtual bool | isFillComplete () const |
true if fillComplete() has been called, else false . More... | |
virtual bool | supportsRowViews () const |
true if row views are supported, else false . More... | |
Extraction methods | |
virtual void | getGlobalRowCopy (global_ordinal_type GlobalRow, nonconst_global_inds_host_view_type &Indices, nonconst_values_host_view_type &Values, size_t &NumEntries) const |
Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage. More... | |
virtual void | getLocalRowCopy (local_ordinal_type LocalRow, nonconst_local_inds_host_view_type &Indices, nonconst_values_host_view_type &Values, size_t &NumEntries) const |
Extract a list of entries in a specified local row of the graph. Put into storage allocated by calling routine. More... | |
virtual void | getGlobalRowView (global_ordinal_type GlobalRow, global_inds_host_view_type &indices, values_host_view_type &values) const |
Extract a const, non-persisting view of global indices in a specified row of the matrix. More... | |
virtual void | getLocalRowView (local_ordinal_type LocalRow, local_inds_host_view_type &indices, values_host_view_type &values) const |
Extract a const, non-persisting view of local indices in a specified row of the matrix. More... | |
virtual void | getLocalDiagCopy (Tpetra::Vector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &diag) const |
Get a copy of the diagonal entries owned by this node, with local row indices. More... | |
Public Member Functions inherited from Ifpack2::Details::RowMatrix< MatrixType > | |
virtual | ~RowMatrix ()=default |
Destructor (virtual for memory safety of derived classes) More... | |
Mathematical methods | |
virtual void | leftScale (const Tpetra::Vector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &x) |
Scales the RowMatrix on the left with the Vector x. More... | |
virtual void | rightScale (const Tpetra::Vector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &x) |
Scales the RowMatrix on the right with the Vector x. More... | |
virtual mag_type | getFrobeniusNorm () const |
Returns the Frobenius norm of the matrix. More... | |
virtual void | apply (const Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &X, Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, scalar_type alpha=Teuchos::ScalarTraits< scalar_type >::one(), scalar_type beta=Teuchos::ScalarTraits< scalar_type >::zero()) const |
Computes the operator-multivector application. More... | |
virtual bool | hasTransposeApply () const |
Whether this operator's apply() method can apply the adjoint (transpose). More... | |
virtual void | importMultiVector (const Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &X, Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &OvX, Tpetra::CombineMode CM=Tpetra::INSERT) |
virtual void | exportMultiVector (const Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &OvX, Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &X, Tpetra::CombineMode CM=Tpetra::ADD) |
std::string | description () const |
void | describe (Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const |
Teuchos::RCP< const crs_matrix_type > | getUnderlyingMatrix () const |
Teuchos::RCP< const crs_matrix_type > | getExtMatrix () const |
Kokkos::View< size_t *, typename OverlappingRowMatrix < MatrixType >::device_type > | getExtHaloStarts () const |
Kokkos::View< size_t *, typename OverlappingRowMatrix < MatrixType >::device_type > ::HostMirror | getExtHaloStartsHost () const |
void | doExtImport () |
Sparse matrix (Tpetra::RowMatrix subclass) with ghost rows.
MatrixType | Tpetra::RowMatrix specialization. |
Ifpack2::OverlappingRowMatrix< MatrixType >::OverlappingRowMatrix | ( | const Teuchos::RCP< const row_matrix_type > & | A, |
const int | overlapLevel | ||
) |
Standard constructor.
A | [in] The input matrix. Currently this class requires that A be a Tpetra::CrsMatrix instance with the same first four template parameters as MatrixType, and with a default fifth template parameter. Furthermore, A must have a nonoverlapping row Map and must be distributed over more than one MPI process. |
overlapLevel | [in] The number of levels of overlap. |
|
default |
Destructor.
|
virtual |
The communicator over which the matrix is distributed.
|
virtual |
The Map that describes the distribution of rows over processes.
|
virtual |
The Map that describes the distribution of columns over processes.
|
virtual |
The Map that describes the domain of this matrix.
The domain is the distribution of valid input vectors of apply().
|
virtual |
The Map that describes the range of this matrix.
The domain is the distribution of valid output vectors of apply().
|
virtual |
This matrix's graph.
|
virtual |
The global number of rows in this matrix.
|
virtual |
The global number of columns in this matrix.
|
virtual |
The number of rows owned by the calling process.
|
virtual |
The number of columns owned by the calling process.
This is the number of columns needed to apply the forward operator on the calling process, that is, the number of elements listed in the column Map on the calling process.
|
virtual |
The index base for global indices for this matrix.
|
virtual |
The global number of entries in this matrix.
|
virtual |
The number of entries in this matrix owned by the calling process.
|
virtual |
The number of entries in the given global row that are owned by the calling process.
globalRow | [in] Global index of the row. |
|
virtual |
The number of entries in the given local row that are owned by the calling process.
globalRow | [in] Local index of the row. |
|
virtual |
The maximum number of entries in any row on any process.
|
virtual |
The maximum number of entries in any row on the calling process.
|
virtual |
The number of degrees of freedom per mesh point.
|
virtual |
Whether this matrix has a column Map.
|
virtual |
Whether this matrix is locally indexed.
|
virtual |
Whether this matrix is globally indexed.
|
virtual |
true
if fillComplete() has been called, else false
.
|
virtual |
true
if row views are supported, else false
.
|
virtual |
Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage.
LocalRow | - (In) Global row number for which indices are desired. |
Indices | - (Out) Global column indices corresponding to values. |
Values | - (Out) Matrix values. |
NumEntries | - (Out) Number of indices. |
Note: A std::runtime_error exception is thrown if either Indices
or Values
is not large enough to hold the data associated with row GlobalRow
. If GlobalRow
does not belong to this node, then Indices
and Values
are unchanged and NumIndices
is returned as Teuchos::OrdinalTraits<size_t>::invalid().
|
virtual |
Extract a list of entries in a specified local row of the graph. Put into storage allocated by calling routine.
LocalRow | - (In) Local row number for which indices are desired. |
Indices | - (Out) Local column indices corresponding to values. |
Values | - (Out) Matrix values. |
NumIndices | - (Out) Number of indices. |
Note: A std::runtime_error exception is thrown if either Indices
or Values
is not large enough to hold the data associated with row LocalRow
. If LocalRow
is not valid for this node, then Indices
and Values
are unchanged and NumIndices
is returned as Teuchos::OrdinalTraits<size_t>::invalid().
|
virtual |
Extract a const, non-persisting view of global indices in a specified row of the matrix.
GlobalRow | - (In) Global row number for which indices are desired. |
Indices | - (Out) Global column indices corresponding to values. |
Values | - (Out) Row values |
isLocallyIndexed() == false
indices.size() == getNumEntriesInGlobalRow(GlobalRow)
Note: If GlobalRow
does not belong to this node, then indices
is set to null.
|
virtual |
Extract a const, non-persisting view of local indices in a specified row of the matrix.
LocalRow | - (In) Local row number for which indices are desired. |
Indices | - (Out) Global column indices corresponding to values. |
Values | - (Out) Row values |
isGloballyIndexed() == false
indices.size() == getNumEntriesInLocalRow(LocalRow)
Note: If LocalRow
does not belong to this node, then indices
is set to null.
|
virtual |
Get a copy of the diagonal entries owned by this node, with local row indices.
Returns a distributed Vector object partitioned according to this matrix's row map, containing the the zero and non-zero diagonals owned by this node.
|
virtual |
Scales the RowMatrix on the left with the Vector x.
This matrix will be scaled such that A(i,j) = x(i)*A(i,j) where i denoes the global row number of A and j denotes the global column number of A.
x | A vector to left scale this matrix. |
|
virtual |
Scales the RowMatrix on the right with the Vector x.
This matrix will be scaled such that A(i,j) = x(j)*A(i,j) where i denoes the global row number of A and j denotes the global column number of A.
x | A vector to right scale this matrix. |
|
virtual |
Returns the Frobenius norm of the matrix.
Computes and returns the Frobenius norm of the matrix, defined as: \( \|A\|_F = \sqrt{\sum_{i,j} \|\a_{ij}\|^2} \)
|
virtual |
Computes the operator-multivector application.
Loosely, performs \(Y = \alpha \cdot A^{\textrm{mode}} \cdot X + \beta \cdot Y\). However, the details of operation vary according to the values of alpha
and beta
. Specifically
beta == 0
, apply() must overwrite Y
, so that any values in Y
(including NaNs) are ignored.alpha == 0
, apply() may short-circuit the operator, so that any values in X
(including NaNs) are ignored.This is analagous to the Multiply function in Ifpack, not the Apply
|
virtual |
Whether this operator's apply() method can apply the adjoint (transpose).