43 #ifndef IFPACK2_LOCALFILTER_DEF_HPP
44 #define IFPACK2_LOCALFILTER_DEF_HPP
46 #include <Ifpack2_LocalFilter_decl.hpp>
47 #include <Tpetra_Map.hpp>
48 #include <Tpetra_MultiVector.hpp>
49 #include <Tpetra_Vector.hpp>
54 # include "Teuchos_DefaultSerialComm.hpp"
60 template<
class MatrixType>
62 LocalFilter<MatrixType>::
63 mapPairsAreFitted (
const row_matrix_type& A)
65 const auto rangeMap = A.getRangeMap();
66 const auto rowMap = A.getRowMap();
67 const bool rangeAndRowFitted = mapPairIsFitted (*rowMap, *rangeMap);
69 const auto domainMap = A.getDomainMap();
70 const auto columnMap = A.getColMap();
71 const bool domainAndColumnFitted = mapPairIsFitted (*columnMap, *domainMap);
78 int localSuccess = rangeAndRowFitted && domainAndColumnFitted;
81 Teuchos::reduceAll<int, int> (*(A.getComm()),
Teuchos::REDUCE_MIN, localSuccess, Teuchos::outArg (globalSuccess));
83 return globalSuccess == 1;
87 template<
class MatrixType>
89 LocalFilter<MatrixType>::
90 mapPairIsFitted (
const map_type& map1,
const map_type& map2)
92 return map1.isLocallyFitted (map2);
96 template<
class MatrixType>
107 #ifdef HAVE_IFPACK2_DEBUG
108 if(! mapPairsAreFitted (*A))
110 std::cout <<
"WARNING: Ifpack2::LocalFilter:\n" <<
111 "A's Map pairs are not fitted to each other on Process "
112 << A_->getRowMap ()->getComm ()->getRank () <<
" of the input matrix's "
114 "This means that LocalFilter may not work with A. "
115 "Please see discussion of Bug 5992.";
117 #endif // HAVE_IFPACK2_DEBUG
120 RCP<const Teuchos::Comm<int> > localComm;
122 localComm =
rcp (
new Teuchos::MpiComm<int> (MPI_COMM_SELF));
151 const size_t numRows = A_->getRangeMap()->getLocalNumElements () / blockSize;
153 const global_ordinal_type indexBase =
static_cast<global_ordinal_type
> (0);
157 Tpetra::GloballyDistributed));
160 localRangeMap_ = localRowMap_;
167 localDomainMap_ = localRangeMap_;
170 const size_t numCols = A_->getDomainMap()->getLocalNumElements () / blockSize;
173 Tpetra::GloballyDistributed));
179 NumEntries_.resize (numRows);
183 MaxNumEntries_ = A_->getLocalMaxNumRowEntries ();
184 MaxNumEntriesA_ = A_->getLocalMaxNumRowEntries ();
187 Kokkos::resize(localIndices_,MaxNumEntries_);
188 Kokkos::resize(localIndicesForGlobalCopy_,MaxNumEntries_);
189 Kokkos::resize(Values_,MaxNumEntries_*blockSize*blockSize);
198 size_t ActualMaxNumEntries = 0;
200 for (
size_t i = 0; i < numRows; ++i) {
202 size_t Nnz, NewNnz = 0;
203 A_->getLocalRowCopy (i, localIndices_, Values_, Nnz);
204 for (
size_t j = 0; j < Nnz; ++j) {
216 if (static_cast<size_t> (localIndices_[j]) < numRows) {
221 if (NewNnz > ActualMaxNumEntries) {
222 ActualMaxNumEntries = NewNnz;
225 NumNonzeros_ += NewNnz;
226 NumEntries_[i] = NewNnz;
229 MaxNumEntries_ = ActualMaxNumEntries;
233 template<
class MatrixType>
238 template<
class MatrixType>
242 return localRowMap_->getComm ();
248 template<
class MatrixType>
249 Teuchos::RCP<
const Tpetra::Map<
typename MatrixType::local_ordinal_type,
250 typename MatrixType::global_ordinal_type,
251 typename MatrixType::node_type> >
258 template<
class MatrixType>
259 Teuchos::RCP<
const Tpetra::Map<
typename MatrixType::local_ordinal_type,
260 typename MatrixType::global_ordinal_type,
261 typename MatrixType::node_type> >
268 template<
class MatrixType>
269 Teuchos::RCP<
const Tpetra::Map<
typename MatrixType::local_ordinal_type,
270 typename MatrixType::global_ordinal_type,
271 typename MatrixType::node_type> >
274 return localDomainMap_;
278 template<
class MatrixType>
279 Teuchos::RCP<
const Tpetra::Map<
typename MatrixType::local_ordinal_type,
280 typename MatrixType::global_ordinal_type,
281 typename MatrixType::node_type> >
284 return localRangeMap_;
288 template<
class MatrixType>
289 Teuchos::RCP<
const Tpetra::RowGraph<
typename MatrixType::local_ordinal_type,
290 typename MatrixType::global_ordinal_type,
291 typename MatrixType::node_type> >
294 if (local_graph_ == Teuchos::null) {
295 local_ordinal_type numRows = this->getLocalNumRows();
297 for(local_ordinal_type i = 0; i < numRows; i++) {
298 entriesPerRow[i] = this->getNumEntriesInLocalRow(i);
305 nonconst_local_inds_host_view_type indices(
"indices",this->getLocalMaxNumRowEntries());
306 nonconst_values_host_view_type values(
"values",this->getLocalMaxNumRowEntries());
307 for(local_ordinal_type i = 0; i < numRows; i++) {
308 size_t numEntries = 0;
309 this->getLocalRowCopy(i, indices, values, numEntries);
310 local_graph_nc->insertLocalIndices (i, numEntries, indices.data());
312 local_graph_nc->fillComplete (this->getDomainMap (), this->getRangeMap ());
313 local_graph_ = Teuchos::rcp_const_cast<
const crs_graph_type> (local_graph_nc);
319 template<
class MatrixType>
322 return static_cast<global_size_t
> (localRangeMap_->getLocalNumElements ());
326 template<
class MatrixType>
329 return static_cast<global_size_t
> (localDomainMap_->getLocalNumElements ());
333 template<
class MatrixType>
336 return static_cast<size_t> (localRangeMap_->getLocalNumElements ());
340 template<
class MatrixType>
343 return static_cast<size_t> (localDomainMap_->getLocalNumElements ());
347 template<
class MatrixType>
348 typename MatrixType::global_ordinal_type
351 return A_->getIndexBase ();
355 template<
class MatrixType>
362 template<
class MatrixType>
368 template<
class MatrixType>
371 return A_->getBlockSize();
374 template<
class MatrixType>
379 const local_ordinal_type localRow = getRowMap ()->getLocalElement (globalRow);
388 return NumEntries_[localRow];
393 template<
class MatrixType>
403 if (getRowMap ()->isNodeLocalElement (localRow)) {
404 return NumEntries_[localRow];
416 template<
class MatrixType>
419 return MaxNumEntries_;
423 template<
class MatrixType>
426 return MaxNumEntries_;
430 template<
class MatrixType>
437 template<
class MatrixType>
440 return A_->isLocallyIndexed ();
444 template<
class MatrixType>
447 return A_->isGloballyIndexed();
451 template<
class MatrixType>
454 return A_->isFillComplete ();
458 template<
class MatrixType>
462 nonconst_global_inds_host_view_type &globalIndices,
463 nonconst_values_host_view_type &values,
464 size_t& numEntries)
const
466 typedef local_ordinal_type LO;
469 const LO localRow = this->getRowMap ()->getLocalElement (globalRow);
482 numEntries = this->getNumEntriesInLocalRow (localRow);
488 this->getLocalRowCopy (localRow, localIndicesForGlobalCopy_, values, numEntries);
490 const map_type& colMap = * (this->getColMap ());
493 const size_type numEnt =
494 std::min (static_cast<size_type> (numEntries),
495 std::min ((size_type)globalIndices.size (), (size_type)values.size ()));
496 for (size_type k = 0; k < numEnt; ++k) {
497 globalIndices[k] = colMap.getGlobalElement (localIndicesForGlobalCopy_[k]);
503 template<
class MatrixType>
507 nonconst_local_inds_host_view_type &Indices,
508 nonconst_values_host_view_type &Values,
509 size_t& NumEntries)
const
511 typedef local_ordinal_type LO;
512 typedef global_ordinal_type GO;
514 if (! A_->getRowMap ()->isNodeLocalElement (LocalRow)) {
520 if (A_->getRowMap()->getComm()->getSize() == 1) {
521 A_->getLocalRowCopy (LocalRow, Indices, Values, NumEntries);
525 const LO blockNumEntr = getBlockSize() * getBlockSize();
527 const size_t numEntInLclRow = NumEntries_[LocalRow];
528 if (static_cast<size_t> (Indices.size ()) < numEntInLclRow ||
529 static_cast<size_t> (Values.size ()) < numEntInLclRow*blockNumEntr) {
534 true, std::runtime_error,
535 "Ifpack2::LocalFilter::getLocalRowCopy: Invalid output array length. "
536 "The output arrays must each have length at least " << numEntInLclRow
537 <<
" for local row " << LocalRow <<
" on Process "
538 << localRowMap_->getComm ()->getRank () <<
".");
540 else if (numEntInLclRow == static_cast<size_t> (0)) {
559 size_t numEntInMat = 0;
560 A_->getLocalRowCopy (LocalRow, localIndices_, Values_ , numEntInMat);
565 const map_type& matrixDomMap = * (A_->getDomainMap ());
566 const map_type& matrixColMap = * (A_->getColMap ());
568 const size_t capacity =
static_cast<size_t> (std::min (Indices.size (),
569 Values.size ()/blockNumEntr));
571 const size_t numRows = localRowMap_->getLocalNumElements ();
572 const bool buggy =
true;
573 for (
size_t j = 0; j < numEntInMat; ++j) {
579 const LO matrixLclCol = localIndices_[j];
580 const GO gblCol = matrixColMap.getGlobalElement (matrixLclCol);
589 if ((
size_t) localIndices_[j] < numRows) {
590 Indices[NumEntries] = localIndices_[j];
591 for (LO k=0;k<blockNumEntr;++k)
592 Values[NumEntries*blockNumEntr + k] = Values_[j*blockNumEntr + k];
596 if (matrixDomMap.isNodeGlobalElement (gblCol)) {
600 if (NumEntries < capacity) {
601 Indices[NumEntries] = matrixLclCol;
602 for (LO k=0;k<blockNumEntr;++k)
603 Values[NumEntries*blockNumEntr + k] = Values_[j*blockNumEntr + k];
612 template<
class MatrixType>
616 global_inds_host_view_type &,
617 values_host_view_type &)
const
620 "Ifpack2::LocalFilter does not implement getGlobalRowView.");
624 template<
class MatrixType>
628 local_inds_host_view_type &,
629 values_host_view_type &)
const
632 "Ifpack2::LocalFilter does not implement getLocalRowView.");
636 template<
class MatrixType>
639 getLocalDiagCopy (Tpetra::Vector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& diag)
const
642 typedef Tpetra::Vector<scalar_type, local_ordinal_type,
643 global_ordinal_type, node_type> vector_type;
646 RCP<vector_type> diagView = diag.offsetViewNonConst (A_->getRowMap (), 0);
647 A_->getLocalDiagCopy (*diagView);
651 template<
class MatrixType>
654 leftScale (
const Tpetra::Vector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& )
657 "Ifpack2::LocalFilter does not implement leftScale.");
661 template<
class MatrixType>
664 rightScale (
const Tpetra::Vector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& )
667 "Ifpack2::LocalFilter does not implement rightScale.");
671 template<
class MatrixType>
674 apply (
const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &X,
675 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &Y,
678 scalar_type beta)
const
680 typedef Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> MV;
682 X.getNumVectors() != Y.getNumVectors(), std::runtime_error,
683 "Ifpack2::LocalFilter::apply: X and Y must have the same number of columns. "
684 "X has " << X.getNumVectors () <<
" columns, but Y has "
685 << Y.getNumVectors () <<
" columns.");
687 #ifdef HAVE_IFPACK2_DEBUG
693 for (
size_t j = 0; j < X.getNumVectors (); ++j) {
694 if (STM::isnaninf (norms[j])) {
699 TEUCHOS_TEST_FOR_EXCEPTION( ! good, std::runtime_error,
"Ifpack2::LocalFilter::apply: The 1-norm of the input X is NaN or Inf.");
701 #endif // HAVE_IFPACK2_DEBUG
704 getBlockSize() > 1, std::runtime_error,
705 "Ifpack2::LocalFilter::apply: Block size greater than zero is not yet supported for "
706 "LocalFilter::apply. Please contact an Ifpack2 developer to request this feature.");
715 applyNonAliased (X_copy, Y, mode, alpha, beta);
717 applyNonAliased (X, Y, mode, alpha, beta);
720 #ifdef HAVE_IFPACK2_DEBUG
726 for (
size_t j = 0; j < Y.getNumVectors (); ++j) {
727 if (STM::isnaninf (norms[j])) {
732 TEUCHOS_TEST_FOR_EXCEPTION( ! good, std::runtime_error,
"Ifpack2::LocalFilter::apply: The 1-norm of the output Y is NaN or Inf.");
734 #endif // HAVE_IFPACK2_DEBUG
737 template<
class MatrixType>
740 applyNonAliased (
const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &X,
741 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &Y,
744 scalar_type beta)
const
750 const scalar_type zero = STS::zero ();
751 const scalar_type one = STS::one ();
756 else if (beta != one) {
760 const size_t NumVectors = Y.getNumVectors ();
761 const size_t numRows = localRowMap_->getLocalNumElements ();
768 const bool constantStride = X.isConstantStride () && Y.isConstantStride ();
769 if (constantStride) {
772 const size_t x_stride = X.getStride();
773 const size_t y_stride = Y.getStride();
774 ArrayRCP<scalar_type> y_rcp = Y.get1dViewNonConst();
775 ArrayRCP<const scalar_type> x_rcp = X.get1dView();
776 ArrayView<scalar_type> y_ptr = y_rcp();
777 ArrayView<const scalar_type> x_ptr = x_rcp();
778 for (
size_t i = 0; i < numRows; ++i) {
781 getLocalRowCopy (i, localIndices_ , Values_ , Nnz);
782 scalar_type* Values =
reinterpret_cast<scalar_type*
>(Values_.data());
784 for (
size_t j = 0; j < Nnz; ++j) {
785 const local_ordinal_type col = localIndices_[j];
786 for (
size_t k = 0; k < NumVectors; ++k) {
787 y_ptr[i + y_stride*k] +=
788 alpha * Values[j] * x_ptr[col + x_stride*k];
793 for (
size_t j = 0; j < Nnz; ++j) {
794 const local_ordinal_type col = localIndices_[j];
795 for (
size_t k = 0; k < NumVectors; ++k) {
796 y_ptr[col + y_stride*k] +=
797 alpha * Values[j] * x_ptr[i + x_stride*k];
802 for (
size_t j = 0; j < Nnz; ++j) {
803 const local_ordinal_type col = localIndices_[j];
804 for (
size_t k = 0; k < NumVectors; ++k) {
805 y_ptr[col + y_stride*k] +=
806 alpha * STS::conjugate (Values[j]) * x_ptr[i + x_stride*k];
815 ArrayRCP<ArrayRCP<const scalar_type> > x_ptr = X.get2dView();
816 ArrayRCP<ArrayRCP<scalar_type> > y_ptr = Y.get2dViewNonConst();
818 for (
size_t i = 0; i < numRows; ++i) {
821 getLocalRowCopy (i, localIndices_ , Values_ , Nnz);
822 scalar_type* Values =
reinterpret_cast<scalar_type*
>(Values_.data());
824 for (
size_t k = 0; k < NumVectors; ++k) {
825 ArrayView<const scalar_type> x_local = (x_ptr())[k]();
826 ArrayView<scalar_type> y_local = (y_ptr())[k]();
827 for (
size_t j = 0; j < Nnz; ++j) {
828 y_local[i] += alpha * Values[j] * x_local[localIndices_[j]];
833 for (
size_t k = 0; k < NumVectors; ++k) {
834 ArrayView<const scalar_type> x_local = (x_ptr())[k]();
835 ArrayView<scalar_type> y_local = (y_ptr())[k]();
836 for (
size_t j = 0; j < Nnz; ++j) {
837 y_local[localIndices_[j]] += alpha * Values[j] * x_local[i];
842 for (
size_t k = 0; k < NumVectors; ++k) {
843 ArrayView<const scalar_type> x_local = (x_ptr())[k]();
844 ArrayView<scalar_type> y_local = (y_ptr())[k]();
845 for (
size_t j = 0; j < Nnz; ++j) {
846 y_local[localIndices_[j]] +=
847 alpha * STS::conjugate (Values[j]) * x_local[i];
855 template<
class MatrixType>
862 template<
class MatrixType>
869 template<
class MatrixType>
871 LocalFilter<MatrixType>::mag_type
874 typedef Kokkos::ArithTraits<scalar_type> STS;
875 typedef Kokkos::ArithTraits<mag_type> STM;
878 const size_type maxNumRowEnt = getLocalMaxNumRowEntries ();
879 const local_ordinal_type blockNumEntr = getBlockSize() * getBlockSize();
880 nonconst_local_inds_host_view_type ind (
"ind",maxNumRowEnt);
881 nonconst_values_host_view_type val (
"val",maxNumRowEnt*blockNumEntr);
882 const size_t numRows =
static_cast<size_t> (localRowMap_->getLocalNumElements ());
885 mag_type sumSquared = STM::zero ();
886 for (
size_t i = 0; i < numRows; ++i) {
887 size_t numEntries = 0;
888 this->getLocalRowCopy (i, ind, val, numEntries);
889 for (size_type k = 0; k < static_cast<size_type> (numEntries*blockNumEntr); ++k) {
890 const mag_type v_k_abs = STS::magnitude (val[k]);
891 sumSquared += v_k_abs * v_k_abs;
894 return STM::squareroot (sumSquared);
897 template<
class MatrixType>
902 std::ostringstream os;
904 os <<
"Ifpack2::LocalFilter: {";
905 os <<
"MatrixType: " << TypeNameTraits<MatrixType>::name ();
906 if (this->getObjectLabel () !=
"") {
907 os <<
", Label: \"" << this->getObjectLabel () <<
"\"";
909 os <<
", Number of rows: " << getGlobalNumRows ()
910 <<
", Number of columns: " << getGlobalNumCols ()
916 template<
class MatrixType>
933 out <<
"Ifpack2::LocalFilter:" << endl;
935 out <<
"MatrixType: " << TypeNameTraits<MatrixType>::name () << endl;
936 if (this->getObjectLabel () !=
"") {
937 out <<
"Label: \"" << this->getObjectLabel () <<
"\"" << endl;
939 out <<
"Number of rows: " << getGlobalNumRows () << endl
940 <<
"Number of columns: " << getGlobalNumCols () << endl
941 <<
"Number of nonzeros: " << NumNonzeros_ << endl;
944 out <<
"Row Map:" << endl;
945 localRowMap_->describe (out, vl);
946 out <<
"Domain Map:" << endl;
947 localDomainMap_->describe (out, vl);
948 out <<
"Range Map:" << endl;
949 localRangeMap_->describe (out, vl);
954 template<
class MatrixType>
955 Teuchos::RCP<
const Tpetra::RowMatrix<
typename MatrixType::scalar_type,
956 typename MatrixType::local_ordinal_type,
957 typename MatrixType::global_ordinal_type,
958 typename MatrixType::node_type> >
967 #define IFPACK2_LOCALFILTER_INSTANT(S,LO,GO,N) \
968 template class Ifpack2::LocalFilter< Tpetra::RowMatrix<S, LO, GO, N> >;
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_LocalFilter_def.hpp:664
virtual Teuchos::RCP< const map_type > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
Definition: Ifpack2_LocalFilter_def.hpp:252
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object to the given output stream.
Definition: Ifpack2_LocalFilter_def.hpp:919
basic_OSTab< char > OSTab
virtual Teuchos::RCP< const map_type > getDomainMap() const
Returns the Map that describes the domain distribution in this matrix.
Definition: Ifpack2_LocalFilter_def.hpp:272
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual size_t getLocalMaxNumRowEntries() const
The maximum number of entries across all rows/columns on this process.
Definition: Ifpack2_LocalFilter_def.hpp:424
virtual void getGlobalRowCopy(global_ordinal_type GlobalRow, nonconst_global_inds_host_view_type &Indices, nonconst_values_host_view_type &Values, size_t &NumEntries) const
Get the entries in the given row, using global indices.
Definition: Ifpack2_LocalFilter_def.hpp:461
virtual bool isFillComplete() const
Returns true if fillComplete() has been called.
Definition: Ifpack2_LocalFilter_def.hpp:452
virtual void getLocalRowCopy(local_ordinal_type LocalRow, nonconst_local_inds_host_view_type &Indices, nonconst_values_host_view_type &Values, size_t &NumEntries) const
Get the entries in the given row, using local indices.
Definition: Ifpack2_LocalFilter_def.hpp:506
virtual Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
Returns the communicator.
Definition: Ifpack2_LocalFilter_def.hpp:240
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_LocalFilter_def.hpp:615
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
Compute Y = beta*Y + alpha*A_local*X.
Definition: Ifpack2_LocalFilter_def.hpp:674
virtual bool hasTransposeApply() const
Whether this operator supports applying the transpose or conjugate transpose.
Definition: Ifpack2_LocalFilter_def.hpp:856
virtual bool isLocallyIndexed() const
Whether the underlying sparse matrix is locally (opposite of globally) indexed.
Definition: Ifpack2_LocalFilter_def.hpp:438
virtual void getLocalDiagCopy(Tpetra::Vector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &diag) const
Get the diagonal entries of the (locally filtered) matrix.
Definition: Ifpack2_LocalFilter_def.hpp:639
virtual size_t getLocalNumEntries() const
Returns the local number of entries in this matrix.
Definition: Ifpack2_LocalFilter_def.hpp:363
Tpetra::Map< local_ordinal_type, global_ordinal_type, node_type > map_type
Type of the Tpetra::Map specialization that this class uses.
Definition: Ifpack2_LocalFilter_decl.hpp:222
virtual mag_type getFrobeniusNorm() const
The Frobenius norm of the (locally filtered) matrix.
Definition: Ifpack2_LocalFilter_def.hpp:872
virtual Teuchos::RCP< const Tpetra::RowGraph< local_ordinal_type, global_ordinal_type, node_type > > getGraph() const
The (locally filtered) matrix's graph.
Definition: Ifpack2_LocalFilter_def.hpp:292
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
virtual bool supportsRowViews() const
Returns true if RowViews are supported.
Definition: Ifpack2_LocalFilter_def.hpp:863
virtual size_t getLocalNumCols() const
The number of columns in the (locally filtered) matrix.
Definition: Ifpack2_LocalFilter_def.hpp:341
virtual size_t getNumEntriesInGlobalRow(global_ordinal_type globalRow) const
The current number of entries on this node in the specified global row.
Definition: Ifpack2_LocalFilter_def.hpp:377
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_LocalFilter_def.hpp:654
virtual global_ordinal_type getIndexBase() const
Returns the index base for global indices for this matrix.
Definition: Ifpack2_LocalFilter_def.hpp:349
virtual bool isGloballyIndexed() const
Whether the underlying sparse matrix is globally (opposite of locally) indexed.
Definition: Ifpack2_LocalFilter_def.hpp:445
virtual size_t getGlobalMaxNumRowEntries() const
The maximum number of entries across all rows/columns on all processes.
Definition: Ifpack2_LocalFilter_def.hpp:417
LocalFilter(const Teuchos::RCP< const row_matrix_type > &A)
Constructor.
Definition: Ifpack2_LocalFilter_def.hpp:98
virtual ~LocalFilter()
Destructor.
Definition: Ifpack2_LocalFilter_def.hpp:234
virtual global_size_t getGlobalNumRows() const
The number of global rows in this matrix.
Definition: Ifpack2_LocalFilter_def.hpp:320
virtual std::string description() const
A one-line description of this object.
Definition: Ifpack2_LocalFilter_def.hpp:899
virtual size_t getLocalNumRows() const
The number of rows owned on the calling process.
Definition: Ifpack2_LocalFilter_def.hpp:334
virtual size_t getNumEntriesInLocalRow(local_ordinal_type localRow) const
The current number of entries on this node in the specified local row.
Definition: Ifpack2_LocalFilter_def.hpp:396
virtual Teuchos::RCP< const map_type > getRangeMap() const
Returns the Map that describes the range distribution in this matrix.
Definition: Ifpack2_LocalFilter_def.hpp:282
virtual global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
Definition: Ifpack2_LocalFilter_def.hpp:356
virtual Teuchos::RCP< const row_matrix_type > getUnderlyingMatrix() const
Return matrix that LocalFilter was built on.
Definition: Ifpack2_LocalFilter_def.hpp:959
Access only local rows and columns of a sparse matrix.
Definition: Ifpack2_LocalFilter_decl.hpp:161
virtual global_size_t getGlobalNumCols() const
The number of global columns in this matrix.
Definition: Ifpack2_LocalFilter_def.hpp:327
virtual Teuchos::RCP< const map_type > getColMap() const
Returns the Map that describes the column distribution in this matrix.
Definition: Ifpack2_LocalFilter_def.hpp:262
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_LocalFilter_def.hpp:627
virtual bool hasColMap() const
Whether this matrix has a well-defined column Map.
Definition: Ifpack2_LocalFilter_def.hpp:431
virtual local_ordinal_type getBlockSize() const
The number of degrees of freedom per mesh point.
Definition: Ifpack2_LocalFilter_def.hpp:369