Ifpack2 Templated Preconditioning Package  Version 1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Ifpack2_LocalSparseTriangularSolver_decl.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 
10 #ifndef IFPACK2_LOCALSPARSETRIANGULARSOLVER_DECL_HPP
11 #define IFPACK2_LOCALSPARSETRIANGULARSOLVER_DECL_HPP
12 
15 #include "Tpetra_CrsMatrix_fwd.hpp"
16 #include "Teuchos_FancyOStream.hpp"
17 #include <type_traits>
18 
19 #include "KokkosSparse_sptrsv.hpp"
20 
21 namespace Ifpack2 {
22 
45 template<class MatrixType>
47  virtual public Ifpack2::Preconditioner<typename MatrixType::scalar_type,
48  typename MatrixType::local_ordinal_type,
49  typename MatrixType::global_ordinal_type,
50  typename MatrixType::node_type>,
51  virtual public Ifpack2::Details::CanChangeMatrix<Tpetra::RowMatrix<typename MatrixType::scalar_type,
52  typename MatrixType::local_ordinal_type,
53  typename MatrixType::global_ordinal_type,
54  typename MatrixType::node_type> >
55 {
56 public:
58  typedef typename MatrixType::scalar_type scalar_type;
60  typedef typename MatrixType::impl_scalar_type impl_scalar_type;
62  typedef typename MatrixType::local_ordinal_type local_ordinal_type;
64  typedef typename MatrixType::global_ordinal_type global_ordinal_type;
66  typedef typename MatrixType::node_type node_type;
67 
69  typedef typename MatrixType::mag_type magnitude_type;
71  typedef Tpetra::Map<local_ordinal_type, global_ordinal_type, node_type> map_type;
73  typedef Tpetra::RowMatrix<scalar_type, local_ordinal_type,
76  typedef Tpetra::CrsMatrix<scalar_type, local_ordinal_type,
78 
79  static_assert (std::is_same<MatrixType, row_matrix_type>::value,
80  "Ifpack2::LocalSparseTriangularSolver: The template parameter "
81  "MatrixType must be a Tpetra::RowMatrix specialization. "
82  "Please don't use Tpetra::CrsMatrix (a subclass of "
83  "Tpetra::RowMatrix) here anymore. The constructor can take "
84  "either a RowMatrix or a CrsMatrix just fine.");
85 
86  // Use the local matrix types
87  using local_matrix_device_type = typename crs_matrix_type::local_matrix_device_type;
88  using local_matrix_graph_device_type = typename local_matrix_device_type::StaticCrsGraphType;
89  using lno_row_view_t = typename local_matrix_graph_device_type::row_map_type;
90  using lno_nonzero_view_t = typename local_matrix_graph_device_type::entries_type;
91  using scalar_nonzero_view_t = typename local_matrix_device_type::values_type;
92  using TemporaryMemorySpace = typename local_matrix_graph_device_type::device_type::memory_space;
93  using PersistentMemorySpace = typename local_matrix_graph_device_type::device_type::memory_space;
94  using HandleExecSpace = typename local_matrix_graph_device_type::device_type::execution_space;
95  using k_handle = typename KokkosKernels::Experimental::KokkosKernelsHandle<typename lno_row_view_t::const_value_type, typename lno_nonzero_view_t::const_value_type, typename scalar_nonzero_view_t::value_type, HandleExecSpace, TemporaryMemorySpace,PersistentMemorySpace >;
96 
125 
136 
141 
153  LocalSparseTriangularSolver (const bool /* unused */, const Teuchos::RCP<Teuchos::FancyOStream>& out);
154 
156  virtual ~LocalSparseTriangularSolver ();
157 
170  void setParameters (const Teuchos::ParameterList& params);
171 
177  void initialize ();
178 
180  inline bool isInitialized () const {
181  return isInitialized_;
182  }
183 
188  void compute ();
189 
191  inline bool isComputed () const {
192  return isComputed_;
193  }
194 
196 
197 
213  void
214  apply (const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
215  Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
219 
222 
225 
234  void
235  applyMat (const Tpetra::MultiVector<scalar_type, local_ordinal_type,
237  Tpetra::MultiVector<scalar_type, local_ordinal_type,
239  Teuchos::ETransp mode = Teuchos::NO_TRANS) const;
240 
243 
246  return A_;
247  }
248 
250  double getComputeFlops () const;
251 
253  double getApplyFlops () const;
254 
256  int getNumInitialize () const;
257 
259  int getNumCompute () const;
260 
262  int getNumApply () const;
263 
265  double getInitializeTime () const;
266 
268  double getComputeTime () const;
269 
271  double getApplyTime () const;
272 
274 
276 
278  std::string description() const;
279 
301  void
303  const Teuchos::EVerbosityLevel verbLevel =
305 
310  virtual void setMatrix (const Teuchos::RCP<const row_matrix_type>& A);
311 
314  void setStreamInfo (const bool& isKokkosKernelsStream, const int& num_streams, const std::vector<HandleExecSpace>& exec_space_instances);
315 
320  void setMatrices (const std::vector< Teuchos::RCP<crs_matrix_type> >& A_crs_v);
321 
323 
324 private:
331  std::vector< Teuchos::RCP<crs_matrix_type> > A_crs_v_;
332 
333  typedef Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type> MV;
334  mutable Teuchos::RCP<MV> X_colMap_;
335  mutable Teuchos::RCP<MV> Y_rowMap_;
336 
337  bool isInitialized_;
338  bool isComputed_;
348  bool isInternallyChanged_;
349  bool reverseStorage_;
350 
351  mutable int numInitialize_;
352  mutable int numCompute_;
353  mutable int numApply_;
354 
355  double initializeTime_;
356  double computeTime_;
357  double applyTime_;
358 
360  class HtsImpl;
361  Teuchos::RCP<HtsImpl> htsImpl_;
362 
364  bool isKokkosKernelsSptrsv_;
366  std::vector< Teuchos::RCP<k_handle> > kh_v_;
367  int num_streams_;
368  bool isKokkosKernelsStream_;
369  bool kh_v_nonnull_;
370  std::vector<HandleExecSpace> exec_space_instances_;
371 
375  std::string uplo_;
378  std::string diag_;
379 
398  void
399  localApply (const MV& X,
400  MV& Y,
401  const Teuchos::ETransp mode,
402  const scalar_type& alpha,
403  const scalar_type& beta) const;
404 
406  void
407  localTriangularSolve (const MV& Y,
408  MV& X,
409  const Teuchos::ETransp mode) const;
410 
411  void initializeState();
412 
413  KokkosSparse::Experimental::SPTRSVAlgorithm kokkosKernelsAlgorithm() const;
414 };
415 
416 } // namespace Ifpack2
417 
418 #endif // IFPACK2_LOCALSPARSETRIANGULARSOLVER_DECL_HPP
Mix-in interface for preconditioners that can change their matrix after construction.
Definition: Ifpack2_Details_CanChangeMatrix.hpp:60
double getApplyFlops() const
Return the number of flops for the application of the preconditioner.
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:756
int getNumCompute() const
Return the number of calls to compute().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1100
std::string description() const
A one-line description of this object.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1135
void compute()
&quot;Numeric&quot; phase of setup
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:652
Teuchos::RCP< const row_matrix_type > getMatrix() const
The original input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:245
MatrixType::global_ordinal_type global_ordinal_type
Type of the global indices of the input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:64
bool isComputed() const
Return true if compute() has been called.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:191
double getApplyTime() const
Return the time spent in apply().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1128
double getComputeTime() const
Return the time spent in compute().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1121
double getInitializeTime() const
Return the time spent in initialize().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1114
double getComputeFlops() const
Return the number of flops in the computation phase.
MatrixType::node_type node_type
Node type of the input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:66
MatrixType::impl_scalar_type impl_scalar_type
Implementation scalar type of the entries of the input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:60
void setMatrices(const std::vector< Teuchos::RCP< crs_matrix_type > > &A_crs_v)
Set this preconditioner&#39;s matrices (used by stream interface of triangular solve).
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1290
void applyMat(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) const
Apply the original input matrix.
Teuchos::RCP< const map_type > getDomainMap() const
The domain of this operator.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1208
&quot;Preconditioner&quot; that solves local sparse triangular systems.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:46
void initialize()
&quot;Symbolic&quot; phase of setup
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:405
Interface for all Ifpack2 preconditioners.
Definition: Ifpack2_Preconditioner.hpp:74
int getNumInitialize() const
Return the number of calls to initialize().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1093
void setStreamInfo(const bool &isKokkosKernelsStream, const int &num_streams, const std::vector< HandleExecSpace > &exec_space_instances)
Set this triangular solver&#39;s stream information.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1279
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:71
virtual void setMatrix(const Teuchos::RCP< const row_matrix_type > &A)
Set this preconditioner&#39;s matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1233
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
This operator&#39;s communicator.
Declaration of interface for preconditioners that can change their matrix after construction.
static const EVerbosityLevel verbLevel_default
int getNumApply() const
Return the number of calls to apply().
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1107
MatrixType::local_ordinal_type local_ordinal_type
Type of the local indices of the input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:62
Teuchos::RCP< const map_type > getRangeMap() const
The range of this operator.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:1221
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:77
void setParameters(const Teuchos::ParameterList &params)
Set this object&#39;s parameters.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:358
MatrixType::scalar_type scalar_type
Scalar type of the entries of the input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:58
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:1170
MatrixType::mag_type magnitude_type
Type of the absolute value (magnitude) of a scalar_type value.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:69
virtual ~LocalSparseTriangularSolver()
Destructor (virtual for memory safety).
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:339
LocalSparseTriangularSolver()
Constructor that takes no input matrix.
Definition: Ifpack2_LocalSparseTriangularSolver_def.hpp:300
bool isInitialized() const
Return true if the preconditioner has been successfully initialized.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:180
Tpetra::RowMatrix< scalar_type, local_ordinal_type, global_ordinal_type, node_type > row_matrix_type
Specialization of Tpetra::RowMatrix used by this class.
Definition: Ifpack2_LocalSparseTriangularSolver_decl.hpp:74