48#ifndef __IFPACK2_FASTILU_BASE_DECL_HPP__
49#define __IFPACK2_FASTILU_BASE_DECL_HPP__
51#include <Tpetra_RowMatrix.hpp>
52#include <Tpetra_CrsMatrix.hpp>
53#include <KokkosCompat_DefaultNode.hpp>
54#include <KokkosSparse_CrsMatrix.hpp>
57#include <shylu_fastutil.hpp>
59#ifdef HAVE_IFPACK2_METIS
70template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
73 Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
81 typedef typename Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::impl_scalar_type
ImplScalar;
83 typedef Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>
TRowMatrix;
85 typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>
TCrsMatrix;
87 typedef Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>
TMultiVec;
89 typedef KokkosSparse::CrsMatrix<Scalar, LocalOrdinal, execution_space>
KCrsMatrix;
91 typedef Kokkos::View<LocalOrdinal *, execution_space>
OrdinalArray;
93 typedef typename Kokkos::View<LocalOrdinal *, execution_space>::HostMirror
OrdinalArrayHost;
96 typedef Kokkos::View< Scalar *, execution_space> ScalarArray;
97 typedef Kokkos::View<const Scalar *, execution_space> ConstScalarArray;
98 #ifdef HAVE_IFPACK2_METIS
99 typedef Kokkos::View<idx_t*, Kokkos::HostSpace> MetisArrayHost;
106 Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >
110 Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >
117 Teuchos::ETransp mode = Teuchos::NO_TRANS,
118 Scalar alpha = Teuchos::ScalarTraits<Scalar>::one(),
119 Scalar beta = Teuchos::ScalarTraits<Scalar>::zero())
const;
147 Teuchos::RCP<const TRowMatrix>
getMatrix()
const;
191 void setMatrix(
const Teuchos::RCP<const TRowMatrix>& A);
194 Teuchos::RCP<const TRowMatrix> mat_;
207 mutable double applyTime_;
214 Params(
const Teuchos::ParameterList& pL, std::string precType);
216 FastILU::SpTRSV sptrsv_algo;
226 static Params getDefaults();
231 #ifdef HAVE_IFPACK2_METIS
232 MetisArrayHost metis_perm_;
233 MetisArrayHost metis_iperm_;
Declaration of interface for preconditioners that can change their matrix after construction.
Mix-in interface for preconditioners that can change their matrix after construction.
Definition Ifpack2_Details_CanChangeMatrix.hpp:93
The base class of the Ifpack2 FastILU wrappers (Filu, Fildl and Fic)
Definition Ifpack2_Details_FastILU_Base_decl.hpp:74
Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >::impl_scalar_type ImplScalar
Kokkos scalar type.
Definition Ifpack2_Details_FastILU_Base_decl.hpp:81
double getInitializeTime() const
Get the time spent in the last initialize() call.
Definition Ifpack2_Details_FastILU_Base_def.hpp:323
virtual void checkLocalILU() const
Verify and print debug information about the underlying ILU preconditioner (only supported if this is...
Definition Ifpack2_Details_FastILU_Base_def.hpp:351
Node::device_type device_type
Kokkos device type.
Definition Ifpack2_Details_FastILU_Base_decl.hpp:77
virtual std::string getSpTrsvType() const =0
Get the name of triangular solve algorithm.
void setMatrix(const Teuchos::RCP< const TRowMatrix > &A)
Definition Ifpack2_Details_FastILU_Base_def.hpp:392
virtual void computeLocalPrec()=0
Get values array from the matrix and then call compute() on the underlying preconditioner.
void compute()
Compute the preconditioner.
Definition Ifpack2_Details_FastILU_Base_def.hpp:262
double getComputeTime() const
Get the time spent in the last compute() call.
Definition Ifpack2_Details_FastILU_Base_def.hpp:330
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TCrsMatrix
Tpetra CRS matrix.
Definition Ifpack2_Details_FastILU_Base_decl.hpp:85
Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TRowMatrix
Tpetra row matrix.
Definition Ifpack2_Details_FastILU_Base_decl.hpp:83
device_type::execution_space execution_space
Kokkos execution space.
Definition Ifpack2_Details_FastILU_Base_decl.hpp:79
double getCopyTime() const
Get the time spent deep copying local 3-array CRS out of the matrix.
Definition Ifpack2_Details_FastILU_Base_def.hpp:344
int getNumApply() const
Get the number of times apply() was called.
Definition Ifpack2_Details_FastILU_Base_def.hpp:316
void initialize()
Initialize the preconditioner.
Definition Ifpack2_Details_FastILU_Base_def.hpp:158
KokkosSparse::CrsMatrix< Scalar, LocalOrdinal, execution_space > KCrsMatrix
Kokkos CRS matrix.
Definition Ifpack2_Details_FastILU_Base_decl.hpp:89
bool isInitialized() const
Whether initialize() has been called since the last time the matrix's structure was changed.
Definition Ifpack2_Details_FastILU_Base_def.hpp:255
virtual std::string getName() const =0
Get the name of the underlying preconditioner ("Filu", "Fildl" or "Fic")
FastILU_Base(Teuchos::RCP< const TRowMatrix > mat_)
Constructor.
Definition Ifpack2_Details_FastILU_Base_def.hpp:61
virtual int getSweeps() const =0
Get the "sweeps" parameter.
void setParameters(const Teuchos::ParameterList &List)
Validate parameters, and set defaults when parameters are not provided.
Definition Ifpack2_Details_FastILU_Base_def.hpp:150
Kokkos::View< LocalOrdinal *, execution_space > OrdinalArray
Array of LocalOrdinal on device.
Definition Ifpack2_Details_FastILU_Base_decl.hpp:91
void apply(const TMultiVec &X, TMultiVec &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const
Apply the preconditioner.
Definition Ifpack2_Details_FastILU_Base_def.hpp:92
double getApplyTime() const
Get the time spent in the last apply() call.
Definition Ifpack2_Details_FastILU_Base_def.hpp:337
Teuchos::RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Get the domain map of the matrix.
Definition Ifpack2_Details_FastILU_Base_def.hpp:77
virtual void checkLocalIC() const
Verify and print debug information about the underlying IC preconditioner.
Definition Ifpack2_Details_FastILU_Base_def.hpp:359
virtual void initLocalPrec()=0
Construct the underlying preconditioner (localPrec_) using given params and then call localPrec_->ini...
int getNumInitialize() const
Get the number of times initialize() was called.
Definition Ifpack2_Details_FastILU_Base_def.hpp:302
Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TMultiVec
Tpetra multivector.
Definition Ifpack2_Details_FastILU_Base_decl.hpp:87
Teuchos::RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Get the range map of the matrix.
Definition Ifpack2_Details_FastILU_Base_def.hpp:85
Teuchos::RCP< const TRowMatrix > getMatrix() const
Get the current matrix.
Definition Ifpack2_Details_FastILU_Base_def.hpp:295
Kokkos::View< LocalOrdinal *, execution_space >::HostMirror OrdinalArrayHost
Array of LocalOrdinal on host.
Definition Ifpack2_Details_FastILU_Base_decl.hpp:93
std::string description() const
Return a brief description of the preconditioner, in YAML format.
Definition Ifpack2_Details_FastILU_Base_def.hpp:366
Kokkos::View< ImplScalar *, execution_space > ImplScalarArray
Array of Scalar on device.
Definition Ifpack2_Details_FastILU_Base_decl.hpp:95
virtual int getNTrisol() const =0
Get the "triangular solve iterations" parameter.
int getNumCompute() const
Get the number of times compute() was called.
Definition Ifpack2_Details_FastILU_Base_def.hpp:309
bool isComputed() const
Whether compute() has been called since the last time the matrix's values or structure were changed.
Definition Ifpack2_Details_FastILU_Base_def.hpp:287
virtual void applyLocalPrec(ImplScalarArray x, ImplScalarArray y) const =0
Apply the local preconditioner with 1-D views of the local parts of X and Y (one vector only)
Interface for all Ifpack2 preconditioners.
Definition Ifpack2_Preconditioner.hpp:108
Ifpack2 implementation details.
Preconditioners and smoothers for Tpetra sparse matrices.
Definition Ifpack2_AdditiveSchwarz_decl.hpp:74