43 #ifndef IFPACK2_CONTAINER_DEF_HPP
44 #define IFPACK2_CONTAINER_DEF_HPP
53 template<
class MatrixType>
58 inputMatrix_ (matrix),
59 inputCrsMatrix_ (Teuchos::rcp_dynamic_cast<const crs_matrix_type>(inputMatrix_)),
60 inputBlockMatrix_ (Teuchos::rcp_dynamic_cast<const block_crs_matrix_type>(inputMatrix_)),
61 pointIndexed_(pointIndexed),
62 IsInitialized_(false),
85 #ifdef HAVE_IFPACK2_DEBUG
93 LO row = blockRows[j];
100 !rowMap.isNodeLocalElement(row),
101 std::invalid_argument,
"Ifpack2::Container: "
102 "On process " << rowMap.getComm()->getRank() <<
" of "
103 << rowMap.getComm()->getSize() <<
", in the given set of local row "
105 "entries is not valid local row index on the calling process: "
112 template<
class MatrixType>
116 template<
class MatrixType>
121 (&blockRows_[blockOffsets_[blockIndex]], blockSizes_[blockIndex]);
124 template<
class MatrixType>
129 LO totalBlockRows = 0;
130 numBlocks_ = partitions.size();
131 blockSizes_.resize(numBlocks_);
132 blockOffsets_.resize(numBlocks_);
134 for(
int i = 0; i < numBlocks_; i++)
136 LO rowsInBlock = partitions[i].size();
137 blockSizes_[i] = rowsInBlock;
138 blockOffsets_[i] = totalBlockRows;
139 totalBlockRows += rowsInBlock;
140 maxBlockSize_ = std::max(maxBlockSize_, rowsInBlock * scalarsPerRow_);
142 blockRows_.resize(totalBlockRows);
145 for(
int i = 0; i < numBlocks_; i++)
147 for(
int j = 0; j < blockSizes_[i]; j++)
149 blockRows_[iter++] = partitions[i][j];
154 template<
class MatrixType>
159 Diag_ =
rcp(
new vector_type(inputMatrix_->getDomainMap()));
160 inputMatrix_->getLocalDiagCopy(*Diag_);
164 template<
class MatrixType>
166 return IsInitialized_;
169 template<
class MatrixType>
174 template<
class MatrixType>
181 template<
class MatrixType>
188 template<
class MatrixType>
195 template<
class MatrixType>
197 SC dampingFactor, LO i)
const
202 template <
class MatrixType>
206 const ISC one = STS::one();
208 size_t numVecs = X.extent(1);
210 for (LO i = 0; i < numBlocks_; i++)
213 if(blockSizes_[i] != 1 || hasBlockCrs_)
215 if(blockSizes_[i] == 0 )
221 LO LRID = blockRows_[blockOffsets_[i]];
223 auto diagView = Diag_->getLocalViewHost(Tpetra::Access::ReadOnly);
224 ISC d = one * (
static_cast<ISC
> (dampingFactor)) / diagView(LRID, 0);
225 for(
size_t nv = 0; nv < numVecs; nv++)
228 Y(LRID, nv) += x * d;
234 template <
class MatrixType>
235 void Container<MatrixType>::DoOverlappingJacobi(ConstHostView X, HostView Y, ConstHostView W, SC dampingFactor,
bool nonsymScaling)
const
239 for(LO i = 0; i < numBlocks_; i++)
242 if(blockSizes_[i] == 0)
244 if(blockSizes_[i] != 1) {
253 HostView tempo(
"", X.extent(0), X.extent(1));
254 size_t numVecs = X.extent(1);
255 LO bOffset = blockOffsets_[i];
256 for (LO ii = 0; ii < blockSizes_[i]; ii++) {
257 LO LRID = blockRows_[bOffset++];
258 for (
size_t jj = 0; jj < numVecs; jj++) tempo(LRID,jj)=X(LRID,jj)/ W(LRID,0);
265 const ISC one = STS::one();
266 size_t numVecs = X.extent(1);
267 LO LRID = blockRows_[blockOffsets_[i]];
269 auto diagView = Diag_->getLocalViewHost(Tpetra::Access::ReadOnly);
270 ISC recip = one / diagView(LRID, 0);
271 ISC wval = W(LRID,0);
272 ISC combo = wval*recip;
273 ISC d = combo*(
static_cast<ISC
> (dampingFactor));
274 for(
size_t nv = 0; nv < numVecs; nv++)
277 Y(LRID, nv) = x * d + Y(LRID, nv);
285 template<
class MatrixType,
typename LocalScalarType>
288 SC dampingFactor, LO i)
const
292 size_t numVecs = X.extent(1);
293 const ISC one = STS::one();
294 if(this->blockSizes_[i] == 0)
296 if(this->hasBlockCrs_ && !this->pointIndexed_)
299 ArrayView<const LO> blockRows = this->getBlockRows(i);
300 const size_t localNumRows = this->blockSizes_[i];
301 using inds_type =
typename block_crs_matrix_type::local_inds_host_view_type;
302 using vals_type =
typename block_crs_matrix_type::values_host_view_type;
303 for(
size_t j = 0; j < localNumRows; j++)
305 LO row = blockRows[j];
308 this->inputBlockMatrix_->getLocalRowView(row, colinds, values);
309 LO numEntries = (LO) colinds.size();
310 for(
size_t m = 0; m < numVecs; m++)
312 for (
int localR = 0; localR < this->bcrsBlockSize_; localR++)
313 Resid(row * this->bcrsBlockSize_ + localR, m) = X(row * this->bcrsBlockSize_ + localR, m);
314 for (LO k = 0; k < numEntries; ++k)
316 const LO col = colinds[k];
317 for(
int localR = 0; localR < this->bcrsBlockSize_; localR++)
319 for(
int localC = 0; localC < this->bcrsBlockSize_; localC++)
321 Resid(row * this->bcrsBlockSize_ + localR, m) -=
322 values[k * this->bcrsBlockSize_ * this->bcrsBlockSize_ + localR + localC * this->bcrsBlockSize_]
323 * Y2(col * this->bcrsBlockSize_ + localC, m); }
336 else if(!this->hasBlockCrs_ && this->blockSizes_[i] == 1)
340 LO LRID = this->blockOffsets_[i];
344 container_exec_space().fence();
345 auto localA = this->inputCrsMatrix_->getLocalMatrixHost();
346 using size_type =
typename crs_matrix_type::local_matrix_host_type::size_type;
347 const auto& rowmap = localA.graph.row_map;
348 const auto& entries = localA.graph.entries;
349 const auto& values = localA.values;
351 auto diagView = this->Diag_->getLocalViewHost(Tpetra::Access::ReadOnly);
352 ISC d = (
static_cast<ISC> (dampingFactor)) / diagView(LRID, 0);
353 for(
size_t m = 0; m < numVecs; m++)
357 for(size_type k = rowmap(LRID); k < rowmap(LRID + 1); k++) {
358 const LO col = entries(k);
359 r -= values(k) * Y2(col, m);
366 else if(!this->inputCrsMatrix_.is_null() &&
367 std::is_same<typename crs_matrix_type::device_type::memory_space, Kokkos::HostSpace>::value)
372 container_exec_space().fence();
373 auto localA = this->inputCrsMatrix_->getLocalMatrixHost();
374 using size_type =
typename crs_matrix_type::local_matrix_host_type::size_type;
375 const auto& rowmap = localA.graph.row_map;
376 const auto& entries = localA.graph.entries;
377 const auto& values = localA.values;
378 ArrayView<const LO> blockRows = this->getBlockRows(i);
379 for(
size_t j = 0; j < size_t(blockRows.size()); j++)
381 const LO row = blockRows[j];
382 for(
size_t m = 0; m < numVecs; m++)
385 for(size_type k = rowmap(row); k < rowmap(row + 1); k++)
387 const LO col = entries(k);
388 r -= values(k) * Y2(col, m);
405 ArrayView<const LO> blockRows = this->getBlockRows(i);
406 for(
size_t j = 0; j < size_t(blockRows.size()); j++)
408 const LO row = blockRows[j];
409 auto rowView = getInputRowView(row);
410 for(
size_t m = 0; m < numVecs; m++)
412 Resid(row, m) = X(row, m);
413 for (
size_t k = 0; k < rowView.size(); ++k)
415 const LO col = rowView.ind(k);
416 Resid(row, m) -= rowView.val(k) * Y2(col, m);
430 template<
class MatrixType>
432 DoGaussSeidel(ConstHostView X, HostView Y, HostView Y2, SC dampingFactor)
const
440 using Teuchos::rcpFromRef;
443 auto numVecs = X.extent(1);
445 HostView Resid(
"", X.extent(0), X.extent(1));
446 for(LO i = 0; i < numBlocks_; i++)
448 DoGSBlock(X, Y, Y2, Resid, dampingFactor, i);
452 auto numMyRows = inputMatrix_->getLocalNumRows();
453 for (
size_t m = 0; m < numVecs; ++m)
455 for (
size_t i = 0; i < numMyRows * bcrsBlockSize_; ++i)
463 template<
class MatrixType>
464 void Container<MatrixType>::
465 DoSGS(ConstHostView X, HostView Y, HostView Y2, SC dampingFactor)
const
475 using Teuchos::rcpFromRef;
476 auto numVecs = X.extent(1);
477 HostView Resid(
"", X.extent(0), X.extent(1));
479 for(LO i = 0; i < numBlocks_; i++)
481 DoGSBlock(X, Y, Y2, Resid, dampingFactor, i);
483 static_assert(std::is_signed<LO>::value,
484 "Local ordinal must be signed (unsigned breaks reverse iteration to 0)");
486 for(LO i = numBlocks_ - 1; i >= 0; --i)
488 DoGSBlock(X, Y, Y2, Resid, dampingFactor, i);
492 auto numMyRows = inputMatrix_->getLocalNumRows();
493 for (
size_t m = 0; m < numVecs; ++m)
495 for (
size_t i = 0; i < numMyRows * bcrsBlockSize_; ++i)
503 template<
class MatrixType>
504 void Container<MatrixType>::
510 blockOffsets_.clear();
511 Diag_ = Teuchos::null;
516 template<
class MatrixType,
class LocalScalarType>
517 ContainerImpl<MatrixType, LocalScalarType>::
522 : Container<MatrixType>(matrix, partitions, pointIndexed) {}
524 template<
class MatrixType,
class LocalScalarType>
528 template<
class MatrixType,
class LocalScalarType>
532 template<
class MatrixType,
class LocalScalarType>
542 template<
class MatrixType,
class LocalScalarType>
546 ConstHostView XView = X.getLocalViewHost(Tpetra::Access::ReadOnly);
547 HostView YView = Y.getLocalViewHost(Tpetra::Access::ReadWrite);
548 this->apply (XView, YView, 0);
551 template<
class MatrixType,
class LocalScalarType>
555 vector_type& W)
const
557 ConstHostView XView = X.getLocalViewHost(Tpetra::Access::ReadOnly);
558 HostView YView = Y.getLocalViewHost(Tpetra::Access::ReadWrite);
559 ConstHostView WView = W.getLocalViewHost(Tpetra::Access::ReadOnly);
560 weightedApply (XView, YView, WView, 0);
563 template<
class MatrixType,
class LocalScalarType>
570 template<
class MatrixType,
class LocalScalarType>
577 const LSC beta)
const
582 template<
class MatrixType,
class LocalScalarType>
583 typename ContainerImpl<MatrixType, LocalScalarType>::LO
589 const map_type& globalRowMap = *(this->inputMatrix_->getRowMap());
590 const map_type& globalColMap = *(this->inputMatrix_->getColMap());
593 if(this->pointIndexed_)
595 rowLID = row / this->bcrsBlockSize_;
596 dofOffset = row % this->bcrsBlockSize_;
598 GO diagGID = globalRowMap.getGlobalElement(rowLID);
600 diagGID == GO_INVALID,
601 std::runtime_error,
"Ifpack2::Container::translateRowToCol: "
602 "On process " << this->inputMatrix_->getRowMap()->getComm()->getRank() <<
603 ", at least one row index in the set of local "
604 "row indices given to the constructor is not a valid local row index in "
605 "the input matrix's row Map on this process. This should be impossible "
606 "because the constructor checks for this case. Here is the complete set "
607 "of invalid local row indices: " << rowLID <<
". "
608 "Please report this bug to the Ifpack2 developers.");
610 LO colLID = globalColMap.getLocalElement(diagGID);
612 colLID == LO_INVALID,
613 std::runtime_error,
"Ifpack2::Container::translateRowToCol: "
614 "On process " << this->inputMatrix_->getRowMap()->getComm()->getRank() <<
", "
615 "at least one row index in the set of row indices given to the constructor "
616 "does not have a corresponding column index in the input matrix's column "
617 "Map. This probably means that the column(s) in question is/are empty on "
618 "this process, which would make the submatrix to extract structurally "
619 "singular. The invalid global column index is " << diagGID <<
".");
621 if(this->pointIndexed_)
622 return colLID * this->bcrsBlockSize_ + dofOffset;
626 template<
class MatrixType,
class LocalScalarType>
649 ! this->IsComputed_, std::runtime_error,
"Ifpack2::Container::apply: "
650 "You must have called the compute() method before you may call apply(). "
651 "You may call the apply() method as many times as you want after calling "
652 "compute() once, but you must have called compute() at least once.");
654 const size_t numVecs = X.extent(1);
685 if(X_localBlocks_.size() == 0 || X.extent(1) != X_local_.extent(1))
688 X_localBlocks_.clear();
689 Y_localBlocks_.clear();
690 X_localBlocks_.reserve(this->numBlocks_);
691 Y_localBlocks_.reserve(this->numBlocks_);
693 X_local_ = HostViewLocal(
"X_local", this->blockRows_.size() * this->scalarsPerRow_, numVecs);
694 Y_local_ = HostViewLocal(
"Y_local", this->blockRows_.size() * this->scalarsPerRow_, numVecs);
698 for(
int i = 0; i < this->numBlocks_; i++)
700 auto blockBounds = std::make_pair(this->blockOffsets_[i] * this->scalarsPerRow_,
701 (this->blockOffsets_[i] + this->blockSizes_[i]) * this->scalarsPerRow_);
702 X_localBlocks_.emplace_back(X_local_, blockBounds, Kokkos::ALL());
703 Y_localBlocks_.emplace_back(Y_local_, blockBounds, Kokkos::ALL());
707 const ArrayView<const LO> blockRows = this->getBlockRows(blockIndex);
709 if(this->scalarsPerRow_ == 1)
710 mvgs.gatherViewToView (X_localBlocks_[blockIndex], X, blockRows);
712 mvgs.gatherViewToViewBlock (X_localBlocks_[blockIndex], X, blockRows, this->scalarsPerRow_);
719 if(this->scalarsPerRow_ == 1)
720 mvgs.gatherViewToView (Y_localBlocks_[blockIndex], Y, blockRows);
722 mvgs.gatherViewToViewBlock (Y_localBlocks_[blockIndex], Y, blockRows, this->scalarsPerRow_);
726 this->solveBlock (X_localBlocks_[blockIndex], Y_localBlocks_[blockIndex], blockIndex, mode,
731 if(this->scalarsPerRow_ == 1)
732 mvgs.scatterViewToView (Y, Y_localBlocks_[blockIndex], blockRows);
734 mvgs.scatterViewToViewBlock (Y, Y_localBlocks_[blockIndex], blockRows, this->scalarsPerRow_);
737 template<
class MatrixType,
class LocalScalarType>
754 using Teuchos::rcp_const_cast;
765 const char prefix[] =
"Ifpack2::Container::weightedApply: ";
767 ! this->IsComputed_, std::runtime_error, prefix <<
"You must have called the "
768 "compute() method before you may call this method. You may call "
769 "weightedApply() as many times as you want after calling compute() once, "
770 "but you must have called compute() at least once first.");
774 this->scalarsPerRow_ > 1, std::logic_error, prefix <<
"Use of block rows isn't allowed "
775 "in overlapping Jacobi (the only method that uses weightedApply");
777 const size_t numVecs = X.extent(1);
780 X.extent(1) != Y.extent(1), std::runtime_error,
781 prefix <<
"X and Y have different numbers of vectors (columns). X has "
782 << X.extent(1) <<
", but Y has " << Y.extent(1) <<
".");
788 const size_t numRows = this->blockSizes_[blockIndex];
813 if(X_localBlocks_.size() == 0 || X.extent(1) != X_local_.extent(1))
816 X_localBlocks_.clear();
817 Y_localBlocks_.clear();
818 X_localBlocks_.reserve(this->numBlocks_);
819 Y_localBlocks_.reserve(this->numBlocks_);
821 X_local_ = HostViewLocal(
"X_local", this->blockRows_.size() * this->scalarsPerRow_, numVecs);
822 Y_local_ = HostViewLocal(
"Y_local", this->blockRows_.size() * this->scalarsPerRow_, numVecs);
826 for(
int i = 0; i < this->numBlocks_; i++)
828 auto blockBounds = std::make_pair(this->blockOffsets_[i] * this->scalarsPerRow_,
829 (this->blockOffsets_[i] + this->blockSizes_[i]) * this->scalarsPerRow_);
830 X_localBlocks_.emplace_back(X_local_, blockBounds, Kokkos::ALL());
831 Y_localBlocks_.emplace_back(Y_local_, blockBounds, Kokkos::ALL());
834 if(
int(weightedApplyScratch_.extent(0)) != 3 * this->maxBlockSize_ ||
835 weightedApplyScratch_.extent(1) != numVecs)
837 weightedApplyScratch_ = HostViewLocal(
"weightedApply scratch", 3 * this->maxBlockSize_, numVecs);
840 ArrayView<const LO> blockRows = this->getBlockRows(blockIndex);
845 mvgs.gatherViewToView (X_localBlocks_[blockIndex], X, blockRows);
851 mvgs.gatherViewToView (Y_localBlocks_[blockIndex], Y, blockRows);
862 auto maxBS = this->maxBlockSize_;
863 auto bs = this->blockSizes_[blockIndex] * this->scalarsPerRow_;
865 HostSubviewLocal D_local(weightedApplyScratch_, std::make_pair(0, bs), std::make_pair(0, 1));
866 mvgs.gatherViewToView (D_local, D, blockRows);
867 HostSubviewLocal X_scaled(weightedApplyScratch_, std::make_pair(maxBS, maxBS + bs), Kokkos::ALL());
868 for(
size_t j = 0; j < numVecs; j++)
869 for(
size_t i = 0; i < numRows; i++)
870 X_scaled(i, j) = X_localBlocks_[blockIndex](i, j) * D_local(i, 0);
872 HostSubviewLocal Y_temp(weightedApplyScratch_, std::make_pair(maxBS * 2, maxBS * 2 + bs), Kokkos::ALL());
874 this->solveBlock (X_scaled, Y_temp, blockIndex, mode, STS::one(), STS::zero());
882 for(
size_t j = 0; j < numVecs; j++)
883 for(
size_t i = 0; i < numRows; i++)
884 Y_localBlocks_[blockIndex](i, j) = b * Y_localBlocks_[blockIndex](i, j) + a * Y_temp(i, j) * D_local(i, 0);
888 mvgs.scatterViewToView (Y, Y_localBlocks_[blockIndex], blockRows);
891 template<
class MatrixType,
class LocalScalarType>
893 typename ContainerImpl<MatrixType, LocalScalarType>::SC,
894 typename ContainerImpl<MatrixType, LocalScalarType>::LO,
895 typename ContainerImpl<MatrixType, LocalScalarType>::GO,
896 typename ContainerImpl<MatrixType, LocalScalarType>::NO>
901 typedef typename MatrixType::nonconst_local_inds_host_view_type nonconst_local_inds_host_view_type;
902 typedef typename MatrixType::nonconst_values_host_view_type nonconst_values_host_view_type;
904 typedef typename MatrixType::local_inds_host_view_type local_inds_host_view_type;
905 typedef typename MatrixType::values_host_view_type values_host_view_type;
906 using IST =
typename row_matrix_type::impl_scalar_type;
908 if(this->hasBlockCrs_)
910 using h_inds_type =
typename block_crs_matrix_type::local_inds_host_view_type;
911 using h_vals_type =
typename block_crs_matrix_type::values_host_view_type;
914 this->inputBlockMatrix_->getLocalRowView(row / this->bcrsBlockSize_, colinds, values);
915 size_t numEntries = colinds.size();
918 h_vals_type subvals = Kokkos::subview(values,std::pair<size_t,size_t>(row % this->bcrsBlockSize_,values.size()));
919 return StridedRowView(subvals, colinds, this->bcrsBlockSize_, numEntries * this->bcrsBlockSize_);
921 else if(!this->inputMatrix_->supportsRowViews())
923 size_t maxEntries = this->inputMatrix_->getLocalMaxNumRowEntries();
926 nonconst_local_inds_host_view_type inds_v(inds.
data(),maxEntries);
927 nonconst_values_host_view_type vals_v(reinterpret_cast<IST*>(vals.
data()),maxEntries);
929 this->inputMatrix_->getLocalRowCopy(row, inds_v, vals_v, numEntries);
936 local_inds_host_view_type colinds;
937 values_host_view_type values;
938 this->inputMatrix_->getLocalRowView(row, colinds, values);
943 template<
class MatrixType,
class LocalScalarType>
947 X_localBlocks_.clear();
948 Y_localBlocks_.clear();
949 X_local_ = HostViewLocal();
950 Y_local_ = HostViewLocal();
957 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
959 StridedRowView(h_vals_type vals_, h_inds_type inds_,
int blockSize_,
size_t nnz_)
960 : vals(vals_), inds(inds_), blockSize(blockSize_), nnz(nnz_)
963 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
966 : vals(), inds(), blockSize(1), nnz(vals_.size())
968 valsCopy.
swap(vals_);
969 indsCopy.
swap(inds_);
972 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
976 #ifdef HAVE_IFPACK2_DEBUG
978 "Out-of-bounds access into Ifpack2::Container::StridedRowView");
985 return vals[i * blockSize];
991 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
992 LocalOrdinal StridedRowView<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
995 #ifdef HAVE_IFPACK2_DEBUG
997 "Out-of-bounds access into Ifpack2::Container::StridedRowView");
1005 return inds[i / blockSize] * blockSize + i % blockSize;
1011 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1012 size_t StridedRowView<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
1021 template <
class MatrixType>
1022 std::ostream& operator<<(std::ostream& os, const Ifpack2::Container<MatrixType>& obj)
1024 return obj.print(os);
1027 #define IFPACK2_CONTAINER_INSTANT(S,LO,GO,N) \
1028 template class Ifpack2::Container<Tpetra::RowMatrix<S, LO, GO, N>>; \
1029 template class Ifpack2::ContainerImpl<Tpetra::RowMatrix<S, LO, GO, N>, S>; \
1030 template class Ifpack2::Details::StridedRowView<S, LO, GO, N>; \
1031 template std::ostream& operator<< <Tpetra::RowMatrix<S, LO, GO, N>>( \
1032 std::ostream& os, const Ifpack2::Container<Tpetra::RowMatrix<S, LO, GO, N>>& obj);
virtual void apply(ConstHostView X, HostView Y, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const
Compute Y := alpha * M^{-1} X + beta*Y.
Definition: Ifpack2_Container_def.hpp:628
void weightedApplyMV(const mv_type &X, mv_type &Y, vector_type &W) const
Wrapper for weightedApply with MVs, used in unit tests (never called by BlockRelaxation) ...
Definition: Ifpack2_Container_def.hpp:553
bool isComputed() const
Whether the container has been successfully computed.
Definition: Ifpack2_Container_def.hpp:170
void DoGSBlock(ConstHostView X, HostView Y, HostView Y2, HostView Resid, SC dampingFactor, LO i) const
Do one step of Gauss-Seidel on block i (used by DoGaussSeidel and DoSGS)
Definition: Ifpack2_Container_def.hpp:286
void setBlockSizes(const Teuchos::Array< Teuchos::Array< LO > > &partitions)
Initialize arrays with information about block sizes.
Definition: Ifpack2_Container_def.hpp:125
Teuchos::Array< LO > blockSizes_
Number of rows in each block.
Definition: Ifpack2_Container_decl.hpp:296
typename Kokkos::ArithTraits< SC >::val_type ISC
Internal representation of Scalar in Kokkos::View.
Definition: Ifpack2_Container_decl.hpp:135
StridedRowView(h_vals_type vals_, h_inds_type inds_, int blockSize_, size_t nnz_)
Constructor for row views (preferred)
Definition: Ifpack2_Container_def.hpp:959
virtual void applyMV(const mv_type &X, mv_type &Y) const
Wrapper for apply with MultiVector.
Definition: Ifpack2_Container_def.hpp:176
bool pointIndexed_
(If hasBlockCrs_) Whether the blocks are described using sub-block row indices instead of full block ...
Definition: Ifpack2_Container_decl.hpp:314
int numBlocks_
The number of blocks (partitions) in the container.
Definition: Ifpack2_Container_decl.hpp:292
The implementation of the numerical features of Container (Jacobi, Gauss-Seidel, SGS). This class allows a custom scalar type (LocalScalarType) to be used for storing blocks and solving the block systems. Hiding this template parameter from the Container interface simplifies the BlockRelaxation and ContainerFactory classes.
Definition: Ifpack2_Container_decl.hpp:343
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual void setParameters(const Teuchos::ParameterList &List)
Set parameters, if any.
Definition: Ifpack2_Container_def.hpp:530
bool isInitialized() const
Whether the container has been successfully initialized.
Definition: Ifpack2_Container_def.hpp:165
LO scalarsPerRow_
Definition: Ifpack2_Container_decl.hpp:317
virtual ~Container()
Destructor.
Definition: Ifpack2_Container_def.hpp:114
bool hasBlockCrs_
Whether the input matrix is a BlockCRS matrix.
Definition: Ifpack2_Container_decl.hpp:310
Details::StridedRowView< SC, LO, GO, NO > getInputRowView(LO row) const
View a row of the input matrix.
Definition: Ifpack2_Container_def.hpp:898
LocalScalarType LSC
The internal representation of LocalScalarType in Kokkos::View.
Definition: Ifpack2_Container_decl.hpp:363
Declaration and definition of the Ifpack2::Details::MultiVectorLocalGatherScatter class...
Container(const Teuchos::RCP< const row_matrix_type > &matrix, const Teuchos::Array< Teuchos::Array< LO > > &partitions, bool pointIndexed)
Constructor.
Definition: Ifpack2_Container_def.hpp:54
int bcrsBlockSize_
If hasBlockCrs_, the number of DOFs per vertex. Otherwise 1.
Definition: Ifpack2_Container_decl.hpp:312
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
#define TEUCHOS_TEST_FOR_EXCEPT_MSG(throw_exception_test, msg)
LO translateRowToCol(LO row)
Definition: Ifpack2_Container_def.hpp:585
Teuchos::RCP< const block_crs_matrix_type > inputBlockMatrix_
The input matrix, dynamic cast to BlockCrsMatrix. May be null.
Definition: Ifpack2_Container_decl.hpp:289
Teuchos::RCP< const row_matrix_type > inputMatrix_
The input matrix to the constructor.
Definition: Ifpack2_Container_decl.hpp:283
virtual void weightedApplyMV(const mv_type &X, mv_type &Y, vector_type &W) const
Wrapper for weightedApply with MultiVector.
Definition: Ifpack2_Container_def.hpp:183
void applyMV(const mv_type &X, mv_type &Y) const
Wrapper for apply with MVs, used in unit tests (never called by BlockRelaxation)
Definition: Ifpack2_Container_def.hpp:544
virtual void applyInverseJacobi(const mv_type &, mv_type &, SC dampingFactor, bool, int) const
Compute Y := (1 - a) Y + a D^{-1} (X - R*Y).
Definition: Ifpack2_Container_def.hpp:534
typename mv_type::dual_view_type::t_host HostView
Definition: Ifpack2_Container_decl.hpp:139
virtual void weightedApply(ConstHostView X, HostView Y, ConstHostView D, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const
Compute Y := alpha * diag(D) * M^{-1} (diag(D) * X) + beta*Y.
Definition: Ifpack2_Container_def.hpp:739
void resize(size_type new_size, const value_type &x=value_type())
virtual ~ContainerImpl()
Destructor.
Definition: Ifpack2_Container_def.hpp:526
bool IsParallel_
Whether the problem is distributed across multiple MPI processes.
Definition: Ifpack2_Container_decl.hpp:302
static std::string getName()
Definition: Ifpack2_Container_def.hpp:565
TypeTo as(const TypeFrom &t)
Interface for creating and solving a set of local linear problems.
Definition: Ifpack2_Container_decl.hpp:112
Teuchos::ArrayView< const LO > getBlockRows(int blockIndex) const
Local indices of the rows of the input matrix that belong to this block.
Definition: Ifpack2_Container_def.hpp:118
Structure for read-only views of general matrix rows.
Definition: Ifpack2_Container_decl.hpp:330
static std::string getName()
Definition: Ifpack2_Container_def.hpp:190
GO NumGlobalRows_
Number of global rows in input matrix.
Definition: Ifpack2_Container_decl.hpp:306
GO NumGlobalNonzeros_
Number of nonzeros in input matrix.
Definition: Ifpack2_Container_decl.hpp:308
LO NumLocalRows_
Number of local rows in input matrix.
Definition: Ifpack2_Container_decl.hpp:304
Implementation detail of Ifpack2::Container subclasses.
Definition: Ifpack2_Details_MultiVectorLocalGatherScatter.hpp:85
virtual void solveBlock(ConstHostSubviewLocal X, HostSubviewLocal Y, int blockIndex, Teuchos::ETransp mode, const LSC alpha, const LSC beta) const
Definition: Ifpack2_Container_def.hpp:572
std::string toString(const T &t)
virtual void DoGSBlock(ConstHostView X, HostView Y, HostView Y2, HostView Resid, SC dampingFactor, LO i) const
Do one step of Gauss-Seidel on block i (used by DoGaussSeidel and DoSGS)
Definition: Ifpack2_Container_def.hpp:196