Public Types | Public Member Functions | Private Member Functions | Private Attributes | List of all members
Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > Class Template Reference

#include <Xpetra_BlockedCrsMatrix_fwd.hpp>

Inheritance diagram for Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >:
Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > Teuchos::Describable Teuchos::LabeledObject

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...
 
global_size_t getGlobalNumDiags () const
 Returns the number of global diagonal entries, based on global row/column index comparisons. More...
 
size_t getNodeNumDiags () const
 Returns the number of local diagonal entries, based on global row/column index comparisons. 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 idices. More...
 
virtual ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm () const
 Get Frobenius norm of the matrix. 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...
 
- 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
 
 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 bool hasTransposeApply () const
 Whether this operator supports applying the transpose or conjugate transpose. More...
 
- Public Member Functions inherited from Teuchos::Describable
void describe (std::ostream &out, const EVerbosityLevel verbLevel=verbLevel_default) const
 
virtual ~Describable ()
 
DescribableStreamManipulatorState describe (const Describable &describable, const EVerbosityLevel verbLevel=Describable::verbLevel_default)
 
std::ostream & operator<< (std::ostream &os, const DescribableStreamManipulatorState &d)
 
void describe (std::ostream &out, const EVerbosityLevel verbLevel=verbLevel_default) const
 
virtual ~Describable ()
 
- Public Member Functions inherited from Teuchos::LabeledObject
 LabeledObject ()
 
virtual ~LabeledObject ()
 
virtual void setObjectLabel (const std::string &objectLabel)
 
virtual std::string getObjectLabel () const
 

Private Member Functions

void CreateDefaultView ()
 

Private Attributes

Teuchos::RCP< const MapExtractordomainmaps_
 
Teuchos::RCP< const MapExtractorrangemaps_
 
Teuchos::RCP< Mapfullrowmap_
 
Teuchos::RCP< Mapfullcolmap_
 
std::vector< Teuchos::RCP< CrsMatrix > > blocks_
 

Constructor/Destructor Methods

 BlockedCrsMatrix (Teuchos::RCP< const MapExtractor > &rangeMaps, Teuchos::RCP< const MapExtractor > &domainMaps, size_t npr, Xpetra::ProfileType pftype=Xpetra::DynamicProfile)
 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, specifying domain and range maps. 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 MapgetDomainMap () const
 Returns the Map associated with the full domain of this operator. This will be null until fillComplete() is called. More...
 
RCP< const MapgetDomainMap (size_t i) const
 Returns the Map associated with the i'th block domain of this operator. This will be null until fillComplete() is called. More...
 
RCP< const MapgetRangeMap () const
 
RCP< const MapgetRangeMap (size_t i) const
 
RCP< const MapExtractorgetRangeMapExtractor ()
 Returns map extractor class for range map. More...
 
RCP< const MapExtractorgetDomainMapExtractor ()
 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...
 
RCP< const CrsGraphgetCrsGraph () const
 Returns the CrsGraph associated with this matrix. More...
 

Block matrix access

virtual size_t Rows () const
 number of row blocks More...
 
virtual size_t Cols () const
 number of column blocks More...
 
Teuchos::RCP< CrsMatrixgetMatrix (size_t r, size_t c) const
 return block (r,c) More...
 
void setMatrix (size_t r, size_t c, Teuchos::RCP< CrsMatrix > mat)
 set matrix block More...
 
Teuchos::RCP< CrsMatrixMerge () const
 merge BlockedCrsMatrix blocks in a CrsMatrix More...
 

helper functions

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

Additional Inherited Members

- Static Public Attributes inherited from Teuchos::Describable
static const EVerbosityLevel verbLevel_default
 
- 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 = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
class Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >

Definition at line 51 of file Xpetra_BlockedCrsMatrix_fwd.hpp.

Member Typedef Documentation

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
typedef Scalar Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::scalar_type

Definition at line 96 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
typedef LocalOrdinal Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::local_ordinal_type

Definition at line 97 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
typedef GlobalOrdinal Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::global_ordinal_type

Definition at line 98 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
typedef Node Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::node_type

Definition at line 99 of file Xpetra_BlockedCrsMatrix.hpp.

Constructor & Destructor Documentation

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::BlockedCrsMatrix ( Teuchos::RCP< const MapExtractor > &  rangeMaps,
Teuchos::RCP< const MapExtractor > &  domainMaps,
size_t  npr,
Xpetra::ProfileType  pftype = Xpetra::DynamicProfile 
)
inline

