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

Sparse matrix that presents a compressed sparse row interface. More...

#include <Tpetra_CrsMatrix_decl.hpp>

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

List of all members.

Public Types

typedef Scalar scalar_type
 The type of the entries of the input and output multivectors.
typedef LocalOrdinal local_ordinal_type
 The local index type.
typedef GlobalOrdinal global_ordinal_type
 The global index type.
typedef Node node_type
 The Kokkos Node type.

Public Member Functions

Constructor/Destructor Methods
 CrsMatrix (const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
 Constructor specifying fixed number of entries for each row.
 CrsMatrix (const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
 Constructor specifying (possibly different) number of entries in each row.
 CrsMatrix (const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
 Constructor specifying column Map and fixed number of entries for each row.
 CrsMatrix (const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
 Constructor specifying column Map and number of entries in each row.
 CrsMatrix (const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node, LocalMatOps > > &graph, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
 Constructor specifying a previously constructed graph.
virtual ~CrsMatrix ()
 Destructor.
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 (const RCP< ParameterList > &params=null)
void fillComplete (const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null)
 Signal that data entry is complete, specifying domain and range maps.
void fillComplete (const RCP< ParameterList > &params=null)
 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
 Number of global elements in the row map of this matrix.
global_size_t getGlobalNumCols () const
 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.
ScalarTraits< Scalar >
::magnitudeType 
getFrobeniusNorm () const
 Returns the Frobenius norm of the 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 indices.
void leftScale (const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
 
void rightScale (const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
 
Advanced templated methods
template<class DomainScalar , class RangeScalar >
void localMultiply (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 localSolve (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.
template<class T >
RCP< CrsMatrix< T,
LocalOrdinal, GlobalOrdinal,
Node > > 
convert () const
 Returns another CrsMatrix with the same entries, but represented as a different scalar type.
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
 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)
 Compare the source and target (this) objects 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 process.
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 Teuchos::ArrayView< const LocalOrdinal > &importLIDs, const Teuchos::ArrayView< const char > &imports, const Teuchos::ArrayView< size_t > &numPacketsPerLID, size_t constantNumPackets, Distributor &distor, CombineMode combineMode)
 Unpack the imported column indices and values, and combine into matrix.
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.
template<class DomainScalar , class RangeScalar >
TPETRA_DEPRECATED void multiply (const MultiVector< DomainScalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< RangeScalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp trans, RangeScalar alpha, RangeScalar beta) const
 Deprecated. Replaced by localMultiply().
template<class DomainScalar , class RangeScalar >
TPETRA_DEPRECATED void solve (const MultiVector< RangeScalar, LocalOrdinal, GlobalOrdinal, Node > &Y, MultiVector< DomainScalar, LocalOrdinal, GlobalOrdinal, Node > &X, Teuchos::ETransp trans) const
 Deprecated. Replaced by localSolve().
TPETRA_DEPRECATED void fillComplete (const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, OptimizeOption os)
 Deprecated. Now takes a ParameterList.
TPETRA_DEPRECATED void fillComplete (OptimizeOption os)
 Deprecated. Now takes a ParameterList.
Public methods for redistributing data
void doImport (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
 Import using an Import object ("forward mode").
void doImport (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
 Import using an Export object ("reverse mode").
void doExport (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
 Export using an Export object ("forward mode").
void doExport (const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
 Export using an Import object ("reverse mode").
Attribute accessor methods
bool isDistributed () const
 Whether this is a globally distributed object.
const Teuchos::RCP< const Map
< LocalOrdinal, GlobalOrdinal,
Node > > & 
getMap () const
 The Map with which this DistObject was constructed.
I/O methods
void print (std::ostream &os) const
 Print this object to the given output stream.

Protected Member Functions

void allocateValues (ELocalGlobal lg, GraphAllocationStatus gas)
 Allocate values (and optionally indices) using the Node.
void sortEntries ()
 Sort the entries of each row by their column indices.
void mergeRedundantEntries ()
 Merge entries in each row with the same column indices.
void clearGlobalConstants ()
 Clear matrix properties that require collectives.
void computeGlobalConstants ()
 Compute matrix properties that require collectives.
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)
 Redistribute data across memory images.
virtual void createViews () const
 Hook for creating a const view.
virtual void createViewsNonConst (Kokkos::ReadWriteOption rwo)
 Hook for creating a nonconst view.
virtual void releaseViews () const
 Hook for releasing views.

Protected Attributes

bool fillComplete_
 Whether the matrix is fill complete.
std::map< GlobalOrdinal, Array
< std::pair< GlobalOrdinal,
Scalar > > > 
nonlocals_
 Nonlocal data added using insertGlobalValues().
Teuchos::RCP< const Map
< LocalOrdinal, GlobalOrdinal,
Node > > 
map_
 The Map over which this object is distributed.

Related Functions

(Note that these are not member functions.)

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node >
Teuchos::RCP< CrsMatrix
< Scalar, LocalOrdinal,
GlobalOrdinal, Node > > 
createCrsMatrix (const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t maxNumEntriesPerRow=0, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
 Non-member function to create an empty CrsMatrix given a row map and a non-zero profile.
template<class CrsMatrixType >
Teuchos::RCP< CrsMatrixType > importAndFillCompleteCrsMatrix (const Teuchos::RCP< const CrsMatrixType > &sourceMatrix, const Import< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > &importer, const Teuchos::RCP< const Map< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
 Nonmember CrsMatrix constructor that fuses Import and fillComplete().
template<class CrsMatrixType >
Teuchos::RCP< CrsMatrixType > exportAndFillCompleteCrsMatrix (const Teuchos::RCP< const CrsMatrixType > &sourceMatrix, const Export< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > &exporter, const Teuchos::RCP< const Map< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
 Nonmember CrsMatrix constructor that fuses Export and fillComplete().
template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node >
void extractBlockDiagonals (const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &matrix, const Teuchos::ArrayView< const LocalOrdinal > &first_points, Teuchos::ArrayRCP< Scalar > &out_diags, Teuchos::ArrayRCP< LocalOrdinal > &out_offsets)
 Extracts the block diagonals from a RowMatrix into contiguous, host storage.
template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node >
void extractBlockDiagonals (const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &matrix, const Tpetra::BlockMap< LocalOrdinal, GlobalOrdinal, Node > &block_map, Teuchos::ArrayRCP< Scalar > &out_diags, Teuchos::ArrayRCP< LocalOrdinal > &out_offsets)
 Extracts the block diagonals from a RowMatrix into contiguous, host storage using a BlockMap.
template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node >
void extractBlockRow (LocalOrdinal localBlockRow, const RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &matrix, const BlockMap< LocalOrdinal, GlobalOrdinal, Node > &block_row_map, const BlockMap< LocalOrdinal, GlobalOrdinal, Node > &block_col_map, ArrayRCP< ArrayRCP< Scalar > > &out_block_entries, ArrayRCP< LocalOrdinal > &out_block_indices)
 Extracts block elements from a RowMatrix into contiguous, host storage.

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 >

Sparse matrix that presents a compressed sparse row interface.

Template Parameters:
ScalarThe type of the numerical entries of the matrix. (You can use real-valued or complex-valued types here, unlike in Epetra, where the scalar type is always double.)
LocalOrdinalThe type of local indices. Same as the LocalOrdinal template parameter of Map objects used by this matrix. (In Epetra, this is just int.)
GlobalOrdinalThe type of global indices. Same as the GlobalOrdinal template parameter of Map objects used by this matrix. (In Epetra, this is just int. One advantage of Tpetra over Epetra is that you can use a 64-bit integer type here if you want to solve big problems.)
NodeA class implementing on-node shared-memory parallel operations. It must implement the Kokkos Node API. The default Node type depends on your Trilinos build options.
LocalMatOpsA local sparse matrix operations class. It must implement the Kokkos CRS Ops API.

This class implements a distributed-memory parallel sparse matrix, and provides sparse matrix-vector multiply (including transpose) and sparse triangular solve operations. It provides access by rows to the elements of the matrix, as if the local data were stored in compressed sparse row format. (Implementations are _not_ required to store the data in this way internally.) This class has an interface like that of Epetra_CrsMatrix, but also allows insertion of data into nonowned rows, much like Epetra_FECrsMatrix.

Local vs. Global

The distinction between local and global indices might confuse new Tpetra users. _Global_ indices represent the rows and columns uniquely over the entire matrix, which may be distributed over multiple processes. _Local_ indices are local to the process that owns them. If global index G is owned by process P, then there is a unique local index L on process P corresponding to G. If the local index L is valid on process P, then there is a unique global index G owned by P corresponding to the pair (L, P). However, multiple processes might own the same global index (an "overlapping Map"), so a global index G might correspond to multiple (L, P) pairs. In summary, local indices on a process correspond to matrix rows or columns owned by that process.

We summarize the different between local and global indices because many of CrsMatrix's methods for adding, modifying, or accessing entries come in versions that take either local or global indices. The matrix itself may store indices either as local or global. You should only use the method version corresponding to the current state of the matrix. For example, getGlobalRowView() returns a view to the indices represented as global; it is incorrect to call this method if the matrix is storing indices as local. Call the isGloballyIndexed() or isLocallyIndexed() methods to find out whether the matrix currently stores indices as local or global.

Method methods that work with global indices only allow operations on indices owned by the calling process. For example, methods that take a global row index expect that row to be owned by the calling process. Access to nonlocal (i.e., not owned by the calling process) 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. It allows you to add data to nonlocal rows. These data are stored locally and communicated to the appropriate node on the next call to globalAssemble() or fillComplete() (the latter calls the former).

Examples:

CrsMatrix_BlockExtraction.cpp, CrsMatrix_NonlocalAfterResume.hpp, LocalMatOpExample.cpp, Tpetra_Power_Method.cpp, and Tpetra_Power_Method_From_File.cpp.

Definition at line 164 of file Tpetra_CrsMatrix_decl.hpp.


Member Typedef Documentation

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

The type of the entries of the input and output multivectors.

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

Definition at line 167 of file Tpetra_CrsMatrix_decl.hpp.

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

The local index type.

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

Definition at line 168 of file Tpetra_CrsMatrix_decl.hpp.

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

The global index type.

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

Definition at line 169 of file Tpetra_CrsMatrix_decl.hpp.

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

The Kokkos Node type.

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

Definition at line 170 of file Tpetra_CrsMatrix_decl.hpp.


Constructor & Destructor Documentation

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::CrsMatrix ( const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  rowMap,
size_t  maxNumEntriesPerRow,
ProfileType  pftype = DynamicProfile,
const Teuchos::RCP< Teuchos::ParameterList > &  params = Teuchos::null 
)

Constructor specifying fixed number of entries for each row.

Parameters:
rowMap[in] Distribution of rows of the matrix.
maxNumEntriesPerRow[in] Maximum number of matrix entries per row. If pftype==DynamicProfile, this is only a hint, and you can set this to zero without affecting correctness. If pftype==StaticProfile, this sets the amount of storage allocated, and you cannot exceed this number of entries in any row.
pftype[in] Whether to allocate storage dynamically (DynamicProfile) or statically (StaticProfile).
params[in/out] Optional list of parameters. If not null, any missing parameters will be filled in with their default values.
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::CrsMatrix ( const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  rowMap,
const ArrayRCP< const size_t > &  NumEntriesPerRowToAlloc,
ProfileType  pftype = DynamicProfile,
const Teuchos::RCP< Teuchos::ParameterList > &  params = Teuchos::null 
)

Constructor specifying (possibly different) number of entries in each row.

Parameters:
rowMap[in] Distribution of rows of the matrix.
NumEntriesPerRowToAlloc[in] Maximum number of matrix entries to allocate for each row. If pftype==DynamicProfile, this is only a hint. If pftype==StaticProfile, this sets the amount of storage allocated, and you cannot exceed the allocated number of entries for any row.
pftype[in] Whether to allocate storage dynamically (DynamicProfile) or statically (StaticProfile).
params[in/out] Optional list of parameters. If not null, any missing parameters will be filled in with their default values.
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::CrsMatrix ( const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  rowMap,
const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  colMap,
size_t  maxNumEntriesPerRow,
ProfileType  pftype = DynamicProfile,
const Teuchos::RCP< Teuchos::ParameterList > &  params = Teuchos::null 
)

Constructor specifying column Map and fixed number of entries for each row.

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

Parameters:
rowMap[in] Distribution of rows of the matrix.
colMap[in] Distribution of columns of the matrix.
maxNumEntriesPerRow[in] Maximum number of matrix entries per row. If pftype==DynamicProfile, this is only a hint, and you can set this to zero without affecting correctness. If pftype==StaticProfile, this sets the amount of storage allocated, and you cannot exceed this number of entries in any row.
pftype[in] Whether to allocate storage dynamically (DynamicProfile) or statically (StaticProfile).
params[in/out] Optional list of parameters. If not null, any missing parameters will be filled in with their default values.
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::CrsMatrix ( const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  rowMap,
const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  colMap,
const ArrayRCP< const size_t > &  NumEntriesPerRowToAlloc,
ProfileType  pftype = DynamicProfile,
const Teuchos::RCP< Teuchos::ParameterList > &  params = Teuchos::null 
)

Constructor specifying column Map and number of entries in each row.

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

Parameters:
rowMap[in] Distribution of rows of the matrix.
colMap[in] Distribution of columns of the matrix.
NumEntriesPerRowToAlloc[in] Maximum number of matrix entries to allocate for each row. If pftype==DynamicProfile, this is only a hint. If pftype==StaticProfile, this sets the amount of storage allocated, and you cannot exceed the allocated number of entries for any row.
pftype[in] Whether to allocate storage dynamically (DynamicProfile) or statically (StaticProfile).
params[in/out] Optional list of parameters. If not null, any missing parameters will be filled in with their default values.
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::CrsMatrix ( const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node, LocalMatOps > > &  graph,
const Teuchos::RCP< Teuchos::ParameterList > &  params = Teuchos::null 
) [explicit]

Constructor specifying a previously constructed graph.

Calling this constructor fixes the graph structure of the sparse matrix. We say in this case that the matrix has a "static graph." If you create a CrsMatrix with this constructor, you are not allowed to insert new entries into the matrix, but you are allowed to change values in the matrix.

The given graph must be fill complete. Note that calling resumeFill() on the graph makes it not fill complete, even if you had previously called fillComplete() on the graph. In that case, you must call fillComplete() on the graph again before invoking this CrsMatrix constructor.

Parameters:
graph[in] The graph structure of the sparse matrix. The graph must be fill complete.
params[in/out] Optional list of parameters. If not null, any missing parameters will be filled in with their default values.
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
virtual Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::~CrsMatrix ( ) [virtual]

Destructor.


Member Function Documentation

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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
isStorageOptimized() == false
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.
If isLocallyIndexed() == true, then the global indices will be translated to local indices via the column map; indices not present in the column map will be discarded.
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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.
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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.
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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.
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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.
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::setAllToScalar ( const Scalar &  alpha)

Set all matrix entries equal to scalarThis.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::scale ( const Scalar &  alpha)

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

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

Communicate non-local contributions to other nodes.

This method only does global assembly if there are nonlocal entries on at least one process. It does an all-reduce to find that out. If not, it returns early, without doing any more communication or work.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::resumeFill ( const RCP< ParameterList > &  params = null)

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
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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,
const RCP< ParameterList > &  params = null 
)

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
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::fillComplete ( const RCP< ParameterList > &  params = null)

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
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
const RCP<const Comm<int> >& Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getComm ( ) const [virtual]
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
RCP<Node> Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getNode ( ) const [virtual]

Returns the underlying node.

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
RCP<const CrsGraph<LocalOrdinal,GlobalOrdinal,Node,LocalMatOps> > Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getCrsGraph ( ) const

Returns the CrsGraph associated with this matrix.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
global_size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getGlobalNumRows ( ) const [virtual]

Number of global elements in the row map of this matrix.

This is <it>not</it> the number of rows in the matrix as a mathematical object. This method returns the global sum of the number of local elements in the row map on each processor, which is the row map's getGlobalNumElements(). Since the row map is not one-to-one in general, that global sum could be different than the number of rows in the matrix. If you want the number of rows in the matrix, ask the range map for its global number of elements, using the following code: global_size_t globalNumRows = getRangeMap()->getGlobalNumElements(); This method retains the behavior of Epetra, which also asks the row map for the global number of rows, rather than asking the range map.

Warning:
Undefined if isFillActive().

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
global_size_t Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getGlobalNumCols ( ) const [virtual]

Number of global columns in the matrix.

Returns the number of entries in the domain map of the matrix.

Warning:
Undefined if isFillActive().

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

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

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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.

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

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
ScalarTraits<Scalar>::magnitudeType Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::getFrobeniusNorm ( ) const [virtual]

Returns the Frobenius norm of the matrix.

Computes and returns the Frobenius norm of the matrix, defined as: $ \|A\|_F = \sqrt{\sum_{i,j} \|\a_{ij}\|^2} $

If the matrix is fill-complete, then the computed value is cached; the cache is cleared whenever resumeFill() is called. Otherwise, the value is computed every time the method is called.

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 indices.

Returns a distributed Vector object partitioned according to this matrix's row map, containing the the zero and non-zero diagonals owned by this node.

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::leftScale ( const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &  x) [virtual]
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::rightScale ( const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &  x) [virtual]
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
template<class DomainScalar , class RangeScalar >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::localMultiply ( 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.

If trans is Teuchos::NO_TRANS, then 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. Otherwise, then X is should be described by the row map of the matrix and Y should be described by the column 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.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
template<class DomainScalar , class RangeScalar >
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::localSolve ( 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/Export operations required to apply a matrix with non-trivial communication needs.

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
template<class T >
RCP<CrsMatrix<T,LocalOrdinal,GlobalOrdinal,Node> > Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::convert ( ) const

Returns another CrsMatrix with the same entries, but represented as a different scalar type.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
std::string Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::description ( ) const [virtual]

A simple one-line description of this object.

Reimplemented from Tpetra::DistObject< char, LocalOrdinal, GlobalOrdinal, Node >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 Tpetra::DistObject< char, LocalOrdinal, GlobalOrdinal, Node >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
bool Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::checkSizes ( const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &  source) [virtual]

Compare the source and target (this) objects for compatibility.

Returns:
True if they are compatible, else false.

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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 process.

Parameters:
source[in] On entry, the source object, from which we are distributing. We distribute to the destination object, which is *this object.
numSameIDs[in] The umber of elements that are the same on the source and destination (this) objects. These elements are owned by the same process in both the source and destination objects. No permutation occurs.
numPermuteIDs[in] The number of elements that are locally permuted between the source and destination objects.
permuteToLIDs[in] List of the elements that are permuted. They are listed by their LID in the destination object.
permuteFromLIDs[in] List of the elements that are permuted. They are listed by their LID in the source object.

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
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:
source[in] Source object for the redistribution.
exportLIDs[in] List of the entries (as local IDs in the source object) we will be sending to other images.
exports[out] On exit, the buffer for data to send.
numPacketsPerLID[out] On exit, numPacketsPerLID[i] contains the number of packets to be exported for exportLIDs[i]. If constantNumPackets is nonzero, you should use that instead, and not rely on numPacketsPerLID[i] being filled.
constantNumPackets[out] 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.
distor[in] The Distributor object we are using.

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::unpackAndCombine ( const Teuchos::ArrayView< const LocalOrdinal > &  importLIDs,
const Teuchos::ArrayView< const char > &  imports,
const Teuchos::ArrayView< size_t > &  numPacketsPerLID,
size_t  constantNumPackets,
Distributor distor,
CombineMode  combineMode 
) [virtual]

Unpack the imported column indices and values, and combine into matrix.

Warning:
The allowed combineMode depends on whether the matrix's graph is static or dynamic. ADD, REPLACE, and ABSMAX are valid for a static graph, but INSERT is not. ADD and INSERT are valid for a dynamic graph; ABSMAX and REPLACE have not yet been implemented (and would require serious changes to matrix assembly in order to implement sensibly).

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

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
TPETRA_DEPRECATED 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.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
TPETRA_DEPRECATED 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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
TPETRA_DEPRECATED 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 >.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
template<class DomainScalar , class RangeScalar >
TPETRA_DEPRECATED 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

Deprecated. Replaced by localMultiply().

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
template<class DomainScalar , class RangeScalar >
TPETRA_DEPRECATED 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

Deprecated. Replaced by localSolve().

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
TPETRA_DEPRECATED 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 
)

Deprecated. Now takes a ParameterList.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
TPETRA_DEPRECATED void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::fillComplete ( OptimizeOption  os)

Deprecated. Now takes a ParameterList.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::allocateValues ( ELocalGlobal  lg,
GraphAllocationStatus  gas 
) [protected]

Allocate values (and optionally indices) using the Node.

Parameters:
gas[in] If GraphNotYetAllocated, allocate the indices of myGraph_ via allocateIndices(lg) before allocating values.
lg[in] Argument passed into myGraph_->allocateIndices(), if applicable.
Precondition:
If the graph (that is, staticGraph_) indices are already allocated, then gas must be GraphAlreadyAllocated. Otherwise, gas must be GraphNotYetAllocated. We only check for this precondition in debug mode.
If the graph indices are not already allocated, then the graph must be owned by the matrix.
template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::sortEntries ( ) [protected]

Sort the entries of each row by their column indices.

This only does anything if the graph isn't already sorted (i.e., ! myGraph_->isSorted ()). This method is called in fillComplete().

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::mergeRedundantEntries ( ) [protected]

Merge entries in each row with the same column indices.

This only does anything if the graph isn't already merged (i.e., ! myGraph_->isMerged ()). This method is called in fillComplete().

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::clearGlobalConstants ( ) [protected]

Clear matrix properties that require collectives.

This clears whatever computeGlobalConstants() (which see) computed, in preparation for changes to the matrix. The current implementation of this method does nothing.

This method is called in resumeFill().

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
void Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::computeGlobalConstants ( ) [protected]

Compute matrix properties that require collectives.

The corresponding Epetra_CrsGraph method computes things like the global number of nonzero entries, that require collectives over the matrix's communicator. The current Tpetra implementation of this method does nothing.

This method is called in fillComplete().

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 using an Import object ("forward mode").

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 Export object ("reverse mode").

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 using an Export object ("forward mode").

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 Import object ("reverse mode").

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

Whether this is a globally distributed object.

For a definition of "globally distributed" (and its opposite, "locally replicated"), see the documentation of Map's isDistributed() method.

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

The Map with which this DistObject was constructed.

Definition at line 172 of file Tpetra_DistObject.hpp.

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

Print this object to the given output stream.

We generally assume that all MPI processes can print to the given stream.

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]

Redistribute data across memory images.

Parameters:
source[in] The source object, to redistribute into the destination object, which is *this object.
CM[in] The combine mode that describes how to combine values that map to the same global ID on the same process.
permuteToLIDs[in] See copyAndPermute().
permuteFromLIDs[in] See copyAndPermute().
remoteLIDs[in] List of entries (as local IDs) in the destination object to receive from other processes.
exportLIDs[in] See packAndPrepare().
distor[in/out] The Distributor object that knows how to redistribute data.
revOp[in] Whether to do a forward or reverse mode redistribution.
virtual void Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::createViews ( ) const [inline, protected, virtual, inherited]

Hook for creating a const view.

doTransfer() calls this on the source object. By default, it does nothing, but the source object can use this as a hint to fetch data from a compute buffer on an off-CPU device (such as a GPU) into host memory.

Definition at line 358 of file Tpetra_DistObject.hpp.

virtual void Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::createViewsNonConst ( Kokkos::ReadWriteOption  rwo) [inline, protected, virtual, inherited]

Hook for creating a nonconst view.

doTransfer() calls this on the destination (*this) object. By default, it does nothing, but the destination object can use this as a hint to fetch data from a compute buffer on an off-CPU device (such as a GPU) into host memory.

Parameters:
rwo[in] Whether to create a write-only or a read-and-write view. For Kokkos Node types where compute buffers live in a separate memory space (e.g., in the device memory of a discrete accelerator like a GPU), a write-only view only requires copying from host memory to the compute buffer, whereas a read-and-write view requires copying both ways (once to read, from the compute buffer to host memory, and once to write, back to the compute buffer).

Definition at line 375 of file Tpetra_DistObject.hpp.

virtual void Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::releaseViews ( ) const [inline, protected, virtual, inherited]

Hook for releasing views.

doTransfer() calls this on both the source and destination objects, once it no longer needs to access that object's data. By default, this method does nothing. Implementations may use this as a hint to free host memory which is a view of a compute buffer, once the host memory view is no longer needed. Some implementations may prefer to mirror compute buffers in host memory; for these implementations, releaseViews() may do nothing.

Definition at line 387 of file Tpetra_DistObject.hpp.


Friends And Related Function Documentation

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node >
Teuchos::RCP< CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > createCrsMatrix ( const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &  map,
size_t  maxNumEntriesPerRow = 0,
const Teuchos::RCP< Teuchos::ParameterList > &  params = Teuchos::null 
) [related]

Non-member function to create an empty CrsMatrix given a row map and a non-zero profile.

Returns:
A dynamically allocated (DynamicProfile) matrix with specified number of nonzeros per row (defaults to zero).

Definition at line 1017 of file Tpetra_CrsMatrix_decl.hpp.

template<class CrsMatrixType >
Teuchos::RCP< CrsMatrixType > importAndFillCompleteCrsMatrix ( const Teuchos::RCP< const CrsMatrixType > &  sourceMatrix,
const Import< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > &  importer,
const Teuchos::RCP< const Map< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > > &  domainMap = Teuchos::null,
const Teuchos::RCP< const Map< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > > &  rangeMap = Teuchos::null,
const Teuchos::RCP< Teuchos::ParameterList > &  params = Teuchos::null 
) [related]

Nonmember CrsMatrix constructor that fuses Import and fillComplete().

Template Parameters:
CrsMatrixTypeA specialization of CrsMatrix.

A common use case is to create an empty destination CrsMatrix, redistribute from a source CrsMatrix (by an Import or Export operation), then call fillComplete() on the destination CrsMatrix. This constructor fuses these three cases, for an Import redistribution.

Fusing redistribution and fillComplete() exposes potential optimizations. For example, it may make constructing the column Map faster, and it may avoid intermediate unoptimized storage in the destination CrsMatrix. These optimizations may improve performance for specialized kernels like sparse matrix-matrix multiply, as well as for redistributing data after doing load balancing.

The resulting matrix is fill complete (in the sense of isFillComplete()) and has optimized storage (in the sense of isStorageOptimized()). By default, its domain Map is the domain Map of the source matrix, and its range Map is the range Map of the source matrix.

Warning:
If the target Map of the Import is a subset of the source Map of the Import, then you cannot use the default range Map. You should instead construct a nonoverlapping version of the target Map and supply that as the nondefault value of the range Map.
Parameters:
sourceMatrix[in] The source matrix from which to import. The source of an Import must have a nonoverlapping distribution.
importer[in] The Import instance containing a precomputed redistribution plan. The source Map of the Import must be the same as the row Map of sourceMatrix.
domainMap[in] Domain Map of the returned matrix. If null, we use the default, which is the domain Map of the source matrix.
rangeMap[in] Range Map of the returned matrix. If null, we use the default, which is the range Map of the source matrix.
params[in/out] Optional list of parameters. If not null, any missing parameters will be filled in with their default values.

Definition at line 1078 of file Tpetra_CrsMatrix_decl.hpp.

template<class CrsMatrixType >
Teuchos::RCP< CrsMatrixType > exportAndFillCompleteCrsMatrix ( const Teuchos::RCP< const CrsMatrixType > &  sourceMatrix,
const Export< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > &  exporter,
const Teuchos::RCP< const Map< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > > &  domainMap = Teuchos::null,
const Teuchos::RCP< const Map< typename CrsMatrixType::local_ordinal_type, typename CrsMatrixType::global_ordinal_type, typename CrsMatrixType::node_type > > &  rangeMap = Teuchos::null,
const Teuchos::RCP< Teuchos::ParameterList > &  params = Teuchos::null 
) [related]

Nonmember CrsMatrix constructor that fuses Export and fillComplete().

Template Parameters:
CrsMatrixTypeA specialization of CrsMatrix.

For justification, see the documentation of importAndFillCompleteCrsMatrix() (which is the Import analog of this function).

The resulting matrix is fill complete (in the sense of isFillComplete()) and has optimized storage (in the sense of isStorageOptimized()). By default, its domain Map is the domain Map of the source matrix, and its range Map is the range Map of the source matrix.

Parameters:
sourceMatrix[in] The source matrix from which to export. Its row Map may be overlapping, since the source of an Export may be overlapping.
exporter[in] The Export instance containing a precomputed redistribution plan. The source Map of the Export must be the same as the row Map of sourceMatrix.
domainMap[in] Domain Map of the returned matrix. If null, we use the default, which is the domain Map of the source matrix.
rangeMap[in] Range Map of the returned matrix. If null, we use the default, which is the range Map of the source matrix.
params[in/out] Optional list of parameters. If not null, any missing parameters will be filled in with their default values.

Definition at line 1153 of file Tpetra_CrsMatrix_decl.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node >
void extractBlockDiagonals ( const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &  matrix,
const Teuchos::ArrayView< const LocalOrdinal > &  first_points,
Teuchos::ArrayRCP< Scalar > &  out_diags,
Teuchos::ArrayRCP< LocalOrdinal > &  out_offsets 
) [related]

Extracts the block diagonals from a RowMatrix into contiguous, host storage.

This method does not initiate any communication, and therefore can be called safely on a single node.

Parameters:
inmatrix - The sparse matrix source for the diagonals.
infirst_points - A list of ordinals, where first_points[b] indicates the first local element in the b-th block.
outout_diags - A reference-counted array, containing the block diagonals in contiguous storage.
outout_offsets - A reference-counted array, indicating the offset to reach each block in out_diags.
Precondition:
- first_points[0] == 0
- elements in first_points are non-decreasing
- matrix.isFillComplete() == true
Postcondition:
- out_offsets.size() == block_sizes.size()
- the non-trivial b-th block is stored contiguously (column-major) in out_diags( out_offsets[b], block_size_b ), where $block\_size\_b = ( first\_points\[b+1\] - first\_points\[b\] )^2 $.
template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node >
void extractBlockDiagonals ( const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &  matrix,
const Tpetra::BlockMap< LocalOrdinal, GlobalOrdinal, Node > &  block_map,
Teuchos::ArrayRCP< Scalar > &  out_diags,
Teuchos::ArrayRCP< LocalOrdinal > &  out_offsets 
) [related]

Extracts the block diagonals from a RowMatrix into contiguous, host storage using a BlockMap.

This method does not initiate any communication, and therefore can be called safely on a single node.

Parameters:
inmatrix - The sparse matrix source for the diagonals.
inblock_map - A BlockMap describing the blocks
outout_diags - A reference-counted array, containing the block diagonals in contiguous storage.
outout_offsets - A reference-counted array, indicating the offset to reach each block in out_diags.
Precondition:
- block_map.getPointMap().isCompatible( matrix.getRowMap() )
- matrix.isFillComplete() == true
Postcondition:
- out_offsets.size() == block_sizes.size()
- the non-trivial b-th block is stored contiguously (column-major) in out_diags( out_offsets[b], block_map.getLocalBlockSize(b) )

Calls extractBlockDiagonals().

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node >
void extractBlockRow ( LocalOrdinal  localBlockRow,
const RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &  matrix,
const BlockMap< LocalOrdinal, GlobalOrdinal, Node > &  block_row_map,
const BlockMap< LocalOrdinal, GlobalOrdinal, Node > &  block_col_map,
ArrayRCP< ArrayRCP< Scalar > > &  out_block_entries,
ArrayRCP< LocalOrdinal > &  out_block_indices 
) [related]

Extracts block elements from a RowMatrix into contiguous, host storage.

This method does not initiate any communication, and therefore can be called safely on a single node.

Parameters:
inblock_row - The block row to be extracted
inblock_row_map - A BlockMap describing the row blocks
inblock_col_map - A BlockMap describing the column blocks
inmatrix - The sparse matrix source for the diagonals.
outblock_entries - The block entries
outblock_indices - The indices for the block entries
Precondition:
- block_row_map.getPointMap().isCompatible( matrix.getRowMap() )
- block_row_map.getPointMap().isCompatible( matrix.getRowMap() )
- block_col_map.getPointMap().isCompatible( matrix.getColMap() )
- block_col_map.getPointMap().isCompatible( matrix.getColMap() )
- matrix.isFillComplete() == true

Member Data Documentation

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

Whether the matrix is fill complete.

Definition at line 984 of file Tpetra_CrsMatrix_decl.hpp.

template<class Scalar, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = Kokkos::DefaultNode::DefaultNodeType, class LocalMatOps = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps>
std::map<GlobalOrdinal, Array<std::pair<GlobalOrdinal,Scalar> > > Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, LocalMatOps >::nonlocals_ [protected]

Nonlocal data added using insertGlobalValues().

These data are cleared by globalAssemble(), once it finishes redistributing them to their owning processes.

Note:
For Epetra developers: Tpetra::CrsMatrix corresponds more to Epetra_FECrsMatrix than to Epetra_CrsMatrix. The insertGlobalValues() method in Tpetra::CrsMatrix, unlike its corresponding method in Epetra_CrsMatrix, allows insertion into rows which are not owned by the calling process. The globalAssemble() method redistributes these to their owning processes.

Definition at line 998 of file Tpetra_CrsMatrix_decl.hpp.

Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > Tpetra::DistObject< char , LocalOrdinal, GlobalOrdinal, Node >::map_ [protected, inherited]

The Map over which this object is distributed.

Definition at line 390 of file Tpetra_DistObject.hpp.


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