49 #ifndef XPETRA_CRSMATRIXWRAP_DEF_HPP
50 #define XPETRA_CRSMATRIXWRAP_DEF_HPP
52 #include <Tpetra_KokkosCompat_DefaultNode.hpp>
54 #include <Teuchos_SerialDenseMatrix.hpp>
55 #include <Teuchos_Hashtable.hpp>
60 #include "Xpetra_MultiVector.hpp"
70 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
72 : finalDefaultView_(false) {
77 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
79 size_t maxNumEntriesPerRow)
80 : finalDefaultView_(false) {
85 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
87 const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc)
88 : finalDefaultView_(false) {
93 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
95 : finalDefaultView_(false) {
103 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
105 : finalDefaultView_(false) {
113 #ifdef HAVE_XPETRA_TPETRA
114 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
116 : finalDefaultView_(false) {
124 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
126 const RCP<const Map> &domainMap,
const RCP<const Map> &rangeMap,
127 const Teuchos::RCP<Teuchos::ParameterList> ¶ms)
128 : finalDefaultView_(false) {
137 #warning "Xpetra Kokkos interface for CrsMatrix is enabled (HAVE_XPETRA_KOKKOS_REFACTOR) but Tpetra is disabled. The Kokkos interface needs Tpetra to be enabled, too."
141 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
143 : finalDefaultView_(matrix->isFillComplete()) {
151 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
153 : finalDefaultView_(false) {
161 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
165 const RCP<ParameterList> ¶mList)
166 : finalDefaultView_(false) {
174 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
177 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
179 matrixData_->insertGlobalValues(globalRow, cols, vals);
182 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
184 matrixData_->insertLocalValues(localRow, cols, vals);
187 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
189 const ArrayView<const GlobalOrdinal> &cols,
190 const ArrayView<const Scalar> &vals) { matrixData_->replaceGlobalValues(globalRow, cols, vals); }
192 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
194 const ArrayView<const LocalOrdinal> &cols,
195 const ArrayView<const Scalar> &vals) { matrixData_->replaceLocalValues(localRow, cols, vals); }
197 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
200 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
202 matrixData_->scale(alpha);
205 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
207 matrixData_->resumeFill(params);
210 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
212 matrixData_->fillComplete(domainMap, rangeMap, params);
218 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
220 matrixData_->fillComplete(params);
226 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
228 return matrixData_->getGlobalNumRows();
231 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
233 return matrixData_->getGlobalNumCols();
236 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
238 return matrixData_->getLocalNumRows();
241 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
243 return matrixData_->getGlobalNumEntries();
246 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
248 return matrixData_->getLocalNumEntries();
251 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
253 return matrixData_->getNumEntriesInLocalRow(localRow);
256 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
258 return matrixData_->getNumEntriesInGlobalRow(globalRow);
261 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
263 return matrixData_->getGlobalMaxNumRowEntries();
266 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
268 return matrixData_->getLocalMaxNumRowEntries();
271 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
273 return matrixData_->isLocallyIndexed();
276 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
278 return matrixData_->isGloballyIndexed();
281 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
283 return matrixData_->isFillComplete();
286 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
288 const ArrayView<LocalOrdinal> &Indices,
289 const ArrayView<Scalar> &Values,
290 size_t &NumEntries)
const {
291 matrixData_->getLocalRowCopy(LocalRow, Indices, Values, NumEntries);
294 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
296 matrixData_->getGlobalRowView(GlobalRow, indices, values);
299 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
301 matrixData_->getLocalRowView(LocalRow, indices, values);
304 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
306 matrixData_->getLocalDiagCopy(diag);
309 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
311 matrixData_->getLocalDiagOffsets(offsets);
314 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
316 matrixData_->getLocalDiagCopy(diag, offsets);
319 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
321 return matrixData_->getFrobeniusNorm();
324 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
326 matrixData_->leftScale(x);
329 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
331 matrixData_->rightScale(x);
334 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
336 return matrixData_->haveGlobalConstants();
339 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
342 Teuchos::ETransp mode,
345 matrixData_->apply(X, Y, mode, alpha, beta);
348 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
351 Teuchos::ETransp mode,
354 bool sumInterfaceValues,
356 const Teuchos::ArrayRCP<LocalOrdinal> ®ionInterfaceLIDs)
const {
357 matrixData_->apply(X, Y, mode, alpha, beta, sumInterfaceValues, regionInterfaceImporter, regionInterfaceLIDs);
360 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
362 return matrixData_->getDomainMap();
365 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
367 return matrixData_->getRangeMap();
370 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
373 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
380 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
382 matrixData_->removeEmptyProcessesInPlace(newMap);
383 this->operatorViewTable_.get(this->GetCurrentViewLabel())->SetRowMap(matrixData_->getRowMap());
384 this->operatorViewTable_.get(this->GetCurrentViewLabel())->SetColMap(matrixData_->getColMap());
387 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
389 return matrixData_->getMap();
392 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
399 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
406 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
413 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
420 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
422 return "Xpetra::CrsMatrixWrap";
425 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
435 matrixData_->describe(out, verbLevel);
440 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
442 Teuchos::LabeledObject::setObjectLabel(objectLabel);
443 matrixData_->setObjectLabel(objectLabel);
446 #ifdef HAVE_XPETRA_TPETRA
447 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
450 return matrixData_->getLocalMatrixHost();
452 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
455 return matrixData_->getLocalMatrixDevice();
459 #warning "Xpetra Kokkos interface for CrsMatrix is enabled (HAVE_XPETRA_KOKKOS_REFACTOR) but Tpetra is disabled. The Kokkos interface needs Tpetra to be enabled, too."
463 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
466 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
469 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
474 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
477 this->defaultViewLabel_ =
"point";
478 this->CreateView(this->GetDefaultViewLabel(), matrixData_->getRowMap(), matrixData_->getColMap());
481 this->currentViewLabel_ = this->GetDefaultViewLabel();
484 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
486 if ((finalDefaultView_ ==
false) && matrixData_->isFillComplete()) {
489 finalDefaultView_ =
true;
494 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
497 Teuchos::Hashtable<viewLabel_t, RCP<MatrixView>> dummy_table;
500 finalDefaultView_ = M->isFillComplete();
508 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
510 return matrixData_->GetStorageBlockSize();
513 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
518 matrixData_->residual(X, B, R);
523 #endif // ifndef XPETRA_CRSMATRIXWRAP_DEF_HPP
void setObjectLabel(const std::string &objectLabel)
RCP< CrsMatrix > getCrsMatrix() const
void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale matrix using the given vector entries.
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs.
virtual local_matrix_type getLocalMatrixDevice() const
global_size_t getGlobalNumRows() const
Returns the number of global rows in this matrix.
size_t getLocalNumEntries() const
Returns the local number of entries in this matrix.
const Teuchos::RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
static RCP< CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Build(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node >> &rowMap)
Constructor for empty matrix (intended use is an import/export target - can't insert entries directly...
const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Returns the Map associated with the domain of this operator. This will be null until fillComplete() i...
bool hasCrsGraph() const
Supports the getCrsGraph() call.
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this.
LocalOrdinal GetStorageBlockSize() const
Returns the block size of the storage mechanism, which is usually 1, except for Tpetra::BlockCrsMatri...
Exception throws to report errors in the internal logical of the program.
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map > &newMap)
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 doExport(const Matrix &dest, const Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
Teuchos::Hashtable< viewLabel_t, RCP< MatrixView > > operatorViewTable_
void updateDefaultView() const
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Get Frobenius norm of the matrix.
void doImport(const Matrix &source, const Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
std::string description() const
Return a simple one-line description of this object.
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
const viewLabel_t & GetCurrentViewLabel() const
CrsMatrix::local_matrix_type local_matrix_type
void getLocalDiagCopy(Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row idices.
CrsMatrixWrap(const RCP< const Map > &rowMap)
Constructor for a dynamic profile matrix (Epetra only)
void residual(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &R) const
Compute a residual R = B - (*this) * X.
RCP< const CrsGraph > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
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.
size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
size_t global_size_t
Global size_t object.
void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs.
void fillComplete(const RCP< const Map > &domainMap, const RCP< const Map > &rangeMap, const RCP< Teuchos::ParameterList > ¶ms=null)
Signal that data entry is complete, specifying domain and range maps.
RCP< CrsMatrix > matrixData_
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
virtual ~CrsMatrixWrap()
Destructor.
const RCP< const Map > & getColMap() const
Returns the Map that describes the column distribution in this matrix. This might be null until fillC...
void resumeFill(const RCP< ParameterList > ¶ms=null)
void replaceCrsMatrix(RCP< CrsMatrix > &M)
Expert only.
bool isFillComplete() const
Returns true if fillComplete() has been called and the matrix is in compute mode. ...
virtual void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalar.
const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Concrete implementation of Xpetra::Matrix.
global_size_t getGlobalNumCols() const
Returns the number of global columns in the matrix.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs.
size_t getLocalNumRows() const
Returns the number of matrix rows owned on the calling node.
CombineMode
Xpetra::Combine Mode enumerable type.
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale matrix using the given vector entries.
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
size_t getLocalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
const viewLabel_t & GetDefaultViewLabel() const
virtual void apply(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, Xpetra::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.
KokkosSparse::CrsMatrix< impl_scalar_type, LocalOrdinal, execution_space, void, typename local_graph_type::size_type > local_matrix_type
The specialization of Kokkos::CrsMatrix that represents the part of the sparse matrix on each MPI pro...
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries in the specified global row.
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.
virtual local_matrix_type::HostMirror getLocalMatrixHost() const
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
Xpetra-specific matrix class.
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 calli...