43 #ifndef IFPACK2_LOCALSPARSETRIANGULARSOLVER_DEF_HPP 44 #define IFPACK2_LOCALSPARSETRIANGULARSOLVER_DEF_HPP 46 #include "Tpetra_CrsMatrix.hpp" 47 #include "Tpetra_Core.hpp" 48 #include "Teuchos_StandardParameterEntryValidators.hpp" 49 #include "Tpetra_Details_determineLocalTriangularStructure.hpp" 51 #ifdef HAVE_IFPACK2_SHYLU_NODEHTS 52 # include "shylu_hts.hpp" 58 struct TrisolverType {
64 static void loadPLTypeOption (Teuchos::Array<std::string>& type_strs, Teuchos::Array<Enum>& type_enums) {
66 type_strs[0] =
"Internal";
69 type_enums[0] = Internal;
75 template<
class MatrixType>
76 class LocalSparseTriangularSolver<MatrixType>::HtsImpl {
78 typedef Tpetra::CrsMatrix<scalar_type, local_ordinal_type, global_ordinal_type, node_type>
crs_matrix_type;
81 #ifdef HAVE_IFPACK2_SHYLU_NODEHTS 82 Timpl_ = Teuchos::null;
83 levelset_block_size_ = 1;
89 #ifdef HAVE_IFPACK2_SHYLU_NODEHTS 90 const char* block_size_s =
"trisolver: block size";
91 if (pl.isParameter(block_size_s)) {
92 TEUCHOS_TEST_FOR_EXCEPT_MSG( ! pl.isType<
int>(block_size_s),
93 "The parameter \"" << block_size_s <<
"\" must be of type int.");
94 levelset_block_size_ = pl.get<
int>(block_size_s);
96 if (levelset_block_size_ < 1)
97 levelset_block_size_ = 1;
106 #ifdef HAVE_IFPACK2_SHYLU_NODEHTS 108 transpose_ = conjugate_ =
false;
115 #ifdef HAVE_IFPACK2_SHYLU_NODEHTS 116 using Teuchos::ArrayRCP;
118 Teuchos::ArrayRCP<const size_t> rowptr;
119 Teuchos::ArrayRCP<const local_ordinal_type> colidx;
120 Teuchos::ArrayRCP<const scalar_type> val;
121 T_in.getAllValues(rowptr, colidx, val);
123 Teuchos::RCP<HtsCrsMatrix> T_hts = Teuchos::rcpWithDealloc(
124 HTST::make_CrsMatrix(rowptr.size() - 1,
125 rowptr.getRawPtr(), colidx.getRawPtr(), val.getRawPtr(),
126 transpose_, conjugate_),
127 HtsCrsMatrixDeleter());
129 if (Teuchos::nonnull(Timpl_)) {
131 HTST::reprocess_numeric(Timpl_.get(), T_hts.get());
134 if (T_in.getCrsGraph().is_null()) {
135 if (Teuchos::nonnull(out))
136 *out <<
"HTS compute failed because T_in.getCrsGraph().is_null().\n";
139 if ( ! T_in.getCrsGraph()->isSorted()) {
140 if (Teuchos::nonnull(out))
141 *out <<
"HTS compute failed because ! T_in.getCrsGraph().isSorted().\n";
144 if ( ! T_in.isStorageOptimized()) {
145 if (Teuchos::nonnull(out))
146 *out <<
"HTS compute failed because ! T_in.isStorageOptimized().\n";
150 typename HTST::PreprocessArgs args;
151 args.T = T_hts.get();
154 args.nthreads = omp_get_max_threads();
158 args.save_for_reprocess =
true;
159 typename HTST::Options opts;
160 opts.levelset_block_size = levelset_block_size_;
161 args.options = &opts;
164 Timpl_ = Teuchos::rcpWithDealloc(HTST::preprocess(args), TImplDeleter());
165 }
catch (
const std::exception& e) {
166 if (Teuchos::nonnull(out))
167 *out <<
"HTS preprocess threw: " << e.what() <<
"\n";
176 #ifdef HAVE_IFPACK2_SHYLU_NODEHTS 177 return Teuchos::nonnull(Timpl_);
184 void localApply (
const MV& X, MV& Y,
185 const Teuchos::ETransp ,
191 #ifdef HAVE_IFPACK2_SHYLU_NODEHTS 192 const auto& X_view = X.template getLocalView<Kokkos::HostSpace>();
193 const auto& Y_view = Y.template getLocalView<Kokkos::HostSpace>();
195 HTST::reset_max_nrhs(Timpl_.get(), X_view.extent(1));
197 HTST::solve_omp(Timpl_.get(),
199 reinterpret_cast<const scalar_type*
>(X_view.data()),
208 #ifdef HAVE_IFPACK2_SHYLU_NODEHTS 209 typedef ::Experimental::HTS<local_ordinal_type, size_t, scalar_type> HTST;
210 typedef typename HTST::Impl TImpl;
211 typedef typename HTST::CrsMatrix HtsCrsMatrix;
213 struct TImplDeleter {
214 void free (TImpl* impl) {
215 HTST::delete_Impl(impl);
219 struct HtsCrsMatrixDeleter {
220 void free (HtsCrsMatrix* T) {
221 HTST::delete_CrsMatrix(T);
225 Teuchos::RCP<TImpl> Timpl_;
226 bool transpose_, conjugate_;
227 int levelset_block_size_;
231 template<
class MatrixType>
239 if (! A.is_null ()) {
240 Teuchos::RCP<const crs_matrix_type> A_crs =
242 TEUCHOS_TEST_FOR_EXCEPTION
243 (A_crs.is_null (), std::invalid_argument,
244 "Ifpack2::LocalSparseTriangularSolver constructor: " 245 "The input matrix A is not a Tpetra::CrsMatrix.");
250 template<
class MatrixType>
253 const Teuchos::RCP<Teuchos::FancyOStream>& out) :
258 if (! out_.is_null ()) {
259 *out_ <<
">>> DEBUG Ifpack2::LocalSparseTriangularSolver constructor" 264 if (! A.is_null ()) {
265 Teuchos::RCP<const crs_matrix_type> A_crs =
267 TEUCHOS_TEST_FOR_EXCEPTION
268 (A_crs.is_null (), std::invalid_argument,
269 "Ifpack2::LocalSparseTriangularSolver constructor: " 270 "The input matrix A is not a Tpetra::CrsMatrix.");
275 template<
class MatrixType>
282 template<
class MatrixType>
288 if (! out_.is_null ()) {
289 *out_ <<
">>> DEBUG Ifpack2::LocalSparseTriangularSolver constructor" 294 template<
class MatrixType>
297 isInitialized_ =
false;
299 reverseStorage_ =
false;
300 isInternallyChanged_ =
false;
304 initializeTime_ = 0.0;
311 template<
class MatrixType>
316 template<
class MatrixType>
322 using Teuchos::ParameterList;
323 using Teuchos::Array;
325 Details::TrisolverType::Enum trisolverType = Details::TrisolverType::Internal;
327 static const char typeName[] =
"trisolver: type";
329 if ( ! pl.isType<std::string>(typeName))
break;
332 Array<std::string> trisolverTypeStrs;
333 Array<Details::TrisolverType::Enum> trisolverTypeEnums;
334 Details::TrisolverType::loadPLTypeOption (trisolverTypeStrs, trisolverTypeEnums);
335 Teuchos::StringToIntegralParameterEntryValidator<Details::TrisolverType::Enum>
336 s2i(trisolverTypeStrs (), trisolverTypeEnums (), typeName,
false);
338 trisolverType = s2i.getIntegralValue(pl.get<std::string>(typeName));
341 if (trisolverType == Details::TrisolverType::HTS) {
342 htsImpl_ = Teuchos::rcp (
new HtsImpl ());
343 htsImpl_->setParameters (pl);
346 if (pl.isParameter(
"trisolver: reverse U"))
347 reverseStorage_ = pl.get<
bool>(
"trisolver: reverse U");
349 TEUCHOS_TEST_FOR_EXCEPTION
350 (reverseStorage_ && trisolverType == Details::TrisolverType::HTS,
351 std::logic_error,
"Ifpack2::LocalSparseTriangularSolver::setParameters: " 352 "You are not allowed to enable both HTS and the \"trisolver: reverse U\" " 353 "options. See GitHub issue #2647.");
356 template<
class MatrixType>
363 using local_matrix_type =
typename crs_matrix_type::local_matrix_type;
365 const char prefix[] =
"Ifpack2::LocalSparseTriangularSolver::initialize: ";
366 if (! out_.is_null ()) {
367 *out_ <<
">>> DEBUG " << prefix << std::endl;
370 TEUCHOS_TEST_FOR_EXCEPTION
371 (A_.is_null (), std::runtime_error, prefix <<
"You must call " 372 "setMatrix() with a nonnull input matrix before you may call " 373 "initialize() or compute().");
374 if (A_crs_.is_null ()) {
376 TEUCHOS_TEST_FOR_EXCEPTION
377 (A_crs.get () ==
nullptr, std::invalid_argument,
378 prefix <<
"The input matrix A is not a Tpetra::CrsMatrix.");
381 auto G = A_crs_->getGraph ();
382 TEUCHOS_TEST_FOR_EXCEPTION
383 (G.is_null (), std::logic_error, prefix <<
"A_ and A_crs_ are nonnull, " 384 "but A_crs_'s RowGraph G is null. " 385 "Please report this bug to the Ifpack2 developers.");
389 TEUCHOS_TEST_FOR_EXCEPTION
390 (! G->isFillComplete (), std::runtime_error,
"If you call this method, " 391 "the matrix's graph must be fill complete. It is not.");
400 if (reverseStorage_ && A_crs_->isUpperTriangularImpl() && htsImpl_.is_null()) {
402 auto Alocal = A_crs_->getLocalMatrix();
403 auto ptr = Alocal.graph.row_map;
404 auto ind = Alocal.graph.entries;
405 auto val = Alocal.values;
407 auto numRows = Alocal.numRows();
408 auto numCols = Alocal.numCols();
409 auto numNnz = Alocal.nnz();
411 typename decltype(ptr)::non_const_type newptr (
"ptr", ptr.extent (0));
412 typename decltype(ind)::non_const_type newind (
"ind", ind.extent (0));
413 decltype(val) newval (
"val", val.extent (0));
416 crs_matrix_type::execution_space::fence();
419 auto A_r = Alocal.row(numRows-1 - row);
421 auto numEnt = A_r.length;
423 newind(rowStart + k) = numCols-1 - A_r.colidx(numEnt-1 - k);
424 newval(rowStart + k) = A_r.value (numEnt-1 - k);
427 newptr(row+1) = rowStart;
429 crs_matrix_type::execution_space::fence();
432 using map_type =
typename crs_matrix_type::map_type;
433 Teuchos::RCP<map_type> newRowMap, newColMap;
436 auto rowMap = A_->getRowMap();
437 auto numElems = rowMap->getNodeNumElements();
438 auto rowElems = rowMap->getNodeElementList();
440 Teuchos::Array<global_ordinal_type> newRowElems(rowElems.size());
441 for (
size_t i = 0; i < numElems; i++)
442 newRowElems[i] = rowElems[numElems-1 - i];
444 newRowMap = Teuchos::rcp(
new map_type(rowMap->getGlobalNumElements(), newRowElems, rowMap->getIndexBase(), rowMap->getComm()));
448 auto colMap = A_->getColMap();
449 auto numElems = colMap->getNodeNumElements();
450 auto colElems = colMap->getNodeElementList();
452 Teuchos::Array<global_ordinal_type> newColElems(colElems.size());
453 for (
size_t i = 0; i < numElems; i++)
454 newColElems[i] = colElems[numElems-1 - i];
456 newColMap = Teuchos::rcp(
new map_type(colMap->getGlobalNumElements(), newColElems, colMap->getIndexBase(), colMap->getComm()));
460 local_matrix_type newLocalMatrix(
"Upermuted", numRows, numCols, numNnz, newval, newptr, newind);
462 A_crs_ = Teuchos::rcp(
new crs_matrix_type(newLocalMatrix, newRowMap, newColMap, A_crs_->getDomainMap(), A_crs_->getRangeMap()));
464 isInternallyChanged_ =
true;
467 if (Teuchos::nonnull (htsImpl_))
469 htsImpl_->initialize (*A_crs_);
470 isInternallyChanged_ =
true;
473 auto lclMatrix = A_crs_->getLocalMatrix ();
474 auto lclRowMap = A_crs_->getRowMap ()->getLocalMap ();
475 auto lclColMap = A_crs_->getColMap ()->getLocalMap ();
476 using Tpetra::Details::determineLocalTriangularStructure;
478 constexpr
bool ignoreMapsForTriangularStructure =
true;
480 determineLocalTriangularStructure (lclMatrix.graph, lclRowMap, lclColMap,
481 ignoreMapsForTriangularStructure);
483 const LO lclNumRows = lclRowMap.getNodeNumElements ();
484 const LO lclNumCols = lclColMap.getNodeNumElements ();
492 this->diag_ = (result.diagCount < lclNumRows) ?
"U" :
"N";
493 this->uplo_ = result.couldBeLowerTriangular ?
"L" :
494 (result.couldBeUpperTriangular ?
"U" :
"N");
496 isInitialized_ =
true;
500 template<
class MatrixType>
505 const char prefix[] =
"Ifpack2::LocalSparseTriangularSolver::compute: ";
506 if (! out_.is_null ()) {
507 *out_ <<
">>> DEBUG " << prefix << std::endl;
510 TEUCHOS_TEST_FOR_EXCEPTION
511 (A_.is_null (), std::runtime_error, prefix <<
"You must call " 512 "setMatrix() with a nonnull input matrix before you may call " 513 "initialize() or compute().");
514 TEUCHOS_TEST_FOR_EXCEPTION
515 (A_crs_.is_null (), std::logic_error, prefix <<
"A_ is nonnull, but " 516 "A_crs_ is null. Please report this bug to the Ifpack2 developers.");
518 TEUCHOS_TEST_FOR_EXCEPTION
519 (! A_crs_->isFillComplete (), std::runtime_error,
"If you call this " 520 "method, the matrix must be fill complete. It is not.");
522 if (! isInitialized_) {
525 TEUCHOS_TEST_FOR_EXCEPTION
526 (! isInitialized_, std::logic_error, prefix <<
"initialize() should have " 527 "been called by this point, but isInitialized_ is false. " 528 "Please report this bug to the Ifpack2 developers.");
530 if (Teuchos::nonnull (htsImpl_))
531 htsImpl_->compute (*A_crs_, out_);
537 template<
class MatrixType>
543 Teuchos::ETransp mode,
549 using Teuchos::rcpFromRef;
551 typedef Teuchos::ScalarTraits<ST> STS;
552 const char prefix[] =
"Ifpack2::LocalSparseTriangularSolver::apply: ";
553 if (! out_.is_null ()) {
554 *out_ <<
">>> DEBUG " << prefix;
555 if (A_crs_.is_null ()) {
556 *out_ <<
"A_crs_ is null!" << std::endl;
559 Teuchos::RCP<const crs_matrix_type> A_crs =
561 const std::string uplo = this->uplo_;
562 const std::string trans = (mode == Teuchos::CONJ_TRANS) ?
"C" :
563 (mode == Teuchos::TRANS ?
"T" :
"N");
564 const std::string diag = this->diag_;
565 *out_ <<
"uplo=\"" << uplo
566 <<
"\", trans=\"" << trans
567 <<
"\", diag=\"" << diag <<
"\"" << std::endl;
571 TEUCHOS_TEST_FOR_EXCEPTION
572 (! isComputed (), std::runtime_error, prefix <<
"If compute() has not yet " 573 "been called, or if you have changed the matrix via setMatrix(), you must " 574 "call compute() before you may call this method.");
577 TEUCHOS_TEST_FOR_EXCEPTION
578 (A_.is_null (), std::logic_error, prefix <<
"A_ is null. " 579 "Please report this bug to the Ifpack2 developers.");
580 TEUCHOS_TEST_FOR_EXCEPTION
581 (A_crs_.is_null (), std::logic_error, prefix <<
"A_crs_ is null. " 582 "Please report this bug to the Ifpack2 developers.");
585 TEUCHOS_TEST_FOR_EXCEPTION
586 (! A_crs_->isFillComplete (), std::runtime_error,
"If you call this " 587 "method, the matrix must be fill complete. It is not. This means that " 588 " you must have called resumeFill() on the matrix before calling apply(). " 589 "This is NOT allowed. Note that this class may use the matrix's data in " 590 "place without copying it. Thus, you cannot change the matrix and expect " 591 "the solver to stay the same. If you have changed the matrix, first call " 592 "fillComplete() on it, then call compute() on this object, before you call" 593 " apply(). You do NOT need to call setMatrix, as long as the matrix " 594 "itself (that is, its address in memory) is the same.");
596 auto G = A_crs_->getGraph ();
597 TEUCHOS_TEST_FOR_EXCEPTION
598 (G.is_null (), std::logic_error, prefix <<
"A_ and A_crs_ are nonnull, " 599 "but A_crs_'s RowGraph G is null. " 600 "Please report this bug to the Ifpack2 developers.");
601 auto importer = G->getImporter ();
602 auto exporter = G->getExporter ();
604 if (! importer.is_null ()) {
605 if (X_colMap_.is_null () || X_colMap_->getNumVectors () != X.getNumVectors ()) {
606 X_colMap_ = rcp (
new MV (importer->getTargetMap (), X.getNumVectors ()));
609 X_colMap_->putScalar (STS::zero ());
614 X_colMap_->doImport (X, *importer, Tpetra::ZERO);
616 RCP<const MV> X_cur = importer.is_null () ? rcpFromRef (X) :
617 Teuchos::rcp_const_cast<
const MV> (X_colMap_);
619 if (! exporter.is_null ()) {
620 if (Y_rowMap_.is_null () || Y_rowMap_->getNumVectors () != Y.getNumVectors ()) {
621 Y_rowMap_ = rcp (
new MV (exporter->getSourceMap (), Y.getNumVectors ()));
624 Y_rowMap_->putScalar (STS::zero ());
626 Y_rowMap_->doExport (Y, *importer, Tpetra::ADD);
628 RCP<MV> Y_cur = exporter.is_null () ? rcpFromRef (Y) : Y_rowMap_;
630 localApply (*X_cur, *Y_cur, mode, alpha, beta);
632 if (! exporter.is_null ()) {
633 Y.putScalar (STS::zero ());
634 Y.doExport (*Y_cur, *exporter, Tpetra::ADD);
640 template<
class MatrixType>
645 const Teuchos::ETransp mode)
const 647 using Teuchos::CONJ_TRANS;
648 using Teuchos::NO_TRANS;
649 using Teuchos::TRANS;
650 typedef Kokkos::HostSpace host_memory_space;
651 using device_type =
typename MV::device_type;
652 using dev_memory_space =
typename device_type::memory_space;
653 const char tfecfFuncName[] =
"localTriangularSolve: ";
655 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
656 (! A_crs_->isFillComplete (), std::runtime_error,
657 "The matrix is not fill complete.");
658 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
659 (! X.isConstantStride () || ! Y.isConstantStride (), std::invalid_argument,
660 "X and Y must be constant stride.");
661 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
662 ( A_crs_->getNodeNumRows() > 0 && this->uplo_ ==
"N", std::runtime_error,
663 "The matrix is neither upper triangular or lower triangular. " 664 "You may only call this method if the matrix is triangular. " 665 "Remember that this is a local (per MPI process) property, and that " 666 "Tpetra only knows how to do a local (per process) triangular solve.");
667 using STS = Teuchos::ScalarTraits<scalar_type>;
668 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
669 (STS::isComplex && mode == TRANS, std::logic_error,
"This method does " 670 "not currently support non-conjugated transposed solve (mode == " 671 "Teuchos::TRANS) for complex scalar types.");
680 const std::string uplo = this->uplo_;
681 const std::string trans = (mode == Teuchos::CONJ_TRANS) ?
"C" :
682 (mode == Teuchos::TRANS ?
"T" :
"N");
683 const std::string diag = this->diag_;
684 auto A_lcl = this->A_crs_->getLocalMatrix ();
691 X.template sync<host_memory_space> ();
692 const_cast<MV&
> (Y).
template sync<host_memory_space> ();
693 X.template modify<host_memory_space> ();
695 if (X.isConstantStride () && Y.isConstantStride ()) {
696 auto X_lcl = X.template getLocalView<host_memory_space> ();
697 auto Y_lcl = Y.template getLocalView<host_memory_space> ();
698 KokkosSparse::trsv (uplo.c_str (), trans.c_str (), diag.c_str (),
699 A_lcl, Y_lcl, X_lcl);
702 const size_t numVecs =
703 std::min (X.getNumVectors (), Y.getNumVectors ());
704 for (
size_t j = 0; j < numVecs; ++j) {
705 auto X_j = X.getVector (j);
706 auto Y_j = X.getVector (j);
707 auto X_lcl = X_j->template getLocalView<host_memory_space> ();
708 auto Y_lcl = Y_j->template getLocalView<host_memory_space> ();
709 KokkosSparse::trsv (uplo.c_str (), trans.c_str (),
710 diag.c_str (), A_lcl, Y_lcl, X_lcl);
714 X.template sync<dev_memory_space> ();
715 const_cast<MV&
> (Y).
template sync<dev_memory_space> ();
718 template<
class MatrixType>
720 LocalSparseTriangularSolver<MatrixType>::
721 localApply (
const MV& X,
723 const Teuchos::ETransp mode,
724 const scalar_type& alpha,
725 const scalar_type& beta)
const 727 if (mode == Teuchos::NO_TRANS && Teuchos::nonnull (htsImpl_) &&
728 htsImpl_->isComputed ()) {
729 htsImpl_->localApply (X, Y, mode, alpha, beta);
734 typedef scalar_type ST;
735 typedef Teuchos::ScalarTraits<ST> STS;
737 if (beta == STS::zero ()) {
738 if (alpha == STS::zero ()) {
739 Y.putScalar (STS::zero ());
742 this->localTriangularSolve (X, Y, mode);
743 if (alpha != STS::one ()) {
749 if (alpha == STS::zero ()) {
753 MV Y_tmp (Y, Teuchos::Copy);
754 this->localTriangularSolve (X, Y_tmp, mode);
755 Y.update (alpha, Y_tmp, beta);
761 template <
class MatrixType>
765 return numInitialize_;
768 template <
class MatrixType>
775 template <
class MatrixType>
782 template <
class MatrixType>
786 return initializeTime_;
789 template<
class MatrixType>
796 template<
class MatrixType>
803 template <
class MatrixType>
808 std::ostringstream os;
813 os <<
"\"Ifpack2::LocalSparseTriangularSolver\": {";
814 if (this->getObjectLabel () !=
"") {
815 os <<
"Label: \"" << this->getObjectLabel () <<
"\", ";
817 os <<
"Initialized: " << (isInitialized () ?
"true" :
"false") <<
", " 818 <<
"Computed: " << (isComputed () ?
"true" :
"false") <<
", ";
821 os <<
"Matrix: null";
824 os <<
"Matrix: not null" 825 <<
", Global matrix dimensions: [" 826 << A_->getGlobalNumRows () <<
", " 827 << A_->getGlobalNumCols () <<
"]";
830 if (Teuchos::nonnull (htsImpl_))
831 os <<
", HTS computed: " << (htsImpl_->isComputed () ?
"true" :
"false");
837 template <
class MatrixType>
840 const Teuchos::EVerbosityLevel verbLevel)
const 845 const Teuchos::EVerbosityLevel vl
846 = (verbLevel == Teuchos::VERB_DEFAULT) ? Teuchos::VERB_LOW : verbLevel;
848 if (vl != Teuchos::VERB_NONE) {
853 auto comm = A_.is_null () ?
854 Tpetra::getDefaultComm () :
859 if (! comm.is_null () && comm->getRank () == 0) {
861 Teuchos::OSTab tab0 (out);
864 out <<
"\"Ifpack2::LocalSparseTriangularSolver\":" << endl;
865 Teuchos::OSTab tab1 (out);
866 out <<
"Scalar: " << Teuchos::TypeNameTraits<scalar_type>::name () << endl
867 <<
"LocalOrdinal: " << Teuchos::TypeNameTraits<local_ordinal_type>::name () << endl
868 <<
"GlobalOrdinal: " << Teuchos::TypeNameTraits<global_ordinal_type>::name () << endl
869 <<
"Node: " << Teuchos::TypeNameTraits<node_type>::name () << endl;
874 template <
class MatrixType>
875 Teuchos::RCP<const typename LocalSparseTriangularSolver<MatrixType>::map_type>
879 TEUCHOS_TEST_FOR_EXCEPTION
880 (A_.is_null (), std::runtime_error,
881 "Ifpack2::LocalSparseTriangularSolver::getDomainMap: " 882 "The matrix is null. Please call setMatrix() with a nonnull input " 883 "before calling this method.");
884 return A_->getDomainMap ();
887 template <
class MatrixType>
888 Teuchos::RCP<const typename LocalSparseTriangularSolver<MatrixType>::map_type>
892 TEUCHOS_TEST_FOR_EXCEPTION
893 (A_.is_null (), std::runtime_error,
894 "Ifpack2::LocalSparseTriangularSolver::getRangeMap: " 895 "The matrix is null. Please call setMatrix() with a nonnull input " 896 "before calling this method.");
897 return A_->getRangeMap ();
900 template<
class MatrixType>
902 setMatrix (
const Teuchos::RCP<const row_matrix_type>& A)
904 const char prefix[] =
"Ifpack2::LocalSparseTriangularSolver::setMatrix: ";
910 if (A.getRawPtr () != A_.getRawPtr () || isInternallyChanged_) {
912 TEUCHOS_TEST_FOR_EXCEPTION
913 (! A.is_null () && A->getComm ()->getSize () == 1 &&
914 A->getNodeNumRows () != A->getNodeNumCols (),
915 std::runtime_error, prefix <<
"If A's communicator only contains one " 916 "process, then A must be square. Instead, you provided a matrix A with " 917 << A->getNodeNumRows () <<
" rows and " << A->getNodeNumCols ()
923 isInitialized_ =
false;
929 A_crs_ = Teuchos::null;
933 Teuchos::RCP<const crs_matrix_type> A_crs =
935 TEUCHOS_TEST_FOR_EXCEPTION
936 (A_crs.is_null (), std::invalid_argument, prefix <<
937 "The input matrix A is not a Tpetra::CrsMatrix.");
942 if (Teuchos::nonnull (htsImpl_))
949 #define IFPACK2_LOCALSPARSETRIANGULARSOLVER_INSTANT(S,LO,GO,N) \ 950 template class Ifpack2::LocalSparseTriangularSolver< Tpetra::RowMatrix<S, LO, GO, N> >; 952 #endif // IFPACK2_LOCALSPARSETRIANGULARSOLVER_DEF_HPP 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
Apply the preconditioner to X, and put the result in Y.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:539
int getNumInitialize() const
Return the number of calls to initialize().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:764
Teuchos::RCP< const map_type > getDomainMap() const
The domain of this operator.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:877
Teuchos::RCP< const map_type > getRangeMap() const
The range of this operator.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:890
int getNumCompute() const
Return the number of calls to compute().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:771
double getApplyTime() const
Return the time spent in apply().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:799
void compute()
"Numeric" phase of setup
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:503
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print this object with given verbosity to the given output stream.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:839
MatrixType::global_ordinal_type global_ordinal_type
Type of the global indices of the input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:93
Ifpack2 implementation details.
MatrixType::node_type node_type
Node type of the input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:95
bool isComputed() const
Return true if compute() has been called.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:209
double getComputeTime() const
Return the time spent in compute().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:792
int getNumApply() const
Return the number of calls to apply().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:778
"Preconditioner" that solves local sparse triangular systems.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:77
void initialize()
"Symbolic" phase of setup
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:359
Tpetra::Map< local_ordinal_type, global_ordinal_type, node_type > map_type
Specialization of Tpetra::Map used by this class.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:100
virtual void setMatrix(const Teuchos::RCP< const row_matrix_type > &A)
Set this preconditioner's matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:902
double getInitializeTime() const
Return the time spent in initialize().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:785
MatrixType::local_ordinal_type local_ordinal_type
Type of the local indices of the input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:91
Tpetra::CrsMatrix< scalar_type, local_ordinal_type, global_ordinal_type, node_type > crs_matrix_type
Specialization of Tpetra::CrsMatrix used by this class.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:106
void setParameters(const Teuchos::ParameterList ¶ms)
Set this object's parameters.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:319
MatrixType::scalar_type scalar_type
Type of the entries of the input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:89
std::string description() const
A one-line description of this object.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:806
Preconditioners and smoothers for Tpetra sparse matrices.
Definition: Ifpack2_AdditiveSchwarz_decl.hpp:72
virtual ~LocalSparseTriangularSolver()
Destructor (virtual for memory safety).
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:313
LocalSparseTriangularSolver()
Constructor that takes no input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:277