Constructor.

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

Definition at line 117 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
virtual Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::~BlockedCrsMatrix ( )
inlinevirtual

Destructor.

Definition at line 249 of file Xpetra_BlockedCrsMatrix.hpp.

Member Function Documentation

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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.

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 280 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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.

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 292 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::removeEmptyProcessesInPlace ( const Teuchos::RCP< const Map > &  newMap)
inlinevirtual
template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 309 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 320 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 327 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 332 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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.

resumeFill() may be called repeatedly.

Postcondition
isFillActive() == true
isFillComplete() == false

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

Definition at line 349 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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, specifying domain and range maps.

Off-node indices are distributed (via globalAssemble()), indices are sorted, redundant indices are eliminated, and global indices are transformed to local indices.

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 366 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 384 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 438 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 454 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 468 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 482 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 494 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 507 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
global_size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getGlobalNumDiags ( ) const
inlinevirtual

Returns the number of global diagonal entries, based on global row/column index comparisons.

Undefined if isFillActive().

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

Definition at line 514 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getNodeNumDiags ( ) const
inlinevirtual

Returns the number of local diagonal entries, based on global row/column index comparisons.

Undefined if isFillActive().

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

Definition at line 521 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 528 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 535 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 543 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 554 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 562 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 584 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 601 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 615 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 idices.

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 623 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 628 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 >.

Definition at line 666 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
RCP<const Map > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getDomainMap ( ) const
inlinevirtual

Returns the Map associated with the full domain of this operator. This will be null until fillComplete() is called.

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

Definition at line 726 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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. This will be null until fillComplete() is called.

Definition at line 730 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
RCP<const Map > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getRangeMap ( ) const
inlinevirtual

Returns the Map associated with the full range of this operator. This will be null until fillComplete() is called.

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

Definition at line 734 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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. This will be null until fillComplete() is called.

Definition at line 738 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
RCP<const MapExtractor> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getRangeMapExtractor ( )
inline

Returns map extractor class for range map.

Definition at line 741 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
RCP<const MapExtractor> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getDomainMapExtractor ( )
inline

Returns map extractor for domain map.

Definition at line 744 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 752 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::doImport ( const Matrix source,
const Import importer,
CombineMode  CM 
)
inlinevirtual
template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::doExport ( const Matrix dest,
const Import importer,
CombineMode  CM 
)
inlinevirtual
template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 767 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 772 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 >.

Definition at line 782 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 >.

Definition at line 785 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
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 809 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
virtual size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::Rows ( ) const
inlinevirtual

number of row blocks

Definition at line 819 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
virtual size_t Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::Cols ( ) const
inlinevirtual

number of column blocks

Definition at line 822 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
Teuchos::RCP<CrsMatrix> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getMatrix ( size_t  r,
size_t  c 
) const
inline

return block (r,c)

Definition at line 825 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::setMatrix ( size_t  r,
size_t  c,
Teuchos::RCP< CrsMatrix mat 
)
inline

set matrix block

Definition at line 828 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
Teuchos::RCP<CrsMatrix> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::Merge ( ) const
inline

merge BlockedCrsMatrix blocks in a CrsMatrix

Definition at line 841 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::Add ( const CrsMatrix A,
const Scalar  scalarA,
CrsMatrix 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 896 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
void Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::CreateDefaultView ( )
inlineprivate

Definition at line 936 of file Xpetra_BlockedCrsMatrix.hpp.

Member Data Documentation

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
Teuchos::RCP<const MapExtractor> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::domainmaps_
private

Definition at line 946 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
Teuchos::RCP<const MapExtractor> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::rangemaps_
private

Definition at line 947 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
Teuchos::RCP<Map> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::fullrowmap_
private

Definition at line 948 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
Teuchos::RCP<Map> Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::fullcolmap_
private

Definition at line 949 of file Xpetra_BlockedCrsMatrix.hpp.

template<class Scalar = Matrix<>::scalar_type, class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type, class GlobalOrdinal = typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type, class Node = typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>
std::vector<Teuchos::RCP<CrsMatrix> > Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::blocks_
private

Definition at line 951 of file Xpetra_BlockedCrsMatrix.hpp.


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