Tpetra Matrix/Vector Services Version of the Day
Protected Types | Protected Member Functions
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps > Class Template Reference

A class for constructing and using sparse compressed matrices with row access. More...

#include <Tpetra_CrsMatrix_decl.hpp>

Inheritance diagram for Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >:
Inheritance graph
[legend]

List of all members.

Public Member Functions

Constructor/Destructor Methods
 CrsMatrix (const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile)
 Constructor specifying the number of non-zeros for all rows.
 CrsMatrix (const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile)
 Constructor specifying the number of non-zeros for each row.
 CrsMatrix (const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile)
 Constructor specifying a column map and the number of non-zeros for all rows.
 CrsMatrix (const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile)
 Constructor specifying a column map and the number of non-zeros for each row.
 CrsMatrix (const RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node, LocalMatOps > > &graph)
 Constructor specifying a pre-constructed graph.
virtual ~CrsMatrix ()
Insertion/Removal Methods
void insertGlobalValues (GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
 Insert matrix entries, using global IDs.
void insertLocalValues (LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
 Insert matrix entries, using local IDs.
void replaceGlobalValues (GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
 Replace matrix entries, using global IDs.
void replaceLocalValues (LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
 Replace matrix entries, using local IDs.
void sumIntoGlobalValues (GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
 Sum into multiple entries, using global IDs.
void sumIntoLocalValues (LocalOrdinal globalRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
 Sum into multiple entries, using local IDs.
void setAllToScalar (const Scalar &alpha)
 Set all matrix entries equal to scalarThis.
void scale (const Scalar &alpha)
 Scale the current values of a matrix, this = alpha*this.
Transformational Methods
void globalAssemble ()
 Communicate non-local contributions to other nodes.
void resumeFill ()
void fillComplete (const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, OptimizeOption os=DoOptimizeStorage)
 Signal that data entry is complete, specifying domain and range maps.
void fillComplete (OptimizeOption os=DoOptimizeStorage)
 Signal that data entry is complete.
Methods implementing RowMatrix
const RCP< const Comm< int > > & getComm () const
 Returns the communicator.
RCP< Node > getNode () const
 Returns the underlying node.
const RCP< const Map
< LocalOrdinal, GlobalOrdinal,
Node > > & 
getRowMap () const
 Returns the Map that describes the row distribution in this matrix.
const RCP< const Map
< LocalOrdinal, GlobalOrdinal,
Node > > & 
getColMap () const
 Returns the Map that describes the column distribution in this matrix.
RCP< const RowGraph
< LocalOrdinal, GlobalOrdinal,
Node > > 
getGraph () const
 Returns the RowGraph associated with this matrix.
RCP< const CrsGraph
< LocalOrdinal, GlobalOrdinal,
Node, LocalMatOps > > 
getCrsGraph () const
 Returns the CrsGraph associated with this matrix.
global_size_t getGlobalNumRows () const
 Returns the number of global rows in this matrix.
global_size_t getGlobalNumCols () const
 Returns the number of global columns in the matrix.
size_t getNodeNumRows () const
 Returns the number of matrix rows owned on the calling node.
size_t getNodeNumCols () const
 Returns the number of columns connected to the locally owned rows of this matrix.
GlobalOrdinal getIndexBase () const
 Returns the index base for global indices for this matrix.
global_size_t getGlobalNumEntries () const
 Returns the global number of entries in this matrix.
size_t getNodeNumEntries () const
 Returns the local number of entries in this matrix.
size_t getNumEntriesInGlobalRow (GlobalOrdinal globalRow) const
 Returns the current number of entries on this node in the specified global row.
size_t getNumEntriesInLocalRow (LocalOrdinal localRow) const
 Returns the current number of entries on this node in the specified local row.
global_size_t getGlobalNumDiags () const
 Returns the number of global diagonal entries, based on global row/column index comparisons.
size_t getNodeNumDiags () const
 Returns the number of local diagonal entries, based on global row/column index comparisons.
size_t getGlobalMaxNumRowEntries () const
 Returns the maximum number of entries across all rows/columns on all nodes.
size_t getNodeMaxNumRowEntries () const
 Returns the maximum number of entries across all rows/columns on this node.
bool hasColMap () const
 Indicates whether the matrix has a well-defined column map.
bool isLowerTriangular () const
 Indicates whether the matrix is lower triangular.
bool isUpperTriangular () const
 Indicates whether the matrix is upper triangular.
bool isLocallyIndexed () const
 If matrix indices are in the local range, this function returns true. Otherwise, this function returns false. */.
bool isGloballyIndexed () const
 If matrix indices are in the global range, this function returns true. Otherwise, this function returns false. */.
bool isFillComplete () const
 Returns true if fillComplete() has been called and the matrix is in compute mode.
bool isFillActive () const
 Returns true if resumeFill() has been called and the matrix is in edit mode.
bool isStorageOptimized () const
 Returns true if storage has been optimized.
ProfileType getProfileType () const
 Returns true if the matrix was allocated with static data structures.
bool isStaticGraph () const
 Indicates that the graph is static, so that new entries cannot be added to this matrix. */.
void getGlobalRowCopy (GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const
 Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage.
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.
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.
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.
void getLocalDiagCopy (Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
 Get a copy of the diagonal entries owned by this node, with local row idices.
Advanced Matrix-vector multiplication and solve methods
template<class DomainScalar , class RangeScalar >
void multiply (const MultiVector< DomainScalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< RangeScalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp trans, RangeScalar alpha, RangeScalar beta) const
 Multiplies this matrix by a MultiVector.
template<class DomainScalar , class RangeScalar >
void solve (const MultiVector< RangeScalar, LocalOrdinal, GlobalOrdinal, Node > &Y, MultiVector< DomainScalar, LocalOrdinal, GlobalOrdinal, Node > &X, Teuchos::ETransp trans) const
 Solves a linear system when the underlying matrix is triangular.
Methods implementing Operator
void apply (const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
 Computes the sparse matrix-multivector multiplication.
bool hasTransposeApply () const
 Indicates whether this operator supports applying the adjoint operator.
const RCP< const Map
< LocalOrdinal, GlobalOrdinal,
Node > > & 
getDomainMap () const
 Returns the Map associated with the domain of this operator. This will be null until fillComplete() is called.
const RCP< const Map
< LocalOrdinal, GlobalOrdinal,
Node > > & 
getRangeMap () const
Overridden from Teuchos::Describable
std::string description () const
 Return a simple one-line description of this object.
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.
Methods implementing Tpetra::DistObject
bool checkSizes (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source)
 Allows the source and target (this) objects to be compared for compatibility.
void copyAndPermute (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, size_t numSameIDs, const ArrayView< const LocalOrdinal > &permuteToLIDs, const ArrayView< const LocalOrdinal > &permuteFromLIDs)
 Perform copies and permutations that are local to this image.
void packAndPrepare (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const ArrayView< const LocalOrdinal > &exportLIDs, Array< char > &exports, const ArrayView< size_t > &numPacketsPerLID, size_t &constantNumPackets, Distributor &distor)
 Perform any packing or preparation required for communication.
void unpackAndCombine (const ArrayView< const LocalOrdinal > &importLIDs, const ArrayView< const char > &imports, const ArrayView< size_t > &numPacketsPerLID, size_t constantNumPackets, Distributor &distor, CombineMode CM)
 Perform any unpacking and combining after communication.
Deprecated routines to be removed at some point in the future.
TPETRA_DEPRECATED void optimizeStorage ()
 Deprecated. Re-allocate the data into contiguous storage.
TPETRA_DEPRECATED void getGlobalRowView (GlobalOrdinal GlobalRow, ArrayRCP< const GlobalOrdinal > &indices, ArrayRCP< const Scalar > &values) const
 Deprecated. Get a persisting const view of the entries in a specified global row of this matrix.
TPETRA_DEPRECATED void getLocalRowView (LocalOrdinal LocalRow, ArrayRCP< const LocalOrdinal > &indices, ArrayRCP< const Scalar > &values) const
 Deprecated. Get a persisting const view of the entries in a specified local row of this matrix.
Import/Export Methods
void doImport (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
 Import.
void doImport (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
 Import (using an Exporter)
void doExport (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
 Export.
void doExport (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
 Export (using an Importer)
Attribute Accessor Methods
bool isDistributed () const
 Accessor for whether or not this is a global object.
const Teuchos::RCP< const Map
< LocalOrdinal, GlobalOrdinal,
Node > > & 
getMap () const
 Access function for the Tpetra::Map this DistObject was constructed with.
I/O methods
void print (std::ostream &os) const
 Print method.

Protected Types

enum  ReverseOption
 Enum indicating whether the transer should be performed in forward or reverse mode. More...

Protected Member Functions

virtual void doTransfer (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, CombineMode CM, size_t numSameIDs, const Teuchos::ArrayView< const LocalOrdinal > &permuteToLIDs, const Teuchos::ArrayView< const LocalOrdinal > &permuteFromLIDs, const Teuchos::ArrayView< const LocalOrdinal > &remoteLIDs, const Teuchos::ArrayView< const LocalOrdinal > &exportLIDs, Distributor &distor, ReverseOption revOp)
 Perform transfer (redistribution) of data across memory images.

Detailed Description

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
class Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >

A class for constructing and using sparse compressed matrices with row access.

Template Parameters:
ScalarThe scalar field describing the numerical entries of the matrix.
LocalOrdinalA ordinal type for lists of local indices. This specifies the LocalOrdinal type for Map objects used by this matrix.
GlobalOrdinalA ordinal type for lists of global indices. This specifies the GlobalOrdinal type for Map objects used by this matrix.
NodeA shared-memory node class, fulfilling the Kokkos Node API
LocalMatOpsA local sparse matrix operations class, fulfiling the Kokkos CRS Ops API. This class allows the construction of sparse matrices with row-access.

Local vs. Global

Matrix entries can be added using either local or global coordinates for the indices. The accessors isGloballyIndexed() and isLocallyIndexed() indicate whether the indices are currently stored as global or local indices. Many of the class methods are divided into global and local versions, which differ only in whether they accept/return indices in the global or local coordinate space. Some of these methods may only be used if the matrix coordinates are in the appropriate coordinates. For example, getGlobalRowView() returns a View to the indices in global coordinates; if the indices are not in global coordinates, then no such View can be created.

The global/local distinction does distinguish between operation on the global/local matrix. Almost all methods operate on the local matrix, i.e., the rows of the matrix associated with the local node, per the distribution specified by the row map. Access to non-local rows requires performing an explicit communication via the import/export capabilities of the CrsMatrix object; see DistObject. However, the method insertGlobalValues() is an exception to this rule, as non-local rows are allowed to be added via the local matrix. These rows are stored in the local matrix and communicated to the appropriate node on the next call to globalAssemble() or fillComplete() (the latter calls the former).

Examples:

LocalMatOpExample.cpp, and Tpetra_Power_Method.cpp.

Definition at line 113 of file Tpetra_CrsMatrix_decl.hpp.


Member Enumeration Documentation

enum Tpetra::DistObject::ReverseOption [protected, inherited]

Enum indicating whether the transer should be performed in forward or reverse mode.

Definition at line 127 of file Tpetra_DistObject.hpp.


Constructor & Destructor Documentation

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::CrsMatrix ( const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  rowMap,
size_t  maxNumEntriesPerRow,
ProfileType  pftype = DynamicProfile 
)

Constructor specifying the number of non-zeros for all rows.

Definition at line 70 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::CrsMatrix ( const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  rowMap,
const ArrayRCP< const size_t > &  NumEntriesPerRowToAlloc,
ProfileType  pftype = DynamicProfile 
)

Constructor specifying the number of non-zeros for each row.

Definition at line 100 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::CrsMatrix ( const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  rowMap,
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  colMap,
size_t  maxNumEntriesPerRow,
ProfileType  pftype = DynamicProfile 
)

Constructor specifying a column map and the number of non-zeros for all rows.

The column map will be used to filter any matrix entries inserted using insertLocalValues() or insertGlobalValues().

Definition at line 130 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::CrsMatrix ( const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  rowMap,
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  colMap,
const ArrayRCP< const size_t > &  NumEntriesPerRowToAlloc,
ProfileType  pftype = DynamicProfile 
)

Constructor specifying a column map and the number of non-zeros for each row.

The column map will be used to filter any matrix entries inserted using insertLocalValues() or insertGlobalValues().

Definition at line 161 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::CrsMatrix ( const RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node, LocalMatOps > > &  graph) [explicit]

Constructor specifying a pre-constructed graph.

Definition at line 192 of file Tpetra_CrsMatrix_def.hpp.


Member Function Documentation

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::insertGlobalValues ( GlobalOrdinal  globalRow,
const ArrayView< const GlobalOrdinal > &  cols,
const ArrayView< const Scalar > &  vals 
)

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.

Definition at line 588 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::insertLocalValues ( LocalOrdinal  localRow,
const ArrayView< const LocalOrdinal > &  cols,
const ArrayView< const Scalar > &  vals 
)

Insert matrix entries, using local IDs.

Precondition:
localRow is a local row belonging to the matrix on this node
isGloballyIndexed() == false
isStorageOptimized() == false
hasColMap() == true
Postcondition:
isLocallyIndexed() == true
Note:
If the matrix row already contains entries 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.

Definition at line 531 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::replaceGlobalValues ( GlobalOrdinal  globalRow,
const ArrayView< const GlobalOrdinal > &  cols,
const ArrayView< const Scalar > &  vals 
)

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.

Definition at line 712 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::replaceLocalValues ( LocalOrdinal  localRow,
const ArrayView< const LocalOrdinal > &  cols,
const ArrayView< const Scalar > &  vals 
)

Replace matrix entries, using local IDs.

All index values must be in the local space.

Definition at line 666 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::sumIntoGlobalValues ( GlobalOrdinal  globalRow,
const ArrayView< const GlobalOrdinal > &  cols,
const ArrayView< const Scalar > &  vals 
)

Sum into multiple 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.

Definition at line 754 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::sumIntoLocalValues ( LocalOrdinal  globalRow,
const ArrayView< const LocalOrdinal > &  cols,
const ArrayView< const Scalar > &  vals 
)

Sum into multiple entries, using local IDs.

All index values must be in the local space.

Precondition:
localRow is a local row belonging to the matrix on this node.

Definition at line 796 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::setAllToScalar ( const Scalar &  alpha)

Set all matrix entries equal to scalarThis.

Definition at line 1052 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::scale ( const Scalar &  alpha)

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

Definition at line 1014 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::globalAssemble ( )

Communicate non-local contributions to other nodes.

Definition at line 1124 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::resumeFill ( )

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

Definition at line 1359 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::fillComplete ( const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  domainMap,
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  rangeMap,
OptimizeOption  os = DoOptimizeStorage 
)

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

Definition at line 1402 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::fillComplete ( OptimizeOption  os = DoOptimizeStorage)

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

Definition at line 1394 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
const RCP< const Comm< int > > & Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getComm ( ) const [virtual]

Returns the communicator.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 223 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
RCP< Node > Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getNode ( ) const [virtual]

Returns the underlying node.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 230 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > & Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getRowMap ( ) const [virtual]

Returns the Map that describes the row distribution in this matrix.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 357 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > & Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getColMap ( ) const [virtual]

Returns the Map that describes the column distribution in this matrix.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 364 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
RCP< const RowGraph< LocalOrdinal, GlobalOrdinal, Node > > Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getGraph ( ) const [virtual]

Returns the RowGraph associated with this matrix.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 385 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node, LocalMatOps > > Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getCrsGraph ( ) const

Returns the CrsGraph associated with this matrix.

Definition at line 393 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
global_size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getGlobalNumRows ( ) const [virtual]

Returns the number of global rows in this matrix.

Undefined if isFillActive().

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 290 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
global_size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getGlobalNumCols ( ) const [virtual]

Returns the number of global columns in the matrix.

Undefined if isFillActive().

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 296 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getNodeNumRows ( ) const [virtual]

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

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 302 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getNodeNumCols ( ) const [virtual]

Returns the number of columns connected to the locally owned rows of this matrix.

Throws std::runtime_error if hasColMap() == false

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 308 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
GlobalOrdinal Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getIndexBase ( ) const [virtual]

Returns the index base for global indices for this matrix.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 350 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
global_size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getGlobalNumEntries ( ) const [virtual]

Returns the global number of entries in this matrix.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 278 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getNodeNumEntries ( ) const [virtual]

Returns the local number of entries in this matrix.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 284 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getNumEntriesInGlobalRow ( GlobalOrdinal  globalRow) const [virtual]

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

Returns OrdinalTraits<size_t>::invalid() if the specified global row does not belong to this matrix.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 326 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getNumEntriesInLocalRow ( LocalOrdinal  localRow) const [virtual]

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 Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 332 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
global_size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getGlobalNumDiags ( ) const [virtual]

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

Undefined if isFillActive().

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 314 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getNodeNumDiags ( ) const [virtual]

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

Undefined if isFillActive().

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 320 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getGlobalMaxNumRowEntries ( ) const [virtual]

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

Undefined if isFillActive().

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 338 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getNodeMaxNumRowEntries ( ) const [virtual]

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

Undefined if isFillActive().

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 344 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::hasColMap ( ) const [virtual]

Indicates whether the matrix has a well-defined column map.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 272 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::isLowerTriangular ( ) const [virtual]

Indicates whether the matrix is lower triangular.

Undefined if isFillActive().

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 400 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::isUpperTriangular ( ) const [virtual]

Indicates whether the matrix is upper triangular.

Undefined if isFillActive().

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 406 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::isLocallyIndexed ( ) const [virtual]

If matrix indices are in the local range, this function returns true. Otherwise, this function returns false. */.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 260 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::isGloballyIndexed ( ) const [virtual]

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

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 266 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::isFillComplete ( ) const [virtual]

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

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 242 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::isFillActive ( ) const

Returns true if resumeFill() has been called and the matrix is in edit mode.

Definition at line 248 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::isStorageOptimized ( ) const

Returns true if storage has been optimized.

Optimized storage means that the allocation of each row is equal to the number of entries. The effect is that a pass through the matrix, i.e., during a mat-vec, requires minimal memory traffic. One limitation of optimized storage is that no new indices can be added to the matrix.

Definition at line 254 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
ProfileType Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getProfileType ( ) const

Returns true if the matrix was allocated with static data structures.

Definition at line 236 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::isStaticGraph ( ) const

Indicates that the graph is static, so that new entries cannot be added to this matrix. */.

Definition at line 412 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getGlobalRowCopy ( GlobalOrdinal  GlobalRow,
const ArrayView< GlobalOrdinal > &  Indices,
const ArrayView< Scalar > &  Values,
size_t &  NumEntries 
) const [virtual]

Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage.

Parameters:
LocalRow- (In) Global row number for which indices are desired.
Indices- (Out) Global column indices corresponding to values.
Values- (Out) Matrix values.
NumEntries- (Out) Number of indices.

Note: A std::runtime_error exception is thrown if either Indices or Values is not large enough to hold the data associated with row GlobalRow. If GlobalRow does not belong to this node, then Indices and Values are unchanged and NumIndices is returned as OrdinalTraits<size_t>::invalid().

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 913 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getLocalRowCopy ( LocalOrdinal  LocalRow,
const ArrayView< LocalOrdinal > &  Indices,
const ArrayView< Scalar > &  Values,
size_t &  NumEntries 
) const [virtual]

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 Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 871 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getGlobalRowView ( GlobalOrdinal  GlobalRow,
ArrayView< const GlobalOrdinal > &  indices,
ArrayView< const Scalar > &  values 
) const [virtual]

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 Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 984 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getLocalRowView ( LocalOrdinal  LocalRow,
ArrayView< const LocalOrdinal > &  indices,
ArrayView< const Scalar > &  values 
) const [virtual]

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 Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 955 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getLocalDiagCopy ( Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &  diag) const [virtual]

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 Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1082 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
template<class DomainScalar , class RangeScalar >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::multiply ( const MultiVector< DomainScalar, LocalOrdinal, GlobalOrdinal, Node > &  X,
MultiVector< RangeScalar, LocalOrdinal, GlobalOrdinal, Node > &  Y,
Teuchos::ETransp  trans,
RangeScalar  alpha,
RangeScalar  beta 
) const

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.

Definition at line 1530 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
template<class DomainScalar , class RangeScalar >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::solve ( const MultiVector< RangeScalar, LocalOrdinal, GlobalOrdinal, Node > &  Y,
MultiVector< DomainScalar, LocalOrdinal, GlobalOrdinal, Node > &  X,
Teuchos::ETransp  trans 
) const

Solves a linear system when the underlying matrix is triangular.

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.

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 solve() not be called directly; instead, use the CrsMatrixSolveOp, as it will handle the import/exprt operations required to apply a matrix with non-trivial communication needs.

Both are required to have constant stride. However, unlike multiply(), it is permissible for &X == &Y. No runtime checking will be performed in a non-debug build.

Definition at line 1564 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::apply ( const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &  X,
MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &  Y,
Teuchos::ETransp  mode = Teuchos::NO_TRANS,
Scalar  alpha = ScalarTraits<Scalar>::one(),
Scalar  beta = ScalarTraits<Scalar>::zero() 
) const [virtual]

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 Tpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1516 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::hasTransposeApply ( ) const [virtual]

Indicates whether this operator supports applying the adjoint operator.

Reimplemented from Tpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 418 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > & Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getDomainMap ( ) const [virtual]

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

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

Definition at line 371 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > & Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getRangeMap ( ) const [virtual]

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

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

Definition at line 378 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
std::string Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::description ( ) const [virtual]

Return a simple one-line description of this object.

Reimplemented from Teuchos::Describable.

Definition at line 1640 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::describe ( Teuchos::FancyOStream out,
const Teuchos::EVerbosityLevel  verbLevel = Teuchos::Describable::verbLevel_default 
) const [virtual]

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

Reimplemented from Teuchos::Describable.

Definition at line 1662 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::checkSizes ( const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &  source) [virtual]

Allows the source and target (this) objects to be compared for compatibility.

Return true if they are compatible, return false if they aren't.

Implements Tpetra::DistObject< char, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1812 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::copyAndPermute ( const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &  source,
size_t  numSameIDs,
const ArrayView< const LocalOrdinal > &  permuteToLIDs,
const ArrayView< const LocalOrdinal > &  permuteFromLIDs 
) [virtual]

Perform copies and permutations that are local to this image.

Parameters:
sourceIn On entry, the DistObject that we are importing from.
numSameIDsIn On entry, the number of elements that are the same on the source and dest objects. (i.e. The element is owned by the same image in both source and dest, and no permutation occurs.)
numPermuteIDsIn On entry, the number of elements that are locally permuted between source and dest objects.
permuteToLIDsIn On entry, contains a list of the elements that are permuted. (Listed by their LID in the destination DistObject.)
permuteFromLIDsIn On entry, contains a list of the elements that are permuted. (Listed by their LID in the source DistObject.)

Implements Tpetra::DistObject< char, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1833 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::packAndPrepare ( const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &  source,
const ArrayView< const LocalOrdinal > &  exportLIDs,
Array< char > &  exports,
const ArrayView< size_t > &  numPacketsPerLID,
size_t &  constantNumPackets,
Distributor distor 
) [virtual]

Perform any packing or preparation required for communication.

Parameters:
sourceIn On entry, the DistObject that we are importing from.
exportLIDsIn On entry, a list of the entries we will be sending to other images. (Listed by their LID in the source DistObject.)
exportsOut On exit, buffer for data we will be sending out.
numPacketsPerLIDOut On exit, numPacketsPerLID[i] contains the number of packets to be exported for exportLIDs[i].
constantNumPacketsOut On exit, 0 if numPacketsPerLID has variable contents (different size for each LID). If nonzero, then it is expected that num-packets-per-LID is constant, and constantNumPackets holds that value.
distorIn On entry, contains the Distributor object we are using.

Implements Tpetra::DistObject< char, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1902 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::unpackAndCombine ( const ArrayView< const LocalOrdinal > &  importLIDs,
const ArrayView< const char > &  imports,
const ArrayView< size_t > &  numPacketsPerLID,
size_t  constantNumPackets,
Distributor distor,
CombineMode  CM 
) [virtual]

Perform any unpacking and combining after communication.

Parameters:
importLIDsIn On entry, a list of the entries we received from other images. (Listed by their LID in the target DistObject.)
importsIn Buffer containing data we received.
numPacketsPerLIDIn numPacketsPerLID[i] contains the number of packets imported for importLIDs[i].
constantNumPacketsIn If nonzero, then numPacketsPerLID is constant (same value in all entries) and constantNumPackets is that value.
distorIn The Distributor object we are using.
CMIn The Tpetra::CombineMode to use when combining the imported entries with existing entries.

Implements Tpetra::DistObject< char, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 1996 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::optimizeStorage ( )

Deprecated. Re-allocate the data into contiguous storage.

This method is deprecated and will be removed in a future version of Tpetra, as the implementation of storage optimization has been below Tpetra to Kokkos.

Currently, the implementation simply calls resumeFill() and then fillComplete(OptimizeStorage). As such, it is required to be called by all nodes that participate in the associated communicator.

Definition at line 2112 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getGlobalRowView ( GlobalOrdinal  GlobalRow,
ArrayRCP< const GlobalOrdinal > &  indices,
ArrayRCP< const Scalar > &  values 
) const [virtual]

Deprecated. Get a persisting const view of the entries in a specified global row of this matrix.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 2060 of file Tpetra_CrsMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class LocalMatOps >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getLocalRowView ( LocalOrdinal  LocalRow,
ArrayRCP< const LocalOrdinal > &  indices,
ArrayRCP< const Scalar > &  values 
) const [virtual]

Deprecated. Get a persisting const view of the entries in a specified local row of this matrix.

Implements Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >.

Definition at line 2087 of file Tpetra_CrsMatrix_def.hpp.

void Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::doImport ( const DistObject< char , LocalOrdinal, GlobalOrdinal, Node > &  source,
const Import< LocalOrdinal, GlobalOrdinal, Node > &  importer,
CombineMode  CM 
) [inherited]

Import.

void Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::doImport ( const DistObject< char , LocalOrdinal, GlobalOrdinal, Node > &  source,
const Export< LocalOrdinal, GlobalOrdinal, Node > &  exporter,
CombineMode  CM 
) [inherited]

Import (using an Exporter)

void Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::doExport ( const DistObject< char , LocalOrdinal, GlobalOrdinal, Node > &  dest,
const Export< LocalOrdinal, GlobalOrdinal, Node > &  exporter,
CombineMode  CM 
) [inherited]

Export.

void Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::doExport ( const DistObject< char , LocalOrdinal, GlobalOrdinal, Node > &  dest,
const Import< LocalOrdinal, GlobalOrdinal, Node > &  importer,
CombineMode  CM 
) [inherited]

Export (using an Importer)

bool Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::isDistributed ( ) const [inline, inherited]

Accessor for whether or not this is a global object.

const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::getMap ( ) const [inline, inherited]

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

Definition at line 111 of file Tpetra_DistObject.hpp.

void Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::print ( std::ostream &  os) const [inherited]

Print method.

virtual void Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::doTransfer ( const DistObject< char , LocalOrdinal, GlobalOrdinal, Node > &  source,
CombineMode  CM,
size_t  numSameIDs,
const Teuchos::ArrayView< const LocalOrdinal > &  permuteToLIDs,
const Teuchos::ArrayView< const LocalOrdinal > &  permuteFromLIDs,
const Teuchos::ArrayView< const LocalOrdinal > &  remoteLIDs,
const Teuchos::ArrayView< const LocalOrdinal > &  exportLIDs,
Distributor distor,
ReverseOption  revOp 
) [protected, virtual, inherited]

Perform transfer (redistribution) of data across memory images.


The documentation for this class was generated from the following files:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines