Ifpack2 Templated Preconditioning Package  Version 1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Ifpack2_Container_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_CONTAINER_DECL_HPP
11 #define IFPACK2_CONTAINER_DECL_HPP
12 
15 
16 #include "Ifpack2_ConfigDefs.hpp"
17 #include "Tpetra_RowMatrix.hpp"
18 #include "Teuchos_Describable.hpp"
19 #include <Tpetra_Map.hpp>
20 #include <Tpetra_BlockCrsMatrix.hpp>
22 #include <iostream>
23 
24 #ifndef DOXYGEN_SHOULD_SKIP_THIS
25 namespace Teuchos {
26 // Forward declaration to avoid include.
27 class ParameterList;
28 } // namespace Teuchos
29 #endif // DOXYGEN_SHOULD_SKIP_THIS
30 
31 namespace Ifpack2 {
32 
78 template <class MatrixType>
80  protected:
81  using scalar_type = typename MatrixType::scalar_type;
82  using local_ordinal_type = typename MatrixType::local_ordinal_type;
83  using global_ordinal_type = typename MatrixType::global_ordinal_type;
84  using node_type = typename MatrixType::node_type;
85  using SC = scalar_type;
86  using LO = local_ordinal_type;
87  using GO = global_ordinal_type;
88  using NO = node_type;
89  using import_type = Tpetra::Import<LO, GO, NO>;
90  using mv_type = Tpetra::MultiVector<SC, LO, GO, NO>;
91  using vector_type = Tpetra::Vector<SC, LO, GO, NO>;
92  using map_type = Tpetra::Map<LO, GO, NO>;
93  using crs_matrix_type = Tpetra::CrsMatrix<SC, LO, GO, NO>;
94  using block_crs_matrix_type = Tpetra::BlockCrsMatrix<SC, LO, GO, NO>;
95  using row_matrix_type = Tpetra::RowMatrix<SC, LO, GO, NO>;
96 
97  static_assert(std::is_same<MatrixType, row_matrix_type>::value,
98  "Ifpack2::Container: Please use MatrixType = Tpetra::RowMatrix.");
99 
101  using ISC = typename Kokkos::ArithTraits<SC>::val_type;
102 
105  using HostView = typename mv_type::dual_view_type::t_host;
106  using ConstHostView = typename HostView::const_type;
107 
108  public:
120  const Teuchos::Array<Teuchos::Array<LO> >& partitions,
121  bool pointIndexed);
122 
124  virtual ~Container();
125 
141  Teuchos::ArrayView<const LO> getBlockRows(int blockIndex) const;
142 
151  virtual void initialize() = 0;
152 
154  void setBlockSizes(const Teuchos::Array<Teuchos::Array<LO> >& partitions);
155 
156  void getMatDiag() const;
157 
159  bool isInitialized() const;
160 
162  bool isComputed() const;
163 
172  virtual void compute() = 0;
173 
174  void DoJacobi(ConstHostView X, HostView Y, SC dampingFactor) const;
175  void DoOverlappingJacobi(ConstHostView X, HostView Y, ConstHostView W, SC dampingFactor, bool nonsymScaling) const;
176  void DoGaussSeidel(ConstHostView X, HostView Y, HostView Y2, SC dampingFactor) const;
177  void DoSGS(ConstHostView X, HostView Y, HostView Y2, SC dampingFactor) const;
178 
180  virtual void setParameters(const Teuchos::ParameterList& List) = 0;
181 
194  virtual void
195  apply(ConstHostView X,
196  HostView Y,
197  int blockIndex,
199  SC alpha = Teuchos::ScalarTraits<SC>::one(),
200  SC beta = Teuchos::ScalarTraits<SC>::zero()) const = 0;
201 
203  virtual void
204  weightedApply(ConstHostView X,
205  HostView Y,
206  ConstHostView D,
207  int blockIndex,
209  SC alpha = Teuchos::ScalarTraits<SC>::one(),
210  SC beta = Teuchos::ScalarTraits<SC>::zero()) const = 0;
211 
214  // The underlying container implements the splitting <tt>A = D + R</tt>. Only
215  // it can have efficient access to D and R, as these are constructed in the
216  // symbolic and numeric phases.
217  //
218  // This is the first performance-portable implementation of a block
219  // relaxation, and it is supported currently only by BlockTriDiContainer.
220  virtual void applyInverseJacobi(const mv_type& /* X */, mv_type& /* Y */,
221  SC dampingFactor,
222  bool /* zeroStartingSolution = false */,
223  int /* numSweeps = 1 */) const = 0;
224 
226  virtual void applyMV(const mv_type& X, mv_type& Y) const;
227 
229  virtual void weightedApplyMV(const mv_type& X,
230  mv_type& Y,
231  vector_type& W) const;
232 
233  virtual void clearBlocks();
234 
236  virtual std::ostream& print(std::ostream& os) const = 0;
237 
240  static std::string getName();
241 
242  protected:
244  virtual void DoGSBlock(ConstHostView X, HostView Y, HostView Y2, HostView Resid,
245  SC dampingFactor, LO i) const;
246 
249 
252 
255 
259  Teuchos::Array<LO> blockRows_; // size: total # of local rows (in all local blocks)
261  Teuchos::Array<LO> blockSizes_; // size: # of blocks
283 
288 
289  LO maxBlockSize_;
290 };
291 
292 namespace Details {
293 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
295 }
296 
305 
306 template <class MatrixType, class LocalScalarType>
307 class ContainerImpl : public Container<MatrixType> {
309 
310  protected:
311  using local_scalar_type = LocalScalarType;
312  using SC = typename Container<MatrixType>::scalar_type;
313  using LO = typename Container<MatrixType>::local_ordinal_type;
314  using GO = typename Container<MatrixType>::global_ordinal_type;
315  using NO = typename Container<MatrixType>::node_type;
317  using typename Container<MatrixType>::import_type;
318  using typename Container<MatrixType>::row_matrix_type;
319  using typename Container<MatrixType>::crs_matrix_type;
320  using typename Container<MatrixType>::block_crs_matrix_type;
321  using typename Container<MatrixType>::mv_type;
322  using typename Container<MatrixType>::vector_type;
323  using typename Container<MatrixType>::map_type;
324  using typename Container<MatrixType>::ISC;
326  using LSC = LocalScalarType;
327  using LISC = typename Kokkos::ArithTraits<LSC>::val_type;
328 
329  using local_mv_type = Tpetra::MultiVector<LSC, LO, GO, NO>;
330 
331  using typename Container<MatrixType>::HostView;
332  using typename Container<MatrixType>::ConstHostView;
333  using HostViewLocal = typename local_mv_type::dual_view_type::t_host;
334  using HostSubviewLocal = Kokkos::View<LISC**, Kokkos::LayoutStride, typename HostViewLocal::memory_space>;
335  using ConstHostSubviewLocal = Kokkos::View<const LISC**, Kokkos::LayoutStride, typename HostViewLocal::memory_space>;
336 
337  static_assert(std::is_same<MatrixType, row_matrix_type>::value,
338  "Ifpack2::Container: Please use MatrixType = Tpetra::RowMatrix.");
340  //
341 
348  LO translateRowToCol(LO row);
349 
352 
357  virtual void extract() = 0;
358 
359  public:
361  const Teuchos::Array<Teuchos::Array<LO> >& partitions,
362  bool pointIndexed);
363 
365  virtual ~ContainerImpl();
366 
375  virtual void initialize() = 0;
376 
385  virtual void compute() = 0;
386 
388  virtual void setParameters(const Teuchos::ParameterList& List);
389 
402  virtual void
403  apply(ConstHostView X,
404  HostView Y,
405  int blockIndex,
407  SC alpha = Teuchos::ScalarTraits<SC>::one(),
408  SC beta = Teuchos::ScalarTraits<SC>::zero()) const;
409 
411  virtual void
412  weightedApply(ConstHostView X,
413  HostView Y,
414  ConstHostView D,
415  int blockIndex,
417  SC alpha = Teuchos::ScalarTraits<SC>::one(),
418  SC beta = Teuchos::ScalarTraits<SC>::zero()) const;
419 
422  // The underlying container implements the splitting <tt>A = D + R</tt>. Only
423  // it can have efficient access to D and R, as these are constructed in the
424  // symbolic and numeric phases.
425  //
426  // This is the first performance-portable implementation of a block
427  // relaxation, and it is supported currently only by BlockTriDiContainer.
428  virtual void applyInverseJacobi(const mv_type& /* X */, mv_type& /* Y */,
429  SC dampingFactor,
430  bool /* zeroStartingSolution = false */,
431  int /* numSweeps = 1 */) const;
432 
434  void applyMV(const mv_type& X, mv_type& Y) const;
435 
437  void weightedApplyMV(const mv_type& X,
438  mv_type& Y,
439  vector_type& W) const;
440 
441  virtual void clearBlocks();
442 
444  virtual std::ostream& print(std::ostream& os) const = 0;
445 
448  static std::string getName();
449 
450  protected:
451  // Do Gauss-Seidel on only block i (this is used by DoGaussSeidel and DoSGS)
452  void DoGSBlock(ConstHostView X, HostView Y, HostView Y2, HostView Resid,
453  SC dampingFactor, LO i) const;
454 
460  virtual void
461  solveBlock(ConstHostSubviewLocal X,
462  HostSubviewLocal Y,
463  int blockIndex,
464  Teuchos::ETransp mode,
465  const LSC alpha,
466  const LSC beta) const;
467 
469  mutable HostViewLocal X_local_; // length: blockRows_.size()
470  mutable HostViewLocal Y_local_; // length: blockRows_.size()
471 
482  mutable HostViewLocal weightedApplyScratch_;
483 
485  mutable std::vector<HostSubviewLocal> X_localBlocks_;
486 
488  mutable std::vector<HostSubviewLocal> Y_localBlocks_;
489 };
490 
491 namespace Details {
497 template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, typename Node>
498 struct StridedRowView {
499  using SC = Scalar;
500  using LO = LocalOrdinal;
501 
502  using block_crs_matrix_type = Tpetra::BlockCrsMatrix<SC, LO, GlobalOrdinal, Node>;
503 
504  using h_inds_type = typename block_crs_matrix_type::local_inds_host_view_type;
505  using h_vals_type = typename block_crs_matrix_type::values_host_view_type;
507  StridedRowView(h_vals_type vals_, h_inds_type inds_, int blockSize_, size_t nnz_);
508 
510  // StridedRowView(const SC* vals_, const LO* inds_, int blockSize_, size_t nnz_);
511 
514 
515  SC val(size_t i) const;
516  LO ind(size_t i) const;
517 
518  size_t size() const;
519 
520  private:
521  h_vals_type vals;
522  h_inds_type inds;
523  int blockSize;
524  size_t nnz;
525  // These arrays are only used if the inputMatrix_ doesn't support row views.
526  Teuchos::Array<SC> valsCopy;
527  Teuchos::Array<LO> indsCopy;
528 };
529 } // namespace Details
530 
531 } // namespace Ifpack2
532 
534 template <class MatrixType>
535 std::ostream& operator<<(std::ostream& os, const Ifpack2::Container<MatrixType>& obj);
536 
537 namespace Teuchos {
538 
543 template <class MatrixType>
544 class TEUCHOSCORE_LIB_DLL_EXPORT TypeNameTraits< ::Ifpack2::Container<MatrixType> > {
545  public:
546  static std::string name() {
547  return std::string("Ifpack2::Container<") +
549  }
550 
551  static std::string concreteName(const ::Ifpack2::Container<MatrixType>&) {
552  return name();
553  }
554 };
555 
556 } // namespace Teuchos
557 
558 #endif // IFPACK2_CONTAINER_HPP
virtual void applyInverseJacobi(const mv_type &, mv_type &, SC dampingFactor, bool, int) const =0
Compute Y := (1 - a) Y + a D^{-1} (X - R*Y).
virtual void apply(ConstHostView X, HostView Y, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const
Compute Y := alpha * M^{-1} X + beta*Y.
Definition: Ifpack2_Container_def.hpp:539
virtual void initialize()=0
Do all set-up operations that only require matrix structure.
void weightedApplyMV(const mv_type &X, mv_type &Y, vector_type &W) const
Wrapper for weightedApply with MVs, used in unit tests (never called by BlockRelaxation) ...
Definition: Ifpack2_Container_def.hpp:464
bool isComputed() const
Whether the container has been successfully computed.
Definition: Ifpack2_Container_def.hpp:128
void DoGSBlock(ConstHostView X, HostView Y, HostView Y2, HostView Resid, SC dampingFactor, LO i) const
Do one step of Gauss-Seidel on block i (used by DoGaussSeidel and DoSGS)
Definition: Ifpack2_Container_def.hpp:231
std::vector< HostSubviewLocal > X_localBlocks_
Views for holding pieces of X corresponding to each block.
Definition: Ifpack2_Container_decl.hpp:485
Teuchos::Array< LO > blockRows_
Local indices of the rows of the input matrix that belong to this block.
Definition: Ifpack2_Container_decl.hpp:259
void setBlockSizes(const Teuchos::Array< Teuchos::Array< LO > > &partitions)
Initialize arrays with information about block sizes.
Definition: Ifpack2_Container_def.hpp:89
bool IsInitialized_
If true, the container has been successfully initialized.
Definition: Ifpack2_Container_decl.hpp:285
Teuchos::Array< LO > blockSizes_
Number of rows in each block.
Definition: Ifpack2_Container_decl.hpp:261
typename Kokkos::ArithTraits< SC >::val_type ISC
Internal representation of Scalar in Kokkos::View.
Definition: Ifpack2_Container_decl.hpp:101
StridedRowView(h_vals_type vals_, h_inds_type inds_, int blockSize_, size_t nnz_)
Constructor for row views (preferred)
Definition: Ifpack2_Container_def.hpp:857
std::vector< HostSubviewLocal > Y_localBlocks_
Views for holding pieces of Y corresponding to each block.
Definition: Ifpack2_Container_decl.hpp:488
virtual void applyMV(const mv_type &X, mv_type &Y) const
Wrapper for apply with MultiVector.
Definition: Ifpack2_Container_def.hpp:134
virtual void initialize()=0
Do all set-up operations that only require matrix structure.
bool pointIndexed_
(If hasBlockCrs_) Whether the blocks are described using sub-block row indices instead of full block ...
Definition: Ifpack2_Container_decl.hpp:279
int numBlocks_
The number of blocks (partitions) in the container.
Definition: Ifpack2_Container_decl.hpp:257
The implementation of the numerical features of Container (Jacobi, Gauss-Seidel, SGS). This class allows a custom scalar type (LocalScalarType) to be used for storing blocks and solving the block systems. Hiding this template parameter from the Container interface simplifies the BlockRelaxation and ContainerFactory classes.
Definition: Ifpack2_Container_decl.hpp:307
virtual void setParameters(const Teuchos::ParameterList &List)
Set parameters, if any.
Definition: Ifpack2_Container_def.hpp:443
bool isInitialized() const
Whether the container has been successfully initialized.
Definition: Ifpack2_Container_def.hpp:123
LO scalarsPerRow_
Definition: Ifpack2_Container_decl.hpp:282
virtual ~Container()
Destructor.
Definition: Ifpack2_Container_def.hpp:80
bool hasBlockCrs_
Whether the input matrix is a BlockCRS matrix.
Definition: Ifpack2_Container_decl.hpp:275
Details::StridedRowView< SC, LO, GO, NO > getInputRowView(LO row) const
View a row of the input matrix.
Definition: Ifpack2_Container_def.hpp:803
HostViewLocal weightedApplyScratch_
Definition: Ifpack2_Container_decl.hpp:482
virtual void compute()=0
Extract the local diagonal blocks and prepare the solver.
InverseType::scalar_type LSC
The internal representation of LocalScalarType in Kokkos::View.
Definition: Ifpack2_Container_decl.hpp:326
Container(const Teuchos::RCP< const row_matrix_type > &matrix, const Teuchos::Array< Teuchos::Array< LO > > &partitions, bool pointIndexed)
Constructor.
Definition: Ifpack2_Container_def.hpp:21
int bcrsBlockSize_
If hasBlockCrs_, the number of DOFs per vertex. Otherwise 1.
Definition: Ifpack2_Container_decl.hpp:277
virtual void apply(ConstHostView X, HostView Y, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const =0
Compute Y := alpha * M^{-1} X + beta*Y.
LO translateRowToCol(LO row)
Definition: Ifpack2_Container_def.hpp:493
HostViewLocal X_local_
Scratch vectors used in apply().
Definition: Ifpack2_Container_decl.hpp:469
Teuchos::RCP< const block_crs_matrix_type > inputBlockMatrix_
The input matrix, dynamic cast to BlockCrsMatrix. May be null.
Definition: Ifpack2_Container_decl.hpp:254
Teuchos::RCP< const row_matrix_type > inputMatrix_
The input matrix to the constructor.
Definition: Ifpack2_Container_decl.hpp:248
virtual void weightedApplyMV(const mv_type &X, mv_type &Y, vector_type &W) const
Wrapper for weightedApply with MultiVector.
Definition: Ifpack2_Container_def.hpp:140
void applyMV(const mv_type &X, mv_type &Y) const
Wrapper for apply with MVs, used in unit tests (never called by BlockRelaxation)
Definition: Ifpack2_Container_def.hpp:456
virtual void applyInverseJacobi(const mv_type &, mv_type &, SC dampingFactor, bool, int) const
Compute Y := (1 - a) Y + a D^{-1} (X - R*Y).
Definition: Ifpack2_Container_def.hpp:447
typename mv_type::dual_view_type::t_host HostView
Definition: Ifpack2_Container_decl.hpp:105
virtual void weightedApply(ConstHostView X, HostView Y, ConstHostView D, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const
Compute Y := alpha * diag(D) * M^{-1} (diag(D) * X) + beta*Y.
Definition: Ifpack2_Container_def.hpp:648
Teuchos::RCP< vector_type > Diag_
Diagonal elements.
Definition: Ifpack2_Container_decl.hpp:265
virtual ~ContainerImpl()
Destructor.
Definition: Ifpack2_Container_def.hpp:439
bool IsParallel_
Whether the problem is distributed across multiple MPI processes.
Definition: Ifpack2_Container_decl.hpp:267
static std::string getName()
Definition: Ifpack2_Container_def.hpp:475
bool IsComputed_
If true, the container has been successfully computed.
Definition: Ifpack2_Container_decl.hpp:287
virtual void weightedApply(ConstHostView X, HostView Y, ConstHostView D, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const =0
Compute Y := alpha * diag(D) * M^{-1} (diag(D) * X) + beta*Y.
virtual void compute()=0
Extract the local diagonal blocks and prepare the solver.
Interface for creating and solving a set of local linear problems.
Definition: Ifpack2_Container_decl.hpp:79
Teuchos::ArrayView< const LO > getBlockRows(int blockIndex) const
Local indices of the rows of the input matrix that belong to this block.
Definition: Ifpack2_Container_def.hpp:84
virtual void extract()=0
Extract the submatrices identified by the local indices set by the constructor. BlockMatrix may be an...
virtual void setParameters(const Teuchos::ParameterList &List)=0
Set parameters, if any.
Structure for read-only views of general matrix rows.
Definition: Ifpack2_Container_decl.hpp:294
virtual std::ostream & print(std::ostream &os) const =0
Print basic information about the container to os.
static std::string getName()
Definition: Ifpack2_Container_def.hpp:146
GO NumGlobalRows_
Number of global rows in input matrix.
Definition: Ifpack2_Container_decl.hpp:271
GO NumGlobalNonzeros_
Number of nonzeros in input matrix.
Definition: Ifpack2_Container_decl.hpp:273
LO NumLocalRows_
Number of local rows in input matrix.
Definition: Ifpack2_Container_decl.hpp:269
Teuchos::Array< LO > blockOffsets_
Starting index in blockRows_ of local row indices for each block.
Definition: Ifpack2_Container_decl.hpp:263
virtual std::ostream & print(std::ostream &os) const =0
Print basic information about the container to os.
static std::string name()
virtual void solveBlock(ConstHostSubviewLocal X, HostSubviewLocal Y, int blockIndex, Teuchos::ETransp mode, const LSC alpha, const LSC beta) const
Definition: Ifpack2_Container_def.hpp:481
Teuchos::RCP< const crs_matrix_type > inputCrsMatrix_
The input matrix, dynamic cast to CrsMatrix. May be null.
Definition: Ifpack2_Container_decl.hpp:251
virtual void DoGSBlock(ConstHostView X, HostView Y, HostView Y2, HostView Resid, SC dampingFactor, LO i) const
Do one step of Gauss-Seidel on block i (used by DoGaussSeidel and DoSGS)
Definition: Ifpack2_Container_def.hpp:151