Ifpack2 Templated Preconditioning Package  Version 1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Ifpack2_Details_Filu_def.hpp
1 // @HEADER
2 // *****************************************************************************
3 // Ifpack2: Templated Object-Oriented Algebraic Preconditioner Package
4 //
5 // Copyright 2009 NTESS and the Ifpack2 contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
11 
12 #ifndef __IFPACK2_FILU_DEF_HPP__
13 #define __IFPACK2_FILU_DEF_HPP__
14 
15 #include "Ifpack2_Details_Filu_decl.hpp"
17 #include "Ifpack2_Details_getCrsMatrix.hpp"
18 #include <Kokkos_Timer.hpp>
19 #include <shylu_fastilu.hpp>
20 
21 namespace Ifpack2 {
22 namespace Details {
23 
24 template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, typename Node, bool BlockCrsEnabled>
27  : FastILU_Base<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A) {}
28 
29 template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, typename Node, bool BlockCrsEnabled>
31  getSweeps() const {
32  return localPrec_->getNFact();
33 }
34 
35 template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, typename Node, bool BlockCrsEnabled>
37  getSpTrsvType() const {
38  return localPrec_->getSpTrsvType();
39 }
40 
41 template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, typename Node, bool BlockCrsEnabled>
43  getNTrisol() const {
44  return localPrec_->getNTrisol();
45 }
46 
47 template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, typename Node, bool BlockCrsEnabled>
49  checkLocalILU() const {
50  localPrec_->checkILU();
51 }
52 
53 template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, typename Node, bool BlockCrsEnabled>
55  checkLocalIC() const {
56  localPrec_->checkIC();
57 }
58 
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_);
65 
66  if (p.blockCrsSize > 1 && !BlockCrsEnabled) {
67  throw std::runtime_error("Must use prec type FAST_ILU_B if you want blockCrs support");
68  }
69 
70  bool skipSortMatrix = !matCrs.is_null() && matCrs->getCrsGraph()->isSorted() &&
71  !p.use_metis;
72  localPrec_ =
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,
75  p.blockCrsSize));
76 
77 #ifdef HAVE_IFPACK2_METIS
78  if (p.use_metis) {
79  localPrec_->setMetisPerm(this->metis_perm_, this->metis_iperm_);
80  }
81 #endif
82 
83  localPrec_->initialize();
84  this->initTime_ = localPrec_->getInitializeTime();
85 }
86 
87 template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, typename Node, bool BlockCrsEnabled>
90  // update values in local prec (until compute(), values aren't needed)
91  localPrec_->setValues(this->localValues_);
92  localPrec_->compute();
93  this->computeTime_ = localPrec_->getComputeTime();
94 }
95 
96 template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, typename Node, bool BlockCrsEnabled>
99  localPrec_->apply(x, y);
100 
101  // since this may be applied to multiple vectors, add to applyTime_ instead of setting it
102  this->applyTime_ += localPrec_->getApplyTime();
103 }
104 
105 template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, typename Node, bool BlockCrsEnabled>
107  getName() const {
108  return "Filu";
109 }
110 
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>;
114 
115 } // namespace Details
116 } // namespace Ifpack2
117 
118 #endif
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 (&quot;Filu&quot;, &quot;Fildl&quot; or &quot;Fic&quot;)
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_-&gt;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...