Xpetra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > Class Template Reference

#include <Xpetra_BlockedCrsMatrix.hpp>

Inheritance diagram for Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >:
Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > Xpetra::ReorderedBlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >

Public Types

typedef Scalar scalar_type
 
typedef LocalOrdinal local_ordinal_type
 
typedef GlobalOrdinal global_ordinal_type
 
typedef Node node_type
 
- Public Types inherited from Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >
typedef Scalar scalar_type
 
typedef LocalOrdinal local_ordinal_type
 
typedef GlobalOrdinal global_ordinal_type
 
typedef Node node_type
 
- Public Types inherited from Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node >
typedef Scalar scalar_type
 The type of the entries of the input and output multivectors. More...
 
typedef LocalOrdinal local_ordinal_type
 The local index type. More...
 
typedef GlobalOrdinal global_ordinal_type
 The global index type. More...
 
typedef Node node_type
 The Kokkos Node type. More...
 

Public Member Functions

global_size_t getGlobalNumRows () const
 Returns the number of global rows. More...
 
global_size_t getGlobalNumCols () const
 Returns the number of global columns in the matrix. More...
 
size_t getNodeNumRows () const
 Returns the number of matrix rows owned on the calling node. More...
 
global_size_t getGlobalNumEntries () const
 Returns the global number of entries in this matrix. More...
 
size_t getNodeNumEntries () const
 Returns the local number of entries in this matrix. More...
 
size_t getNumEntriesInLocalRow (LocalOrdinal localRow) const
 Returns the current number of entries on this node in the specified local row. More...
 
size_t getNumEntriesInGlobalRow (GlobalOrdinal globalRow) const
 Returns the current number of entries in the specified (locally owned) global row. More...
 
size_t getGlobalMaxNumRowEntries () const
 Returns the maximum number of entries across all rows/columns on all nodes. More...
 
size_t getNodeMaxNumRowEntries () const
 Returns the maximum number of entries across all rows/columns on this node. More...
 
bool isLocallyIndexed () const
 If matrix indices of all matrix blocks are in the local range, this function returns true. Otherwise, this function returns false. More...
 
bool isGloballyIndexed () const
 If matrix indices are in the global range, this function returns true. Otherwise, this function returns false. More...
 
bool isFillComplete () const
 Returns true if fillComplete() has been called and the matrix is in compute mode. More...
 
