42 #ifndef BELOS_PSEUDO_BLOCK_GMRES_SOLMGR_HPP 43 #define BELOS_PSEUDO_BLOCK_GMRES_SOLMGR_HPP 59 #ifdef HAVE_BELOS_TSQR 61 #endif // HAVE_BELOS_TSQR 65 #include "Teuchos_BLAS.hpp" 66 #ifdef BELOS_TEUCHOS_TIME_MONITOR 67 #include "Teuchos_TimeMonitor.hpp" 125 template<
class ScalarType,
class MV,
class OP>
131 typedef Teuchos::ScalarTraits<ScalarType> SCT;
132 typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
133 typedef Teuchos::ScalarTraits<MagnitudeType> MT;
259 const Teuchos::RCP<Teuchos::ParameterList> &pl );
283 Teuchos::Array<Teuchos::RCP<Teuchos::Time> >
getTimers()
const {
284 return Teuchos::tuple(timerSolve_);
374 void setParameters (
const Teuchos::RCP<Teuchos::ParameterList> ¶ms);
442 describe (Teuchos::FancyOStream& out,
443 const Teuchos::EVerbosityLevel verbLevel =
444 Teuchos::Describable::verbLevel_default)
const;
467 bool checkStatusTest();
470 Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > problem_;
473 Teuchos::RCP<OutputManager<ScalarType> > printer_;
474 Teuchos::RCP<std::ostream> outputStream_;
477 Teuchos::RCP<StatusTest<ScalarType,MV,OP> > userConvStatusTest_;
478 Teuchos::RCP<StatusTest<ScalarType,MV,OP> > debugStatusTest_;
479 Teuchos::RCP<StatusTest<ScalarType,MV,OP> > sTest_;
480 Teuchos::RCP<StatusTestMaxIters<ScalarType,MV,OP> > maxIterTest_;
481 Teuchos::RCP<StatusTest<ScalarType,MV,OP> > convTest_;
482 Teuchos::RCP<StatusTestResNorm<ScalarType,MV,OP> > impConvTest_, expConvTest_;
483 Teuchos::RCP<StatusTestOutput<ScalarType,MV,OP> > outputTest_;
485 Teuchos::RCP<std::map<std::string, Teuchos::RCP<StatusTest<ScalarType, MV, OP> > > > taggedTests_;
488 Teuchos::RCP<MatOrthoManager<ScalarType,MV,OP> > ortho_;
491 Teuchos::RCP<Teuchos::ParameterList> params_;
494 static const MagnitudeType convtol_default_;
495 static const MagnitudeType orthoKappa_default_;
496 static const int maxRestarts_default_;
497 static const int maxIters_default_;
498 static const bool showMaxResNormOnly_default_;
499 static const int blockSize_default_;
500 static const int numBlocks_default_;
501 static const int verbosity_default_;
502 static const int outputStyle_default_;
503 static const int outputFreq_default_;
504 static const int defQuorum_default_;
505 static const std::string impResScale_default_;
506 static const std::string expResScale_default_;
507 static const MagnitudeType resScaleFactor_default_;
508 static const std::string label_default_;
509 static const std::string orthoType_default_;
510 static const Teuchos::RCP<std::ostream> outputStream_default_;
513 MagnitudeType convtol_, orthoKappa_, achievedTol_;
514 int maxRestarts_, maxIters_, numIters_;
515 int blockSize_, numBlocks_, verbosity_, outputStyle_, outputFreq_, defQuorum_;
516 bool showMaxResNormOnly_;
517 std::string orthoType_;
518 std::string impResScale_, expResScale_;
519 MagnitudeType resScaleFactor_;
523 Teuchos::RCP<Teuchos::Time> timerSolve_;
526 bool isSet_, isSTSet_, expResTest_;
532 template<
class ScalarType,
class MV,
class OP>
533 const typename PseudoBlockGmresSolMgr<ScalarType,MV,OP>::MagnitudeType PseudoBlockGmresSolMgr<ScalarType,MV,OP>::convtol_default_ = 1e-8;
535 template<
class ScalarType,
class MV,
class OP>
536 const typename PseudoBlockGmresSolMgr<ScalarType,MV,OP>::MagnitudeType PseudoBlockGmresSolMgr<ScalarType,MV,OP>::orthoKappa_default_ = -1.0;
538 template<
class ScalarType,
class MV,
class OP>
539 const int PseudoBlockGmresSolMgr<ScalarType,MV,OP>::maxRestarts_default_ = 20;
541 template<
class ScalarType,
class MV,
class OP>
542 const int PseudoBlockGmresSolMgr<ScalarType,MV,OP>::maxIters_default_ = 1000;
544 template<
class ScalarType,
class MV,
class OP>
545 const bool PseudoBlockGmresSolMgr<ScalarType,MV,OP>::showMaxResNormOnly_default_ =
false;
547 template<
class ScalarType,
class MV,
class OP>
548 const int PseudoBlockGmresSolMgr<ScalarType,MV,OP>::blockSize_default_ = 1;
550 template<
class ScalarType,
class MV,
class OP>
551 const int PseudoBlockGmresSolMgr<ScalarType,MV,OP>::numBlocks_default_ = 300;
553 template<
class ScalarType,
class MV,
class OP>
554 const int PseudoBlockGmresSolMgr<ScalarType,MV,OP>::verbosity_default_ =
Belos::Errors;
556 template<
class ScalarType,
class MV,
class OP>
557 const int PseudoBlockGmresSolMgr<ScalarType,MV,OP>::outputStyle_default_ =
Belos::General;
559 template<
class ScalarType,
class MV,
class OP>
560 const int PseudoBlockGmresSolMgr<ScalarType,MV,OP>::outputFreq_default_ = -1;
562 template<
class ScalarType,
class MV,
class OP>
563 const int PseudoBlockGmresSolMgr<ScalarType,MV,OP>::defQuorum_default_ = 1;
565 template<
class ScalarType,
class MV,
class OP>
566 const std::string PseudoBlockGmresSolMgr<ScalarType,MV,OP>::impResScale_default_ =
"Norm of Preconditioned Initial Residual";
568 template<
class ScalarType,
class MV,
class OP>
569 const std::string PseudoBlockGmresSolMgr<ScalarType,MV,OP>::expResScale_default_ =
"Norm of Initial Residual";
571 template<
class ScalarType,
class MV,
class OP>
572 const typename PseudoBlockGmresSolMgr<ScalarType,MV,OP>::MagnitudeType PseudoBlockGmresSolMgr<ScalarType,MV,OP>::resScaleFactor_default_ = 1.0;
574 template<
class ScalarType,
class MV,
class OP>
575 const std::string PseudoBlockGmresSolMgr<ScalarType,MV,OP>::label_default_ =
"Belos";
577 template<
class ScalarType,
class MV,
class OP>
578 const std::string PseudoBlockGmresSolMgr<ScalarType,MV,OP>::orthoType_default_ =
"DGKS";
580 template<
class ScalarType,
class MV,
class OP>
581 const Teuchos::RCP<std::ostream> PseudoBlockGmresSolMgr<ScalarType,MV,OP>::outputStream_default_ = Teuchos::rcp(&std::cout,
false);
584 template<
class ScalarType,
class MV,
class OP>
586 outputStream_(outputStream_default_),
587 taggedTests_(Teuchos::null),
588 convtol_(convtol_default_),
589 orthoKappa_(orthoKappa_default_),
590 achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
591 maxRestarts_(maxRestarts_default_),
592 maxIters_(maxIters_default_),
594 blockSize_(blockSize_default_),
595 numBlocks_(numBlocks_default_),
596 verbosity_(verbosity_default_),
597 outputStyle_(outputStyle_default_),
598 outputFreq_(outputFreq_default_),
599 defQuorum_(defQuorum_default_),
600 showMaxResNormOnly_(showMaxResNormOnly_default_),
601 orthoType_(orthoType_default_),
602 impResScale_(impResScale_default_),
603 expResScale_(expResScale_default_),
604 resScaleFactor_(resScaleFactor_default_),
605 label_(label_default_),
613 template<
class ScalarType,
class MV,
class OP>
616 const Teuchos::RCP<Teuchos::ParameterList> &pl) :
618 outputStream_(outputStream_default_),
619 taggedTests_(Teuchos::null),
620 convtol_(convtol_default_),
621 orthoKappa_(orthoKappa_default_),
622 achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
623 maxRestarts_(maxRestarts_default_),
624 maxIters_(maxIters_default_),
626 blockSize_(blockSize_default_),
627 numBlocks_(numBlocks_default_),
628 verbosity_(verbosity_default_),
629 outputStyle_(outputStyle_default_),
630 outputFreq_(outputFreq_default_),
631 defQuorum_(defQuorum_default_),
632 showMaxResNormOnly_(showMaxResNormOnly_default_),
633 orthoType_(orthoType_default_),
634 impResScale_(impResScale_default_),
635 expResScale_(expResScale_default_),
636 resScaleFactor_(resScaleFactor_default_),
637 label_(label_default_),
643 TEUCHOS_TEST_FOR_EXCEPTION(problem_ == Teuchos::null, std::invalid_argument,
"Problem not given to solver manager.");
652 template<
class ScalarType,
class MV,
class OP>
657 using Teuchos::ParameterList;
658 using Teuchos::parameterList;
660 using Teuchos::rcp_dynamic_cast;
663 if (params_ == Teuchos::null) {
664 params_ = parameterList (*getValidParameters ());
671 params->validateParameters (*getValidParameters (), 0);
675 if (params->isParameter (
"Maximum Restarts")) {
676 maxRestarts_ = params->get (
"Maximum Restarts", maxRestarts_default_);
679 params_->set (
"Maximum Restarts", maxRestarts_);
683 if (params->isParameter (
"Maximum Iterations")) {
684 maxIters_ = params->get (
"Maximum Iterations", maxIters_default_);
687 params_->set (
"Maximum Iterations", maxIters_);
688 if (! maxIterTest_.is_null ()) {
689 maxIterTest_->setMaxIters (maxIters_);
694 if (params->isParameter (
"Block Size")) {
695 blockSize_ = params->get (
"Block Size", blockSize_default_);
696 TEUCHOS_TEST_FOR_EXCEPTION(
697 blockSize_ <= 0, std::invalid_argument,
698 "Belos::PseudoBlockGmresSolMgr::setParameters: " 699 "The \"Block Size\" parameter must be strictly positive, " 700 "but you specified a value of " << blockSize_ <<
".");
703 params_->set (
"Block Size", blockSize_);
707 if (params->isParameter (
"Num Blocks")) {
708 numBlocks_ = params->get (
"Num Blocks", numBlocks_default_);
709 TEUCHOS_TEST_FOR_EXCEPTION(
710 numBlocks_ <= 0, std::invalid_argument,
711 "Belos::PseudoBlockGmresSolMgr::setParameters: " 712 "The \"Num Blocks\" parameter must be strictly positive, " 713 "but you specified a value of " << numBlocks_ <<
".");
716 params_->set (
"Num Blocks", numBlocks_);
720 if (params->isParameter (
"Timer Label")) {
721 const std::string tempLabel = params->get (
"Timer Label", label_default_);
724 if (tempLabel != label_) {
726 params_->set (
"Timer Label", label_);
727 const std::string solveLabel =
728 label_ +
": PseudoBlockGmresSolMgr total solve time";
729 #ifdef BELOS_TEUCHOS_TIME_MONITOR 730 timerSolve_ = Teuchos::TimeMonitor::getNewCounter (solveLabel);
731 #endif // BELOS_TEUCHOS_TIME_MONITOR 732 if (ortho_ != Teuchos::null) {
733 ortho_->setLabel( label_ );
739 if (params->isParameter (
"Orthogonalization")) {
740 std::string tempOrthoType = params->get (
"Orthogonalization", orthoType_default_);
741 #ifdef HAVE_BELOS_TSQR 742 TEUCHOS_TEST_FOR_EXCEPTION(
743 tempOrthoType !=
"DGKS" && tempOrthoType !=
"ICGS" &&
744 tempOrthoType !=
"IMGS" && tempOrthoType !=
"TSQR",
745 std::invalid_argument,
746 "Belos::PseudoBlockGmresSolMgr::setParameters: " 747 "The \"Orthogonalization\" parameter must be one of \"DGKS\", \"ICGS\", " 748 "\"IMGS\", or \"TSQR\".");
750 TEUCHOS_TEST_FOR_EXCEPTION(
751 tempOrthoType !=
"DGKS" && tempOrthoType !=
"ICGS" &&
752 tempOrthoType !=
"IMGS",
753 std::invalid_argument,
754 "Belos::PseudoBlockGmresSolMgr::setParameters: " 755 "The \"Orthogonalization\" parameter must be one of \"DGKS\", \"ICGS\", " 757 #endif // HAVE_BELOS_TSQR 759 if (tempOrthoType != orthoType_) {
760 orthoType_ = tempOrthoType;
762 if (orthoType_ ==
"DGKS") {
764 if (orthoKappa_ <= 0) {
765 ortho_ = rcp (
new ortho_type (label_));
768 ortho_ = rcp (
new ortho_type (label_));
769 rcp_dynamic_cast<ortho_type> (ortho_)->setDepTol (orthoKappa_);
772 else if (orthoType_ ==
"ICGS") {
774 ortho_ = rcp (
new ortho_type (label_));
776 else if (orthoType_ ==
"IMGS") {
778 ortho_ = rcp (
new ortho_type (label_));
780 #ifdef HAVE_BELOS_TSQR 781 else if (orthoType_ ==
"TSQR") {
783 ortho_ = rcp (
new ortho_type (label_));
785 #endif // HAVE_BELOS_TSQR 790 if (params->isParameter (
"Orthogonalization Constant")) {
791 orthoKappa_ = params->get (
"Orthogonalization Constant", orthoKappa_default_);
794 params_->set (
"Orthogonalization Constant", orthoKappa_);
795 if (orthoType_ ==
"DGKS") {
796 if (orthoKappa_ > 0 && ! ortho_.is_null ()) {
798 rcp_dynamic_cast<ortho_type> (ortho_)->setDepTol (orthoKappa_);
804 if (params->isParameter (
"Verbosity")) {
805 if (Teuchos::isParameterType<int> (*params,
"Verbosity")) {
806 verbosity_ = params->get (
"Verbosity", verbosity_default_);
808 verbosity_ = (int) Teuchos::getParameter<Belos::MsgType> (*params,
"Verbosity");
812 params_->set (
"Verbosity", verbosity_);
813 if (! printer_.is_null ()) {
814 printer_->setVerbosity (verbosity_);
819 if (params->isParameter (
"Output Style")) {
820 if (Teuchos::isParameterType<int> (*params,
"Output Style")) {
821 outputStyle_ = params->get (
"Output Style", outputStyle_default_);
823 outputStyle_ = (int) Teuchos::getParameter<Belos::OutputType> (*params,
"Output Style");
827 params_->set (
"Output Style", verbosity_);
828 if (! outputTest_.is_null ()) {
833 if (params->isSublist (
"User Status Tests")) {
834 Teuchos::ParameterList userStatusTestsList = params->sublist(
"User Status Tests",
true);
835 if ( userStatusTestsList.numParams() > 0 ) {
836 std::string userCombo_string = params->get<std::string>(
"User Status Tests Combo Type",
"SEQ");
838 setUserConvStatusTest( testFactory->buildStatusTests(userStatusTestsList), testFactory->stringToComboType(userCombo_string) );
839 taggedTests_ = testFactory->getTaggedTests();
845 if (params->isParameter (
"Output Stream")) {
846 outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> > (*params,
"Output Stream");
849 params_->set(
"Output Stream", outputStream_);
850 if (! printer_.is_null ()) {
851 printer_->setOStream (outputStream_);
857 if (params->isParameter (
"Output Frequency")) {
858 outputFreq_ = params->get (
"Output Frequency", outputFreq_default_);
862 params_->set (
"Output Frequency", outputFreq_);
863 if (! outputTest_.is_null ()) {
864 outputTest_->setOutputFrequency (outputFreq_);
869 if (printer_.is_null ()) {
876 if (params->isParameter (
"Convergence Tolerance")) {
877 convtol_ = params->get (
"Convergence Tolerance", convtol_default_);
880 params_->set (
"Convergence Tolerance", convtol_);
881 if (! impConvTest_.is_null ()) {
882 impConvTest_->setTolerance (convtol_);
884 if (! expConvTest_.is_null ()) {
885 expConvTest_->setTolerance (convtol_);
890 bool userDefinedResidualScalingUpdated =
false;
891 if (params->isParameter (
"User Defined Residual Scaling")) {
892 const MagnitudeType tempResScaleFactor =
893 Teuchos::getParameter<MagnitudeType> (*params,
"User Defined Residual Scaling");
896 if (resScaleFactor_ != tempResScaleFactor) {
897 resScaleFactor_ = tempResScaleFactor;
898 userDefinedResidualScalingUpdated =
true;
901 if(userDefinedResidualScalingUpdated)
903 if (! params->isParameter (
"Implicit Residual Scaling") && ! impConvTest_.is_null ()) {
905 if(impResScale_ ==
"User Provided")
908 catch (std::exception& e) {
913 if (! params->isParameter (
"Explicit Residual Scaling") && ! expConvTest_.is_null ()) {
915 if(expResScale_ ==
"User Provided")
918 catch (std::exception& e) {
927 if (params->isParameter (
"Implicit Residual Scaling")) {
928 const std::string tempImpResScale =
929 Teuchos::getParameter<std::string> (*params,
"Implicit Residual Scaling");
932 if (impResScale_ != tempImpResScale) {
934 impResScale_ = tempImpResScale;
937 params_->set (
"Implicit Residual Scaling", impResScale_);
938 if (! impConvTest_.is_null ()) {
940 if(impResScale_ ==
"User Provided")
941 impConvTest_->defineScaleForm (impResScaleType,
Belos::TwoNorm, resScaleFactor_);
945 catch (std::exception& e) {
951 else if (userDefinedResidualScalingUpdated) {
954 if (! impConvTest_.is_null ()) {
956 if(impResScale_ ==
"User Provided")
957 impConvTest_->defineScaleForm (impResScaleType,
Belos::TwoNorm, resScaleFactor_);
959 catch (std::exception& e) {
967 if (params->isParameter (
"Explicit Residual Scaling")) {
968 const std::string tempExpResScale =
969 Teuchos::getParameter<std::string> (*params,
"Explicit Residual Scaling");
972 if (expResScale_ != tempExpResScale) {
974 expResScale_ = tempExpResScale;
977 params_->set (
"Explicit Residual Scaling", expResScale_);
978 if (! expConvTest_.is_null ()) {
980 if(expResScale_ ==
"User Provided")
981 expConvTest_->defineScaleForm (expResScaleType,
Belos::TwoNorm, resScaleFactor_);
985 catch (std::exception& e) {
991 else if (userDefinedResidualScalingUpdated) {
994 if (! expConvTest_.is_null ()) {
996 if(expResScale_ ==
"User Provided")
997 expConvTest_->defineScaleForm (expResScaleType,
Belos::TwoNorm, resScaleFactor_);
999 catch (std::exception& e) {
1008 if (params->isParameter (
"Show Maximum Residual Norm Only")) {
1009 showMaxResNormOnly_ =
1010 Teuchos::getParameter<bool> (*params,
"Show Maximum Residual Norm Only");
1013 params_->set (
"Show Maximum Residual Norm Only", showMaxResNormOnly_);
1014 if (! impConvTest_.is_null ()) {
1015 impConvTest_->setShowMaxResNormOnly (showMaxResNormOnly_);
1017 if (! expConvTest_.is_null ()) {
1018 expConvTest_->setShowMaxResNormOnly (showMaxResNormOnly_);
1025 if (params->isParameter(
"Deflation Quorum")) {
1026 defQuorum_ = params->get(
"Deflation Quorum", defQuorum_);
1027 TEUCHOS_TEST_FOR_EXCEPTION(
1028 defQuorum_ > blockSize_, std::invalid_argument,
1029 "Belos::PseudoBlockGmresSolMgr::setParameters: " 1030 "The \"Deflation Quorum\" parameter (= " << defQuorum_ <<
") must not be " 1031 "larger than \"Block Size\" (= " << blockSize_ <<
").");
1032 params_->set (
"Deflation Quorum", defQuorum_);
1033 if (! impConvTest_.is_null ()) {
1034 impConvTest_->setQuorum (defQuorum_);
1036 if (! expConvTest_.is_null ()) {
1037 expConvTest_->setQuorum (defQuorum_);
1042 if (ortho_.is_null ()) {
1043 if (orthoType_ ==
"DGKS") {
1045 if (orthoKappa_ <= 0) {
1046 ortho_ = rcp (
new ortho_type (label_));
1049 ortho_ = rcp (
new ortho_type (label_));
1050 rcp_dynamic_cast<ortho_type> (ortho_)->setDepTol (orthoKappa_);
1053 else if (orthoType_ ==
"ICGS") {
1055 ortho_ = rcp (
new ortho_type (label_));
1057 else if (orthoType_ ==
"IMGS") {
1059 ortho_ = rcp (
new ortho_type (label_));
1061 #ifdef HAVE_BELOS_TSQR 1062 else if (orthoType_ ==
"TSQR") {
1064 ortho_ = rcp (
new ortho_type (label_));
1066 #endif // HAVE_BELOS_TSQR 1068 #ifdef HAVE_BELOS_TSQR 1069 TEUCHOS_TEST_FOR_EXCEPTION(
1070 orthoType_ !=
"ICGS" && orthoType_ !=
"DGKS" &&
1071 orthoType_ !=
"IMGS" && orthoType_ !=
"TSQR",
1073 "Belos::PseudoBlockGmresSolMgr::setParameters(): " 1074 "Invalid orthogonalization type \"" << orthoType_ <<
"\".");
1076 TEUCHOS_TEST_FOR_EXCEPTION(
1077 orthoType_ !=
"ICGS" && orthoType_ !=
"DGKS" &&
1078 orthoType_ !=
"IMGS",
1080 "Belos::PseudoBlockGmresSolMgr::setParameters(): " 1081 "Invalid orthogonalization type \"" << orthoType_ <<
"\".");
1082 #endif // HAVE_BELOS_TSQR 1087 if (timerSolve_ == Teuchos::null) {
1088 std::string solveLabel = label_ +
": PseudoBlockGmresSolMgr total solve time";
1089 #ifdef BELOS_TEUCHOS_TIME_MONITOR 1090 timerSolve_ = Teuchos::TimeMonitor::getNewCounter (solveLabel);
1099 template<
class ScalarType,
class MV,
class OP>
1106 userConvStatusTest_ = userConvStatusTest;
1107 comboType_ = comboType;
1110 template<
class ScalarType,
class MV,
class OP>
1116 debugStatusTest_ = debugStatusTest;
1121 template<
class ScalarType,
class MV,
class OP>
1122 Teuchos::RCP<const Teuchos::ParameterList>
1125 static Teuchos::RCP<const Teuchos::ParameterList> validPL;
1126 if (is_null(validPL)) {
1127 Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
1129 pl= Teuchos::rcp(
new Teuchos::ParameterList() );
1130 pl->set(
"Convergence Tolerance", convtol_default_,
1131 "The relative residual tolerance that needs to be achieved by the\n" 1132 "iterative solver in order for the linear system to be declared converged.");
1133 pl->set(
"Maximum Restarts", maxRestarts_default_,
1134 "The maximum number of restarts allowed for each\n" 1135 "set of RHS solved.");
1136 pl->set(
"Maximum Iterations", maxIters_default_,
1137 "The maximum number of block iterations allowed for each\n" 1138 "set of RHS solved.");
1139 pl->set(
"Num Blocks", numBlocks_default_,
1140 "The maximum number of vectors allowed in the Krylov subspace\n" 1141 "for each set of RHS solved.");
1142 pl->set(
"Block Size", blockSize_default_,
1143 "The number of RHS solved simultaneously.");
1144 pl->set(
"Verbosity", verbosity_default_,
1145 "What type(s) of solver information should be outputted\n" 1146 "to the output stream.");
1147 pl->set(
"Output Style", outputStyle_default_,
1148 "What style is used for the solver information outputted\n" 1149 "to the output stream.");
1150 pl->set(
"Output Frequency", outputFreq_default_,
1151 "How often convergence information should be outputted\n" 1152 "to the output stream.");
1153 pl->set(
"Deflation Quorum", defQuorum_default_,
1154 "The number of linear systems that need to converge before\n" 1155 "they are deflated. This number should be <= block size.");
1156 pl->set(
"Output Stream", outputStream_default_,
1157 "A reference-counted pointer to the output stream where all\n" 1158 "solver output is sent.");
1159 pl->set(
"Show Maximum Residual Norm Only", showMaxResNormOnly_default_,
1160 "When convergence information is printed, only show the maximum\n" 1161 "relative residual norm when the block size is greater than one.");
1162 pl->set(
"Implicit Residual Scaling", impResScale_default_,
1163 "The type of scaling used in the implicit residual convergence test.");
1164 pl->set(
"Explicit Residual Scaling", expResScale_default_,
1165 "The type of scaling used in the explicit residual convergence test.");
1166 pl->set(
"Timer Label", label_default_,
1167 "The string to use as a prefix for the timer labels.");
1169 pl->set(
"Orthogonalization", orthoType_default_,
1170 "The type of orthogonalization to use: DGKS, ICGS, IMGS.");
1171 pl->set(
"Orthogonalization Constant",orthoKappa_default_,
1172 "The constant used by DGKS orthogonalization to determine\n" 1173 "whether another step of classical Gram-Schmidt is necessary.");
1174 pl->sublist(
"User Status Tests");
1175 pl->set(
"User Status Tests Combo Type",
"SEQ",
1176 "Type of logical combination operation of user-defined\n" 1177 "and/or solver-specific status tests.");
1184 template<
class ScalarType,
class MV,
class OP>
1196 if ( !Teuchos::is_null(problem_->getLeftPrec()) ) {
1203 Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
1204 Teuchos::rcp(
new StatusTestGenResNorm_t( convtol_, defQuorum_ ) );
1205 if(impResScale_ ==
"User Provided")
1210 tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
1211 impConvTest_ = tmpImpConvTest;
1214 Teuchos::RCP<StatusTestGenResNorm_t> tmpExpConvTest =
1215 Teuchos::rcp(
new StatusTestGenResNorm_t( convtol_, defQuorum_ ) );
1216 tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit,
Belos::TwoNorm );
1217 if(expResScale_ ==
"User Provided")
1221 tmpExpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
1222 expConvTest_ = tmpExpConvTest;
1225 convTest_ = Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
1231 Teuchos::RCP<StatusTestImpResNorm_t> tmpImpConvTest =
1232 Teuchos::rcp(
new StatusTestImpResNorm_t( convtol_, defQuorum_ ) );
1233 if(impResScale_ ==
"User Provided")
1237 tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
1238 impConvTest_ = tmpImpConvTest;
1241 expConvTest_ = impConvTest_;
1242 convTest_ = impConvTest_;
1245 if (nonnull(debugStatusTest_) ) {
1247 convTest_ = Teuchos::rcp(
1248 new StatusTestCombo_t( StatusTestCombo_t::AND, convTest_, debugStatusTest_ ) );
1251 if (nonnull(userConvStatusTest_) ) {
1253 convTest_ = Teuchos::rcp(
1254 new StatusTestCombo_t( comboType_, convTest_, userConvStatusTest_ ) );
1261 sTest_ = Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
1269 StatusTestOutputFactory<ScalarType,MV,OP> stoFactory( outputStyle_, taggedTests_ );
1273 std::string solverDesc =
" Pseudo Block Gmres ";
1274 outputTest_->setSolverDesc( solverDesc );
1285 template<
class ScalarType,
class MV,
class OP>
1291 if (!isSet_) { setParameters( params_ ); }
1293 Teuchos::BLAS<int,ScalarType> blas;
1296 "Belos::PseudoBlockGmresSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
1299 if (!isSTSet_ || (!expResTest_ && !Teuchos::is_null(problem_->getLeftPrec())) ) {
1301 "Belos::BlockGmresSolMgr::solve(): Linear problem and requested status tests are incompatible.");
1306 int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
1307 int numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
1309 std::vector<int> currIdx( numCurrRHS );
1310 blockSize_ = numCurrRHS;
1311 for (
int i=0; i<numCurrRHS; ++i)
1312 { currIdx[i] = startPtr+i; }
1315 problem_->setLSIndex( currIdx );
1319 Teuchos::ParameterList plist;
1320 plist.set(
"Num Blocks",numBlocks_);
1323 outputTest_->reset();
1324 loaDetected_ =
false;
1327 bool isConverged =
true;
1332 Teuchos::RCP<PseudoBlockGmresIter<ScalarType,MV,OP> > block_gmres_iter
1337 #ifdef BELOS_TEUCHOS_TIME_MONITOR 1338 Teuchos::TimeMonitor slvtimer(*timerSolve_);
1341 while ( numRHS2Solve > 0 ) {
1344 std::vector<int> convRHSIdx;
1345 std::vector<int> currRHSIdx( currIdx );
1346 currRHSIdx.resize(numCurrRHS);
1349 block_gmres_iter->setNumBlocks( numBlocks_ );
1352 block_gmres_iter->resetNumIters();
1355 outputTest_->resetNumCalls();
1361 std::vector<int> index(1);
1362 Teuchos::RCP<MV> tmpV, R_0 = MVT::CloneCopy( *(problem_->getInitPrecResVec()), currIdx );
1363 newState.
V.resize( blockSize_ );
1364 newState.
Z.resize( blockSize_ );
1365 for (
int i=0; i<blockSize_; ++i) {
1367 tmpV = MVT::CloneViewNonConst( *R_0, index );
1370 Teuchos::RCP<Teuchos::SerialDenseVector<int,ScalarType> > tmpZ
1371 = Teuchos::rcp(
new Teuchos::SerialDenseVector<int,ScalarType>( 1 ));
1374 int rank = ortho_->normalize( *tmpV, tmpZ );
1376 "Belos::PseudoBlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors.");
1378 newState.
V[i] = tmpV;
1379 newState.
Z[i] = tmpZ;
1383 block_gmres_iter->initialize(newState);
1384 int numRestarts = 0;
1390 block_gmres_iter->iterate();
1397 if ( convTest_->getStatus() ==
Passed ) {
1399 if ( expConvTest_->getLOADetected() ) {
1409 loaDetected_ =
true;
1411 "Belos::PseudoBlockGmresSolMgr::solve(): Warning! Solver has experienced a loss of accuracy!" << std::endl;
1412 isConverged =
false;
1416 std::vector<int> convIdx = expConvTest_->convIndices();
1420 if (convIdx.size() == currRHSIdx.size())
1427 problem_->setCurrLS();
1431 int curDim = oldState.
curDim;
1436 std::vector<int> oldRHSIdx( currRHSIdx );
1437 std::vector<int> defRHSIdx;
1438 for (
unsigned int i=0; i<currRHSIdx.size(); ++i) {
1440 for (
unsigned int j=0; j<convIdx.size(); ++j) {
1441 if (currRHSIdx[i] == convIdx[j]) {
1447 defRHSIdx.push_back( i );
1450 defState.
V.push_back( Teuchos::rcp_const_cast<MV>( oldState.
V[i] ) );
1451 defState.
Z.push_back( Teuchos::rcp_const_cast<Teuchos::SerialDenseVector<int,ScalarType> >( oldState.
Z[i] ) );
1452 defState.
H.push_back( Teuchos::rcp_const_cast<Teuchos::SerialDenseMatrix<int,ScalarType> >( oldState.
H[i] ) );
1453 defState.
sn.push_back( Teuchos::rcp_const_cast<Teuchos::SerialDenseVector<int,ScalarType> >( oldState.
sn[i] ) );
1454 defState.
cs.push_back( Teuchos::rcp_const_cast<Teuchos::SerialDenseVector<int,MagnitudeType> >(oldState.
cs[i] ) );
1455 currRHSIdx[have] = currRHSIdx[i];
1459 defRHSIdx.resize(currRHSIdx.size()-have);
1460 currRHSIdx.resize(have);
1464 Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1465 Teuchos::RCP<MV> defUpdate = MVT::CloneViewNonConst( *update, defRHSIdx );
1468 problem_->setLSIndex( convIdx );
1471 problem_->updateSolution( defUpdate,
true );
1475 problem_->setLSIndex( currRHSIdx );
1478 defState.
curDim = curDim;
1481 block_gmres_iter->initialize(defState);
1488 else if ( maxIterTest_->getStatus() ==
Passed ) {
1490 isConverged =
false;
1498 else if ( block_gmres_iter->getCurSubspaceDim() == block_gmres_iter->getMaxSubspaceDim() ) {
1500 if ( numRestarts >= maxRestarts_ ) {
1501 isConverged =
false;
1506 printer_->stream(
Debug) <<
" Performing restart number " << numRestarts <<
" of " << maxRestarts_ << std::endl << std::endl;
1509 Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1510 problem_->updateSolution( update,
true );
1517 newstate.
V.resize(currRHSIdx.size());
1518 newstate.
Z.resize(currRHSIdx.size());
1522 R_0 = MVT::Clone( *(problem_->getInitPrecResVec()), currRHSIdx.size() );
1523 problem_->computeCurrPrecResVec( &*R_0 );
1524 for (
unsigned int i=0; i<currRHSIdx.size(); ++i) {
1527 tmpV = MVT::CloneViewNonConst( *R_0, index );
1530 Teuchos::RCP<Teuchos::SerialDenseVector<int,ScalarType> > tmpZ
1531 = Teuchos::rcp(
new Teuchos::SerialDenseVector<int,ScalarType>( 1 ));
1534 int rank = ortho_->normalize( *tmpV, tmpZ );
1536 "Belos::PseudoBlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors after the restart.");
1538 newstate.
V[i] = tmpV;
1539 newstate.
Z[i] = tmpZ;
1544 block_gmres_iter->initialize(newstate);
1556 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
1557 "Belos::PseudoBlockGmresSolMgr::solve(): Invalid return from PseudoBlockGmresIter::iterate().");
1563 block_gmres_iter->updateLSQR( block_gmres_iter->getCurSubspaceDim() );
1566 sTest_->checkStatus( &*block_gmres_iter );
1567 if (convTest_->getStatus() !=
Passed)
1568 isConverged =
false;
1571 catch (
const std::exception &e) {
1572 printer_->stream(
Errors) <<
"Error! Caught std::exception in PseudoBlockGmresIter::iterate() at iteration " 1573 << block_gmres_iter->getNumIters() << std::endl
1574 << e.what() << std::endl;
1581 if (nonnull(userConvStatusTest_)) {
1583 Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1584 problem_->updateSolution( update,
true );
1586 else if (nonnull(expConvTest_->getSolution())) {
1588 Teuchos::RCP<MV> newX = expConvTest_->getSolution();
1589 Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
1590 MVT::MvAddMv( 0.0, *newX, 1.0, *newX, *curX );
1594 Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1595 problem_->updateSolution( update,
true );
1599 problem_->setCurrLS();
1602 startPtr += numCurrRHS;
1603 numRHS2Solve -= numCurrRHS;
1604 if ( numRHS2Solve > 0 ) {
1605 numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
1607 blockSize_ = numCurrRHS;
1608 currIdx.resize( numCurrRHS );
1609 for (
int i=0; i<numCurrRHS; ++i)
1610 { currIdx[i] = startPtr+i; }
1613 if (defQuorum_ > blockSize_) {
1614 if (impConvTest_ != Teuchos::null)
1615 impConvTest_->setQuorum( blockSize_ );
1616 if (expConvTest_ != Teuchos::null)
1617 expConvTest_->setQuorum( blockSize_ );
1621 problem_->setLSIndex( currIdx );
1624 currIdx.resize( numRHS2Solve );
1635 #ifdef BELOS_TEUCHOS_TIME_MONITOR 1640 Teuchos::TimeMonitor::summarize( printer_->stream(
TimingDetails) );
1644 numIters_ = maxIterTest_->getNumIters();
1657 const std::vector<MagnitudeType>* pTestValues = NULL;
1659 pTestValues = expConvTest_->getTestValue();
1660 if (pTestValues == NULL || pTestValues->size() < 1) {
1661 pTestValues = impConvTest_->getTestValue();
1666 pTestValues = impConvTest_->getTestValue();
1668 TEUCHOS_TEST_FOR_EXCEPTION(pTestValues == NULL, std::logic_error,
1669 "Belos::PseudoBlockGmresSolMgr::solve(): The implicit convergence test's " 1670 "getTestValue() method returned NULL. Please report this bug to the " 1671 "Belos developers.");
1672 TEUCHOS_TEST_FOR_EXCEPTION(pTestValues->size() < 1, std::logic_error,
1673 "Belos::PseudoBlockGmresSolMgr::solve(): The implicit convergence test's " 1674 "getTestValue() method returned a vector of length zero. Please report " 1675 "this bug to the Belos developers.");
1680 achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
1683 if (!isConverged || loaDetected_) {
1690 template<
class ScalarType,
class MV,
class OP>
1693 std::ostringstream out;
1695 out <<
"\"Belos::PseudoBlockGmresSolMgr\": {";
1696 if (this->getObjectLabel () !=
"") {
1697 out <<
"Label: " << this->getObjectLabel () <<
", ";
1699 out <<
"Num Blocks: " << numBlocks_
1700 <<
", Maximum Iterations: " << maxIters_
1701 <<
", Maximum Restarts: " << maxRestarts_
1702 <<
", Convergence Tolerance: " << convtol_
1708 template<
class ScalarType,
class MV,
class OP>
1712 const Teuchos::EVerbosityLevel verbLevel)
const 1714 using Teuchos::TypeNameTraits;
1715 using Teuchos::VERB_DEFAULT;
1716 using Teuchos::VERB_NONE;
1717 using Teuchos::VERB_LOW;
1724 const Teuchos::EVerbosityLevel vl =
1725 (verbLevel == VERB_DEFAULT) ? VERB_LOW : verbLevel;
1727 if (vl != VERB_NONE) {
1728 Teuchos::OSTab tab0 (out);
1730 out <<
"\"Belos::PseudoBlockGmresSolMgr\":" << endl;
1731 Teuchos::OSTab tab1 (out);
1732 out <<
"Template parameters:" << endl;
1734 Teuchos::OSTab tab2 (out);
1735 out <<
"ScalarType: " << TypeNameTraits<ScalarType>::name () << endl
1736 <<
"MV: " << TypeNameTraits<MV>::name () << endl
1737 <<
"OP: " << TypeNameTraits<OP>::name () << endl;
1739 if (this->getObjectLabel () !=
"") {
1740 out <<
"Label: " << this->getObjectLabel () << endl;
1742 out <<
"Num Blocks: " << numBlocks_ << endl
1743 <<
"Maximum Iterations: " << maxIters_ << endl
1744 <<
"Maximum Restarts: " << maxRestarts_ << endl
1745 <<
"Convergence Tolerance: " << convtol_ << endl;
ScaleType convertStringToScaleType(const std::string &scaleType)
Convert the given string to its ScaleType enum value.
std::vector< Teuchos::RCP< const Teuchos::SerialDenseMatrix< int, ScalarType > > > H
The current Hessenberg matrix.
Collection of types and exceptions used within the Belos solvers.
ComboType
The test can be either the AND of all the component tests, or the OR of all the component tests...
virtual ~PseudoBlockGmresSolMgr()
Destructor.
Belos's basic output manager for sending information of select verbosity levels to the appropriate ou...
MagnitudeType achievedTol() const
Tolerance achieved by the last solve() invocation.
bool isLOADetected() const
Whether a "loss of accuracy" was detected during the last solve().
Class which manages the output and verbosity of the Belos solvers.
PseudoBlockGmresSolMgr()
Empty constructor.
void setParameters(const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Set the parameters the solver manager should use to solve the linear problem.
PseudoBlockGmresSolMgrLinearProblemFailure(const std::string &what_arg)
ScaleType
The type of scaling to use on the residual norm value.
Convergence test using the implicit residual norm(s), with an explicit residual norm(s) check for los...
PseudoBlockGmresSolMgrOrthoFailure(const std::string &what_arg)
Teuchos::Array< Teuchos::RCP< Teuchos::Time > > getTimers() const
Return the timers for this object.
std::vector< Teuchos::RCP< const MV > > V
The current Krylov basis.
An implementation of StatusTestResNorm using a family of residual norms.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
A list of valid default parameters for this solver.
int getNumIters() const
Iteration count for the most recent call to solve().
A pure virtual class for defining the status tests for the Belos iterative solvers.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with the given verbosity level to a FancyOStream.
Structure to contain pointers to PseudoBlockGmresIter state variables.
A factory class for generating StatusTestOutput objects.
Iterated Modified Gram-Schmidt (IMGS) implementation of the Belos::OrthoManager class.
virtual void setUserConvStatusTest(const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &userConvStatusTest, const typename StatusTestCombo< ScalarType, MV, OP >::ComboType &comboType=StatusTestCombo< ScalarType, MV, OP >::SEQ)
Set a custom status test.
Traits class which defines basic operations on multivectors.
Classical Gram-Schmidt (with DGKS correction) implementation of the Belos::OrthoManager class...
A Belos::StatusTest class for specifying a maximum number of iterations.
ResetType
How to reset the solver.
Pure virtual base class which describes the basic interface for a solver manager. ...
PseudoBlockGmresSolMgrLinearProblemFailure is thrown when the linear problem is not setup (i...
std::vector< Teuchos::RCP< const Teuchos::SerialDenseVector< int, MagnitudeType > > > cs
A linear system to solve, and its associated information.
Class which describes the linear problem to be solved by the iterative solver.
std::string description() const
Return a one-line description of this object.
Interface to standard and "pseudoblock" GMRES.
ReturnType
Whether the Belos solve converged for all linear systems.
The Belos::SolverManager is a templated virtual base class that defines the basic interface that any ...
Iterated Classical Gram-Schmidt (ICGS) implementation of the Belos::OrthoManager class.
An implementation of the Belos::MatOrthoManager that performs orthogonalization using multiple steps ...
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const
The current parameters for this solver.
An implementation of the Belos::MatOrthoManager that performs orthogonalization using multiple steps ...
void reset(const ResetType type)
Performs a reset of the solver manager specified by the ResetType. This informs the solver manager th...
std::vector< Teuchos::RCP< const Teuchos::SerialDenseVector< int, ScalarType > > > Z
The current right-hand side of the least squares system RY = Z.
Belos concrete class for performing the pseudo-block GMRES iteration.
virtual void setDebugStatusTest(const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &debugStatusTest)
Set a debug status test.
An implementation of the Belos::MatOrthoManager that performs orthogonalization using (potentially) m...
A class for extending the status testing capabilities of Belos via logical combinations.
Orthogonalization manager based on Tall Skinny QR (TSQR)
Class which defines basic traits for the operator type.
int curDim
The current dimension of the reduction.
void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem)
Set the linear problem to solve.
Parent class to all Belos exceptions.
std::vector< Teuchos::RCP< const Teuchos::SerialDenseVector< int, ScalarType > > > sn
The current Given's rotation coefficients.
const LinearProblem< ScalarType, MV, OP > & getProblem() const
Return a reference to the linear problem being solved by this solver manager.
ReturnType solve()
This method performs possibly repeated calls to the underlying linear solver's iterate() routine unti...
MatOrthoManager subclass using TSQR or DGKS.
Belos header file which uses auto-configuration information to include necessary C++ headers...
This class implements the pseudo-block GMRES iteration, where a block Krylov subspace is constructed ...
Factory to build a set of status tests from a parameter list.
PseudoBlockGmresIterOrthoFailure is thrown when the orthogonalization manager is unable to generate o...
PseudoBlockGmresSolMgrOrthoFailure is thrown when the orthogonalization manager is unable to generate...