Xpetra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Xpetra_CrsMatrixWrap_decl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Xpetra: A linear algebra interface package
4 //
5 // Copyright 2012 NTESS and the Xpetra contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 // WARNING: This code is experimental. Backwards compatibility should not be expected.
11 
12 #ifndef XPETRA_CRSMATRIXWRAP_DECL_HPP
13 #define XPETRA_CRSMATRIXWRAP_DECL_HPP
14 
15 #include <Tpetra_KokkosCompat_DefaultNode.hpp>
16 
17 #include "Xpetra_ConfigDefs.hpp"
18 #include "Xpetra_Exceptions.hpp"
19 
21 #include "Xpetra_CrsGraph.hpp"
22 #include "Xpetra_CrsMatrix.hpp"
24 
25 #include "Xpetra_Matrix.hpp"
27 #include "Xpetra_TpetraBlockCrsMatrix.hpp"
29 #include "Xpetra_TpetraCrsMatrix.hpp"
30 
31 #include <Teuchos_SerialDenseMatrix.hpp>
32 #include <Teuchos_Hashtable.hpp>
33 
38 namespace Xpetra {
39 
40 typedef std::string viewLabel_t;
41 
46 template <class Scalar,
47  class LocalOrdinal,
48  class GlobalOrdinal,
49  class Node = Tpetra::KokkosClassic::DefaultNode::DefaultNodeType>
50 class CrsMatrixWrap : public Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> {
55 #ifdef HAVE_XPETRA_TPETRA
57 #endif
60 #ifdef HAVE_XPETRA_TPETRA
62 #endif
63 
64  public:
66 
67 
69  CrsMatrixWrap(const RCP<const Map> &rowMap);
70 
72  CrsMatrixWrap(const RCP<const Map> &rowMap,
73  size_t maxNumEntriesPerRow);
74 
76  CrsMatrixWrap(const RCP<const Map> &rowMap,
77  const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc);
78 
80  CrsMatrixWrap(const RCP<const Map> &rowMap, const RCP<const Map> &colMap, size_t maxNumEntriesPerRow);
81 
83  CrsMatrixWrap(const RCP<const Map> &rowMap, const RCP<const Map> &colMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc);
84 
85 #ifdef HAVE_XPETRA_TPETRA
86  CrsMatrixWrap(const RCP<const Map> &rowMap, const RCP<const Map> &colMap, const local_matrix_type &lclMatrix, const Teuchos::RCP<Teuchos::ParameterList> &params = null);
88 
90  CrsMatrixWrap(const local_matrix_type &lclMatrix, const RCP<const Map> &rowMap, const RCP<const Map> &colMap,
91  const RCP<const Map> &domainMap = Teuchos::null, const RCP<const Map> &rangeMap = Teuchos::null,
92  const Teuchos::RCP<Teuchos::ParameterList> &params = null);
93 #else
94 #ifdef __GNUC__
95 #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."
96 #endif
97 #endif
98 
99  CrsMatrixWrap(RCP<CrsMatrix> matrix);
100 
101  CrsMatrixWrap(const RCP<const CrsGraph> &graph, const RCP<ParameterList> &paramList = Teuchos::null);
102 
103  CrsMatrixWrap(const RCP<const CrsGraph> &graph,
105  const RCP<ParameterList> &paramList = Teuchos::null);
106 
108  virtual ~CrsMatrixWrap();
109 
111 
113 
114 
116 
127  void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView<const GlobalOrdinal> &cols, const ArrayView<const Scalar> &vals);
128 
130 
137  void insertLocalValues(LocalOrdinal localRow, const ArrayView<const LocalOrdinal> &cols, const ArrayView<const Scalar> &vals);
138 
140 
145  void replaceGlobalValues(GlobalOrdinal globalRow,
146  const ArrayView<const GlobalOrdinal> &cols,
147  const ArrayView<const Scalar> &vals);
148 
150 
153  void replaceLocalValues(LocalOrdinal localRow,
154  const ArrayView<const LocalOrdinal> &cols,
155  const ArrayView<const Scalar> &vals);
156 
158  virtual void setAllToScalar(const Scalar &alpha);
159 
161  void scale(const Scalar &alpha);
162 
164 
166 
167 
176  void resumeFill(const RCP<ParameterList> &params = null);
177 
189  void fillComplete(const RCP<const Map> &domainMap, const RCP<const Map> &rangeMap, const RCP<Teuchos::ParameterList> &params = null);
190 
204  // TODO : Get ride of "Tpetra"::OptimizeOption
205  void fillComplete(const RCP<ParameterList> &params = null);
206 
208 
210 
213 
215 
218 
220  size_t getLocalNumRows() const;
221 
224 
226  size_t getLocalNumEntries() const;
227 
229 
230  size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const;
231 
233 
234  size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const;
235 
237 
239  size_t getGlobalMaxNumRowEntries() const;
240 
242 
244  size_t getLocalMaxNumRowEntries() const;
245 
247  bool isLocallyIndexed() const;
248 
250  bool isGloballyIndexed() const;
251 
253  bool isFillComplete() const;
254 
256 
268  void getLocalRowCopy(LocalOrdinal LocalRow,
269  const ArrayView<LocalOrdinal> &Indices,
270  const ArrayView<Scalar> &Values,
271  size_t &NumEntries) const;
272 
274 
283  void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView<const GlobalOrdinal> &indices, ArrayView<const Scalar> &values) const;
284 
286 
295  void getLocalRowView(LocalOrdinal LocalRow, ArrayView<const LocalOrdinal> &indices, ArrayView<const Scalar> &values) const;
296 
298 
301 
303  void getLocalDiagOffsets(Teuchos::ArrayRCP<size_t> &offsets) const;
304 
306  void getLocalDiagCopy(Xpetra::Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &diag, const Teuchos::ArrayView<const size_t> &offsets) const;
307 
309  typename ScalarTraits<Scalar>::magnitudeType getFrobeniusNorm() const;
310 
313 
316 
318  bool haveGlobalConstants() const;
319 
321 
323 
324 
326 
335  // TODO virtual=0 // TODO: Add default parameters ?
336  // void multiply(const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> & X, MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &Y, Teuchos::ETransp trans, Scalar alpha, Scalar beta) const {
337  // matrixData_->multiply(X, Y, trans, alpha, beta);
338  // }
339 
341 
343 
344 
346 
351  Teuchos::ETransp mode = Teuchos::NO_TRANS,
352  Scalar alpha = ScalarTraits<Scalar>::one(),
353  Scalar beta = ScalarTraits<Scalar>::zero()) const;
354 
358  Teuchos::ETransp mode,
359  Scalar alpha,
360  Scalar beta,
361  bool sumInterfaceValues,
362  const RCP<Import<LocalOrdinal, GlobalOrdinal, Node>> &regionInterfaceImporter,
363  const Teuchos::ArrayRCP<LocalOrdinal> &regionInterfaceLIDs) const;
364 
367  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> getDomainMap() const;
368 
371  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> getRangeMap() const;
372 
375  const RCP<const Map> &getColMap() const;
376 
378  const RCP<const Map> &getColMap(viewLabel_t viewLabel) const;
379 
380  void removeEmptyProcessesInPlace(const Teuchos::RCP<const Map> &newMap);
381 
383 
385  //{@
386 
388  const Teuchos::RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> getMap() const;
389 
391  void doImport(const Matrix &source,
393 
395  void doExport(const Matrix &dest,
397 
399  void doImport(const Matrix &source,
401 
403  void doExport(const Matrix &dest,
405 
406  // @}
407 
409 
410 
412  std::string description() const;
413 
415  void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel = Teuchos::Describable::verbLevel_default) const;
416 
418 
419  void setObjectLabel(const std::string &objectLabel);
421 
422 #ifdef HAVE_XPETRA_TPETRA
423  virtual local_matrix_type getLocalMatrixDevice() const;
424  virtual typename local_matrix_type::HostMirror getLocalMatrixHost() const;
425 #else
426 #ifdef __GNUC__
427 #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."
428 #endif
429 #endif
430 
431  // JG: Added:
432 
433  bool hasCrsGraph() const;
434 
436  RCP<const CrsGraph> getCrsGraph() const;
437 
438  RCP<CrsMatrix> getCrsMatrix() const;
439 
441  LocalOrdinal GetStorageBlockSize() const;
442 
447 
449  void replaceCrsMatrix(RCP<CrsMatrix> &M);
450 
452  private:
453  // Default view is created after fillComplete()
454  // Because ColMap might not be available before fillComplete().
455  void CreateDefaultView();
456 
457  // The colMap can be <tt>null</tt> until fillComplete() is called. The default view of the Matrix have to be updated when fillComplete() is called.
458  // If CrsMatrix::fillComplete() have been used instead of CrsMatrixWrap::fillComplete(), the default view is updated when getColMap() is called.
459  void updateDefaultView() const;
460  // The boolean finalDefaultView_ keep track of the status of the default view (= already updated or not)
461  // See also CrsMatrixWrap::updateDefaultView()
462  mutable bool finalDefaultView_;
463 
464  // The underlying matrix object
465  RCP<CrsMatrix> matrixData_;
466 
467 }; // class CrsMatrixWrap
468 
469 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
470 Teuchos::RCP<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
472  return Teuchos::rcp_dynamic_cast<Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(matrix, true)->getCrsMatrix();
473 }
474 
475 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
476 Teuchos::RCP<const Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
478  return Teuchos::rcp_dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(matrix, true)->getCrsMatrix();
479 }
480 
481 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
484  return *dynamic_cast<Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node> &>(matrix).getCrsMatrix();
485 }
486 
487 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
490  return *dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node> &>(matrix).getCrsMatrix();
491 }
492 
493 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
494 Teuchos::RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
495 toXpetra(const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> &A) {
496  auto xpTpCrs = Teuchos::rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A));
497  auto xpCrs = Teuchos::rcp_dynamic_cast<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(xpTpCrs, true);
499 }
500 
501 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
502 Teuchos::RCP<const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
503 toXpetra(Teuchos::RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> &A) {
504  auto xpTpCrs = Teuchos::rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A));
505  auto xpCrs = Teuchos::rcp_dynamic_cast<const Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(xpTpCrs, true);
507 }
508 
509 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
510 Teuchos::RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
511 toXpetra(Teuchos::RCP<Tpetra::FECrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> &A) {
512  auto xpTpCrs = Teuchos::rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A));
513  auto xpCrs = Teuchos::rcp_dynamic_cast<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(xpTpCrs, true);
515 }
516 
517 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
518 Teuchos::RCP<const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
519 toXpetra(Teuchos::RCP<const Tpetra::FECrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> &A) {
520  auto xpTpCrs = Teuchos::rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A));
521  auto xpCrs = Teuchos::rcp_dynamic_cast<const Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(xpTpCrs, true);
523 }
524 
525 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
526 Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
528  return toTpetra(toCrsMatrix(A));
529 }
530 
531 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
532 Teuchos::RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
534  return toTpetra(toCrsMatrix(A));
535 }
536 
537 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
538 Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> &
540  return toTpetra(toCrsMatrix(A));
541 }
542 
543 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
544 Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> &
546  return toTpetra(toCrsMatrix(A));
547 }
548 
549 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
550 Teuchos::RCP<Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
552  return Teuchos::rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(A, true)->getTpetra_BlockCrsMatrix();
553 }
554 
555 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
556 Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
558  return Teuchos::rcp_dynamic_cast<Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(A, true)->getTpetra_BlockCrsMatrixNonConst();
559 }
560 
561 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
562 Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> &
564  return *(dynamic_cast<Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> &>(A).getTpetra_BlockCrsMatrixNonConst());
565 }
566 
567 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
568 const Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> &
570  return *(dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> &>(A).getTpetra_BlockCrsMatrix());
571 }
572 
573 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
574 Teuchos::RCP<Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
576  return toTpetraBlock(toCrsMatrix(A));
577 }
578 
579 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
580 Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
582  return toTpetraBlock(toCrsMatrix(A));
583 }
584 
585 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
586 Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> &
588  return toTpetraBlock(toCrsMatrix(A));
589 }
590 
591 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
592 const Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> &
594  return toTpetraBlock(toCrsMatrix(A));
595 }
596 
597 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
598 Teuchos::RCP<Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
600  auto crsMat = toCrsMatrix(mat);
601  auto tmp_Crs = Teuchos::rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(crsMat);
602  if (!tmp_Crs.is_null()) {
603  return tmp_Crs->getTpetra_CrsMatrixNonConst();
604  } else {
605  auto tmp_BlockCrs = Teuchos::rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(crsMat);
606  if (tmp_BlockCrs.is_null())
607  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
608  return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
609  }
610 }
611 
612 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
613 Teuchos::RCP<const Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
615  auto crsMat = toCrsMatrix(mat);
616  auto tmp_Crs = Teuchos::rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(crsMat);
617  if (!tmp_Crs.is_null()) {
618  return tmp_Crs->getTpetra_CrsMatrix();
619  } else {
620  auto tmp_BlockCrs = Teuchos::rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(crsMat);
621  if (tmp_BlockCrs.is_null())
622  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
623  return tmp_BlockCrs->getTpetra_BlockCrsMatrix();
624  }
625 }
626 
627 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
628 Teuchos::RCP<Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
630  auto mat = Teuchos::rcp_dynamic_cast<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(A);
631  auto rmat = Teuchos::rcp_dynamic_cast<Xpetra::TpetraRowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(A);
632  if (!mat.is_null()) {
633  return toTpetraRowMatrix(mat);
634  } else if (!rmat.is_null()) {
635  return rmat->getTpetra_RowMatrixNonConst();
636  } else {
637  auto tpOp = Teuchos::rcp_dynamic_cast<Xpetra::TpetraOperator<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(A, true);
638  auto tOp = tpOp->getOperator();
639  auto tRow = Teuchos::rcp_dynamic_cast<Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tOp, true);
640  return tRow;
641  }
642 }
643 
644 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
645 Teuchos::RCP<const Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
647  auto mat = Teuchos::rcp_dynamic_cast<const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(A);
648  auto rmat = Teuchos::rcp_dynamic_cast<const Xpetra::TpetraRowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(A);
649  if (!mat.is_null()) {
650  return toTpetraRowMatrix(mat);
651  } else if (!rmat.is_null()) {
652  return rmat->getTpetra_RowMatrix();
653  } else {
654  auto tpOp = Teuchos::rcp_dynamic_cast<const Xpetra::TpetraOperator<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(A, true);
655  auto tOp = tpOp->getOperator();
656  auto tRow = Teuchos::rcp_dynamic_cast<const Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tOp, true);
657  return tRow;
658  }
659 }
660 
661 } // namespace Xpetra
662 
663 #define XPETRA_CRSMATRIXWRAP_SHORT
664 #endif // XPETRA_CRSMATRIXWRAP_DECL_HPP
665 
666 // NOTE: if CrsMatrix and VbrMatrix share a common interface for fillComplete() etc, I can move some stuff in Xpetra_Matrix.hpp
667 // TODO: getUnderlyingMatrix() method
void setObjectLabel(const std::string &objectLabel)
RCP< CrsMatrix > getCrsMatrix() const
Teuchos::RCP< Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > toTpetraRowMatrix(const Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >> &mat)
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...
Xpetra::TpetraCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraCrsMatrix
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::HostMirror getLocalMatrixHost() const
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.
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...
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map > &newMap)
Xpetra::MatrixView< Scalar, LocalOrdinal, GlobalOrdinal, Node > MatrixView
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.
Teuchos::RCP< Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > toCrsMatrix(const Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >> &matrix)
void doExport(const Matrix &dest, const Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
virtual RCP< Tpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getOperator()
Gets the operator out.
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.
Exception indicating invalid cast attempted.
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.
Xpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > CrsGraph
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.
RCP< const Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > toTpetra(const RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph)
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 > &params=null)
Signal that data entry is complete, specifying domain and range maps.
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 > &params=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. ...
RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrix() const
Get the underlying Tpetra matrix.
virtual void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalar.
Teuchos::RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > toTpetraBlock(const Teuchos::RCP< Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >> &A)
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrixNonConst() const
Get the underlying Tpetra matrix.
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.
Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > Map
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.
Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrix
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.
Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > 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.
Xpetra::CrsMatrixFactory< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrixFactory
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...
RCP< const CrsGraph< int, GlobalOrdinal, Node > > toXpetra(const Epetra_CrsGraph &g)