10 #ifndef IFPACK2_DETAILS_CHEBYSHEV_DECL_HPP
11 #define IFPACK2_DETAILS_CHEBYSHEV_DECL_HPP
20 #include "Ifpack2_ConfigDefs.hpp"
22 #include "Teuchos_Describable.hpp"
23 #include "Tpetra_CrsMatrix.hpp"
28 #ifndef DOXYGEN_SHOULD_SKIP_THIS
29 template <
class TpetraOperatorType>
30 class ChebyshevKernel;
31 #endif // DOXYGEN_SHOULD_SKIP_THIS
74 template <
class ScalarType,
class MV>
81 typedef ScalarType
ST;
87 typedef Tpetra::Operator<
typename MV::scalar_type,
88 typename MV::local_ordinal_type,
89 typename MV::global_ordinal_type,
90 typename MV::node_type>
93 typedef Tpetra::RowMatrix<
typename MV::scalar_type,
94 typename MV::local_ordinal_type,
95 typename MV::global_ordinal_type,
96 typename MV::node_type>
99 typedef Tpetra::Vector<
typename MV::scalar_type,
100 typename MV::local_ordinal_type,
101 typename MV::global_ordinal_type,
102 typename MV::node_type>
105 typedef Tpetra::Map<
typename MV::local_ordinal_type,
106 typename MV::global_ordinal_type,
107 typename MV::node_type>
241 void setZeroStartingSolution(
bool zeroStartingSolution) { zeroStartingSolution_ = zeroStartingSolution; }
296 ST getLambdaMaxForApply()
const;
312 void print(std::ostream& out);
355 typedef Kokkos::View<size_t*, typename MV::node_type::device_type> offsets_type;
362 offsets_type diagOffsets_;
371 bool savedDiagOffsets_;
388 ST computedLambdaMax_;
395 ST computedLambdaMin_;
403 ST lambdaMaxForApply_;
414 ST lambdaMinForApply_;
417 ST eigRatioForApply_;
460 bool eigKeepVectors_;
463 std::string eigenAnalysisType_;
470 int eigNormalizationFreq_;
473 bool zeroStartingSolution_;
481 bool assumeMatrixUnchanged_;
484 std::string chebyshevAlgorithm_;
487 bool computeMaxResNorm_;
490 bool computeSpectralRadius_;
494 bool ckUseNativeSpMV_;
499 bool preAllocateTempVector_;
514 void checkConstructorInput()
const;
517 void checkInputMatrix()
const;
545 firstIterationWithZeroStartingSolution(MV& W,
546 const ScalarType& alpha,
553 computeResidual(MV& R,
const MV& B,
const op_type& A,
const MV& X,
561 static void solve(MV& Z,
const V& D_inv,
const MV& R);
568 static void solve(MV& Z,
const ST alpha,
const V& D_inv,
const MV& R);
578 makeInverseDiagonal(
const row_matrix_type& A,
const bool useDiagOffsets =
false)
const;
614 textbookApplyImpl(
const op_type& A,
621 const V& D_inv)
const;
645 fourthKindApplyImpl(
const op_type& A,
675 ifpackApplyImpl(
const op_type& A,
696 ST cgMethodWithInitGuess(
const op_type& A,
const V& D_inv,
const int numIters,
V& x);
707 ST cgMethod(
const op_type& A,
const V& D_inv,
const int numIters);
710 static MT maxNormInf(
const MV& X);
748 const ST zero = Teuchos::as<ST> (0);
749 const ST one = Teuchos::as<ST> (1);
750 const ST two = Teuchos::as<ST> (2);
752 MV pAux (B.getMap (), B.getNumVectors ());
753 MV dk (B.getMap (), B.getNumVectors ());
754 MV R (B.getMap (), B.getNumVectors ());
756 ST beta = Teuchos::as<ST> (1.1) * lambdaMax;
757 ST alpha = lambdaMax / eigRatio;
759 ST delta = (beta - alpha) / two;
760 ST theta = (beta + alpha) / two;
761 ST s1 = theta / delta;
772 if (lambdaMin == one && lambdaMin == lambdaMax) {
781 if (! zeroStartingSolution_) {
785 R.update (-one, pAux, one);
786 dk.elementWiseMultiply (one/theta, D_inv, R, zero);
787 X.update (one, dk, one);
789 dk.elementWiseMultiply (one/theta, D_inv, B, zero);
793 ST rhokp1, dtemp1, dtemp2;
794 for (
int k = 0; k < numIters-1; ++k) {
796 rhokp1 = one / (two*s1 - rhok);
797 dtemp1 = rhokp1*rhok;
798 dtemp2 = two*rhokp1/delta;
802 R.update (-one, pAux, one);
804 dk.elementWiseMultiply (dtemp2, D_inv, B, dtemp1);
805 X.update (one, dk, one);
815 #endif // IFPACK2_DETAILS_CHEBYSHEV_DECL_HPP
void setParameters(Teuchos::ParameterList &plist)
Set (or reset) parameters.
Definition: Ifpack2_Details_Chebyshev_def.hpp:340
MT apply(const MV &B, MV &X)
Solve Ax=b for x with Chebyshev iteration with left diagonal scaling.
Definition: Ifpack2_Details_Chebyshev_def.hpp:999
Tpetra::Operator< typename MV::scalar_type, typename MV::local_ordinal_type, typename MV::global_ordinal_type, typename MV::node_type > op_type
Specialization of Tpetra::Operator.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:91
scalar_type magnitudeType
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print a description of the Chebyshev solver to out.
Definition: Ifpack2_Details_Chebyshev_def.hpp:1661
void setMatrix(const Teuchos::RCP< const row_matrix_type > &A)
Set the matrix.
Definition: Ifpack2_Details_Chebyshev_def.hpp:755
Teuchos::RCP< const row_matrix_type > getMatrix() const
Get the matrix given to the constructor.
Definition: Ifpack2_Details_Chebyshev_def.hpp:1597
Chebyshev(Teuchos::RCP< const row_matrix_type > A)
Definition: Ifpack2_Details_Chebyshev_def.hpp:275
Tpetra::Map< typename MV::local_ordinal_type, typename MV::global_ordinal_type, typename MV::node_type > map_type
Specialization of Tpetra::Map.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:108
void print(std::ostream &out)
Print instance data to the given output stream.
Definition: Ifpack2_Details_Chebyshev_def.hpp:1049
STS::magnitudeType MT
The type of the absolute value of a ScalarType.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:85
ScalarType ST
The type of entries in the matrix and vectors.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:81
Tpetra::RowMatrix< typename MV::scalar_type, typename MV::local_ordinal_type, typename MV::global_ordinal_type, typename MV::node_type > row_matrix_type
Specialization of Tpetra::RowMatrix.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:97
Left-scaled Chebyshev iteration.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:75
static const EVerbosityLevel verbLevel_default
std::string description() const
A single-line description of the Chebyshev solver.
Definition: Ifpack2_Details_Chebyshev_def.hpp:1641
Teuchos::ScalarTraits< ScalarType > STS
Traits class for ST.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:83
bool hasTransposeApply() const
Whether it's possible to apply the transpose of this operator.
Definition: Ifpack2_Details_Chebyshev_def.hpp:1603
void compute()
(Re)compute the left scaling D_inv, and estimate min and max eigenvalues of D_inv * A...
Definition: Ifpack2_Details_Chebyshev_def.hpp:791
Tpetra::Vector< typename MV::scalar_type, typename MV::local_ordinal_type, typename MV::global_ordinal_type, typename MV::node_type > V
Specialization of Tpetra::Vector.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:103