46 #ifndef IFPACK2_MDF_DECL_HPP
47 #define IFPACK2_MDF_DECL_HPP
51 #include "Tpetra_CrsMatrix_decl.hpp"
54 #include "Ifpack2_LocalSparseTriangularSolver_decl.hpp"
55 #include "KokkosSparse_mdf.hpp"
57 #include <type_traits>
82 template<
class MatrixType>
85 typename MatrixType::local_ordinal_type,
86 typename MatrixType::global_ordinal_type,
87 typename MatrixType::node_type>,
89 typename MatrixType::local_ordinal_type,
90 typename MatrixType::global_ordinal_type,
91 typename MatrixType::node_type> >
122 static_assert(std::is_same<MatrixType, row_matrix_type>::value,
"Ifpack2::MDF: The template parameter MatrixType must be a Tpetra::RowMatrix specialization. Please don't use Tpetra::CrsMatrix (a subclass of Tpetra::RowMatrix) here anymore.");
133 template <
class NewMatrixType>
friend class MDF;
135 typedef typename crs_matrix_type::global_inds_host_view_type global_inds_host_view_type;
136 typedef typename crs_matrix_type::local_inds_host_view_type local_inds_host_view_type;
137 typedef typename crs_matrix_type::values_host_view_type values_host_view_type;
140 typedef typename crs_matrix_type::nonconst_global_inds_host_view_type nonconst_global_inds_host_view_type;
141 typedef typename crs_matrix_type::nonconst_local_inds_host_view_type nonconst_local_inds_host_view_type;
142 typedef typename crs_matrix_type::nonconst_values_host_view_type nonconst_values_host_view_type;
149 typedef typename crs_matrix_type::local_matrix_device_type local_matrix_device_type;
150 typedef typename local_matrix_device_type::StaticCrsGraphType::row_map_type lno_row_view_t;
151 typedef typename local_matrix_device_type::StaticCrsGraphType::entries_type lno_nonzero_view_t;
152 typedef typename local_matrix_device_type::values_type scalar_nonzero_view_t;
153 typedef typename local_matrix_device_type::StaticCrsGraphType::device_type::memory_space TemporaryMemorySpace;
154 typedef typename local_matrix_device_type::StaticCrsGraphType::device_type::memory_space PersistentMemorySpace;
155 typedef typename local_matrix_device_type::StaticCrsGraphType::device_type::execution_space HandleExecSpace;
156 typedef typename KokkosKernels::Experimental::KokkosKernelsHandle
157 <
typename lno_row_view_t::const_value_type,
typename lno_nonzero_view_t::const_value_type,
typename scalar_nonzero_view_t::value_type,
158 HandleExecSpace, TemporaryMemorySpace,PersistentMemorySpace > kk_handle_type;
181 virtual ~MDF () =
default;
206 return isInitialized_;
215 return numInitialize_;
228 return initializeTime_;
319 apply (
const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
320 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
329 void apply_impl (
const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
330 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
357 multiply (
const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
358 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
361 using MDF_handle_device_type = KokkosSparse::Experimental::MDF_handle<local_matrix_device_type>;
373 true, std::logic_error,
"Ifpack2::MDF::SetOverlapMode: "
374 "MDF no longer implements overlap on its own. "
375 "Use MDF with AdditiveSchwarz if you want overlap.");
380 return getL ().getGlobalNumEntries () +
getU ().getGlobalNumEntries ();
399 typedef Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> MV;
403 void allocateSolvers ();
404 void allocatePermutations(
bool force =
false);
417 typedef Tpetra::Vector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> vec_type;
429 lno_row_view_t A_local_rowmap_;
430 lno_nonzero_view_t A_local_entries_;
431 scalar_nonzero_view_t A_local_values_;
459 mutable int numApply_;
461 double initializeTime_;
463 mutable double applyTime_;
Tpetra::CrsMatrix< scalar_type, local_ordinal_type, global_ordinal_type, node_type > crs_matrix_type
Tpetra::CrsMatrix specialization used by this class for representing L and U.
Definition: Ifpack2_MDF_decl.hpp:122
virtual ~MDF()=default
Destructor (declared virtual for memory safety).
Mix-in interface for preconditioners that can change their matrix after construction.
Definition: Ifpack2_Details_CanChangeMatrix.hpp:93
Declaration and definition of IlukGraph.
Ifpack2::ScalingType enumerable type.
void compute()
Compute the (numeric) incomplete factorization.
Definition: Ifpack2_MDF_def.hpp:564
MatrixType::local_ordinal_type local_ordinal_type
The type of local indices in the input MatrixType.
Definition: Ifpack2_MDF_decl.hpp:98
MDF (incomplete LU factorization with minimum discarded fill reordering) of a Tpetra sparse matrix...
Definition: Ifpack2_MDF_decl.hpp:83
bool isComputed() const
Whether compute() has been called on this object.
Definition: Ifpack2_MDF_decl.hpp:209
permutations_type permutations_
The computed permuations from MDF factorization.
Definition: Ifpack2_MDF_decl.hpp:443
Teuchos::RCP< const row_matrix_type > A_
The (original) input matrix for which to compute ILU(k).
Definition: Ifpack2_MDF_decl.hpp:420
permutations_type & getReversePermutations() const
Return the reverse permutations of the MDF factorization.
Definition: Ifpack2_MDF_def.hpp:286
crs_matrix_type::impl_scalar_type impl_scalar_type
Scalar type stored in Kokkos::Views (CrsMatrix and MultiVector)
Definition: Ifpack2_MDF_decl.hpp:131
int getNumInitialize() const
Number of successful initialize() calls for this object.
Definition: Ifpack2_MDF_decl.hpp:214
MatrixType::node_type node_type
The Node type used by the input MatrixType.
Definition: Ifpack2_MDF_decl.hpp:104
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
permutations_type reversePermutations_
The reverse permuations from MDF factorization.
Definition: Ifpack2_MDF_decl.hpp:446
std::string description() const
A one-line description of this object.
Definition: Ifpack2_MDF_def.hpp:765
Teuchos::RCP< const Tpetra::Map< local_ordinal_type, global_ordinal_type, node_type > > getRangeMap() const
Returns the Tpetra::Map object associated with the range of this operator.
Definition: Ifpack2_MDF_def.hpp:348
int getLevelOfFill() const
Get level of fill (the "k" in ILU(k)).
Definition: Ifpack2_MDF_decl.hpp:368
MatrixType::global_ordinal_type global_ordinal_type
The type of global indices in the input MatrixType.
Definition: Ifpack2_MDF_decl.hpp:101
Teuchos::RCP< const crs_matrix_type > getCrsMatrix() const
Return the input matrix A as a Tpetra::CrsMatrix, if possible; else throws.
Definition: Ifpack2_MDF_def.hpp:428
size_t getNodeSmootherComplexity() const
Get a rough estimate of cost per iteration.
Definition: Ifpack2_MDF_def.hpp:312
node_type::device_type device_type
The Kokkos device type of the input MatrixType.
Definition: Ifpack2_MDF_decl.hpp:110
Tpetra::CombineMode getOverlapMode()
Get overlap mode type.
Definition: Ifpack2_MDF_decl.hpp:371
Teuchos::RCP< crs_matrix_type > U_
The U (upper triangular) factor of ILU(k).
Definition: Ifpack2_MDF_decl.hpp:438
permutations_type & getPermutations() const
Return the permutations of the MDF factorization.
Definition: Ifpack2_MDF_def.hpp:274
Teuchos::RCP< LocalSparseTriangularSolver< row_matrix_type > > L_solver_
Sparse triangular solver for L.
Definition: Ifpack2_MDF_decl.hpp:436
Teuchos::ScalarTraits< scalar_type >::magnitudeType magnitude_type
The type of the magnitude (absolute value) of a matrix entry.
Definition: Ifpack2_MDF_decl.hpp:107
Teuchos::RCP< crs_matrix_type > L_
The L (lower triangular) factor of ILU(k).
Definition: Ifpack2_MDF_decl.hpp:434
int getNumApply() const
Number of successful apply() calls for this object.
Definition: Ifpack2_MDF_decl.hpp:222
node_type::execution_space execution_space
The Kokkos execution space of the input MatrixType.
Definition: Ifpack2_MDF_decl.hpp:113
Interface for all Ifpack2 preconditioners.
Definition: Ifpack2_Preconditioner.hpp:107
virtual void setMatrix(const Teuchos::RCP< const row_matrix_type > &A)
Change the matrix to be preconditioned.
Definition: Ifpack2_MDF_def.hpp:224
void setParameters(const Teuchos::ParameterList ¶ms)
Definition: Ifpack2_MDF_def.hpp:365
double getApplyTime() const
Total time in seconds taken by all successful apply() calls for this object.
Definition: Ifpack2_MDF_decl.hpp:235
Declaration of interface for preconditioners that can change their matrix after construction.
double getInitializeTime() const
Total time in seconds taken by all successful initialize() calls for this object. ...
Definition: Ifpack2_MDF_decl.hpp:227
const crs_matrix_type & getL() const
Return the L factor of the MDF factorization.
Definition: Ifpack2_MDF_def.hpp:261
int getNumCompute() const
Number of successful compute() calls for this object.
Definition: Ifpack2_MDF_decl.hpp:218
Teuchos::RCP< const Tpetra::Map< local_ordinal_type, global_ordinal_type, node_type > > getDomainMap() const
Returns the Tpetra::Map object associated with the domain of this operator.
Definition: Ifpack2_MDF_def.hpp:329
Tpetra::RowMatrix< scalar_type, local_ordinal_type, global_ordinal_type, node_type > row_matrix_type
Tpetra::RowMatrix specialization used by this class.
Definition: Ifpack2_MDF_decl.hpp:119
MatrixType::scalar_type scalar_type
The type of the entries of the input MatrixType.
Definition: Ifpack2_MDF_decl.hpp:95
Teuchos::RCP< const row_matrix_type > A_local_
The matrix whos numbers are used to to compute ILU(k). The graph may be computed using a crs_matrix_t...
Definition: Ifpack2_MDF_decl.hpp:428
Teuchos::RCP< const row_matrix_type > getMatrix() const
Get the input matrix.
Definition: Ifpack2_MDF_def.hpp:421
const crs_matrix_type & getU() const
Return the U factor of the MDF factorization.
Definition: Ifpack2_MDF_def.hpp:299
Teuchos::RCP< LocalSparseTriangularSolver< row_matrix_type > > U_solver_
Sparse triangular solver for U.
Definition: Ifpack2_MDF_decl.hpp:440
void initialize()
Initialize by computing the symbolic incomplete factorization.
Definition: Ifpack2_MDF_def.hpp:468
bool isInitialized() const
Whether initialize() has been called on this object.
Definition: Ifpack2_MDF_decl.hpp:205
double getComputeTime() const
Total time in seconds taken by all successful compute() calls for this object.
Definition: Ifpack2_MDF_decl.hpp:231
Tpetra::global_size_t getGlobalNumEntries() const
Returns the number of nonzero entries in the global graph.
Definition: Ifpack2_MDF_decl.hpp:379
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 (inverse of the) incomplete factorization to X, resulting in Y.
Definition: Ifpack2_MDF_def.hpp:696