43 #ifndef IFPACK2_SPARSECONTAINER_DEF_HPP
44 #define IFPACK2_SPARSECONTAINER_DEF_HPP
51 #include "Teuchos_DefaultSerialComm.hpp"
58 template<
class MatrixType,
class InverseType>
64 ContainerImpl<MatrixType, InverseScalar> (matrix, partitions, pointIndexed),
66 localComm_ (Teuchos::
rcp (new Teuchos::MpiComm<int> (MPI_COMM_SELF)))
68 localComm_ (Teuchos::
rcp (new Teuchos::SerialComm<int> ()))
73 template<
class MatrixType,
class InverseType>
78 template<
class MatrixType,
class InverseType>
85 template<
class MatrixType,
class InverseType>
99 diagBlocks_.assign(this->numBlocks_, Teuchos::null);
100 Inverses_.assign(this->numBlocks_, Teuchos::null);
103 this->extractGraph();
106 for(
int i = 0; i < this->numBlocks_; i++)
108 Inverses_[i]->setParameters(List_);
109 Inverses_[i]->initialize ();
112 this->IsInitialized_ =
true;
113 this->IsComputed_ =
false;
117 template<
class MatrixType,
class InverseType>
124 this->IsComputed_ =
false;
125 if (!this->isInitialized ()) {
130 this->extractValues();
133 for(
int i = 0; i < this->numBlocks_; i++) {
134 Inverses_[i]->compute ();
137 this->IsComputed_ =
true;
141 template<
class MatrixType,
class InverseType>
144 for(
auto inv : Inverses_)
152 template<
class MatrixType,
class InverseType>
159 InverseScalar beta)
const
162 Inverses_[blockIndex]->getDomainMap()->getLocalNumElements() != X.getLocalLength(),
163 std::logic_error,
"Ifpack2::SparseContainer::apply: Inverse_ "
164 "operator and X have incompatible dimensions (" <<
165 Inverses_[blockIndex]->getDomainMap()->getLocalNumElements() <<
" resp. "
166 << X.getLocalLength() <<
"). Please report this bug to "
167 "the Ifpack2 developers.");
169 Inverses_[blockIndex]->getRangeMap()->getLocalNumElements() != Y.getLocalLength(),
170 std::logic_error,
"Ifpack2::SparseContainer::apply: Inverse_ "
171 "operator and Y have incompatible dimensions (" <<
172 Inverses_[blockIndex]->getRangeMap()->getLocalNumElements() <<
" resp. "
173 << Y.getLocalLength() <<
"). Please report this bug to "
174 "the Ifpack2 developers.");
175 Inverses_[blockIndex]->apply(X, Y, mode, alpha, beta);
178 template<
class MatrixType,
class InverseType>
205 size_t numVecs = X.extent(1);
208 ! this->IsComputed_, std::runtime_error,
"Ifpack2::SparseContainer::apply: "
209 "You must have called the compute() method before you may call apply(). "
210 "You may call the apply() method as many times as you want after calling "
211 "compute() once, but you must have called compute() at least once.");
213 X.extent(1) != Y.extent(1), std::runtime_error,
214 "Ifpack2::SparseContainer::apply: X and Y have different numbers of "
215 "vectors. X has " << X.extent(1)
216 <<
", but Y has " << Y.extent(1) <<
".");
222 const LO numRows = this->blockSizes_[blockIndex];
251 for(LO i = 0; i < this->numBlocks_; i++)
252 invX.emplace_back(Inverses_[i]->getDomainMap(), numVecs);
253 for(LO i = 0; i < this->numBlocks_; i++)
254 invY.emplace_back(Inverses_[i]->getDomainMap(), numVecs);
256 inverse_mv_type& X_local = invX[blockIndex];
258 X_local.getLocalLength() != size_t(numRows * this->scalarsPerRow_), std::logic_error,
259 "Ifpack2::SparseContainer::apply: "
260 "X_local has length " << X_local.getLocalLength() <<
", which does "
261 "not match numRows = " << numRows * this->scalarsPerRow_ <<
". Please report this bug to "
262 "the Ifpack2 developers.");
263 const ArrayView<const LO> blockRows = this->getBlockRows(blockIndex);
264 if(this->scalarsPerRow_ == 1)
265 mvgs.gatherMVtoView(X_local, X, blockRows);
267 mvgs.gatherMVtoViewBlock(X_local, X, blockRows, this->scalarsPerRow_);
274 inverse_mv_type& Y_local = invY[blockIndex];
276 Y_local.getLocalLength () != size_t(numRows * this->scalarsPerRow_), std::logic_error,
277 "Ifpack2::SparseContainer::apply: "
278 "Y_local has length " << Y_local.getLocalLength () <<
", which does "
279 "not match numRows = " << numRows * this->scalarsPerRow_ <<
". Please report this bug to "
280 "the Ifpack2 developers.");
282 if(this->scalarsPerRow_ == 1)
283 mvgs.gatherMVtoView(Y_local, Y, blockRows);
285 mvgs.gatherMVtoViewBlock(Y_local, Y, blockRows, this->scalarsPerRow_);
289 this->solveBlockMV(X_local, Y_local, blockIndex, mode,
290 InverseScalar(alpha), InverseScalar(beta));
295 if(this->scalarsPerRow_ == 1)
296 mvgs.scatterMVtoView(Y, Y_local, blockRows);
298 mvgs.scatterMVtoViewBlock(Y, Y_local, blockRows, this->scalarsPerRow_);
302 template<
class MatrixType,
class InverseType>
328 typedef Tpetra::Vector<InverseScalar, InverseLocalOrdinal, InverseGlobalOrdinal, InverseNode> inverse_vector_type;
331 const size_t numVecs = X.extent(1);
334 ! this->IsComputed_, std::runtime_error,
"Ifpack2::SparseContainer::"
335 "weightedApply: You must have called the compute() method before you may "
336 "call apply(). You may call the apply() method as many times as you want "
337 "after calling compute() once, but you must have called compute() at least "
340 X.extent(1) != Y.extent(1), std::runtime_error,
341 "Ifpack2::SparseContainer::weightedApply: X and Y have different numbers "
342 "of vectors. X has " << X.extent(1) <<
", but Y has "
343 << Y.extent(1) <<
".");
347 this->scalarsPerRow_ > 1, std::logic_error,
"Ifpack2::SparseContainer::weightedApply: "
348 "Use of block rows isn't allowed in overlapping Jacobi (the only method that uses weightedApply");
378 const LO numRows = this->blockSizes_[blockIndex];
382 for(LO i = 0; i < this->numBlocks_; i++)
383 invX.emplace_back(Inverses_[i]->getDomainMap(), numVecs);
384 for(LO i = 0; i < this->numBlocks_; i++)
385 invY.emplace_back(Inverses_[i]->getDomainMap(), numVecs);
387 inverse_mv_type& X_local = invX[blockIndex];
388 const ArrayView<const LO> blockRows = this->getBlockRows(blockIndex);
389 mvgs.gatherMVtoView(X_local, X, blockRows);
396 inverse_mv_type Y_local = invY[blockIndex];
397 TEUCHOS_TEST_FOR_EXCEPTION(
398 Y_local.getLocalLength() != size_t(numRows), std::logic_error,
399 "Ifpack2::SparseContainer::weightedApply: "
400 "Y_local has length " << X_local.getLocalLength() <<
", which does "
401 "not match numRows = " << numRows <<
". Please report this bug to "
402 "the Ifpack2 developers.");
403 mvgs.gatherMVtoView(Y_local, Y, blockRows);
415 inverse_vector_type D_local(Inverses_[blockIndex]->getDomainMap());
416 TEUCHOS_TEST_FOR_EXCEPTION(
417 D_local.getLocalLength() != size_t(this->blockSizes_[blockIndex]), std::logic_error,
418 "Ifpack2::SparseContainer::weightedApply: "
419 "D_local has length " << X_local.getLocalLength () <<
", which does "
420 "not match numRows = " << this->blockSizes_[blockIndex] <<
". Please report this bug to "
421 "the Ifpack2 developers.");
422 mvgs.gatherMVtoView(D_local, D, blockRows);
423 inverse_mv_type X_scaled(Inverses_[blockIndex]->getDomainMap(), numVecs);
424 X_scaled.elementWiseMultiply(STS::one(), D_local, X_local, STS::zero());
431 inverse_mv_type* Y_temp;
432 if (InverseScalar(beta) == STS::zero ()) {
435 Y_temp =
new inverse_mv_type(Inverses_[blockIndex]->getRangeMap(), numVecs);
438 Inverses_[blockIndex]->apply(X_scaled, *Y_temp, mode);
444 Y_local.elementWiseMultiply(alpha, D_local, *Y_temp, beta);
445 if(Y_temp != &Y_local)
449 mvgs.scatterMVtoView(Y, Y_local, blockRows);
453 template<
class MatrixType,
class InverseType>
463 template<
class MatrixType,
class InverseType>
466 std::ostringstream oss;
467 oss <<
"\"Ifpack2::SparseContainer\": {";
468 if (this->isInitialized()) {
469 if (this->isComputed()) {
470 oss <<
"status = initialized, computed";
473 oss <<
"status = initialized, not computed";
477 oss <<
"status = not initialized, not computed";
479 for(
int i = 0; i < this->numBlocks_; i++)
481 oss <<
", Block Inverse " << i <<
": {";
482 oss << Inverses_[i]->description();
490 template<
class MatrixType,
class InverseType>
495 os <<
"================================================================================" << endl;
496 os <<
"Ifpack2::SparseContainer" << endl;
497 for(
int i = 0; i < this->numBlocks_; i++)
499 os <<
"Block " << i <<
" rows: = " << this->blockSizes_[i] << endl;
501 os <<
"isInitialized() = " << this->IsInitialized_ << endl;
502 os <<
"isComputed() = " << this->IsComputed_ << endl;
503 os <<
"================================================================================" << endl;
508 template<
class MatrixType,
class InverseType>
524 if(this->scalarsPerRow_ > 1)
526 Array<LO> colToBlockOffset(this->inputBlockMatrix_->getLocalNumCols(), INVALID);
527 for(
int i = 0; i < this->numBlocks_; i++)
530 LO blockStart = this->blockOffsets_[i];
531 LO blockSize = this->blockSizes_[i];
532 LO blockPointSize = this->blockSizes_[i] * this->scalarsPerRow_;
533 LO blockEnd = blockStart + blockSize;
534 ArrayView<const LO> blockRows = this->getBlockRows(i);
538 for(LO j = 0; j < blockSize; j++)
540 LO localCol = this->translateRowToCol(blockRows[j]);
541 colToBlockOffset[localCol] = blockStart + j;
545 Array<size_t> rowEntryCounts(blockPointSize, 0);
548 using inds_type =
typename block_crs_matrix_type::local_inds_host_view_type;
549 using vals_type =
typename block_crs_matrix_type::values_host_view_type;
550 for(LO blockRow = 0; blockRow < blockSize; blockRow++)
555 LO inputRow = this->blockRows_[blockStart + blockRow];
556 this->inputBlockMatrix_->getLocalRowView(inputRow, indices, values);
557 LO numEntries = (LO) indices.size();
558 for(LO br = 0; br < this->bcrsBlockSize_; br++)
560 for(LO k = 0; k < numEntries; k++)
562 LO colOffset = colToBlockOffset[indices[k]];
563 if(blockStart <= colOffset && colOffset < blockEnd)
565 rowEntryCounts[blockRow * this->bcrsBlockSize_ + br] += this->bcrsBlockSize_;
571 RCP<InverseMap> tempMap(
new InverseMap(blockPointSize, 0, this->localComm_));
572 diagBlocks_[i] =
rcp(
new InverseCrs(tempMap, rowEntryCounts));
573 Inverses_[i] =
rcp(
new InverseType(diagBlocks_[i]));
575 for(LO blockRow = 0; blockRow < blockSize; blockRow++)
580 LO inputRow = this->blockRows_[blockStart + blockRow];
581 this->inputBlockMatrix_->getLocalRowView(inputRow, indices, values);
582 LO numEntries = (LO) indices.size();
583 for(LO br = 0; br < this->bcrsBlockSize_; br++)
585 indicesToInsert.
clear();
586 valuesToInsert.
clear();
587 for(LO k = 0; k < numEntries; k++)
589 LO colOffset = colToBlockOffset[indices[k]];
590 if(blockStart <= colOffset && colOffset < blockEnd)
592 LO blockCol = colOffset - blockStart;
594 for(LO bc = 0; bc < this->bcrsBlockSize_; bc++)
596 indicesToInsert.
push_back(blockCol * this->bcrsBlockSize_ + bc);
597 valuesToInsert.
push_back(values[k * this->bcrsBlockSize_ * this->bcrsBlockSize_ + bc * this->bcrsBlockSize_ + br]);
601 InverseGlobalOrdinal rowToInsert = blockRow * this->bcrsBlockSize_ + br;
602 if(indicesToInsert.
size())
603 diagBlocks_[i]->insertGlobalValues(rowToInsert, indicesToInsert(), valuesToInsert());
606 diagBlocks_[i]->fillComplete();
613 Array<LO> colToBlockOffset(this->inputMatrix_->getLocalNumCols() * this->bcrsBlockSize_, INVALID);
614 for(
int i = 0; i < this->numBlocks_; i++)
617 LO blockStart = this->blockOffsets_[i];
618 LO blockSize = this->blockSizes_[i];
619 LO blockEnd = blockStart + blockSize;
620 ArrayView<const LO> blockRows = this->getBlockRows(i);
624 for(LO j = 0; j < blockSize; j++)
627 LO localCol = this->translateRowToCol(blockRows[j]);
628 colToBlockOffset[localCol] = blockStart + j;
631 for(LO j = 0; j < blockSize; j++)
633 rowEntryCounts[j] = this->getInputRowView(this->blockRows_[blockStart + j]).size();
635 RCP<InverseMap> tempMap(
new InverseMap(blockSize, 0, this->localComm_));
636 diagBlocks_[i] =
rcp(
new InverseCrs(tempMap, rowEntryCounts));
637 Inverses_[i] =
rcp(
new InverseType(diagBlocks_[i]));
638 for(LO blockRow = 0; blockRow < blockSize; blockRow++)
640 valuesToInsert.
clear();
641 indicesToInsert.
clear();
643 LO inputSplitRow = this->blockRows_[blockStart + blockRow];
644 auto rowView = this->getInputRowView(inputSplitRow);
645 for(
size_t k = 0; k < rowView.size(); k++)
647 LO colOffset = colToBlockOffset[rowView.ind(k)];
648 if(blockStart <= colOffset && colOffset < blockEnd)
650 LO blockCol = colOffset - blockStart;
652 valuesToInsert.
push_back(rowView.val(k));
655 if(indicesToInsert.
size())
656 diagBlocks_[i]->insertGlobalValues(blockRow, indicesToInsert(), valuesToInsert());
658 diagBlocks_[i]->fillComplete ();
664 template<
class MatrixType,
class InverseType>
665 void SparseContainer<MatrixType,InverseType>::
683 if(this->scalarsPerRow_ > 1)
685 Array<LO> colToBlockOffset(this->inputBlockMatrix_->getLocalNumCols(), INVALID);
686 for(
int i = 0; i < this->numBlocks_; i++)
689 LO blockStart = this->blockOffsets_[i];
690 LO blockSize = this->blockSizes_[i];
691 LO blockPointSize = this->blockSizes_[i] * this->scalarsPerRow_;
692 LO blockEnd = blockStart + blockSize;
693 ArrayView<const LO> blockRows = this->getBlockRows(i);
697 for(LO j = 0; j < blockSize; j++)
699 LO localCol = this->translateRowToCol(blockRows[j]);
700 colToBlockOffset[localCol] = blockStart + j;
704 Array<size_t> rowEntryCounts(blockPointSize, 0);
707 using inds_type =
typename block_crs_matrix_type::local_inds_host_view_type;
708 for(LO blockRow = 0; blockRow < blockSize; blockRow++)
712 LO inputRow = this->blockRows_[blockStart + blockRow];
713 this->inputBlockMatrix_->getGraph()->getLocalRowView(inputRow, indices);
714 LO numEntries = (LO) indices.size();
715 for(LO br = 0; br < this->bcrsBlockSize_; br++)
717 for(LO k = 0; k < numEntries; k++)
719 LO colOffset = colToBlockOffset[indices[k]];
720 if(blockStart <= colOffset && colOffset < blockEnd)
722 rowEntryCounts[blockRow * this->bcrsBlockSize_ + br] += this->bcrsBlockSize_;
728 RCP<InverseMap> tempMap(
new InverseMap(blockPointSize, 0, this->localComm_));
729 auto diagGraph =
rcp(
new InverseGraph(tempMap, rowEntryCounts));
731 for(LO blockRow = 0; blockRow < blockSize; blockRow++)
735 LO inputRow = this->blockRows_[blockStart + blockRow];
736 this->inputBlockMatrix_->getGraph()->getLocalRowView(inputRow, indices);
737 LO numEntries = (LO) indices.size();
738 for(LO br = 0; br < this->bcrsBlockSize_; br++)
740 indicesToInsert.
clear();
741 for(LO k = 0; k < numEntries; k++)
743 LO colOffset = colToBlockOffset[indices[k]];
744 if(blockStart <= colOffset && colOffset < blockEnd)
746 LO blockCol = colOffset - blockStart;
748 for(LO bc = 0; bc < this->bcrsBlockSize_; bc++)
750 indicesToInsert.
push_back(blockCol * this->bcrsBlockSize_ + bc);
754 InverseGlobalOrdinal rowToInsert = blockRow * this->bcrsBlockSize_ + br;
755 if(indicesToInsert.
size())
756 diagGraph->insertGlobalIndices(rowToInsert, indicesToInsert());
759 diagGraph->fillComplete();
762 diagBlocks_[i] =
rcp(
new InverseCrs(diagGraph));
763 Inverses_[i] =
rcp(
new InverseType(diagBlocks_[i]));
770 Array<LO> colToBlockOffset(this->inputMatrix_->getLocalNumCols() * this->bcrsBlockSize_, INVALID);
771 for(
int i = 0; i < this->numBlocks_; i++)
774 LO blockStart = this->blockOffsets_[i];
775 LO blockSize = this->blockSizes_[i];
776 LO blockEnd = blockStart + blockSize;
777 ArrayView<const LO> blockRows = this->getBlockRows(i);
781 for(LO j = 0; j < blockSize; j++)
784 LO localCol = this->translateRowToCol(blockRows[j]);
785 colToBlockOffset[localCol] = blockStart + j;
788 for(LO j = 0; j < blockSize; j++)
790 rowEntryCounts[j] = this->getInputRowView(this->blockRows_[blockStart + j]).size();
792 RCP<InverseMap> tempMap(
new InverseMap(blockSize, 0, this->localComm_));
793 auto diagGraph =
rcp(
new InverseGraph(tempMap, rowEntryCounts));
794 for(LO blockRow = 0; blockRow < blockSize; blockRow++)
796 indicesToInsert.
clear();
798 LO inputSplitRow = this->blockRows_[blockStart + blockRow];
799 auto rowView = this->getInputRowView(inputSplitRow);
800 for(
size_t k = 0; k < rowView.size(); k++)
802 LO colOffset = colToBlockOffset[rowView.ind(k)];
803 if(blockStart <= colOffset && colOffset < blockEnd)
805 LO blockCol = colOffset - blockStart;
809 if(indicesToInsert.
size())
810 diagGraph->insertGlobalIndices(blockRow, indicesToInsert());
812 diagGraph->fillComplete();
815 diagBlocks_[i] =
rcp(
new InverseCrs(diagGraph));
816 Inverses_[i] =
rcp(
new InverseType(diagBlocks_[i]));
822 template<
class MatrixType,
class InverseType>
823 void SparseContainer<MatrixType,InverseType>::
842 if(this->scalarsPerRow_ > 1)
844 Array<LO> colToBlockOffset(this->inputBlockMatrix_->getLocalNumCols(), INVALID);
845 for(
int i = 0; i < this->numBlocks_; i++)
848 LO blockStart = this->blockOffsets_[i];
849 LO blockSize = this->blockSizes_[i];
850 LO blockEnd = blockStart + blockSize;
851 ArrayView<const LO> blockRows = this->getBlockRows(i);
855 for(LO j = 0; j < blockSize; j++)
857 LO localCol = this->translateRowToCol(blockRows[j]);
858 colToBlockOffset[localCol] = blockStart + j;
860 using inds_type =
typename block_crs_matrix_type::local_inds_host_view_type;
861 using vals_type =
typename block_crs_matrix_type::values_host_view_type;
863 diagBlocks_[i]->resumeFill();
864 for(LO blockRow = 0; blockRow < blockSize; blockRow++)
869 LO inputRow = this->blockRows_[blockStart + blockRow];
870 this->inputBlockMatrix_->getLocalRowView(inputRow, indices, values);
871 LO numEntries = (LO) indices.size();
872 for(LO br = 0; br < this->bcrsBlockSize_; br++)
874 indicesToInsert.
clear();
875 valuesToInsert.
clear();
876 for(LO k = 0; k < numEntries; k++)
878 LO colOffset = colToBlockOffset[indices[k]];
879 if(blockStart <= colOffset && colOffset < blockEnd)
881 LO blockCol = colOffset - blockStart;
883 for(LO bc = 0; bc < this->bcrsBlockSize_; bc++)
885 indicesToInsert.
push_back(blockCol * this->bcrsBlockSize_ + bc);
886 valuesToInsert.
push_back(values[k * this->bcrsBlockSize_ * this->bcrsBlockSize_ + bc * this->bcrsBlockSize_ + br]);
890 InverseGlobalOrdinal rowToInsert = blockRow * this->bcrsBlockSize_ + br;
891 if(indicesToInsert.
size())
892 diagBlocks_[i]->replaceGlobalValues(rowToInsert, indicesToInsert(), valuesToInsert());
895 diagBlocks_[i]->fillComplete();
902 Array<LO> colToBlockOffset(this->inputMatrix_->getLocalNumCols() * this->bcrsBlockSize_, INVALID);
903 for(
int i = 0; i < this->numBlocks_; i++)
906 LO blockStart = this->blockOffsets_[i];
907 LO blockSize = this->blockSizes_[i];
908 LO blockEnd = blockStart + blockSize;
909 ArrayView<const LO> blockRows = this->getBlockRows(i);
913 for(LO j = 0; j < blockSize; j++)
916 LO localCol = this->translateRowToCol(blockRows[j]);
917 colToBlockOffset[localCol] = blockStart + j;
919 diagBlocks_[i]->resumeFill();
920 for(LO blockRow = 0; blockRow < blockSize; blockRow++)
922 valuesToInsert.
clear();
923 indicesToInsert.
clear();
925 LO inputSplitRow = this->blockRows_[blockStart + blockRow];
926 auto rowView = this->getInputRowView(inputSplitRow);
927 for(
size_t k = 0; k < rowView.size(); k++)
929 LO colOffset = colToBlockOffset[rowView.ind(k)];
930 if(blockStart <= colOffset && colOffset < blockEnd)
932 LO blockCol = colOffset - blockStart;
934 valuesToInsert.
push_back(rowView.val(k));
937 if(indicesToInsert.
size())
938 diagBlocks_[i]->replaceGlobalValues(blockRow, indicesToInsert, valuesToInsert);
940 diagBlocks_[i]->fillComplete ();
945 template<
typename MatrixType,
typename InverseType>
949 #ifdef HAVE_IFPACK2_AMESOS2
951 if(std::is_same<InverseType, ILUTInverse>::value)
955 else if(std::is_same<InverseType, AmesosInverse>::value)
957 return "SparseAmesos";
961 throw std::logic_error(
"InverseType for SparseContainer must be Ifpack2::ILUT or Details::Amesos2Wrapper");
966 constexpr
bool inverseTypeIsILUT = std::is_same<InverseType, ILUTInverse>::value;
968 "InverseType for SparseContainer must be Ifpack2::ILUT<ROW>");
976 #include "Ifpack2_ILUT.hpp"
977 #ifdef HAVE_IFPACK2_AMESOS2
978 #include "Ifpack2_Details_Amesos2Wrapper.hpp"
985 #ifdef HAVE_IFPACK2_AMESOS2
986 # define IFPACK2_SPARSECONTAINER_INSTANT(S,LO,GO,N) \
987 template class Ifpack2::SparseContainer< Tpetra::RowMatrix<S, LO, GO, N>, \
988 Ifpack2::ILUT<Tpetra::RowMatrix<S,LO,GO,N> > >; \
989 template class Ifpack2::SparseContainer< Tpetra::RowMatrix<S, LO, GO, N>, \
990 Ifpack2::Details::Amesos2Wrapper<Tpetra::RowMatrix<S,LO,GO,N> > >;
992 # define IFPACK2_SPARSECONTAINER_INSTANT(S,LO,GO,N) \
993 template class Ifpack2::SparseContainer< Tpetra::RowMatrix<S,LO,GO,N>, \
994 Ifpack2::ILUT<Tpetra::RowMatrix<S, LO, GO, N> > >;
996 #endif // IFPACK2_SPARSECONTAINER_DEF_HPP
static std::string getName()
Get the name of this container type for Details::constructContainer()
Definition: Ifpack2_SparseContainer_def.hpp:946
SparseContainer(const Teuchos::RCP< const row_matrix_type > &matrix, const Teuchos::Array< Teuchos::Array< LO > > &partitions, const Teuchos::RCP< const import_type > &importer, bool pointIndexed)
Constructor.
Definition: Ifpack2_SparseContainer_def.hpp:60
Store and solve a local sparse linear problem.
Definition: Ifpack2_SparseContainer_decl.hpp:134
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to the given FancyOStream.
Definition: Ifpack2_SparseContainer_def.hpp:491
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 initialize()
Do all set-up operations that only require matrix structure.
Definition: Ifpack2_SparseContainer_def.hpp:86
ILUT (incomplete LU factorization with threshold) of a Tpetra sparse matrix.
Definition: Ifpack2_ILUT_decl.hpp:93
virtual std::ostream & print(std::ostream &os) const
Print information about this object to the given output stream.
Definition: Ifpack2_SparseContainer_def.hpp:454
virtual std::string description() const
A one-line description of this object.
Definition: Ifpack2_SparseContainer_def.hpp:464
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Wrapper class for direct solvers in Amesos2.
Definition: Ifpack2_Details_Amesos2Wrapper_decl.hpp:102
typename mv_type::dual_view_type::t_host HostView
Definition: Ifpack2_Container_decl.hpp:139
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_SparseContainer_def.hpp:180
basic_FancyOStream & setOutputToRootOnly(const int rootRank)
virtual void weightedApply(ConstHostView X, HostView Y, ConstHostView W, 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_SparseContainer_def.hpp:304
void push_back(const value_type &x)
virtual ~SparseContainer()
Destructor (declared virtual for memory safety of derived classes).
Definition: Ifpack2_SparseContainer_def.hpp:75
TypeTo as(const TypeFrom &t)
virtual void setParameters(const Teuchos::ParameterList &List)
Set all necessary parameters.
Definition: Ifpack2_SparseContainer_def.hpp:79
virtual void compute()
Initialize and compute all blocks.
Definition: Ifpack2_SparseContainer_def.hpp:118
Ifpack2::SparseContainer class declaration.
void clearBlocks()
Definition: Ifpack2_SparseContainer_def.hpp:142
Implementation detail of Ifpack2::Container subclasses.
Definition: Ifpack2_Details_MultiVectorLocalGatherScatter.hpp:85