43 #ifndef IFPACK2_OVERLAPPINGROWMATRIX_DEF_HPP
44 #define IFPACK2_OVERLAPPINGROWMATRIX_DEF_HPP
48 #include <Ifpack2_Details_OverlappingRowGraph.hpp>
49 #include <Tpetra_CrsMatrix.hpp>
50 #include <Tpetra_Import.hpp>
51 #include "Tpetra_Map.hpp"
52 #include <Teuchos_CommHelpers.hpp>
53 #include <unordered_set>
57 template<
class MatrixType>
60 const int overlapLevel) :
61 A_ (Teuchos::rcp_dynamic_cast<const crs_matrix_type> (A, true)),
62 OverlapLevel_ (overlapLevel)
67 using Teuchos::outArg;
69 using Teuchos::reduceAll;
70 typedef Tpetra::global_size_t GST;
71 typedef Tpetra::CrsGraph<local_ordinal_type,
72 global_ordinal_type, node_type> crs_graph_type;
74 OverlapLevel_ <= 0, std::runtime_error,
75 "Ifpack2::OverlappingRowMatrix: OverlapLevel must be > 0.");
77 (A_.
is_null (), std::runtime_error,
78 "Ifpack2::OverlappingRowMatrix: The input matrix must be a "
79 "Tpetra::CrsMatrix with the same scalar_type, local_ordinal_type, "
80 "global_ordinal_type, and device_type typedefs as MatrixType.");
82 A_->getComm()->getSize() == 1, std::runtime_error,
83 "Ifpack2::OverlappingRowMatrix: Matrix must be "
84 "distributed over more than one MPI process.");
87 RCP<const crs_graph_type> A_crsGraph = A_->getCrsGraph ();
88 const size_t numMyRowsA = A_->getLocalNumRows ();
89 const global_ordinal_type global_invalid =
93 Array<global_ordinal_type> ExtElements;
97 std::unordered_set<global_ordinal_type> ExtElementSet;
99 RCP<crs_graph_type> TmpGraph;
100 RCP<import_type> TmpImporter;
101 RCP<const map_type> RowMap, ColMap;
102 Kokkos::resize(ExtHaloStarts_, OverlapLevel_+1);
103 ExtHaloStarts_h = Kokkos::create_mirror_view(ExtHaloStarts_);
106 for (
int overlap = 0 ; overlap < OverlapLevel_ ; ++overlap) {
107 ExtHaloStarts_h(overlap) = (size_t) ExtElements.size();
111 RowMap = A_->getRowMap ();
112 ColMap = A_->getColMap ();
115 RowMap = TmpGraph->getRowMap ();
116 ColMap = TmpGraph->getColMap ();
119 const size_t size = ColMap->getLocalNumElements () - RowMap->getLocalNumElements ();
120 Array<global_ordinal_type> mylist (size);
124 for (local_ordinal_type i = 0 ; (size_t) i < ColMap->getLocalNumElements() ; ++i) {
125 const global_ordinal_type GID = ColMap->getGlobalElement (i);
126 if (A_->getRowMap ()->getLocalElement (GID) == global_invalid) {
129 if(ExtElementSet.insert(GID).second)
131 ExtElements.push_back (GID);
140 if (overlap + 1 < OverlapLevel_) {
142 TmpMap =
rcp (
new map_type (global_invalid, mylist (0, count),
146 TmpGraph =
rcp (
new crs_graph_type (TmpMap, 0));
147 TmpImporter =
rcp (
new import_type (A_->getRowMap (), TmpMap));
150 TmpGraph->doImport (*A_crsGraph, *TmpImporter, Tpetra::INSERT);
151 TmpGraph->fillComplete (A_->getDomainMap (), TmpMap);
154 ExtHaloStarts_h[OverlapLevel_] = (size_t) ExtElements.size();
155 Kokkos::deep_copy(ExtHaloStarts_,ExtHaloStarts_h);
159 Array<global_ordinal_type> mylist (numMyRowsA + ExtElements.size ());
160 for (local_ordinal_type i = 0; (size_t)i < numMyRowsA; ++i) {
161 mylist[i] = A_->getRowMap ()->getGlobalElement (i);
163 for (local_ordinal_type i = 0; i < ExtElements.size (); ++i) {
164 mylist[i + numMyRowsA] = ExtElements[i];
168 RowMap_ =
rcp (
new map_type (global_invalid, mylist (),
171 Importer_ =
rcp (
new import_type (A_->getRowMap (), RowMap_));
176 ExtMap_ =
rcp (
new map_type (global_invalid, ExtElements (),
179 ExtImporter_ =
rcp (
new import_type (A_->getRowMap (), ExtMap_));
182 ExtMatrix_ =
rcp (
new crs_matrix_type (ExtMap_, ColMap_, 0));
183 ExtMatrix_->doImport (*A_, *ExtImporter_, Tpetra::INSERT);
184 ExtMatrix_->fillComplete (A_->getDomainMap (), RowMap_);
188 const size_t numMyRowsB = ExtMatrix_->getLocalNumRows ();
190 GST NumMyNonzeros_tmp = A_->getLocalNumEntries () + ExtMatrix_->getLocalNumEntries ();
191 GST NumMyRows_tmp = numMyRowsA + numMyRowsB;
193 GST inArray[2], outArray[2];
194 inArray[0] = NumMyNonzeros_tmp;
195 inArray[1] = NumMyRows_tmp;
198 reduceAll<int, GST> (* (A_->getComm ()), REDUCE_SUM, 2, inArray, outArray);
199 NumGlobalNonzeros_ = outArray[0];
200 NumGlobalRows_ = outArray[1];
207 MaxNumEntries_ = A_->getLocalMaxNumRowEntries ();
209 MaxNumEntries_ = ExtMatrix_->getLocalMaxNumRowEntries ();
214 RCP<row_graph_impl_type> graph =
215 rcp (
new row_graph_impl_type (A_->getGraph (),
216 ExtMatrix_->getGraph (),
225 graph_ = Teuchos::rcp_const_cast<
const row_graph_type>
226 (Teuchos::rcp_implicit_cast<row_graph_type> (graph));
228 Kokkos::resize(Indices_,MaxNumEntries_);
229 Kokkos::resize(Values_,MaxNumEntries_);
233 template<
class MatrixType>
237 return A_->getComm ();
243 template<
class MatrixType>
252 template<
class MatrixType>
261 template<
class MatrixType>
277 template<
class MatrixType>
285 template<
class MatrixType>
293 template<
class MatrixType>
296 return NumGlobalRows_;
300 template<
class MatrixType>
303 return NumGlobalRows_;
307 template<
class MatrixType>
310 return A_->getLocalNumRows () + ExtMatrix_->getLocalNumRows ();
314 template<
class MatrixType>
317 return this->getLocalNumRows ();
321 template<
class MatrixType>
322 typename MatrixType::global_ordinal_type
325 return A_->getIndexBase();
329 template<
class MatrixType>
332 return NumGlobalNonzeros_;
336 template<
class MatrixType>
339 return A_->getLocalNumEntries () + ExtMatrix_->getLocalNumEntries ();
343 template<
class MatrixType>
348 const local_ordinal_type localRow = RowMap_->getLocalElement (globalRow);
352 return getNumEntriesInLocalRow (localRow);
357 template<
class MatrixType>
363 const size_t numMyRowsA = A_->getLocalNumRows ();
364 if (as<size_t> (localRow) < numMyRowsA) {
365 return A_->getNumEntriesInLocalRow (localRow);
367 return ExtMatrix_->getNumEntriesInLocalRow (as<local_ordinal_type> (localRow - numMyRowsA));
372 template<
class MatrixType>
375 return std::max<size_t>(A_->getGlobalMaxNumRowEntries(), ExtMatrix_->getGlobalMaxNumRowEntries());
379 template<
class MatrixType>
382 return MaxNumEntries_;
385 template<
class MatrixType>
388 return A_->getBlockSize();
391 template<
class MatrixType>
398 template<
class MatrixType>
405 template<
class MatrixType>
412 template<
class MatrixType>
419 template<
class MatrixType>
423 nonconst_global_inds_host_view_type &Indices,
424 nonconst_values_host_view_type &Values,
425 size_t& NumEntries)
const
427 throw std::runtime_error(
"Ifpack2::OverlappingRowMatrix::getGlobalRowCopy() not supported.");
430 template<
class MatrixType>
434 nonconst_local_inds_host_view_type &Indices,
435 nonconst_values_host_view_type &Values,
436 size_t& NumEntries)
const
439 const size_t numMyRowsA = A_->getLocalNumRows ();
440 if (as<size_t> (LocalRow) < numMyRowsA) {
441 A_->getLocalRowCopy (LocalRow, Indices, Values, NumEntries);
443 ExtMatrix_->getLocalRowCopy (LocalRow - as<local_ordinal_type> (numMyRowsA),
444 Indices, Values, NumEntries);
448 template<
class MatrixType>
452 global_inds_host_view_type &indices,
453 values_host_view_type &values)
const {
454 const local_ordinal_type LocalRow = RowMap_->getLocalElement (GlobalRow);
456 indices = global_inds_host_view_type();
457 values = values_host_view_type();
459 if (Teuchos::as<size_t> (LocalRow) < A_->getLocalNumRows ()) {
460 A_->getGlobalRowView (GlobalRow, indices, values);
462 ExtMatrix_->getGlobalRowView (GlobalRow, indices, values);
467 template<
class MatrixType>
471 local_inds_host_view_type & indices,
472 values_host_view_type & values)
const {
474 const size_t numMyRowsA = A_->getLocalNumRows ();
475 if (as<size_t> (LocalRow) < numMyRowsA) {
476 A_->getLocalRowView (LocalRow, indices, values);
478 ExtMatrix_->getLocalRowView (LocalRow - as<local_ordinal_type> (numMyRowsA),
484 template<
class MatrixType>
487 getLocalDiagCopy (Tpetra::Vector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& diag)
const
492 vector_type baseDiag(A_->getRowMap());
493 A_->getLocalDiagCopy(baseDiag);
494 Array<scalar_type> baseDiagVals(baseDiag.getLocalLength());
495 baseDiag.get1dCopy(baseDiagVals());
497 vector_type extDiag(ExtMatrix_->getRowMap());
498 ExtMatrix_->getLocalDiagCopy(extDiag);
499 Array<scalar_type> extDiagVals(extDiag.getLocalLength());
500 extDiag.get1dCopy(extDiagVals());
503 if (allDiagVals.
size() != baseDiagVals.size() + extDiagVals.size()) {
504 std::ostringstream errStr;
505 errStr <<
"Ifpack2::OverlappingRowMatrix::getLocalDiagCopy : Mismatch in diagonal lengths, "
506 << allDiagVals.
size() <<
" != " << baseDiagVals.size() <<
"+" << extDiagVals.size();
507 throw std::runtime_error(errStr.str());
509 for (Teuchos::Ordinal i=0; i<baseDiagVals.size(); ++i)
510 allDiagVals[i] = baseDiagVals[i];
511 Teuchos_Ordinal offset=baseDiagVals.
size();
512 for (Teuchos::Ordinal i=0; i<extDiagVals.size(); ++i)
513 allDiagVals[i+offset] = extDiagVals[i];
517 template<
class MatrixType>
520 leftScale (
const Tpetra::Vector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& )
522 throw std::runtime_error(
"Ifpack2::OverlappingRowMatrix does not support leftScale.");
526 template<
class MatrixType>
529 rightScale (
const Tpetra::Vector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& )
531 throw std::runtime_error(
"Ifpack2::OverlappingRowMatrix does not support leftScale.");
535 template<
class MatrixType>
536 typename OverlappingRowMatrix<MatrixType>::mag_type
539 throw std::runtime_error(
"Ifpack2::OverlappingRowMatrix does not support getFrobeniusNorm.");
543 template<
class MatrixType>
546 apply (
const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &X,
547 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &Y,
550 scalar_type beta)
const
552 using MV = Tpetra::MultiVector<scalar_type, local_ordinal_type,
553 global_ordinal_type, node_type>;
555 (X.getNumVectors() != Y.getNumVectors(), std::runtime_error,
556 "Ifpack2::OverlappingRowMatrix::apply: X.getNumVectors() = "
557 << X.getNumVectors() <<
" != Y.getNumVectors() = " << Y.getNumVectors()
560 bool aliases = X.aliases(Y);
563 this->apply (X_copy, Y, mode, alpha, beta);
567 const auto& rowMap0 = * (A_->getRowMap ());
568 const auto& colMap0 = * (A_->getColMap ());
571 A_->localApply (X_0, Y_0, mode, alpha, beta);
573 const auto& rowMap1 = * (ExtMatrix_->getRowMap ());
574 const auto& colMap1 = * (ExtMatrix_->getColMap ());
576 MV Y_1 (Y, mode ==
Teuchos::NO_TRANS ? rowMap1 : colMap1, A_->getLocalNumRows ());
577 ExtMatrix_->localApply (X_1, Y_1, mode, alpha, beta);
581 template<
class MatrixType>
584 importMultiVector (
const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &X,
585 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &OvX,
586 Tpetra::CombineMode CM)
588 OvX.doImport (X, *Importer_, CM);
592 template<
class MatrixType>
594 OverlappingRowMatrix<MatrixType>::
595 exportMultiVector (
const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &OvX,
596 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &X,
597 Tpetra::CombineMode CM)
599 X.doExport (OvX, *Importer_, CM);
603 template<
class MatrixType>
610 template<
class MatrixType>
616 template<
class MatrixType>
619 std::ostringstream oss;
620 if (isFillComplete()) {
621 oss <<
"{ isFillComplete: true"
622 <<
", global rows: " << getGlobalNumRows()
623 <<
", global columns: " << getGlobalNumCols()
624 <<
", global entries: " << getGlobalNumEntries()
628 oss <<
"{ isFillComplete: false"
629 <<
", global rows: " << getGlobalNumRows()
635 template<
class MatrixType>
656 RCP<const Teuchos::Comm<int> > comm = this->getComm();
657 const int myRank = comm->getRank();
658 const int numProcs = comm->getSize();
660 for (
size_t dec=10; dec<getGlobalNumRows(); dec *= 10) {
663 width = std::max<size_t> (width, as<size_t> (11)) + 2;
674 out << this->description() << std::endl;
683 out << endl <<
"Row map:" << endl;
685 getRowMap()->describe(out,vl);
687 if (getColMap() != null) {
688 if (getColMap() == getRowMap()) {
690 out << endl <<
"Column map is row map.";
695 out << endl <<
"Column map:" << endl;
697 getColMap()->describe(out,vl);
700 if (getDomainMap() != null) {
701 if (getDomainMap() == getRowMap()) {
703 out << endl <<
"Domain map is row map.";
706 else if (getDomainMap() == getColMap()) {
708 out << endl <<
"Domain map is column map.";
713 out << endl <<
"Domain map:" << endl;
715 getDomainMap()->describe(out,vl);
718 if (getRangeMap() != null) {
719 if (getRangeMap() == getDomainMap()) {
721 out << endl <<
"Range map is domain map." << endl;
724 else if (getRangeMap() == getRowMap()) {
726 out << endl <<
"Range map is row map." << endl;
731 out << endl <<
"Range map: " << endl;
733 getRangeMap()->describe(out,vl);
742 for (
int curRank = 0; curRank < numProcs; ++curRank) {
743 if (myRank == curRank) {
744 out <<
"Process rank: " << curRank << std::endl;
745 out <<
" Number of entries: " << getLocalNumEntries() << std::endl;
746 out <<
" Max number of entries per row: " << getLocalMaxNumRowEntries() << std::endl;
755 for (
int curRank = 0; curRank < numProcs; ++curRank) {
756 if (myRank == curRank) {
757 out << std::setw(width) <<
"Proc Rank"
758 << std::setw(width) <<
"Global Row"
759 << std::setw(width) <<
"Num Entries";
761 out << std::setw(width) <<
"(Index,Value)";
764 for (
size_t r = 0; r < getLocalNumRows (); ++r) {
765 const size_t nE = getNumEntriesInLocalRow(r);
766 typename MatrixType::global_ordinal_type gid = getRowMap()->getGlobalElement(r);
767 out << std::setw(width) << myRank
768 << std::setw(width) << gid
769 << std::setw(width) << nE;
771 if (isGloballyIndexed()) {
772 global_inds_host_view_type rowinds;
773 values_host_view_type rowvals;
774 getGlobalRowView (gid, rowinds, rowvals);
775 for (
size_t j = 0; j < nE; ++j) {
776 out <<
" (" << rowinds[j]
777 <<
", " << rowvals[j]
781 else if (isLocallyIndexed()) {
782 local_inds_host_view_type rowinds;
783 values_host_view_type rowvals;
784 getLocalRowView (r, rowinds, rowvals);
785 for (
size_t j=0; j < nE; ++j) {
786 out <<
" (" << getColMap()->getGlobalElement(rowinds[j])
787 <<
", " << rowvals[j]
802 out <<
"===========\nlocal matrix\n=================" << std::endl;
806 out <<
"===========\nend of local matrix\n=================" << std::endl;
809 out <<
"=================\nghost matrix\n=================" << std::endl;
813 out <<
"===========\nend of ghost matrix\n=================" << std::endl;
818 template<
class MatrixType>
820 OverlappingRowMatrix<MatrixType>::getUnderlyingMatrix()
const
825 template<
class MatrixType>
827 OverlappingRowMatrix<MatrixType>::getExtMatrix()
const
832 template<
class MatrixType>
833 Kokkos::View<size_t*, typename OverlappingRowMatrix<MatrixType>::device_type>
834 OverlappingRowMatrix<MatrixType>::getExtHaloStarts()
const
836 return ExtHaloStarts_;
839 template<
class MatrixType>
840 typename Kokkos::View<size_t*, typename OverlappingRowMatrix<MatrixType>::device_type>::HostMirror
841 OverlappingRowMatrix<MatrixType>::getExtHaloStartsHost()
const
843 return ExtHaloStarts_h;
846 template<
class MatrixType>
847 void OverlappingRowMatrix<MatrixType>::doExtImport()
852 ExtMatrix_ =
rcp (
new crs_matrix_type (ExtMap_, ColMap_, 0));
853 ExtMatrix_->doImport (*A_, *ExtImporter_, Tpetra::INSERT);
854 ExtMatrix_->fillComplete (A_->getDomainMap (), RowMap_);
859 #define IFPACK2_OVERLAPPINGROWMATRIX_INSTANT(S,LO,GO,N) \
860 template class Ifpack2::OverlappingRowMatrix< Tpetra::RowMatrix<S, LO, GO, N> >;
862 #endif // IFPACK2_OVERLAPPINGROWMATRIX_DEF_HPP
virtual bool isGloballyIndexed() const
Whether this matrix is globally indexed.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:406
virtual mag_type getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:537
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual bool isFillComplete() const
true if fillComplete() has been called, else false.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:413
virtual void getLocalRowCopy(local_ordinal_type LocalRow, nonconst_local_inds_host_view_type &Indices, nonconst_values_host_view_type &Values, size_t &NumEntries) const
Extract a list of entries in a specified local row of the graph. Put into storage allocated by callin...
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:433
virtual void leftScale(const Tpetra::Vector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &x)
Scales the RowMatrix on the left with the Vector x.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:520
virtual size_t getLocalMaxNumRowEntries() const
The maximum number of entries in any row on the calling process.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:380
virtual size_t getLocalNumRows() const
The number of rows owned by the calling process.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:308
virtual Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
The communicator over which the matrix is distributed.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:235
virtual Teuchos::RCP< const Tpetra::Map< local_ordinal_type, global_ordinal_type, node_type > > getRowMap() const
The Map that describes the distribution of rows over processes.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:245
virtual global_size_t getGlobalNumRows() const
The global number of rows in this matrix.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:294
virtual global_size_t getGlobalNumCols() const
The global number of columns in this matrix.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:301
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
virtual bool isLocallyIndexed() const
Whether this matrix is locally indexed.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:399
virtual void getGlobalRowCopy(global_ordinal_type GlobalRow, nonconst_global_inds_host_view_type &Indices, nonconst_values_host_view_type &Values, size_t &NumEntries) const
Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage...
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:422
virtual local_ordinal_type getBlockSize() const
The number of degrees of freedom per mesh point.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:386
virtual Teuchos::RCP< const Tpetra::Map< local_ordinal_type, global_ordinal_type, node_type > > getRangeMap() const
The Map that describes the range of this matrix.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:279
virtual size_t getLocalNumCols() const
The number of columns owned by the calling process.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:315
virtual Teuchos::RCP< const Tpetra::Map< local_ordinal_type, global_ordinal_type, node_type > > getDomainMap() const
The Map that describes the domain of this matrix.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:263
Sparse graph (Tpetra::RowGraph subclass) with ghost rows.
Definition: Ifpack2_Details_OverlappingRowGraph_decl.hpp:65
virtual void getLocalRowView(local_ordinal_type LocalRow, local_inds_host_view_type &indices, values_host_view_type &values) const
Extract a const, non-persisting view of local indices in a specified row of the matrix.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:470
basic_FancyOStream & setOutputToRootOnly(const int rootRank)
virtual size_t getNumEntriesInGlobalRow(global_ordinal_type globalRow) const
The number of entries in the given global row that are owned by the calling process.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:346
virtual global_size_t getGlobalNumEntries() const
The global number of entries in this matrix.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:330
virtual bool hasTransposeApply() const
Whether this operator's apply() method can apply the adjoint (transpose).
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:604
virtual global_ordinal_type getIndexBase() const
The index base for global indices for this matrix.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:323
virtual bool hasColMap() const
Whether this matrix has a column Map.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:392
virtual void getGlobalRowView(global_ordinal_type GlobalRow, global_inds_host_view_type &indices, values_host_view_type &values) const
Extract a const, non-persisting view of global indices in a specified row of the matrix.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:451
virtual size_t getNumEntriesInLocalRow(local_ordinal_type localRow) const
The number of entries in the given local row that are owned by the calling process.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:360
virtual void rightScale(const Tpetra::Vector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &x)
Scales the RowMatrix on the right with the Vector x.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:529
virtual Teuchos::RCP< const Tpetra::RowGraph< local_ordinal_type, global_ordinal_type, node_type > > getGraph() const
This matrix's graph.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:287
TypeTo as(const TypeFrom &t)
virtual void getLocalDiagCopy(Tpetra::Vector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &diag) const
Get a copy of the diagonal entries owned by this node, with local row indices.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:487
Sparse matrix (Tpetra::RowMatrix subclass) with ghost rows.
Definition: Ifpack2_OverlappingRowMatrix_decl.hpp:58
OverlappingRowMatrix(const Teuchos::RCP< const row_matrix_type > &A, const int overlapLevel)
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:59
virtual void apply(const Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &X, Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, scalar_type alpha=Teuchos::ScalarTraits< scalar_type >::one(), scalar_type beta=Teuchos::ScalarTraits< scalar_type >::zero()) const
Computes the operator-multivector application.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:546
virtual size_t getLocalNumEntries() const
The number of entries in this matrix owned by the calling process.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:337
virtual Teuchos::RCP< const Tpetra::Map< local_ordinal_type, global_ordinal_type, node_type > > getColMap() const
The Map that describes the distribution of columns over processes.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:254
virtual bool supportsRowViews() const
true if row views are supported, else false.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:611
virtual size_t getGlobalMaxNumRowEntries() const
The maximum number of entries in any row on any process.
Definition: Ifpack2_OverlappingRowMatrix_def.hpp:373