virtual void getLocalRowCopy (LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const
 Extract a list of entries in a specified local row of the matrix. Put into storage allocated by calling routine. More...
 
void getGlobalRowView (GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values) const
 Extract a const, non-persisting view of global indices in a specified row of the matrix. More...
 
void getLocalRowView (LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values) const
 Extract a const, non-persisting view of local indices in a specified row of the matrix. More...
 
void getLocalDiagCopy (Vector &diag) const
 Get a copy of the diagonal entries owned by this node, with local row indices. More...
 
void leftScale (const Vector &x)
 Left scale matrix using the given vector entries. More...
 
void rightScale (const Vector &x)
 Right scale matrix using the given vector entries. More...
 
virtual ScalarTraits< Scalar >
::magnitudeType 
getFrobeniusNorm () const
 Get Frobenius norm of the matrix. More...
 
virtual bool haveGlobalConstants () const
 Returns true if globalConstants have been computed; false otherwise. More...
 
virtual void bgs_apply (const MultiVector &X, MultiVector &Y, size_t row, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
 Special multiplication routine (for BGS/Jacobi smoother) More...
 
const Teuchos::RCP< const MapgetMap () const
 Implements DistObject interface. More...
 
void doImport (const Matrix &source, const Import &importer, CombineMode CM)
 Import. More...
 
void doExport (const Matrix &dest, const Import &importer, CombineMode CM)
 Export. More...
 
void doImport (const Matrix &source, const Export &exporter, CombineMode CM)
 Import (using an Exporter). More...
 
void doExport (const Matrix &dest, const Export &exporter, CombineMode CM)
 Export (using an Importer). More...
 
bool hasCrsGraph () const
 Supports the getCrsGraph() call. More...
 
RCP< const CrsGraphgetCrsGraph () const
 Returns the CrsGraph associated with this matrix. More...
 
void residual (const MultiVector &X, const MultiVector &B, MultiVector &R) const
 Compute a residual R = B - (*this) * X. More...
 
- Public Member Functions inherited from Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >
void SetFixedBlockSize (LocalOrdinal blksize, GlobalOrdinal offset=0)
 
LocalOrdinal GetFixedBlockSize () const
 
virtual void SetMaxEigenvalueEstimate (Scalar const &sigma)
 
virtual Scalar GetMaxEigenvalueEstimate () const
 
virtual void residual (const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &R) const =0
 Compute a residual R = B - (*this) * X. More...
 
 Matrix ()
 
virtual ~Matrix ()
 Destructor. More...
 
void CreateView (viewLabel_t viewLabel, const RCP< const Map > &rowMap, const RCP< const Map > &colMap)
 
void CreateView (const viewLabel_t viewLabel, const RCP< const Matrix > &A, bool transposeA=false, const RCP< const Matrix > &B=Teuchos::null, bool transposeB=false)
 
void PrintViews (Teuchos::FancyOStream &out) const
 Print all of the views associated with the Matrix. More...
 
void RemoveView (const viewLabel_t viewLabel)
 
const viewLabel_t SwitchToView (const viewLabel_t viewLabel)
 
bool IsView (const viewLabel_t viewLabel) const
 
const viewLabel_t SwitchToDefaultView ()
 
const viewLabel_tGetDefaultViewLabel () const
 
const viewLabel_tGetCurrentViewLabel () const
 
virtual const RCP< const Map > & getRowMap () const
 Returns the Map that describes the row distribution in this matrix. More...
 
virtual const RCP< const Map > & getRowMap (viewLabel_t viewLabel) const
 Returns the Map that describes the row distribution in this matrix. More...
 
virtual const RCP< const Map > & getColMap () const
 Returns the Map that describes the column distribution in this matrix. This might be null until fillComplete() is called. More...
 
virtual const RCP< const Map > & getColMap (viewLabel_t viewLabel) const
 Returns the Map that describes the column distribution in this matrix. More...
 
- Public Member Functions inherited from Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node >
virtual ~Operator ()
 
virtual void removeEmptyProcessesInPlace (const RCP< const Map > &)
 
virtual bool hasTransposeApply () const
 Whether this operator supports applying the transpose or conjugate transpose. More...
 

Private Member Functions

void CreateDefaultView ()
 

Private Attributes

bool is_diagonal_
 
Teuchos::RCP< const MapExtractordomainmaps_
 
Teuchos::RCP< const MapExtractorrangemaps_
 
std::vector< Teuchos::RCP
< Matrix > > 
blocks_
 
bool bRangeThyraMode_
 boolean flag, which is true, if BlockedCrsMatrix has been created using Thyra-style numbering for sub blocks, i.e. all GIDs of submaps are contiguous and start from 0. More...
 
bool bDomainThyraMode_
 boolean flag, which is true, if BlockedCrsMatrix has been created using Thyra-style numbering for sub blocks, i.e. all GIDs of submaps are contiguous and start from 0. More...
 

Constructor/Destructor Methods

 BlockedCrsMatrix (const Teuchos::RCP< const BlockedMap > &rangeMaps, const Teuchos::RCP< const BlockedMap > &domainMaps, size_t npr)
 Constructor. More...
 
 BlockedCrsMatrix (Teuchos::RCP< const MapExtractor > &rangeMaps, Teuchos::RCP< const MapExtractor > &domainMaps, size_t npr)
 Constructor. More...
 
virtual ~BlockedCrsMatrix ()
 Destructor. More...
 

Insertion/Removal Methods

void insertGlobalValues (GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
 Insert matrix entries, using global IDs. More...
 
void insertLocalValues (LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
 Insert matrix entries, using local IDs. More...
 
void removeEmptyProcessesInPlace (const Teuchos::RCP< const Map > &newMap)
 
void replaceGlobalValues (GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
 Replace matrix entries, using global IDs. More...
 
void replaceLocalValues (LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
 Replace matrix entries, using local IDs. More...
 
virtual void setAllToScalar (const Scalar &alpha)
 Set all matrix entries equal to scalar. More...
 
void scale (const Scalar &alpha)
 Scale the current values of a matrix, this = alpha*this. More...
 

Transformational Methods

void resumeFill (const RCP< ParameterList > &params=null)
 
void fillComplete (const RCP< const Map > &domainMap, const RCP< const Map > &rangeMap, const RCP< ParameterList > &params=null)
 Signal that data entry is complete. More...
 
void fillComplete (const RCP< ParameterList > &params=null)
 Signal that data entry is complete. More...
 

Methods implementing Matrix

Multiplies this matrix by a MultiVector.

X is required to be post-imported, i.e., described by the column map of the matrix. Y is required to be pre-exported, i.e., described by the row map of the matrix.

Both are required to have constant stride, and they are not permitted to ocupy overlapping space. No runtime checking will be performed in a non-debug build.

This method is templated on the scalar type of MultiVector objects, allowing this method to be applied to MultiVector objects of arbitrary type. However, it is recommended that multiply() not be called directly; instead, use the CrsMatrixMultiplyOp, as it will handle the import/exprt operations required to apply a matrix with non-trivial communication needs.

If beta is equal to zero, the operation will enjoy overwrite semantics (Y will be overwritten with the result of the multiplication). Otherwise, the result of the multiplication will be accumulated into Y.

virtual void apply (const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
 Computes the sparse matrix-multivector multiplication. More...
 
RCP< const MapgetFullDomainMap () const
 Returns the Map associated with the full domain of this operator. More...
 
RCP< const BlockedMapgetBlockedDomainMap () const
 Returns the BlockedMap associated with the domain of this operator. More...
 
RCP< const MapgetDomainMap () const
 Returns the Map associated with the domain of this operator. More...
 
RCP< const MapgetDomainMap (size_t i) const
 Returns the Map associated with the i'th block domain of this operator. More...
 
RCP< const MapgetDomainMap (size_t i, bool bThyraMode) const
 Returns the Map associated with the i'th block domain of this operator. More...
 
RCP< const MapgetFullRangeMap () const
 Returns the Map associated with the full range of this operator. More...
 
RCP< const BlockedMapgetBlockedRangeMap () const
 Returns the BlockedMap associated with the range of this operator. More...
 
RCP< const MapgetRangeMap () const
 Returns the Map associated with the range of this operator. More...
 
RCP< const MapgetRangeMap (size_t i) const
 Returns the Map associated with the i'th block range of this operator. More...
 
RCP< const MapgetRangeMap (size_t i, bool bThyraMode) const
 Returns the Map associated with the i'th block range of this operator. More...
 
RCP< const MapExtractorgetRangeMapExtractor () const
 Returns map extractor class for range map. More...
 
RCP< const MapExtractorgetDomainMapExtractor () const
 Returns map extractor for domain map. More...
 

Overridden from Teuchos::Describable

std::string description () const
 Return a simple one-line description of this object. More...
 
void describe (Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
 Print the object with some verbosity level to an FancyOStream object. More...
 

Overridden from Teuchos::LabeledObject

void setObjectLabel (const std::string &objectLabel)
 

Block matrix access

virtual bool isDiagonal () const
 
virtual size_t Rows () const
 number of row blocks More...
 
virtual size_t Cols () const
 number of column blocks More...
 
Teuchos::RCP< MatrixgetCrsMatrix () const
 return unwrap 1x1 blocked operators More...
 
Teuchos::RCP< MatrixgetInnermostCrsMatrix ()
 helper routine recursively returns the first inner-most non-null matrix block from a (nested) blocked operator More...
 
Teuchos::RCP< MatrixgetMatrix (size_t r, size_t c) const
 return block (r,c) More...
 
void setMatrix (size_t r, size_t c, Teuchos::RCP< Matrix > mat)
 set matrix block More...
 
Teuchos::RCP< MatrixMerge () const
 merge BlockedCrsMatrix blocks in a CrsMatrix More...
 

helper functions

void Add (const Matrix &A, const Scalar scalarA, Matrix &B, const Scalar scalarB) const
 Add a Xpetra::CrsMatrix to another: B = B*scalarB + A*scalarA. More...
 

Additional Inherited Members

- Protected Attributes inherited from Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >
Teuchos::Hashtable
< viewLabel_t, RCP< MatrixView > > 
operatorViewTable_
 
viewLabel_t defaultViewLabel_
 
viewLabel_t currentViewLabel_
 

Detailed Description

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
class Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >

Definition at line 96 of file Xpetra_BlockedCrsMatrix.hpp.

Member Typedef Documentation

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
typedef Scalar Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::scalar_type

Definition at line 99 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
typedef LocalOrdinal Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::local_ordinal_type

Definition at line 100 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
typedef GlobalOrdinal Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::global_ordinal_type

Definition at line 101 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
typedef Node Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::node_type

Definition at line 102 of file Xpetra_BlockedCrsMatrix.hpp.

Constructor & Destructor Documentation

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::BlockedCrsMatrix ( const Teuchos::RCP< const BlockedMap > &  rangeMaps,
const Teuchos::RCP< const BlockedMap > &  domainMaps,
size_t  npr 
)
inline

Constructor.

Parameters
rangeMapsrange maps for all blocks
domainMapsdomain maps for all blocks
nprextimated number of entries per row in each block(!)

Definition at line 119 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::BlockedCrsMatrix ( Teuchos::RCP< const MapExtractor > &  rangeMaps,
Teuchos::RCP< const MapExtractor > &  domainMaps,
size_t  npr 
)
inline

Constructor.

Parameters
rangeMapsrange maps for all blocks
domainMapsdomain maps for all blocks
nprextimated number of entries per row in each block(!)
Note
This constructor will be deprecated. Please use the constructor which takes BlockedMap objects instead.

Definition at line 147 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
virtual Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::~BlockedCrsMatrix ( )
inlinevirtual

Destructor.

Definition at line 305 of file Xpetra_BlockedCrsMatrix.hpp.

Member Function Documentation

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::insertGlobalValues ( GlobalOrdinal  globalRow,
const ArrayView< const GlobalOrdinal > &  cols,
const ArrayView< const Scalar > &  vals 
)
inlinevirtual

Insert matrix entries, using global IDs.

Note: this routine throws for Rows() > 1 and/or Cols() > 1

All index values must be in the global space.

Precondition
globalRow exists as an ID in the global row map
isLocallyIndexed() == false
isStorageOptimized() == false
Postcondition
isGloballyIndexed() == true
Note
If globalRow does not belong to the matrix on this node, then it will be communicated to the appropriate node when globalAssemble() is called (which will, at the latest, occur during the next call to fillComplete().) Otherwise, the entries will be inserted in the local matrix.
If the matrix row already contains values at the indices corresponding to values in cols, then the new values will be summed with the old values; this may happen at insertion or during the next call to fillComplete().
If hasColMap() == true, only (cols[i],vals[i]) where cols[i] belongs to the column map on this node will be inserted into the matrix.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 339 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::insertLocalValues ( LocalOrdinal  localRow,
const ArrayView< const LocalOrdinal > &  cols,
const ArrayView< const Scalar > &  vals 
)
inlinevirtual

Insert matrix entries, using local IDs.

Note: this routine throws if Rows() > 1 and/or Cols() > 1

All index values must be in the local space.

Precondition
localRow exists as an ID in the local row map
isGloballyIndexed() == false
isStorageOptimized() == false
Postcondition
isLocallyIndexed() == true

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 359 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::removeEmptyProcessesInPlace ( const Teuchos::RCP< const Map > &  newMap)
inline

Definition at line 368 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::replaceGlobalValues ( GlobalOrdinal  globalRow,
const ArrayView< const GlobalOrdinal > &  cols,
const ArrayView< const Scalar > &  vals 
)
inlinevirtual

Replace matrix entries, using global IDs.

All index values must be in the global space.

Precondition
globalRow is a global row belonging to the matrix on this node.
Note
If (globalRow,cols[i]) corresponds to an entry that is duplicated in this matrix row (likely because it was inserted more than once and fillComplete() has not been called in the interim), the behavior of this function is not defined.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 386 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::replaceLocalValues ( LocalOrdinal  localRow,
const ArrayView< const LocalOrdinal > &  cols,
const ArrayView< const Scalar > &  vals 
)
inlinevirtual

Replace matrix entries, using local IDs.

All index values must be in the local space. Note that if a value is not already present for the specified location in the matrix, the input value will be ignored silently.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 402 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
virtual void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::setAllToScalar ( const Scalar &  alpha)
inlinevirtual

Set all matrix entries equal to scalar.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 414 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::scale ( const Scalar &  alpha)
inlinevirtual

Scale the current values of a matrix, this = alpha*this.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 426 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::resumeFill ( const RCP< ParameterList > &  params = null)
inlinevirtual

Resume fill operations. After calling fillComplete(), resumeFill() must be called before initiating any changes to the matrix.

For BlockedCrsMatrix objects we call the routine iteratively for all sub-blocks.

resumeFill() may be called repeatedly.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 449 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::fillComplete ( const RCP< const Map > &  domainMap,
const RCP< const Map > &  rangeMap,
const RCP< ParameterList > &  params = null 
)
inlinevirtual

Signal that data entry is complete.

Note: for blocked operators the specified domain and range maps have no meaning. We just call fillComplete for all underlying blocks

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 465 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::fillComplete ( const RCP< ParameterList > &  params = null)
inlinevirtual

Signal that data entry is complete.

Off-node entries are distributed (via globalAssemble()), repeated entries are summed, and global indices are transformed to local indices.

Note
This method calls fillComplete( getRowMap(), getRowMap(), os ).
Precondition
isFillActive() == true
isFillComplete()() == false
Postcondition
isFillActive() == false
isFillComplete() == true
if os == DoOptimizeStorage, then isStorageOptimized() == true

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 488 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
global_size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getGlobalNumRows ( ) const
inlinevirtual

Returns the number of global rows.

Undefined if isFillActive().

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 548 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
global_size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getGlobalNumCols ( ) const
inlinevirtual

Returns the number of global columns in the matrix.

Undefined if isFillActive().

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 565 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getNodeNumRows ( ) const
inlinevirtual

Returns the number of matrix rows owned on the calling node.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 580 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
global_size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getGlobalNumEntries ( ) const
inlinevirtual

Returns the global number of entries in this matrix.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 595 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getNodeNumEntries ( ) const
inlinevirtual

Returns the local number of entries in this matrix.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 608 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getNumEntriesInLocalRow ( LocalOrdinal  localRow) const
inlinevirtual

Returns the current number of entries on this node in the specified local row.

Returns OrdinalTraits<size_t>::invalid() if the specified local row is not valid for this matrix.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 622 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getNumEntriesInGlobalRow ( GlobalOrdinal  globalRow) const
inlinevirtual

Returns the current number of entries in the specified (locally owned) global row.

Returns OrdinalTraits<size_t>::invalid() if the specified local row is not valid for this matrix.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 636 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getGlobalMaxNumRowEntries ( ) const
inlinevirtual

Returns the maximum number of entries across all rows/columns on all nodes.

Undefined if isFillActive().

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 649 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getNodeMaxNumRowEntries ( ) const
inlinevirtual

Returns the maximum number of entries across all rows/columns on this node.

Undefined if isFillActive().

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 669 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
bool Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::isLocallyIndexed ( ) const
inlinevirtual

If matrix indices of all matrix blocks are in the local range, this function returns true. Otherwise, this function returns false.

if false, then this does not automatically mean that all blocks are globally indexed. The user has to make sure, that all matrix blocks are indexed in the same way (locally or globally). Otherwise the block matrix is not valid...

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 690 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
bool Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::isGloballyIndexed ( ) const
inlinevirtual

If matrix indices are in the global range, this function returns true. Otherwise, this function returns false.

if false, then this does not automatically mean that all blocks are locally indexed. The user has to make sure, that all matrix blocks are indexed in the same way (locally or globally). Otherwise the block matrix is not valid...

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 702 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
bool Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::isFillComplete ( ) const
inlinevirtual

Returns true if fillComplete() has been called and the matrix is in compute mode.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 711 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
virtual void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getLocalRowCopy ( LocalOrdinal  LocalRow,
const ArrayView< LocalOrdinal > &  Indices,
const ArrayView< Scalar > &  Values,
size_t &  NumEntries 
) const
inlinevirtual

Extract a list of entries in a specified local row of the matrix. Put into storage allocated by calling routine.

Parameters
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 OrdinalTraits<size_t>::invalid().

Precondition
isLocallyIndexed()==true or hasColMap() == true

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 734 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getGlobalRowView ( GlobalOrdinal  GlobalRow,
ArrayView< const GlobalOrdinal > &  indices,
ArrayView< const Scalar > &  values 
) const
inlinevirtual

Extract a const, non-persisting view of global indices in a specified row of the matrix.

Parameters
GlobalRow- (In) Global row number for which indices are desired.
Indices- (Out) Global column indices corresponding to values.
Values- (Out) Row values
Precondition
isLocallyIndexed() == false
Postcondition
indices.size() == getNumEntriesInGlobalRow(GlobalRow)

Note: If GlobalRow does not belong to this node, then indices is set to null.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 756 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getLocalRowView ( LocalOrdinal  LocalRow,
ArrayView< const LocalOrdinal > &  indices,
ArrayView< const Scalar > &  values 
) const
inlinevirtual

Extract a const, non-persisting view of local indices in a specified row of the matrix.

Parameters
LocalRow- (In) Local row number for which indices are desired.
Indices- (Out) Global column indices corresponding to values.
Values- (Out) Row values
Precondition
isGloballyIndexed() == false
Postcondition
indices.size() == getNumEntriesInLocalRow(LocalRow)

Note: If LocalRow does not belong to this node, then indices is set to null.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 775 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getLocalDiagCopy ( Vector diag) const
inlinevirtual

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.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 794 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::leftScale ( const Vector x)
inlinevirtual

Left scale matrix using the given vector entries.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 829 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::rightScale ( const Vector x)
inlinevirtual

Right scale matrix using the given vector entries.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 867 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
virtual ScalarTraits<Scalar>::magnitudeType Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getFrobeniusNorm ( ) const
inlinevirtual

Get Frobenius norm of the matrix.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 906 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
virtual bool Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::haveGlobalConstants ( ) const
inlinevirtual

Returns true if globalConstants have been computed; false otherwise.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 922 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
virtual void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::apply ( const MultiVector X,
MultiVector Y,
Teuchos::ETransp  mode = Teuchos::NO_TRANS,
Scalar  alpha = ScalarTraits<Scalar>::one(),
Scalar  beta = ScalarTraits<Scalar>::zero() 
) const
inlinevirtual

Computes the sparse matrix-multivector multiplication.

Performs \(Y = \alpha A^{\textrm{mode}} X + \beta Y\), with one special exceptions:

  • if beta == 0, apply() overwrites Y, so that any values in Y (including NaNs) are ignored.

Implements Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Reimplemented in Xpetra::ReorderedBlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 959 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const Map > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getFullDomainMap ( ) const
inline

Returns the Map associated with the full domain of this operator.

Definition at line 1052 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const BlockedMap > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getBlockedDomainMap ( ) const
inline

Returns the BlockedMap associated with the domain of this operator.

Definition at line 1055 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const Map > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getDomainMap ( ) const
inlinevirtual

Returns the Map associated with the domain of this operator.

Implements Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1058 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const Map > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getDomainMap ( size_t  i) const
inline

Returns the Map associated with the i'th block domain of this operator.

Definition at line 1061 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const Map > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getDomainMap ( size_t  i,
bool  bThyraMode 
) const
inline

Returns the Map associated with the i'th block domain of this operator.

Definition at line 1064 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const Map > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getFullRangeMap ( ) const
inline

Returns the Map associated with the full range of this operator.

Definition at line 1067 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const BlockedMap > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getBlockedRangeMap ( ) const
inline

Returns the BlockedMap associated with the range of this operator.

Definition at line 1070 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const Map > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getRangeMap ( ) const
inlinevirtual

Returns the Map associated with the range of this operator.

Implements Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1073 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const Map > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getRangeMap ( size_t  i) const
inline

Returns the Map associated with the i'th block range of this operator.

Definition at line 1076 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const Map > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getRangeMap ( size_t  i,
bool  bThyraMode 
) const
inline

Returns the Map associated with the i'th block range of this operator.

Definition at line 1079 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const MapExtractor> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getRangeMapExtractor ( ) const
inline

Returns map extractor class for range map.

Definition at line 1082 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const MapExtractor> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getDomainMapExtractor ( ) const
inline

Returns map extractor for domain map.

Definition at line 1085 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
virtual void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::bgs_apply ( const MultiVector X,
MultiVector Y,
size_t  row,
Teuchos::ETransp  mode = Teuchos::NO_TRANS,
Scalar  alpha = ScalarTraits<Scalar>::one(),
Scalar  beta = ScalarTraits<Scalar>::zero() 
) const
inlinevirtual

Special multiplication routine (for BGS/Jacobi smoother)

Computes the sparse matrix-multivector multiplication (plus linear combination with input/result vector)

Performs \(Y = \alpha A^{\textrm{mode}} X + \beta Y\), with one special exception:

  • if beta == 0, apply() overwrites Y, so that any values in Y (including NaNs) are ignored.
  • calculates result only for blocked row "row"
  • useful for BGS/Jacobi smoother in MueLu: there we have to calculate the residual for the current block row we can skip the MatVec calls in all other block rows
Parameters
XVector to be multiplied by matrix (input)
Yresult vector
rowIndex of block row to be treated
modeTranspose mode
alphascaling factor for result of matrix-vector product
betascaling factor for linear combination with result vector

Definition at line 1100 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
const Teuchos::RCP< const Map > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getMap ( ) const
inlinevirtual

Implements DistObject interface.

Access function for the Tpetra::Map this DistObject was constructed with.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1172 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::doImport ( const Matrix source,
const Import importer,
CombineMode  CM 
)
inlinevirtual
template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::doExport ( const Matrix dest,
const Import importer,
CombineMode  CM 
)
inlinevirtual
template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::doImport ( const Matrix source,
const Export exporter,
CombineMode  CM 
)
inlinevirtual

Import (using an Exporter).

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1201 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::doExport ( const Matrix dest,
const Export exporter,
CombineMode  CM 
)
inlinevirtual

Export (using an Importer).

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1211 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
std::string Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::description ( ) const
inlinevirtual

Return a simple one-line description of this object.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Reimplemented in Xpetra::ReorderedBlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1226 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::describe ( Teuchos::FancyOStream &  out,
const Teuchos::EVerbosityLevel  verbLevel = Teuchos::Describable::verbLevel_default 
) const
inlinevirtual

Print the object with some verbosity level to an FancyOStream object.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Reimplemented in Xpetra::ReorderedBlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1229 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::setObjectLabel ( const std::string &  objectLabel)
inlinevirtual
template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
bool Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::hasCrsGraph ( ) const
inlinevirtual
template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
RCP<const CrsGraph> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getCrsGraph ( ) const
inlinevirtual

Returns the CrsGraph associated with this matrix.

Implements Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1280 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
virtual bool Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::isDiagonal ( ) const
inlinevirtual

Definition at line 1293 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
virtual size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::Rows ( ) const
inlinevirtual

number of row blocks

Definition at line 1296 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
virtual size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::Cols ( ) const
inlinevirtual

number of column blocks

Definition at line 1299 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
Teuchos::RCP<Matrix> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getCrsMatrix ( ) const
inline

return unwrap 1x1 blocked operators

Definition at line 1302 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
Teuchos::RCP<Matrix> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getInnermostCrsMatrix ( )
inline

helper routine recursively returns the first inner-most non-null matrix block from a (nested) blocked operator

Definition at line 1314 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
Teuchos::RCP<Matrix> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getMatrix ( size_t  r,
size_t  c 
) const
inline

return block (r,c)

Definition at line 1332 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::setMatrix ( size_t  r,
size_t  c,
Teuchos::RCP< Matrix mat 
)
inline

set matrix block

Definition at line 1346 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
Teuchos::RCP<Matrix> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::Merge ( ) const
inline

merge BlockedCrsMatrix blocks in a CrsMatrix

Definition at line 1360 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::residual ( const MultiVector X,
const MultiVector B,
MultiVector R 
) const
inlinevirtual

Compute a residual R = B - (*this) * X.

Implements Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1505 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::Add ( const Matrix A,
const Scalar  scalarA,
Matrix B,
const Scalar  scalarB 
) const
inlineprivate

Add a Xpetra::CrsMatrix to another: B = B*scalarB + A*scalarA.

Note, that this routine works only correctly if A only has entries which are empty (zero) in B. We use the insertGlobalValues routine for inserting the new values from A in B. The sumIntoGlobalValues routine is not implemented in Xpetra (and would not extend the graph of B for new entries). Here we need something to catch the exceptions of a future implementation of sumIntoGlobalValues that then adds the remaining new entries with insertGlobal Values.

This routine is private and used only by Merge. Since the blocks in BlockedCrsMatrix are seperated, this routine works for merging a BlockedCrsMatrix.

Definition at line 1529 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::CreateDefaultView ( )
inlineprivate

Definition at line 1576 of file Xpetra_BlockedCrsMatrix.hpp.

Member Data Documentation

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
bool Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::is_diagonal_
private

Definition at line 1588 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
Teuchos::RCP<const MapExtractor> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::domainmaps_
private

Definition at line 1589 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
Teuchos::RCP<const MapExtractor> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::rangemaps_
private

Definition at line 1590 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
std::vector<Teuchos::RCP<Matrix> > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::blocks_
private

Definition at line 1592 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
bool Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::bRangeThyraMode_
private

boolean flag, which is true, if BlockedCrsMatrix has been created using Thyra-style numbering for sub blocks, i.e. all GIDs of submaps are contiguous and start from 0.

Definition at line 1596 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
bool Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::bDomainThyraMode_
private

boolean flag, which is true, if BlockedCrsMatrix has been created using Thyra-style numbering for sub blocks, i.e. all GIDs of submaps are contiguous and start from 0.

Definition at line 1597 of file Xpetra_BlockedCrsMatrix.hpp.


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