12 #ifndef __IFPACK2_FILU_DEF_HPP__
13 #define __IFPACK2_FILU_DEF_HPP__
15 #include "Ifpack2_Details_Filu_decl.hpp"
17 #include "Ifpack2_Details_getCrsMatrix.hpp"
18 #include <Kokkos_Timer.hpp>
19 #include <shylu_fastilu.hpp>
24 template <
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal,
typename Node,
bool BlockCrsEnabled>
27 :
FastILU_Base<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A) {}
29 template <
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal,
typename Node,
bool BlockCrsEnabled>
32 return localPrec_->getNFact();
35 template <
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal,
typename Node,
bool BlockCrsEnabled>
38 return localPrec_->getSpTrsvType();
41 template <
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal,
typename Node,
bool BlockCrsEnabled>
44 return localPrec_->getNTrisol();
47 template <
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal,
typename Node,
bool BlockCrsEnabled>
50 localPrec_->checkILU();
53 template <
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal,
typename Node,
bool BlockCrsEnabled>
56 localPrec_->checkIC();
59 template <
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal,
typename Node,
bool BlockCrsEnabled>
62 auto nRows = this->mat_->getLocalNumRows();
63 auto& p = this->params_;
64 auto matCrs = Ifpack2::Details::getCrsMatrix(this->mat_);
66 if (p.blockCrsSize > 1 && !BlockCrsEnabled) {
67 throw std::runtime_error(
"Must use prec type FAST_ILU_B if you want blockCrs support");
70 bool skipSortMatrix = !matCrs.is_null() && matCrs->getCrsGraph()->isSorted() &&
73 Teuchos::rcp(
new LocalFILU(skipSortMatrix, this->localRowPtrs_, this->localColInds_, this->localValues_, nRows, p.sptrsv_algo,
74 p.nFact, p.nTrisol, p.level, p.omega, p.shift, p.guessFlag ? 1 : 0, p.blockSizeILU, p.blockSize,
77 #ifdef HAVE_IFPACK2_METIS
79 localPrec_->setMetisPerm(this->metis_perm_, this->metis_iperm_);
83 localPrec_->initialize();
84 this->initTime_ = localPrec_->getInitializeTime();
87 template <
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal,
typename Node,
bool BlockCrsEnabled>
91 localPrec_->setValues(this->localValues_);
92 localPrec_->compute();
93 this->computeTime_ = localPrec_->getComputeTime();
96 template <
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal,
typename Node,
bool BlockCrsEnabled>
99 localPrec_->apply(x, y);
102 this->applyTime_ += localPrec_->getApplyTime();
105 template <
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal,
typename Node,
bool BlockCrsEnabled>
111 #define IFPACK2_DETAILS_FILU_INSTANT(S, L, G, N) \
112 template class Ifpack2::Details::Filu<S, L, G, N, false>; \
113 template class Ifpack2::Details::Filu<S, L, G, N, true>;
std::string getSpTrsvType() const
Get the name of triangular solve algorithm.
Definition: Ifpack2_Details_Filu_def.hpp:37
int getSweeps() const
Get the sweeps ("nFact") from localPrec_.
Definition: Ifpack2_Details_Filu_def.hpp:31
void computeLocalPrec()
Get values array from the matrix and then call compute() on the underlying preconditioner.
Definition: Ifpack2_Details_Filu_def.hpp:89
std::string getName() const
Get the name of the underlying preconditioner ("Filu", "Fildl" or "Fic")
Definition: Ifpack2_Details_Filu_def.hpp:107
void checkLocalILU() const
Verify and print debug info about the internal ILU preconditioner.
Definition: Ifpack2_Details_Filu_def.hpp:49
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void initLocalPrec()
Construct the underlying preconditioner (localPrec_) using given params and then call localPrec_->ini...
Definition: Ifpack2_Details_Filu_def.hpp:61
int getNTrisol() const
Get the number of triangular solves ("nTrisol") from localPrec_.
Definition: Ifpack2_Details_Filu_def.hpp:43
The base class of the Ifpack2 FastILU wrappers (Filu, Fildl and Fic)
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:36
void applyLocalPrec(ImplScalarArray x, ImplScalarArray y) const
Apply the local preconditioner with 1-D views of the local parts of X and Y (one vector only) ...
Definition: Ifpack2_Details_Filu_def.hpp:98
Filu(Teuchos::RCP< const TRowMatrix > mat_)
Constructor.
Definition: Ifpack2_Details_Filu_def.hpp:26
void checkLocalIC() const
Verify and print debug info about the internal IC preconditioner.
Definition: Ifpack2_Details_Filu_def.hpp:55
Kokkos::View< ImplScalar *, execution_space > ImplScalarArray
Array of Scalar on device.
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:59
Provides functions for retrieving local CRS arrays (row pointers, column indices, and values) from Tp...