43 #ifndef IFPACK2_REORDERFILTER_DEF_HPP
44 #define IFPACK2_REORDERFILTER_DEF_HPP
45 #include "Ifpack2_ReorderFilter_decl.hpp"
48 #include "Tpetra_ConfigDefs.hpp"
49 #include "Tpetra_RowMatrix.hpp"
50 #include "Tpetra_Map.hpp"
51 #include "Tpetra_MultiVector.hpp"
52 #include "Tpetra_Vector.hpp"
56 template<
class MatrixType>
63 reverseperm_ (reverseperm)
66 A_.
is_null (), std::invalid_argument,
67 "Ifpack2::ReorderFilter: The input matrix is null.");
71 A_->getComm()->getSize() != 1, std::invalid_argument,
72 "Ifpack2::ReorderFilter: This class may only be used if the input matrix's "
73 "communicator has one process. This class is an implementation detail of "
74 "Ifpack2::AdditiveSchwarz, and it is not meant to be used otherwise.");
77 A_->getNodeNumRows () != A_->getGlobalNumRows (),
78 std::invalid_argument,
79 "Ifpack2::ReorderFilter: The input matrix is not square.");
82 Indices_.
resize (A_->getNodeMaxNumRowEntries ());
83 Values_.
resize (A_->getNodeMaxNumRowEntries ());
87 template<
class MatrixType>
91 template<
class MatrixType>
98 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
99 template<
class MatrixType>
104 return Teuchos::null;
106 #endif // TPETRA_ENABLE_DEPRECATED_CODE
109 template<
class MatrixType>
114 A_.is_null (), std::runtime_error,
"Ifpack2::ReorderFilter::"
115 "getRowMap: The matrix A is null, so there is no row Map.");
117 return A_->getRowMap ();
121 template<
class MatrixType>
126 A_.is_null (), std::runtime_error,
"Ifpack2::ReorderFilter::"
127 "getColMap: The matrix A is null, so there is no column Map.");
129 return A_->getColMap();
133 template<
class MatrixType>
138 A_.is_null (), std::runtime_error,
"Ifpack2::ReorderFilter::"
139 "getDomainMap: The matrix A is null, so there is no domain Map.");
141 return A_->getDomainMap();
145 template<
class MatrixType>
150 A_.is_null (), std::runtime_error,
"Ifpack2::ReorderFilter::"
151 "getRangeMap: The matrix A is null, so there is no range Map.");
153 return A_->getRangeMap();
157 template<
class MatrixType>
158 Teuchos::RCP<
const Tpetra::RowGraph<
typename MatrixType::local_ordinal_type,
159 typename MatrixType::global_ordinal_type,
160 typename MatrixType::node_type> >
163 throw std::runtime_error(
"Ifpack2::ReorderFilter: does not support getGraph.");
167 template<
class MatrixType>
170 return A_->getGlobalNumRows();
174 template<
class MatrixType>
177 return A_->getGlobalNumCols();
181 template<
class MatrixType>
184 return A_->getNodeNumRows();
188 template<
class MatrixType>
191 return A_->getNodeNumCols();
195 template<
class MatrixType>
198 return A_->getIndexBase();
202 template<
class MatrixType>
205 return A_->getGlobalNumEntries();
209 template<
class MatrixType>
212 return A_->getNodeNumEntries();
216 template<
class MatrixType>
220 if (A_.is_null () || A_->getRowMap ().is_null ()) {
224 const local_ordinal_type lclRow =
225 A_->getRowMap ()->getLocalElement (globalRow);
228 return static_cast<size_t> (0);
230 const local_ordinal_type origLclRow = reverseperm_[lclRow];
231 return A_->getNumEntriesInLocalRow (origLclRow);
237 template<
class MatrixType>
243 if (A_->getRowMap ()->isNodeLocalElement (localRow)) {
245 const local_ordinal_type localReorderedRow = reverseperm_[localRow];
246 return A_->getNumEntriesInLocalRow (localReorderedRow);
249 return static_cast<size_t> (0);
254 template<
class MatrixType>
257 return A_->getGlobalMaxNumRowEntries();
261 template<
class MatrixType>
264 return A_->getNodeMaxNumRowEntries();
268 template<
class MatrixType>
275 template<
class MatrixType>
278 return A_->isLocallyIndexed();
282 template<
class MatrixType>
285 return A_->isGloballyIndexed();
289 template<
class MatrixType>
292 return A_->isFillComplete();
296 template<
class MatrixType>
301 size_t& numEntries)
const
305 using Teuchos::av_reinterpret_cast;
306 typedef local_ordinal_type LO;
307 typedef global_ordinal_type GO;
310 const map_type& rowMap = * (A_->getRowMap ());
311 const local_ordinal_type localRow = rowMap.getLocalElement (globalRow);
313 localRow == OTLO::invalid (), std::invalid_argument,
"Ifpack2::Reorder"
314 "Filter::getGlobalRowCopy: The given global row index " << globalRow
315 <<
" is not owned by the calling process with rank "
316 << rowMap.getComm ()->getRank () <<
".");
318 if (
sizeof (GO) ==
sizeof (LO)) {
320 ArrayView<LO> localInd = av_reinterpret_cast<LO> (globalInd);
321 this->getLocalRowCopy (localRow, localInd, val, numEntries);
324 for (
size_t k = 0; k < numEntries; ++k) {
325 globalInd[k] = rowMap.getGlobalElement (localInd[k]);
331 numEntries = this->getNumEntriesInLocalRow (localRow);
332 Array<LO> localInd (numEntries);
333 this->getLocalRowCopy (localRow, localInd, val, numEntries);
336 for (
size_t k = 0; k < numEntries; ++k) {
337 globalInd[k] = rowMap.getGlobalElement (localInd[k]);
343 template<
class MatrixType>
348 size_t &NumEntries)
const
351 ! A_->getRowMap ()->isNodeLocalElement (LocalRow),
352 std::invalid_argument,
353 "Ifpack2::ReorderFilter::getLocalRowCopy: The given local row index "
354 << LocalRow <<
" is not a valid local row index on the calling process "
355 "with rank " << A_->getRowMap ()->getComm ()->getRank () <<
".");
359 const local_ordinal_type origLclRow = reverseperm_[LocalRow];
360 const size_t numEntries = A_->getNumEntriesInLocalRow (origLclRow);
363 static_cast<size_t> (Indices.
size ()) < numEntries ||
364 static_cast<size_t> (Values.
size ()) < numEntries,
365 std::invalid_argument,
366 "Ifpack2::ReorderFilter::getLocalRowCopy: The given array views are not "
367 "long enough to store all the data in the given row " << LocalRow
368 <<
". Indices.size() = " << Indices.
size () <<
", Values.size() = "
369 << Values.
size () <<
", but the (original) row has " << numEntries
372 A_->getLocalRowCopy (origLclRow, Indices, Values, NumEntries);
377 for (
size_t i = 0; i < NumEntries; ++i) {
378 Indices[i] = perm_[Indices[i]];
383 template<
class MatrixType>
389 throw std::runtime_error(
"Ifpack2::ReorderFilter: does not support getGlobalRowView.");
393 template<
class MatrixType>
399 throw std::runtime_error(
"Ifpack2::ReorderFilter: does not support getLocalRowView.");
403 template<
class MatrixType>
405 getLocalDiagCopy (Tpetra::Vector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &diag)
const
408 return A_->getLocalDiagCopy(diag);
412 template<
class MatrixType>
415 throw std::runtime_error(
"Ifpack2::ReorderFilter does not support leftScale.");
419 template<
class MatrixType>
422 throw std::runtime_error(
"Ifpack2::ReorderFilter does not support rightScale.");
426 template<
class MatrixType>
428 apply (
const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &X,
429 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &Y,
432 scalar_type beta)
const
437 alpha != STS::one () || beta != STS::zero (), std::logic_error,
438 "Ifpack2::ReorderFilter::apply is only implemented for alpha = 1 and "
439 "beta = 0. You set alpha = " << alpha <<
" and beta = " << beta <<
".");
444 X.getNumVectors() != Y.getNumVectors(), std::runtime_error,
445 "Ifpack2::ReorderFilter::apply: X.getNumVectors() != Y.getNumVectors().");
447 const scalar_type zero = STS::zero ();
452 const size_t NumVectors = Y.getNumVectors ();
454 for (
size_t i = 0; i < A_->getNodeNumRows (); ++i) {
457 getLocalRowCopy (i, Indices_ (), Values_ (), Nnz);
459 for (
size_t j = 0; j < Nnz; ++j) {
460 for (
size_t k = 0; k < NumVectors; ++k) {
461 y_ptr[k][i] += Values_[j] * x_ptr[k][Indices_[j]];
466 for (
size_t j = 0; j < Nnz; ++j) {
467 for (
size_t k = 0; k < NumVectors; ++k) {
468 y_ptr[k][Indices_[j]] += Values_[j] * x_ptr[k][i];
473 for (
size_t j = 0; j < Nnz; ++j) {
474 for (
size_t k = 0; k < NumVectors; ++k) {
475 y_ptr[k][Indices_[j]] += STS::conjugate(Values_[j]) * x_ptr[k][i];
483 template<
class MatrixType>
490 template<
class MatrixType>
497 template<
class MatrixType>
501 return A_->getFrobeniusNorm ();
505 template<
class MatrixType>
508 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &reorderedY)
const
510 this->
template permuteOriginalToReorderedTempl<scalar_type,scalar_type>(originalX, reorderedY);
514 template<
class MatrixType>
515 template<
class DomainScalar,
class RangeScalar>
517 Tpetra::MultiVector<RangeScalar,local_ordinal_type,global_ordinal_type,node_type> &reorderedY)
const
520 "Ifpack2::ReorderFilter::permuteOriginalToReordered ERROR: X.getNumVectors() != Y.getNumVectors().");
525 for(
size_t k=0; k < originalX.getNumVectors(); k++)
526 for(local_ordinal_type i=0; (size_t)i< originalX.getLocalLength(); i++)
527 y_ptr[k][perm_[i]] = (RangeScalar)x_ptr[k][i];
531 template<
class MatrixType>
533 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> &originalY)
const
535 this->
template permuteReorderedToOriginalTempl<scalar_type,scalar_type>(reorderedX, originalY);
539 template<
class MatrixType>
540 template<
class DomainScalar,
class RangeScalar>
543 Tpetra::MultiVector<RangeScalar,local_ordinal_type,global_ordinal_type,node_type> &originalY)
const
546 reorderedX.getNumVectors() != originalY.getNumVectors(),
548 "Ifpack2::ReorderFilter::permuteReorderedToOriginal: "
549 "X.getNumVectors() != Y.getNumVectors().");
551 #ifdef HAVE_IFPACK2_DEBUG
554 typedef typename STS::magnitudeType magnitude_type;
557 reorderedX.norm2 (norms ());
560 j < reorderedX.getNumVectors (); ++j) {
561 if (STM::isnaninf (norms[j])) {
567 ! good, std::runtime_error,
"Ifpack2::ReorderFilter::"
568 "permuteReorderedToOriginalTempl: The 2-norm of the input reorderedX is "
571 #endif // HAVE_IFPACK2_DEBUG
576 for (
size_t k = 0; k < reorderedX.getNumVectors (); ++k) {
577 for (local_ordinal_type i = 0; (size_t)i < reorderedX.getLocalLength (); ++i) {
578 y_ptr[k][reverseperm_[i]] = (RangeScalar) x_ptr[k][i];
582 #ifdef HAVE_IFPACK2_DEBUG
585 typedef typename STS::magnitudeType magnitude_type;
588 originalY.norm2 (norms ());
591 j < originalY.getNumVectors (); ++j) {
592 if (STM::isnaninf (norms[j])) {
598 ! good, std::runtime_error,
"Ifpack2::ReorderFilter::"
599 "permuteReorderedToOriginalTempl: The 2-norm of the output originalY is "
602 #endif // HAVE_IFPACK2_DEBUG
607 #define IFPACK2_REORDERFILTER_INSTANT(S,LO,GO,N) \
608 template class Ifpack2::ReorderFilter< Tpetra::RowMatrix<S, LO, GO, N> >;
virtual void getGlobalRowView(global_ordinal_type GlobalRow, Teuchos::ArrayView< const global_ordinal_type > &indices, Teuchos::ArrayView< const scalar_type > &values) const
Extract a const, non-persisting view of global indices in a specified row of the matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:385
Wraps a Tpetra::RowMatrix in a filter that reorders local rows and columns.
Definition: Ifpack2_ReorderFilter_decl.hpp:69
virtual mag_type getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:498
virtual Teuchos::RCP< const map_type > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:111
virtual bool supportsRowViews() const
Returns true if RowViews are supported.
Definition: Ifpack2_ReorderFilter_def.hpp:491
virtual void getGlobalRowCopy(global_ordinal_type GlobalRow, const Teuchos::ArrayView< global_ordinal_type > &Indices, const Teuchos::ArrayView< scalar_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_ReorderFilter_def.hpp:298
virtual size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
Definition: Ifpack2_ReorderFilter_def.hpp:255
virtual Teuchos::RCP< const map_type > getDomainMap() const
Returns the Map that describes the domain distribution in this matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:135
virtual size_t getNodeMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
Definition: Ifpack2_ReorderFilter_def.hpp:262
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual global_size_t getGlobalNumRows() const
Returns the number of global rows in this matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:168
virtual Teuchos::RCP< const map_type > getColMap() const
Returns the Map that describes the column distribution in this matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:123
ReorderFilter(const Teuchos::RCP< const row_matrix_type > &A, const Teuchos::ArrayRCP< local_ordinal_type > &perm, const Teuchos::ArrayRCP< local_ordinal_type > &reverseperm)
Constructor.
Definition: Ifpack2_ReorderFilter_def.hpp:58
virtual Teuchos::RCP< const map_type > getRangeMap() const
Returns the Map that describes the range distribution in this matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:147
virtual size_t getNodeNumEntries() const
Returns the local number of entries in this matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:210
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
, where Op(A) is either A, , or .
Definition: Ifpack2_ReorderFilter_def.hpp:428
virtual size_t getNodeNumRows() const
Returns the number of rows owned on the calling node.
Definition: Ifpack2_ReorderFilter_def.hpp:182
virtual size_t getNodeNumCols() const
Returns the number of columns needed to apply the forward operator on this node, i.e., the number of elements listed in the column map.
Definition: Ifpack2_ReorderFilter_def.hpp:189
virtual bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
Definition: Ifpack2_ReorderFilter_def.hpp:283
virtual Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
The matrix's communicator.
Definition: Ifpack2_ReorderFilter_def.hpp:92
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_ReorderFilter_def.hpp:405
void resize(size_type new_size, const value_type &x=value_type())
virtual size_t getNumEntriesInLocalRow(local_ordinal_type localRow) const
The current number of entries in this matrix, stored on the calling process, in the row whose local i...
Definition: Ifpack2_ReorderFilter_def.hpp:239
virtual bool hasColMap() const
Indicates whether this matrix has a well-defined column map.
Definition: Ifpack2_ReorderFilter_def.hpp:269
virtual size_t getNumEntriesInGlobalRow(global_ordinal_type globalRow) const
The current number of entries in this matrix, stored on the calling process, in the row whose global ...
Definition: Ifpack2_ReorderFilter_def.hpp:218
virtual void getLocalRowCopy(local_ordinal_type DropRow, const Teuchos::ArrayView< local_ordinal_type > &Indices, const Teuchos::ArrayView< scalar_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_ReorderFilter_def.hpp:345
virtual global_size_t getGlobalNumCols() const
Returns the number of global columns in this matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:175
virtual bool hasTransposeApply() const
Whether apply() can apply the transpose or conjugate transpose.
Definition: Ifpack2_ReorderFilter_def.hpp:484
virtual global_ordinal_type getIndexBase() const
Returns the index base for global indices for this matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:196
virtual void getLocalRowView(local_ordinal_type LocalRow, Teuchos::ArrayView< const local_ordinal_type > &indices, Teuchos::ArrayView< const scalar_type > &values) const
Extract a const, non-persisting view of local indices in a specified row of the matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:395
virtual global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:203
virtual bool isFillComplete() const
Returns true if fillComplete() has been called.
Definition: Ifpack2_ReorderFilter_def.hpp:290
virtual void permuteReorderedToOriginal(const Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &reorderedX, Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &originalY) const
Permute multivector: reordered-to-original.
Definition: Ifpack2_ReorderFilter_def.hpp:532
virtual Teuchos::RCP< const Tpetra::RowGraph< local_ordinal_type, global_ordinal_type, node_type > > getGraph() const
Returns the RowGraph associated with this matrix.
Definition: Ifpack2_ReorderFilter_def.hpp:161
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_ReorderFilter_def.hpp:413
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_ReorderFilter_def.hpp:420
virtual void permuteOriginalToReordered(const Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &originalX, Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &reorderedY) const
Permute multivector: original-to-reordered.
Definition: Ifpack2_ReorderFilter_def.hpp:507
virtual ~ReorderFilter()
Destructor.
Definition: Ifpack2_ReorderFilter_def.hpp:88
virtual bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
Definition: Ifpack2_ReorderFilter_def.hpp:276