10 #ifndef IFPACK2_RELAXATION_DECL_HPP
11 #define IFPACK2_RELAXATION_DECL_HPP
15 #include "Ifpack2_Parameters.hpp"
16 #include "Tpetra_Vector.hpp"
18 #include "Tpetra_CrsMatrix.hpp"
19 #include "Tpetra_BlockCrsMatrix.hpp"
20 #include <type_traits>
21 #include <KokkosKernels_Handle.hpp>
22 #include "Ifpack2_Details_GaussSeidel.hpp"
23 #include "Ifpack2_Details_InverseDiagonalKernel.hpp"
25 #ifndef DOXYGEN_SHOULD_SKIP_THIS
29 template <
class TpetraOperatorType>
30 class ScaledDampedResidual;
40 #endif // DOXYGEN_SHOULD_SKIP_THIS
205 template <
class MatrixType>
207 typename MatrixType::scalar_type,
208 typename MatrixType::local_ordinal_type,
209 typename MatrixType::global_ordinal_type,
210 typename MatrixType::node_type>,
212 Tpetra::RowMatrix<typename MatrixType::scalar_type,
213 typename MatrixType::local_ordinal_type,
214 typename MatrixType::global_ordinal_type,
215 typename MatrixType::node_type> > {
250 static_assert(std::is_same<MatrixType, row_matrix_type>::value,
"Ifpack2::Relaxation: Please use MatrixType = Tpetra::RowMatrix. This saves build times, library sizes, and executable sizes. Don't worry, this class still works with CrsMatrix and BlockCrsMatrix; those are both subclasses of RowMatrix.");
375 bool supportsZeroStartingSolution() {
return true; }
388 return isInitialized_;
448 apply(
const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
449 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& Y,
569 #if KOKKOS_VERSION >= 40799
570 typedef typename KokkosKernels::ArithTraits<scalar_type>::val_type impl_scalar_type;
572 typedef typename Kokkos::ArithTraits<scalar_type>::val_type impl_scalar_type;
590 block_crs_matrix_type;
593 block_multivector_type;
595 typedef Tpetra::Map<local_ordinal_type, global_ordinal_type, node_type> map_type;
596 typedef Tpetra::Import<local_ordinal_type, global_ordinal_type, node_type> import_type;
598 typedef typename crs_matrix_type::nonconst_local_inds_host_view_type nonconst_local_inds_host_view_type;
599 typedef typename crs_matrix_type::nonconst_values_host_view_type nonconst_values_host_view_type;
607 typedef typename crs_matrix_type::local_matrix_device_type local_matrix_device_type;
608 typedef typename local_matrix_device_type::StaticCrsGraphType::row_map_type lno_row_view_t;
609 typedef typename local_matrix_device_type::StaticCrsGraphType::entries_type lno_nonzero_view_t;
610 typedef typename local_matrix_device_type::values_type scalar_nonzero_view_t;
611 typedef typename local_matrix_device_type::StaticCrsGraphType::device_type TemporaryWorkSpace;
612 typedef typename local_matrix_device_type::StaticCrsGraphType::device_type PersistentWorkSpace;
613 typedef typename local_matrix_device_type::StaticCrsGraphType::execution_space MyExecSpace;
614 typedef typename KokkosKernels::Experimental::KokkosKernelsHandle<
typename lno_row_view_t::const_value_type,
local_ordinal_type,
typename scalar_nonzero_view_t::value_type,
615 MyExecSpace, TemporaryWorkSpace, PersistentWorkSpace>
616 mt_kernel_handle_type;
623 using SerialGaussSeidel = Ifpack2::Details::GaussSeidel<scalar_type, local_ordinal_type, global_ordinal_type, node_type>;
647 void ApplyInverseRichardson(
648 const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
649 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& Y)
const;
652 void ApplyInverseJacobi(
653 const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
654 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& Y)
const;
657 void ApplyInverseJacobi_BlockCrsMatrix(
658 const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
659 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& Y)
const;
662 void ApplyInverseMTGS_CrsMatrix(
663 const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& B,
664 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
665 Tpetra::ESweepDirection direction)
const;
668 void ApplyInverseSerialGS(
669 const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
670 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& Y,
671 Tpetra::ESweepDirection direction)
const;
674 void ApplyInverseSerialGS_CrsMatrix(
const crs_matrix_type& A,
675 const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& B,
676 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
677 Tpetra::ESweepDirection direction)
const;
680 void ApplyInverseSerialGS_RowMatrix(
681 const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
682 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& Y,
683 Tpetra::ESweepDirection direction)
const;
686 void ApplyInverseSerialGS_BlockCrsMatrix(
687 const block_crs_matrix_type& A,
688 const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
689 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& Y,
690 Tpetra::ESweepDirection direction);
693 void ApplyInverseMTSGS_CrsMatrix(
694 const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
695 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& Y)
const;
698 const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& B,
699 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
700 const Tpetra::ESweepDirection direction)
const;
702 void computeBlockCrs();
705 void updateCachedMultiVector(
const Teuchos::RCP<
const Tpetra::Map<local_ordinal_type, global_ordinal_type, node_type> >& map,
size_t numVecs)
const;
731 typedef Kokkos::View<
typename block_crs_matrix_type::impl_scalar_type***,
732 typename block_crs_matrix_type::device_type>
734 typedef Kokkos::View<
typename block_crs_matrix_type::impl_scalar_type***,
735 typename block_crs_matrix_type::device_type,
736 Kokkos::MemoryUnmanaged>
737 unmanaged_block_diag_type;
753 block_diag_type blockDiag_;
760 Details::RelaxationType PrecType_ = Ifpack2::Details::JACOBI;
766 bool ZeroStartingSolution_ =
true;
768 bool DoBackwardGS_ =
false;
770 bool DoL1Method_ =
false;
776 bool fixTinyDiagEntries_ =
false;
778 bool checkDiagEntries_ =
false;
780 int clusterSize_ = 1;
782 int longRowThreshold_ = 0;
784 KokkosGraph::ColoringAlgorithm mtColoringAlgorithm_ = KokkosGraph::COLORING_DEFAULT;
787 int NumOuterSweeps_ = 1;
789 int NumInnerSweeps_ = 1;
791 bool InnerSpTrsv_ =
false;
795 bool CompactForm_ =
false;
798 bool is_matrix_structurally_symmetric_ =
false;
801 bool ifpack2_dump_matrix_ =
false;
804 bool isInitialized_ =
false;
806 bool IsComputed_ =
false;
808 int NumInitialize_ = 0;
812 bool TimerForApply_ =
true;
814 mutable int NumApply_ = 0;
816 double InitializeTime_ = 0.0;
818 double ComputeTime_ = 0.0;
820 mutable double ApplyTime_ = 0.0;
822 double ComputeFlops_ = 0.0;
824 mutable double ApplyFlops_ = 0.0;
831 size_t globalNumSmallDiagEntries_ = 0;
833 size_t globalNumZeroDiagEntries_ = 0;
835 size_t globalNumNegDiagEntries_ = 0;
847 Kokkos::View<size_t*, typename node_type::device_type> diagOffsets_;
854 bool savedDiagOffsets_ =
false;
856 bool hasBlockCrsMatrix_ =
false;
866 #endif // IFPACK2_RELAXATION_DECL_HPP
bool hasTransposeApply() const
Whether apply() and applyMat() let you apply the transpose or conjugate transpose.
Definition: Ifpack2_Relaxation_def.hpp:488
Mix-in interface for preconditioners that can change their matrix after construction.
Definition: Ifpack2_Details_CanChangeMatrix.hpp:60
Tpetra::Operator< scalar_type, local_ordinal_type, global_ordinal_type, node_type > op_type
Tpetra::Operator specialization used by this class.
Definition: Ifpack2_Relaxation_decl.hpp:248
double getComputeFlops() const
Total number of floating-point operations over all calls to compute().
Definition: Ifpack2_Relaxation_def.hpp:523
double getApplyFlops() const
Total number of floating-point operations over all calls to apply().
Definition: Ifpack2_Relaxation_def.hpp:528
void setZeroStartingSolution(bool zeroStartingSolution)
Set this preconditioner's parameters.
Definition: Ifpack2_Relaxation_decl.hpp:377
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object's attributes to the given output stream.
Definition: Ifpack2_Relaxation_def.hpp:2194
Relaxation(const Teuchos::RCP< const row_matrix_type > &A)
Constructor.
Definition: Ifpack2_Relaxation_def.hpp:186
void compute()
Compute the preconditioner ("numeric setup");.
Definition: Ifpack2_Relaxation_def.hpp:964
size_t getNodeSmootherComplexity() const
Get a rough estimate of cost per iteration.
Definition: Ifpack2_Relaxation_def.hpp:533
MatrixType::node_type::device_type device_type
The Kokkos device type used by the input MatrixType.
Definition: Ifpack2_Relaxation_decl.hpp:233
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
The communicator over which the matrix and vectors are distributed.
Definition: Ifpack2_Relaxation_def.hpp:446
virtual void setMatrix(const Teuchos::RCP< const row_matrix_type > &A)
Change the matrix to be preconditioned.
Definition: Ifpack2_Relaxation_def.hpp:166
int getNumApply() const
Total number of calls to apply().
Definition: Ifpack2_Relaxation_def.hpp:503
double getComputeTime() const
Total time in seconds spent in all calls to compute().
Definition: Ifpack2_Relaxation_def.hpp:513
MatrixType::node_type node_type
The Node type used by the input MatrixType.
Definition: Ifpack2_Relaxation_decl.hpp:230
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_Relaxation_def.hpp:465
void setParameters(const Teuchos::ParameterList ¶ms)
Set the relaxation / preconditioner parameters.
Definition: Ifpack2_Relaxation_def.hpp:438
std::string description() const
A simple one-line description of this object.
Definition: Ifpack2_Relaxation_def.hpp:2093
int getNumInitialize() const
Total number of calls to initialize().
Definition: Ifpack2_Relaxation_def.hpp:493
Interface for all Ifpack2 preconditioners.
Definition: Ifpack2_Preconditioner.hpp:74
bool isComputed() const
Return true if compute() has been called.
Definition: Ifpack2_Relaxation_decl.hpp:395
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_Relaxation_def.hpp:478
MatrixType::global_ordinal_type global_ordinal_type
The type of global indices in the input MatrixType.
Definition: Ifpack2_Relaxation_decl.hpp:227
double getInitializeTime() const
Total time in seconds spent in all calls to initialize().
Definition: Ifpack2_Relaxation_def.hpp:508
Declaration of interface for preconditioners that can change their matrix after construction.
static const EVerbosityLevel verbLevel_default
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 input matrix to X, returning the result in Y.
Definition: Ifpack2_Relaxation_def.hpp:657
Teuchos::RCP< const row_matrix_type > getMatrix() const
The matrix to be preconditioned.
Definition: Ifpack2_Relaxation_def.hpp:457
MatrixType::local_ordinal_type local_ordinal_type
The type of local indices in the input MatrixType.
Definition: Ifpack2_Relaxation_decl.hpp:224
virtual ~Relaxation()=default
Destructor.
MatrixType::scalar_type scalar_type
The type of the entries of the input MatrixType.
Definition: Ifpack2_Relaxation_decl.hpp:221
static scalar_type zero()
int getNumCompute() const
Total number of calls to compute().
Definition: Ifpack2_Relaxation_def.hpp:498
bool isInitialized() const
Returns true if the preconditioner has been successfully initialized.
Definition: Ifpack2_Relaxation_decl.hpp:387
Relaxation preconditioners for Tpetra::RowMatrix and Tpetra::CrsMatrix sparse matrices.
Definition: Ifpack2_Relaxation_decl.hpp:206
double getApplyTime() const
Total time in seconds spent in all calls to apply().
Definition: Ifpack2_Relaxation_def.hpp:518
void initialize()
Initialize the preconditioner ("symbolic setup").
Definition: Ifpack2_Relaxation_def.hpp:678
Tpetra::RowMatrix< scalar_type, local_ordinal_type, global_ordinal_type, node_type > row_matrix_type
Tpetra::RowMatrix specialization used by this class.
Definition: Ifpack2_Relaxation_decl.hpp:241
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, returning the result in Y.
Definition: Ifpack2_Relaxation_def.hpp:545
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Return a list of all the parameters that this class accepts.
Definition: Ifpack2_Relaxation_def.hpp:195
Teuchos::ScalarTraits< scalar_type >::magnitudeType magnitude_type
The type of the magnitude (absolute value) of a matrix entry.
Definition: Ifpack2_Relaxation_decl.hpp:236