44 #ifndef IFPACK2_DETAILS_CHEBYSHEV_DEF_HPP 45 #define IFPACK2_DETAILS_CHEBYSHEV_DEF_HPP 55 #include "Kokkos_ArithTraits.hpp" 56 #include "Teuchos_FancyOStream.hpp" 57 #include "Teuchos_oblackholestream.hpp" 67 const char computeBeforeApplyReminder[] =
68 "This means one of the following:\n" 69 " - you have not yet called compute() on this instance, or \n" 70 " - you didn't call compute() after calling setParameters().\n\n" 71 "After creating an Ifpack2::Chebyshev instance,\n" 72 "you must _always_ call compute() at least once before calling apply().\n" 73 "After calling compute() once, you do not need to call it again,\n" 74 "unless the matrix has changed or you have changed parameters\n" 75 "(by calling setParameters()).";
81 template<
class XV,
class SizeType =
typename XV::
size_type>
82 struct V_ReciprocalThresholdSelfFunctor
84 typedef typename XV::execution_space execution_space;
85 typedef typename XV::non_const_value_type value_type;
86 typedef SizeType size_type;
87 typedef Kokkos::Details::ArithTraits<value_type> KAT;
88 typedef typename KAT::mag_type mag_type;
91 const value_type minVal_;
92 const mag_type minValMag_;
94 V_ReciprocalThresholdSelfFunctor (
const XV& X,
95 const value_type& min_val) :
98 minValMag_ (KAT::abs (min_val))
101 KOKKOS_INLINE_FUNCTION
102 void operator() (
const size_type& i)
const 104 const mag_type X_i_abs = KAT::abs (X_[i]);
106 if (X_i_abs < minValMag_) {
110 X_[i] = KAT::one () / X_[i];
115 template<
class XV,
class SizeType =
typename XV::
size_type>
116 struct LocalReciprocalThreshold {
118 compute (
const XV& X,
119 const typename XV::non_const_value_type& minVal)
121 typedef typename XV::execution_space execution_space;
122 Kokkos::RangePolicy<execution_space, SizeType> policy (0, X.extent (0));
123 V_ReciprocalThresholdSelfFunctor<XV, SizeType> op (X, minVal);
124 Kokkos::parallel_for (policy, op);
128 template <
class TpetraVectorType,
129 const bool classic = TpetraVectorType::node_type::classic>
130 struct GlobalReciprocalThreshold {};
132 template <
class TpetraVectorType>
133 struct GlobalReciprocalThreshold<TpetraVectorType, true> {
135 compute (TpetraVectorType& V,
136 const typename TpetraVectorType::scalar_type& min_val)
138 typedef typename TpetraVectorType::scalar_type scalar_type;
139 typedef typename TpetraVectorType::mag_type mag_type;
140 typedef Kokkos::Details::ArithTraits<scalar_type> STS;
142 const scalar_type ONE = STS::one ();
143 const mag_type min_val_abs = STS::abs (min_val);
145 Teuchos::ArrayRCP<scalar_type> V_0 = V.getDataNonConst (0);
146 const size_t lclNumRows = V.getLocalLength ();
148 for (
size_t i = 0; i < lclNumRows; ++i) {
149 const scalar_type V_0i = V_0[i];
150 if (STS::abs (V_0i) < min_val_abs) {
159 template <
class TpetraVectorType>
160 struct GlobalReciprocalThreshold<TpetraVectorType, false> {
162 compute (TpetraVectorType& X,
163 const typename TpetraVectorType::scalar_type& minVal)
165 typedef typename TpetraVectorType::impl_scalar_type value_type;
166 typedef typename TpetraVectorType::device_type::memory_space memory_space;
168 X.template sync<memory_space> ();
169 X.template modify<memory_space> ();
171 const value_type minValS =
static_cast<value_type
> (minVal);
172 auto X_0 = Kokkos::subview (X.template getLocalView<memory_space> (),
174 LocalReciprocalThreshold<decltype (X_0) >::compute (X_0, minValS);
179 template <
typename S,
typename L,
typename G,
typename N>
181 reciprocal_threshold (Tpetra::Vector<S,L,G,N>& V,
const S& minVal)
183 GlobalReciprocalThreshold<Tpetra::Vector<S,L,G,N> >::compute (V, minVal);
191 template<
class OneDViewType,
192 class LocalOrdinal =
typename OneDViewType::size_type>
193 class PositivizeVector {
194 static_assert (Kokkos::Impl::is_view<OneDViewType>::value,
195 "OneDViewType must be a 1-D Kokkos::View.");
196 static_assert (static_cast<int> (OneDViewType::rank) == 1,
197 "This functor only works with a 1-D View.");
198 static_assert (std::is_integral<LocalOrdinal>::value,
199 "The second template parameter, LocalOrdinal, " 200 "must be an integer type.");
202 PositivizeVector (
const OneDViewType& x) : x_ (x) {}
204 KOKKOS_INLINE_FUNCTION
void 205 operator () (
const LocalOrdinal& i)
const 207 typedef typename OneDViewType::non_const_value_type IST;
208 typedef Kokkos::Details::ArithTraits<IST> STS;
209 typedef Kokkos::Details::ArithTraits<typename STS::mag_type> STM;
211 if (STS::real (x_(i)) < STM::zero ()) {
222 template<
class ScalarType,
class MV>
223 void Chebyshev<ScalarType, MV>::checkInputMatrix ()
const 225 TEUCHOS_TEST_FOR_EXCEPTION(
226 ! A_.is_null () && A_->getGlobalNumRows () != A_->getGlobalNumCols (),
227 std::invalid_argument,
228 "Ifpack2::Chebyshev: The input matrix A must be square. " 229 "A has " << A_->getGlobalNumRows () <<
" rows and " 230 << A_->getGlobalNumCols () <<
" columns.");
234 if (debug_ && ! A_.is_null ()) {
235 Teuchos::RCP<const map_type> domainMap = A_->getDomainMap ();
236 Teuchos::RCP<const map_type> rangeMap = A_->getRangeMap ();
241 TEUCHOS_TEST_FOR_EXCEPTION(
242 ! domainMap->isSameAs (*rangeMap), std::invalid_argument,
243 "Ifpack2::Chebyshev: The domain Map and range Map of the matrix must be " 244 "the same (in the sense of isSameAs())." << std::endl <<
"We only check " 245 "for this in debug mode.");
250 template<
class ScalarType,
class MV>
252 Chebyshev<ScalarType, MV>::
253 checkConstructorInput ()
const 258 if (STS::isComplex) {
259 TEUCHOS_TEST_FOR_EXCEPTION
260 (
true, std::logic_error,
"Ifpack2::Chebyshev: This class' implementation " 261 "of Chebyshev iteration only works for real-valued, symmetric positive " 262 "definite matrices. However, you instantiated this class for ScalarType" 263 " = " << Teuchos::TypeNameTraits<ScalarType>::name () <<
", which is a " 264 "complex-valued type. While this may be algorithmically correct if all " 265 "of the complex numbers in the matrix have zero imaginary part, we " 266 "forbid using complex ScalarType altogether in order to remind you of " 267 "the limitations of our implementation (and of the algorithm itself).");
274 template<
class ScalarType,
class MV>
278 savedDiagOffsets_ (false),
279 computedLambdaMax_ (
STS::nan ()),
280 computedLambdaMin_ (
STS::nan ()),
281 lambdaMaxForApply_ (
STS::nan ()),
282 lambdaMinForApply_ (
STS::nan ()),
283 eigRatioForApply_ (
STS::nan ()),
284 userLambdaMax_ (
STS::nan ()),
285 userLambdaMin_ (
STS::nan ()),
287 minDiagVal_ (
STS::eps ()),
290 zeroStartingSolution_ (true),
291 assumeMatrixUnchanged_ (false),
292 textbookAlgorithm_ (false),
293 computeMaxResNorm_ (false),
296 checkConstructorInput ();
299 template<
class ScalarType,
class MV>
301 Chebyshev (Teuchos::RCP<const row_matrix_type> A, Teuchos::ParameterList& params) :
303 savedDiagOffsets_ (false),
304 computedLambdaMax_ (
STS::nan ()),
305 computedLambdaMin_ (
STS::nan ()),
306 lambdaMaxForApply_ (
STS::nan ()),
307 boostFactor_ (static_cast<
MT> (1.1)),
308 lambdaMinForApply_ (
STS::nan ()),
309 eigRatioForApply_ (
STS::nan ()),
310 userLambdaMax_ (
STS::nan ()),
311 userLambdaMin_ (
STS::nan ()),
313 minDiagVal_ (
STS::eps ()),
316 zeroStartingSolution_ (true),
317 assumeMatrixUnchanged_ (false),
318 textbookAlgorithm_ (false),
319 computeMaxResNorm_ (false),
322 checkConstructorInput ();
323 setParameters (params);
326 template<
class ScalarType,
class MV>
333 using Teuchos::rcp_const_cast;
345 const ST defaultLambdaMax = STS::nan ();
346 const ST defaultLambdaMin = STS::nan ();
353 const ST defaultEigRatio = Teuchos::as<ST> (30);
354 const MT defaultBoostFactor =
static_cast<MT> (1.1);
355 const ST defaultMinDiagVal = STS::eps ();
356 const int defaultNumIters = 1;
357 const int defaultEigMaxIters = 10;
358 const bool defaultZeroStartingSolution =
true;
359 const bool defaultAssumeMatrixUnchanged =
false;
360 const bool defaultTextbookAlgorithm =
false;
361 const bool defaultComputeMaxResNorm =
false;
362 const bool defaultDebug =
false;
368 RCP<const V> userInvDiagCopy;
369 ST lambdaMax = defaultLambdaMax;
370 ST lambdaMin = defaultLambdaMin;
371 ST eigRatio = defaultEigRatio;
372 MT boostFactor = defaultBoostFactor;
373 ST minDiagVal = defaultMinDiagVal;
374 int numIters = defaultNumIters;
375 int eigMaxIters = defaultEigMaxIters;
376 bool zeroStartingSolution = defaultZeroStartingSolution;
377 bool assumeMatrixUnchanged = defaultAssumeMatrixUnchanged;
378 bool textbookAlgorithm = defaultTextbookAlgorithm;
379 bool computeMaxResNorm = defaultComputeMaxResNorm;
380 bool debug = defaultDebug;
387 if (plist.isParameter (
"debug")) {
389 debug = plist.get<
bool> (
"debug");
391 catch (Teuchos::Exceptions::InvalidParameterType&) {}
394 int debugInt = plist.get<
bool> (
"debug");
395 debug = debugInt != 0;
397 catch (Teuchos::Exceptions::InvalidParameterType&) {}
409 if (plist.isParameter (
"chebyshev: operator inv diagonal")) {
411 RCP<const V> userInvDiag;
414 const V* rawUserInvDiag =
415 plist.get<
const V*> (
"chebyshev: operator inv diagonal");
417 userInvDiag = rcp (rawUserInvDiag,
false);
418 }
catch (Teuchos::Exceptions::InvalidParameterType&) {
420 if (userInvDiag.is_null ()) {
422 V* rawUserInvDiag = plist.get<
V*> (
"chebyshev: operator inv diagonal");
424 userInvDiag = rcp (const_cast<const V*> (rawUserInvDiag),
false);
425 }
catch (Teuchos::Exceptions::InvalidParameterType&) {
428 if (userInvDiag.is_null ()) {
431 plist.get<RCP<const V> > (
"chebyshev: operator inv diagonal");
432 }
catch (Teuchos::Exceptions::InvalidParameterType&) {
435 if (userInvDiag.is_null ()) {
437 RCP<V> userInvDiagNonConst =
438 plist.get<RCP<V> > (
"chebyshev: operator inv diagonal");
439 userInvDiag = rcp_const_cast<
const V> (userInvDiagNonConst);
440 }
catch (Teuchos::Exceptions::InvalidParameterType&) {
443 if (userInvDiag.is_null ()) {
452 const V& userInvDiagRef =
453 plist.get<
const V> (
"chebyshev: operator inv diagonal");
454 userInvDiagCopy = rcp (
new V (userInvDiagRef, Teuchos::Copy));
456 userInvDiag = userInvDiagCopy;
457 }
catch (Teuchos::Exceptions::InvalidParameterType&) {
460 TEUCHOS_TEST_FOR_EXCEPTION(
461 true, std::runtime_error,
462 "Ifpack2::Chebyshev::setParameters: \"chebyshev: operator inv diagonal\" " 463 "plist.get<const V> does not compile around return held == other_held " 464 "in Teuchos::any in Visual Studio. Can't fix it now, so throwing " 465 "in case someone builds there.");
468 if (userInvDiag.is_null ()) {
477 V& userInvDiagNonConstRef =
478 plist.get<
V> (
"chebyshev: operator inv diagonal");
479 const V& userInvDiagRef =
const_cast<const V&
> (userInvDiagNonConstRef);
480 userInvDiagCopy = rcp (
new V (userInvDiagRef, Teuchos::Copy));
482 userInvDiag = userInvDiagCopy;
483 }
catch (Teuchos::Exceptions::InvalidParameterType&) {
486 TEUCHOS_TEST_FOR_EXCEPTION(
487 true, std::runtime_error,
488 "Ifpack2::Chebyshev::setParameters: \"chebyshev: operator inv diagonal\" " 489 "plist.get<V> does not compile around return held == other_held " 490 "in Teuchos::any in Visual Studio. Can't fix it now, so throwing " 491 "in case someone builds there.");
502 if (! userInvDiag.is_null () && userInvDiagCopy.is_null ()) {
503 userInvDiagCopy = rcp (
new V (*userInvDiag, Teuchos::Copy));
515 if (plist.isParameter (
"chebyshev: max eigenvalue")) {
516 if (plist.isType<
double>(
"chebyshev: max eigenvalue"))
517 lambdaMax = plist.get<
double> (
"chebyshev: max eigenvalue");
519 lambdaMax = plist.get<
ST> (
"chebyshev: max eigenvalue");
520 TEUCHOS_TEST_FOR_EXCEPTION(
521 STS::isnaninf (lambdaMax), std::invalid_argument,
522 "Ifpack2::Chebyshev::setParameters: \"chebyshev: max eigenvalue\" " 523 "parameter is NaN or Inf. This parameter is optional, but if you " 524 "choose to supply it, it must have a finite value.");
526 if (plist.isParameter (
"chebyshev: min eigenvalue")) {
527 if (plist.isType<
double>(
"chebyshev: min eigenvalue"))
528 lambdaMin = plist.get<
double> (
"chebyshev: min eigenvalue");
530 lambdaMin = plist.get<
ST> (
"chebyshev: min eigenvalue");
531 TEUCHOS_TEST_FOR_EXCEPTION(
532 STS::isnaninf (lambdaMin), std::invalid_argument,
533 "Ifpack2::Chebyshev::setParameters: \"chebyshev: min eigenvalue\" " 534 "parameter is NaN or Inf. This parameter is optional, but if you " 535 "choose to supply it, it must have a finite value.");
539 if (plist.isParameter (
"smoother: Chebyshev alpha")) {
540 if (plist.isType<
double>(
"smoother: Chebyshev alpha"))
541 eigRatio = plist.get<
double> (
"smoother: Chebyshev alpha");
543 eigRatio = plist.get<
ST> (
"smoother: Chebyshev alpha");
546 eigRatio = plist.get (
"chebyshev: ratio eigenvalue", eigRatio);
547 TEUCHOS_TEST_FOR_EXCEPTION(
548 STS::isnaninf (eigRatio), std::invalid_argument,
549 "Ifpack2::Chebyshev::setParameters: \"chebyshev: ratio eigenvalue\" " 550 "parameter (also called \"smoother: Chebyshev alpha\") is NaN or Inf. " 551 "This parameter is optional, but if you choose to supply it, it must have " 558 TEUCHOS_TEST_FOR_EXCEPTION(
559 STS::real (eigRatio) < STS::real (STS::one ()),
560 std::invalid_argument,
561 "Ifpack2::Chebyshev::setParameters: \"chebyshev: ratio eigenvalue\"" 562 "parameter (also called \"smoother: Chebyshev alpha\") must be >= 1, " 563 "but you supplied the value " << eigRatio <<
".");
568 const char paramName[] =
"chebyshev: boost factor";
570 if (plist.isParameter (paramName)) {
571 if (plist.isType<
MT> (paramName)) {
572 boostFactor = plist.get<
MT> (paramName);
574 else if (! std::is_same<double, MT>::value &&
575 plist.isType<
double> (paramName)) {
576 const double dblBF = plist.get<
double> (paramName);
577 boostFactor =
static_cast<MT> (dblBF);
580 TEUCHOS_TEST_FOR_EXCEPTION
581 (
true, std::invalid_argument,
582 "Ifpack2::Chebyshev::setParameters: \"chebyshev: boost factor\"" 583 "parameter must have type magnitude_type (MT) or double.");
593 plist.set (paramName, defaultBoostFactor);
595 TEUCHOS_TEST_FOR_EXCEPTION
596 (boostFactor < Teuchos::ScalarTraits<MT>::one (), std::invalid_argument,
597 "Ifpack2::Chebyshev::setParameters: \"" << paramName <<
"\" parameter " 598 "must be >= 1, but you supplied the value " << boostFactor <<
".");
602 minDiagVal = plist.get (
"chebyshev: min diagonal value", minDiagVal);
603 TEUCHOS_TEST_FOR_EXCEPTION(
604 STS::isnaninf (minDiagVal), std::invalid_argument,
605 "Ifpack2::Chebyshev::setParameters: \"chebyshev: min diagonal value\" " 606 "parameter is NaN or Inf. This parameter is optional, but if you choose " 607 "to supply it, it must have a finite value.");
610 if (plist.isParameter (
"smoother: sweeps")) {
611 numIters = plist.get<
int> (
"smoother: sweeps");
613 if (plist.isParameter (
"relaxation: sweeps")) {
614 numIters = plist.get<
int> (
"relaxation: sweeps");
616 numIters = plist.get (
"chebyshev: degree", numIters);
617 TEUCHOS_TEST_FOR_EXCEPTION(
618 numIters < 0, std::invalid_argument,
619 "Ifpack2::Chebyshev::setParameters: \"chebyshev: degree\" parameter (also " 620 "called \"smoother: sweeps\" or \"relaxation: sweeps\") must be a " 621 "nonnegative integer. You gave a value of " << numIters <<
".");
624 if (plist.isParameter (
"eigen-analysis: iterations")) {
625 eigMaxIters = plist.get<
int> (
"eigen-analysis: iterations");
627 eigMaxIters = plist.get (
"chebyshev: eigenvalue max iterations", eigMaxIters);
628 TEUCHOS_TEST_FOR_EXCEPTION(
629 eigMaxIters < 0, std::invalid_argument,
630 "Ifpack2::Chebyshev::setParameters: \"chebyshev: eigenvalue max iterations" 631 "\" parameter (also called \"eigen-analysis: iterations\") must be a " 632 "nonnegative integer. You gave a value of " << eigMaxIters <<
".");
634 zeroStartingSolution = plist.get (
"chebyshev: zero starting solution",
635 zeroStartingSolution);
636 assumeMatrixUnchanged = plist.get (
"chebyshev: assume matrix does not change",
637 assumeMatrixUnchanged);
641 if (plist.isParameter (
"chebyshev: textbook algorithm")) {
642 textbookAlgorithm = plist.get<
bool> (
"chebyshev: textbook algorithm");
644 if (plist.isParameter (
"chebyshev: compute max residual norm")) {
645 computeMaxResNorm = plist.get<
bool> (
"chebyshev: compute max residual norm");
651 TEUCHOS_TEST_FOR_EXCEPTION(
652 plist.isParameter (
"chebyshev: use block mode") &&
653 ! plist.get<
bool> (
"chebyshev: use block mode"),
654 std::invalid_argument,
655 "Ifpack2::Chebyshev requires that if you supply the Ifpack parameter " 656 "\"chebyshev: use block mode\", it must be set to false. Ifpack2's " 657 "Chebyshev does not implement Ifpack's block mode.");
658 TEUCHOS_TEST_FOR_EXCEPTION(
659 plist.isParameter (
"chebyshev: solve normal equations") &&
660 ! plist.get<
bool> (
"chebyshev: solve normal equations"),
661 std::invalid_argument,
662 "Ifpack2::Chebyshev does not and will never implement the Ifpack " 663 "parameter \"chebyshev: solve normal equations\". If you want to solve " 664 "the normal equations, construct a Tpetra::Operator that implements " 665 "A^* A, and use Chebyshev to solve A^* A x = A^* b.");
673 std::string eigenAnalysisType (
"power-method");
674 if (plist.isParameter (
"eigen-analysis: type")) {
675 eigenAnalysisType = plist.get<std::string> (
"eigen-analysis: type");
676 TEUCHOS_TEST_FOR_EXCEPTION(
677 eigenAnalysisType !=
"power-method" &&
678 eigenAnalysisType !=
"power method",
679 std::invalid_argument,
680 "Ifpack2::Chebyshev: This class supports the ML parameter \"eigen-" 681 "analysis: type\" for backwards compatibility. However, Ifpack2 " 682 "currently only supports the \"power-method\" option for this " 683 "parameter. This imitates Ifpack, which only implements the power " 684 "method for eigenanalysis.");
688 userInvDiag_ = userInvDiagCopy;
689 userLambdaMax_ = lambdaMax;
690 userLambdaMin_ = lambdaMin;
691 userEigRatio_ = eigRatio;
692 boostFactor_ =
static_cast<MT> (boostFactor);
693 minDiagVal_ = minDiagVal;
694 numIters_ = numIters;
695 eigMaxIters_ = eigMaxIters;
696 zeroStartingSolution_ = zeroStartingSolution;
697 assumeMatrixUnchanged_ = assumeMatrixUnchanged;
698 textbookAlgorithm_ = textbookAlgorithm;
699 computeMaxResNorm_ = computeMaxResNorm;
705 if (A_.is_null () || A_->getComm ().is_null ()) {
711 myRank = A_->getComm ()->getRank ();
715 out_ = Teuchos::getFancyOStream (Teuchos::rcpFromRef (std::cerr));
718 Teuchos::RCP<Teuchos::oblackholestream> blackHole (
new Teuchos::oblackholestream ());
719 out_ = Teuchos::getFancyOStream (blackHole);
724 out_ = Teuchos::null;
729 template<
class ScalarType,
class MV>
734 diagOffsets_ = offsets_type ();
735 savedDiagOffsets_ =
false;
738 computedLambdaMax_ = STS::nan ();
739 computedLambdaMin_ = STS::nan ();
743 template<
class ScalarType,
class MV>
747 if (A.getRawPtr () != A_.getRawPtr ()) {
748 if (! assumeMatrixUnchanged_) {
759 if (A.is_null () || A->getComm ().is_null ()) {
765 myRank = A->getComm ()->getRank ();
769 out_ = Teuchos::getFancyOStream (Teuchos::rcpFromRef (std::cerr));
772 Teuchos::RCP<Teuchos::oblackholestream> blackHole (
new Teuchos::oblackholestream ());
773 out_ = Teuchos::getFancyOStream (blackHole);
778 out_ = Teuchos::null;
784 template<
class ScalarType,
class MV>
793 typedef Tpetra::CrsMatrix<
typename MV::scalar_type,
794 typename MV::local_ordinal_type,
795 typename MV::global_ordinal_type,
796 typename MV::node_type> crs_matrix_type;
798 TEUCHOS_TEST_FOR_EXCEPTION(
799 A_.is_null (), std::runtime_error,
"Ifpack2::Chebyshev::compute: The input " 800 "matrix A is null. Please call setMatrix() with a nonnull input matrix " 801 "before calling this method.");
816 if (userInvDiag_.is_null ()) {
817 Teuchos::RCP<const crs_matrix_type> A_crsMat =
818 Teuchos::rcp_dynamic_cast<
const crs_matrix_type> (A_);
821 if (! A_crsMat.is_null () && A_crsMat->isStaticGraph ()) {
823 const size_t lclNumRows = A_crsMat->getNodeNumRows ();
824 if (diagOffsets_.extent (0) < lclNumRows) {
825 diagOffsets_ = offsets_type ();
826 diagOffsets_ = offsets_type (
"offsets", lclNumRows);
828 A_crsMat->getCrsGraph ()->getLocalDiagOffsets (diagOffsets_);
829 savedDiagOffsets_ =
true;
830 D_ = makeInverseDiagonal (*A_,
true);
833 D_ = makeInverseDiagonal (*A_);
836 else if (! assumeMatrixUnchanged_) {
837 if (! A_crsMat.is_null () && A_crsMat->isStaticGraph ()) {
840 if (! savedDiagOffsets_) {
841 const size_t lclNumRows = A_crsMat->getNodeNumRows ();
842 if (diagOffsets_.extent (0) < lclNumRows) {
843 diagOffsets_ = offsets_type ();
844 diagOffsets_ = offsets_type (
"offsets", lclNumRows);
846 A_crsMat->getCrsGraph ()->getLocalDiagOffsets (diagOffsets_);
847 savedDiagOffsets_ =
true;
850 D_ = makeInverseDiagonal (*A_,
true);
853 D_ = makeInverseDiagonal (*A_);
858 D_ = makeRangeMapVectorConst (userInvDiag_);
862 const bool computedEigenvalueEstimates =
863 STS::isnaninf (computedLambdaMax_) || STS::isnaninf (computedLambdaMin_);
872 if (! assumeMatrixUnchanged_ ||
873 (! computedEigenvalueEstimates && STS::isnaninf (userLambdaMax_))) {
874 const ST computedLambdaMax = powerMethod (*A_, *D_, eigMaxIters_);
875 TEUCHOS_TEST_FOR_EXCEPTION(
876 STS::isnaninf (computedLambdaMax),
878 "Ifpack2::Chebyshev::compute: Estimation of the max eigenvalue " 879 "of D^{-1} A failed, by producing Inf or NaN. This probably means that " 880 "the matrix contains Inf or NaN values, or that it is badly scaled.");
881 TEUCHOS_TEST_FOR_EXCEPTION(
882 STS::isnaninf (userEigRatio_),
884 "Ifpack2::Chebyshev::compute: userEigRatio_ is Inf or NaN." 885 << endl <<
"This should be impossible." << endl <<
886 "Please report this bug to the Ifpack2 developers.");
892 const ST computedLambdaMin = computedLambdaMax / userEigRatio_;
895 computedLambdaMax_ = computedLambdaMax;
896 computedLambdaMin_ = computedLambdaMin;
898 TEUCHOS_TEST_FOR_EXCEPTION(
899 STS::isnaninf (userLambdaMax_) && STS::isnaninf (computedLambdaMax_),
901 "Ifpack2::Chebyshev::compute: " << endl <<
902 "Both userLambdaMax_ and computedLambdaMax_ are Inf or NaN." 903 << endl <<
"This should be impossible." << endl <<
904 "Please report this bug to the Ifpack2 developers.");
912 if (STS::isnaninf (userLambdaMax_)) {
913 lambdaMaxForApply_ = computedLambdaMax_;
915 lambdaMaxForApply_ = userLambdaMax_;
928 lambdaMinForApply_ = lambdaMaxForApply_ / userEigRatio_;
929 eigRatioForApply_ = userEigRatio_;
931 if (! textbookAlgorithm_) {
934 const ST one = Teuchos::as<ST> (1);
937 if (STS::magnitude (lambdaMaxForApply_ - one) < Teuchos::as<MT> (1.0e-6)) {
938 lambdaMinForApply_ = one;
939 lambdaMaxForApply_ = lambdaMinForApply_;
940 eigRatioForApply_ = one;
946 template<
class ScalarType,
class MV>
950 return lambdaMaxForApply_;
954 template<
class ScalarType,
class MV>
958 const char prefix[] =
"Ifpack2::Chebyshev::apply: ";
961 *out_ <<
"apply: " << std::endl;
963 TEUCHOS_TEST_FOR_EXCEPTION
964 (A_.is_null (), std::runtime_error, prefix <<
"The input matrix A is null. " 965 " Please call setMatrix() with a nonnull input matrix, and then call " 966 "compute(), before calling this method.");
967 TEUCHOS_TEST_FOR_EXCEPTION
968 (STS::isnaninf (lambdaMaxForApply_), std::runtime_error,
969 prefix <<
"There is no estimate for the max eigenvalue." 970 << std::endl << computeBeforeApplyReminder);
971 TEUCHOS_TEST_FOR_EXCEPTION
972 (STS::isnaninf (lambdaMinForApply_), std::runtime_error,
973 prefix <<
"There is no estimate for the min eigenvalue." 974 << std::endl << computeBeforeApplyReminder);
975 TEUCHOS_TEST_FOR_EXCEPTION
976 (STS::isnaninf (eigRatioForApply_), std::runtime_error,
977 prefix <<
"There is no estimate for the ratio of the max " 978 "eigenvalue to the min eigenvalue." 979 << std::endl << computeBeforeApplyReminder);
980 TEUCHOS_TEST_FOR_EXCEPTION
981 (D_.is_null (), std::runtime_error, prefix <<
"The vector of inverse " 982 "diagonal entries of the matrix has not yet been computed." 983 << std::endl << computeBeforeApplyReminder);
985 if (textbookAlgorithm_) {
986 textbookApplyImpl (*A_, B, X, numIters_, lambdaMaxForApply_,
987 lambdaMinForApply_, eigRatioForApply_, *D_);
990 ifpackApplyImpl (*A_, B, X, numIters_, lambdaMaxForApply_,
991 lambdaMinForApply_, eigRatioForApply_, *D_);
994 if (computeMaxResNorm_ && B.getNumVectors () > 0) {
995 MV R (B.getMap (), B.getNumVectors ());
996 computeResidual (R, B, *A_, X);
997 Teuchos::Array<MT> norms (B.getNumVectors ());
999 return *std::max_element (norms.begin (), norms.end ());
1002 return Teuchos::ScalarTraits<MT>::zero ();
1006 template<
class ScalarType,
class MV>
1010 this->describe (* (Teuchos::getFancyOStream (Teuchos::rcpFromRef (out))),
1011 Teuchos::VERB_MEDIUM);
1014 template<
class ScalarType,
class MV>
1018 const Teuchos::ETransp mode)
1020 Tpetra::deep_copy(R, B);
1021 A.apply (X, R, mode, -STS::one(), STS::one());
1024 template<
class ScalarType,
class MV>
1027 solve (MV& Z,
const V& D_inv,
const MV& R) {
1028 Z.elementWiseMultiply (STS::one(), D_inv, R, STS::zero());
1031 template<
class ScalarType,
class MV>
1033 Chebyshev<ScalarType, MV>::
1034 solve (MV& Z,
const ST alpha,
const V& D_inv,
const MV& R) {
1035 Z.elementWiseMultiply (alpha, D_inv, R, STS::zero());
1038 template<
class ScalarType,
class MV>
1039 Teuchos::RCP<const typename Chebyshev<ScalarType, MV>::V>
1040 Chebyshev<ScalarType, MV>::
1041 makeInverseDiagonal (
const row_matrix_type& A,
const bool useDiagOffsets)
const 1044 using Teuchos::rcpFromRef;
1045 using Teuchos::rcp_dynamic_cast;
1047 RCP<V> D_rowMap (
new V (A.getGraph ()->getRowMap ()));
1048 if (useDiagOffsets) {
1052 typedef Tpetra::CrsMatrix<
typename MV::scalar_type,
1053 typename MV::local_ordinal_type,
1054 typename MV::global_ordinal_type,
1055 typename MV::node_type> crs_matrix_type;
1056 RCP<const crs_matrix_type> A_crsMat =
1057 rcp_dynamic_cast<
const crs_matrix_type> (rcpFromRef (A));
1058 if (! A_crsMat.is_null ()) {
1059 TEUCHOS_TEST_FOR_EXCEPTION(
1060 ! savedDiagOffsets_, std::logic_error,
1061 "Ifpack2::Details::Chebyshev::makeInverseDiagonal: " 1062 "It is not allowed to call this method with useDiagOffsets=true, " 1063 "if you have not previously saved offsets of diagonal entries. " 1064 "This situation should never arise if this class is used properly. " 1065 "Please report this bug to the Ifpack2 developers.");
1066 A_crsMat->getLocalDiagCopy (*D_rowMap, diagOffsets_);
1072 A.getLocalDiagCopy (*D_rowMap);
1074 RCP<V> D_rangeMap = makeRangeMapVector (D_rowMap);
1080 D_rangeMap->template sync<Kokkos::HostSpace> ();
1081 auto D_lcl = D_rangeMap->template getLocalView<Kokkos::HostSpace> ();
1082 auto D_lcl_1d = Kokkos::subview (D_lcl, Kokkos::ALL (), 0);
1084 typedef typename MV::impl_scalar_type IST;
1085 typedef typename MV::local_ordinal_type LO;
1086 typedef Kokkos::Details::ArithTraits<IST> STS;
1087 typedef Kokkos::Details::ArithTraits<typename STS::mag_type> STM;
1089 const LO lclNumRows =
static_cast<LO
> (D_rangeMap->getLocalLength ());
1090 bool foundNonpositiveValue =
false;
1091 for (LO i = 0; i < lclNumRows; ++i) {
1092 if (STS::real (D_lcl_1d(i)) <= STM::zero ()) {
1093 foundNonpositiveValue =
true;
1098 using Teuchos::outArg;
1099 using Teuchos::REDUCE_MIN;
1100 using Teuchos::reduceAll;
1102 const int lclSuccess = foundNonpositiveValue ? 0 : 1;
1103 int gblSuccess = lclSuccess;
1104 if (! D_rangeMap->getMap ().is_null () && D_rangeMap->getMap ()->getComm ().is_null ()) {
1105 const Teuchos::Comm<int>& comm = * (D_rangeMap->getMap ()->getComm ());
1106 reduceAll<int, int> (comm, REDUCE_MIN, lclSuccess, outArg (gblSuccess));
1108 if (gblSuccess == 1) {
1109 *out_ <<
"makeInverseDiagonal: The matrix's diagonal entries all have " 1110 "positive real part (this is good for Chebyshev)." << std::endl;
1113 *out_ <<
"makeInverseDiagonal: The matrix's diagonal has at least one " 1114 "entry with nonpositive real part, on at least one process in the " 1115 "matrix's communicator. This is bad for Chebyshev." << std::endl;
1121 reciprocal_threshold (*D_rangeMap, minDiagVal_);
1122 return Teuchos::rcp_const_cast<
const V> (D_rangeMap);
1126 template<
class ScalarType,
class MV>
1127 Teuchos::RCP<const typename Chebyshev<ScalarType, MV>::V>
1128 Chebyshev<ScalarType, MV>::
1129 makeRangeMapVectorConst (
const Teuchos::RCP<const V>& D)
const 1133 typedef Tpetra::Export<
typename MV::local_ordinal_type,
1134 typename MV::global_ordinal_type,
1135 typename MV::node_type> export_type;
1139 TEUCHOS_TEST_FOR_EXCEPTION(
1140 A_.is_null (), std::logic_error,
"Ifpack2::Details::Chebyshev::" 1141 "makeRangeMapVector: The input matrix A is null. Please call setMatrix() " 1142 "with a nonnull input matrix before calling this method. This is probably " 1143 "a bug in Ifpack2; please report this bug to the Ifpack2 developers.");
1144 TEUCHOS_TEST_FOR_EXCEPTION(
1145 D.is_null (), std::logic_error,
"Ifpack2::Details::Chebyshev::" 1146 "makeRangeMapVector: The input Vector D is null. This is probably " 1147 "a bug in Ifpack2; please report this bug to the Ifpack2 developers.");
1149 RCP<const map_type> sourceMap = D->getMap ();
1150 RCP<const map_type> rangeMap = A_->getRangeMap ();
1151 RCP<const map_type> rowMap = A_->getRowMap ();
1153 if (rangeMap->isSameAs (*sourceMap)) {
1159 RCP<const export_type> exporter;
1163 if (sourceMap->isSameAs (*rowMap)) {
1165 exporter = A_->getGraph ()->getExporter ();
1168 exporter = rcp (
new export_type (sourceMap, rangeMap));
1171 if (exporter.is_null ()) {
1175 RCP<V> D_out = rcp (
new V (*D, Teuchos::Copy));
1176 D_out->doExport (*D, *exporter, Tpetra::ADD);
1177 return Teuchos::rcp_const_cast<
const V> (D_out);
1183 template<
class ScalarType,
class MV>
1184 Teuchos::RCP<typename Chebyshev<ScalarType, MV>::V>
1185 Chebyshev<ScalarType, MV>::
1186 makeRangeMapVector (
const Teuchos::RCP<V>& D)
const 1188 using Teuchos::rcp_const_cast;
1189 return rcp_const_cast<V> (makeRangeMapVectorConst (rcp_const_cast<V> (D)));
1193 template<
class ScalarType,
class MV>
1195 Chebyshev<ScalarType, MV>::
1196 textbookApplyImpl (
const op_type& A,
1203 const V& D_inv)
const 1206 const ST myLambdaMin = lambdaMax / eigRatio;
1208 const ST zero = Teuchos::as<ST> (0);
1209 const ST one = Teuchos::as<ST> (1);
1210 const ST two = Teuchos::as<ST> (2);
1211 const ST d = (lambdaMax + myLambdaMin) / two;
1212 const ST c = (lambdaMax - myLambdaMin) / two;
1214 if (zeroStartingSolution_ && numIters > 0) {
1218 MV R (B.getMap (), B.getNumVectors (),
false);
1219 MV P (B.getMap (), B.getNumVectors (),
false);
1220 MV Z (B.getMap (), B.getNumVectors (),
false);
1222 for (
int i = 0; i < numIters; ++i) {
1223 computeResidual (R, B, A, X);
1224 solve (Z, D_inv, R);
1232 beta = alpha * (c/two) * (c/two);
1233 alpha = one / (d - beta);
1234 P.update (one, Z, beta);
1236 X.update (alpha, P, one);
1243 template<
class ScalarType,
class MV>
1244 typename Chebyshev<ScalarType, MV>::MT
1245 Chebyshev<ScalarType, MV>::maxNormInf (
const MV& X) {
1246 Teuchos::Array<MT> norms (X.getNumVectors ());
1247 X.normInf (norms());
1248 return *std::max_element (norms.begin (), norms.end ());
1251 template<
class ScalarType,
class MV>
1253 Chebyshev<ScalarType, MV>::
1254 ifpackApplyImpl (
const op_type& A,
1264 #ifdef HAVE_IFPACK2_DEBUG 1265 const bool debug = debug_;
1267 const bool debug =
false;
1271 *out_ <<
" \\|B\\|_{\\infty} = " << maxNormInf (B) << endl;
1272 *out_ <<
" \\|X\\|_{\\infty} = " << maxNormInf (X) << endl;
1275 if (numIters <= 0) {
1278 const ST one =
static_cast<ST
> (1.0);
1279 const ST two =
static_cast<ST
> (2.0);
1282 if (lambdaMin == one && lambdaMax == lambdaMin) {
1283 solve (X, D_inv, B);
1288 const ST alpha = lambdaMax / eigRatio;
1289 const ST beta = boostFactor_ * lambdaMax;
1290 const ST delta = two / (beta - alpha);
1291 const ST theta = (beta + alpha) / two;
1292 const ST s1 = theta * delta;
1295 *out_ <<
" alpha = " << alpha << endl
1296 <<
" beta = " << beta << endl
1297 <<
" delta = " << delta << endl
1298 <<
" theta = " << theta << endl
1299 <<
" s1 = " << s1 << endl;
1303 Teuchos::RCP<MV> V_ptr, W_ptr;
1304 makeTempMultiVectors (V_ptr, W_ptr, B);
1314 *out_ <<
" Iteration " << 1 <<
":" << endl
1315 <<
" - \\|D\\|_{\\infty} = " << D_->normInf () << endl;
1319 if (! zeroStartingSolution_) {
1320 computeResidual (V1, B, A, X);
1322 *out_ <<
" - \\|B - A*X\\|_{\\infty} = " << maxNormInf (V1) << endl;
1325 solve (W, one/theta, D_inv, V1);
1327 *out_ <<
" - \\|W\\|_{\\infty} = " << maxNormInf (W) << endl;
1329 X.update (one, W, one);
1332 solve (W, one/theta, D_inv, B);
1334 *out_ <<
" - \\|W\\|_{\\infty} = " << maxNormInf (W) << endl;
1336 Tpetra::deep_copy (X, W);
1339 *out_ <<
" - \\|X\\|_{\\infty} = " << maxNormInf (X) << endl;
1344 ST rhokp1, dtemp1, dtemp2;
1345 for (
int deg = 1; deg < numIters; ++deg) {
1347 *out_ <<
" Iteration " << deg+1 <<
":" << endl
1348 <<
" - \\|D\\|_{\\infty} = " << D_->normInf () << endl
1349 <<
" - \\|B\\|_{\\infty} = " << maxNormInf (B) << endl
1350 <<
" - \\|A\\|_{\\text{frob}} = " << A_->getFrobeniusNorm () << endl
1351 <<
" - rhok = " << rhok << endl;
1352 V1.putScalar (STS::zero ());
1355 computeResidual (V1, B, A, X);
1357 *out_ <<
" - \\|B - A*X\\|_{\\infty} = " << maxNormInf (V1) << endl;
1360 rhokp1 = one / (two * s1 - rhok);
1361 dtemp1 = rhokp1 * rhok;
1362 dtemp2 = two * rhokp1 * delta;
1366 *out_ <<
" - dtemp1 = " << dtemp1 << endl
1367 <<
" - dtemp2 = " << dtemp2 << endl;
1371 W.elementWiseMultiply (dtemp2, D_inv, V1, dtemp1);
1372 X.update (one, W, one);
1375 *out_ <<
" - \\|W\\|_{\\infty} = " << maxNormInf (W) << endl;
1376 *out_ <<
" - \\|X\\|_{\\infty} = " << maxNormInf (X) << endl;
1381 template<
class ScalarType,
class MV>
1382 typename Chebyshev<ScalarType, MV>::ST
1383 Chebyshev<ScalarType, MV>::
1384 powerMethodWithInitGuess (
const op_type& A,
1391 *out_ <<
" powerMethodWithInitGuess:" << endl;
1394 const ST zero =
static_cast<ST
> (0.0);
1395 const ST one =
static_cast<ST
> (1.0);
1396 ST lambdaMax = zero;
1397 ST RQ_top, RQ_bottom, norm;
1399 V y (A.getRangeMap ());
1401 TEUCHOS_TEST_FOR_EXCEPTION
1402 (norm == zero, std::runtime_error,
1403 "Ifpack2::Chebyshev::powerMethodWithInitGuess: " 1404 "The initial guess has zero norm. This could be either because Tpetra::" 1405 "Vector::randomize() filled the vector with zeros (if that was used to " 1406 "compute the initial guess), or because the norm2 method has a bug. The " 1407 "first is not impossible, but unlikely.");
1410 *out_ <<
" Original norm1(x): " << x.norm1 () <<
", norm2(x): " << norm << endl;
1413 x.scale (one / norm);
1416 *out_ <<
" norm1(x.scale(one/norm)): " << x.norm1 () << endl;
1419 for (
int iter = 0; iter < numIters; ++iter) {
1421 *out_ <<
" Iteration " << iter << endl;
1424 solve (y, D_inv, y);
1426 RQ_bottom = x.dot (x);
1428 *out_ <<
" RQ_top: " << RQ_top <<
", RQ_bottom: " << RQ_bottom << endl;
1430 lambdaMax = RQ_top / RQ_bottom;
1434 *out_ <<
" norm is zero; returning zero" << endl;
1438 x.update (one / norm, y, zero);
1441 *out_ <<
" lambdaMax: " << lambdaMax << endl;
1446 template<
class ScalarType,
class MV>
1448 Chebyshev<ScalarType, MV>::
1449 computeInitialGuessForPowerMethod (V& x,
const bool nonnegativeRealParts)
const 1451 typedef typename MV::device_type::execution_space dev_execution_space;
1452 typedef typename MV::device_type::memory_space dev_memory_space;
1453 typedef typename MV::local_ordinal_type LO;
1457 if (nonnegativeRealParts) {
1458 x.template modify<dev_memory_space> ();
1459 auto x_lcl = x.template getLocalView<dev_memory_space> ();
1460 auto x_lcl_1d = Kokkos::subview (x_lcl, Kokkos::ALL (), 0);
1462 const LO lclNumRows =
static_cast<LO
> (x.getLocalLength ());
1463 Kokkos::RangePolicy<dev_execution_space, LO> range (0, lclNumRows);
1464 PositivizeVector<decltype (x_lcl_1d), LO> functor (x_lcl_1d);
1465 Kokkos::parallel_for (range, functor);
1469 template<
class ScalarType,
class MV>
1470 typename Chebyshev<ScalarType, MV>::ST
1471 Chebyshev<ScalarType, MV>::
1472 powerMethod (
const op_type& A,
const V& D_inv,
const int numIters)
1476 *out_ <<
"powerMethod:" << endl;
1479 const ST zero =
static_cast<ST
> (0.0);
1480 V x (A.getDomainMap ());
1484 computeInitialGuessForPowerMethod (x,
false);
1486 ST lambdaMax = powerMethodWithInitGuess (A, D_inv, numIters, x);
1495 if (STS::real (lambdaMax) < STS::real (zero)) {
1497 *out_ <<
"real(lambdaMax) = " << STS::real (lambdaMax) <<
" < 0; try " 1498 "again with a different random initial guess" << endl;
1507 computeInitialGuessForPowerMethod (x,
true);
1508 lambdaMax = powerMethodWithInitGuess (A, D_inv, numIters, x);
1513 template<
class ScalarType,
class MV>
1514 Teuchos::RCP<const typename Chebyshev<ScalarType, MV>::row_matrix_type>
1519 template<
class ScalarType,
class MV>
1527 template<
class ScalarType,
class MV>
1531 Teuchos::RCP<MV>& W,
1537 if (V_.is_null () || V_->getNumVectors () != X.getNumVectors ()) {
1538 V_ = Teuchos::rcp (
new MV (X.getMap (), X.getNumVectors (),
false));
1541 if (W_.is_null () || W_->getNumVectors () != X.getNumVectors ()) {
1542 W_ = Teuchos::rcp (
new MV (X.getMap (), X.getNumVectors (),
true));
1548 template<
class ScalarType,
class MV>
1550 Chebyshev<ScalarType, MV>::
1551 description ()
const {
1552 std::ostringstream oss;
1555 oss <<
"\"Ifpack2::Details::Chebyshev\":" 1557 <<
"degree: " << numIters_
1558 <<
", lambdaMax: " << lambdaMaxForApply_
1559 <<
", alpha: " << eigRatioForApply_
1560 <<
", lambdaMin: " << lambdaMinForApply_
1561 <<
", boost factor: " << boostFactor_
1566 template<
class ScalarType,
class MV>
1570 const Teuchos::EVerbosityLevel verbLevel)
const 1572 using Teuchos::TypeNameTraits;
1575 const Teuchos::EVerbosityLevel vl =
1576 (verbLevel == Teuchos::VERB_DEFAULT) ? Teuchos::VERB_LOW : verbLevel;
1577 if (vl == Teuchos::VERB_NONE) {
1587 Teuchos::OSTab tab0 (out);
1593 if (A_.is_null () || A_->getComm ().is_null ()) {
1597 myRank = A_->getComm ()->getRank ();
1602 out <<
"\"Ifpack2::Details::Chebyshev\":" << endl;
1604 Teuchos::OSTab tab1 (out);
1606 if (vl == Teuchos::VERB_LOW) {
1608 out <<
"degree: " << numIters_ << endl
1609 <<
"lambdaMax: " << lambdaMaxForApply_ << endl
1610 <<
"alpha: " << eigRatioForApply_ << endl
1611 <<
"lambdaMin: " << lambdaMinForApply_ << endl
1612 <<
"boost factor: " << boostFactor_ << endl;
1620 out <<
"Template parameters:" << endl;
1622 Teuchos::OSTab tab2 (out);
1623 out <<
"ScalarType: " << TypeNameTraits<ScalarType>::name () << endl
1624 <<
"MV: " << TypeNameTraits<MV>::name () << endl;
1630 out <<
"Computed parameters:" << endl;
1636 Teuchos::OSTab tab2 (out);
1642 if (D_.is_null ()) {
1644 out <<
"unset" << endl;
1647 else if (vl <= Teuchos::VERB_HIGH) {
1649 out <<
"set" << endl;
1658 D_->describe (out, vl);
1663 out <<
"V_: " << (V_.is_null () ?
"unset" :
"set") << endl
1664 <<
"W_: " << (W_.is_null () ?
"unset" :
"set") << endl
1665 <<
"computedLambdaMax_: " << computedLambdaMax_ << endl
1666 <<
"computedLambdaMin_: " << computedLambdaMin_ << endl
1667 <<
"lambdaMaxForApply_: " << lambdaMaxForApply_ << endl
1668 <<
"lambdaMinForApply_: " << lambdaMinForApply_ << endl
1669 <<
"eigRatioForApply_: " << eigRatioForApply_ << endl;
1674 out <<
"User parameters:" << endl;
1679 Teuchos::OSTab tab2 (out);
1680 out <<
"userInvDiag_: ";
1681 if (userInvDiag_.is_null ()) {
1682 out <<
"unset" << endl;
1684 else if (vl <= Teuchos::VERB_HIGH) {
1685 out <<
"set" << endl;
1691 userInvDiag_->describe (out, vl);
1694 out <<
"userLambdaMax_: " << userLambdaMax_ << endl
1695 <<
"userLambdaMin_: " << userLambdaMin_ << endl
1696 <<
"userEigRatio_: " << userEigRatio_ << endl
1697 <<
"numIters_: " << numIters_ << endl
1698 <<
"eigMaxIters_: " << eigMaxIters_ << endl
1699 <<
"zeroStartingSolution_: " << zeroStartingSolution_ << endl
1700 <<
"assumeMatrixUnchanged_: " << assumeMatrixUnchanged_ << endl
1701 <<
"textbookAlgorithm_: " << textbookAlgorithm_ << endl
1702 <<
"computeMaxResNorm_: " << computeMaxResNorm_ << endl;
1710 #define IFPACK2_DETAILS_CHEBYSHEV_INSTANT(S,LO,GO,N) \ 1711 template class Ifpack2::Details::Chebyshev< S, Tpetra::MultiVector<S, LO, GO, N> >; 1713 #endif // IFPACK2_DETAILS_CHEBYSHEV_DEF_HPP Chebyshev(Teuchos::RCP< const row_matrix_type > A)
Definition: Ifpack2_Details_Chebyshev_def.hpp:276
Ifpack2 implementation details.
Diagonally scaled Chebyshev iteration for Tpetra sparse matrices.
Definition: Ifpack2_Chebyshev_decl.hpp:199
STS::magnitudeType MT
The type of the absolute value of a ScalarType.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:111
scalar_type ST
The type of entries in the matrix and vectors.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:107
Left-scaled Chebyshev iteration.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:101
Declaration of Chebyshev implementation.
Definition: Ifpack2_Container.hpp:774
Teuchos::ScalarTraits< scalar_type > STS
Traits class for ST.
Definition: Ifpack2_Details_Chebyshev_decl.hpp:109
Preconditioners and smoothers for Tpetra sparse matrices.
Definition: Ifpack2_AdditiveSchwarz_decl.hpp:72
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:126