42 #ifndef BELOS_PSEUDO_BLOCK_GMRES_SOLMGR_HPP
43 #define BELOS_PSEUDO_BLOCK_GMRES_SOLMGR_HPP
60 #ifdef BELOS_TEUCHOS_TIME_MONITOR
119 template<
class ScalarType,
class MV,
class OP>
284 return Teuchos::tuple(timerSolve_);
464 bool checkStatusTest();
491 static constexpr
int maxRestarts_default_ = 20;
492 static constexpr
int maxIters_default_ = 1000;
493 static constexpr
bool showMaxResNormOnly_default_ =
false;
494 static constexpr
int blockSize_default_ = 1;
495 static constexpr
int numBlocks_default_ = 300;
498 static constexpr
int outputFreq_default_ = -1;
499 static constexpr
int defQuorum_default_ = 1;
500 static constexpr
const char * impResScale_default_ =
"Norm of Preconditioned Initial Residual";
501 static constexpr
const char * expResScale_default_ =
"Norm of Initial Residual";
502 static constexpr
const char * label_default_ =
"Belos";
503 static constexpr
const char * orthoType_default_ =
"ICGS";
506 MagnitudeType convtol_, orthoKappa_, achievedTol_;
507 int maxRestarts_, maxIters_, numIters_;
508 int blockSize_, numBlocks_, verbosity_, outputStyle_, outputFreq_, defQuorum_;
509 bool showMaxResNormOnly_;
510 std::string orthoType_;
511 std::string impResScale_, expResScale_;
512 MagnitudeType resScaleFactor_;
519 bool isSet_, isSTSet_, expResTest_;
525 template<
class ScalarType,
class MV,
class OP>
527 outputStream_(Teuchos::rcpFromRef(std::cout)),
528 taggedTests_(Teuchos::null),
531 achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
532 maxRestarts_(maxRestarts_default_),
533 maxIters_(maxIters_default_),
535 blockSize_(blockSize_default_),
536 numBlocks_(numBlocks_default_),
537 verbosity_(verbosity_default_),
538 outputStyle_(outputStyle_default_),
539 outputFreq_(outputFreq_default_),
540 defQuorum_(defQuorum_default_),
541 showMaxResNormOnly_(showMaxResNormOnly_default_),
542 orthoType_(orthoType_default_),
543 impResScale_(impResScale_default_),
544 expResScale_(expResScale_default_),
546 label_(label_default_),
554 template<
class ScalarType,
class MV,
class OP>
559 outputStream_(Teuchos::rcpFromRef(std::cout)),
560 taggedTests_(Teuchos::null),
563 achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
564 maxRestarts_(maxRestarts_default_),
565 maxIters_(maxIters_default_),
567 blockSize_(blockSize_default_),
568 numBlocks_(numBlocks_default_),
569 verbosity_(verbosity_default_),
570 outputStyle_(outputStyle_default_),
571 outputFreq_(outputFreq_default_),
572 defQuorum_(defQuorum_default_),
573 showMaxResNormOnly_(showMaxResNormOnly_default_),
574 orthoType_(orthoType_default_),
575 impResScale_(impResScale_default_),
576 expResScale_(expResScale_default_),
578 label_(label_default_),
593 template<
class ScalarType,
class MV,
class OP>
599 using Teuchos::parameterList;
601 using Teuchos::rcp_dynamic_cast;
604 if (params_ == Teuchos::null) {
605 params_ = parameterList (*getValidParameters ());
617 maxRestarts_ = params->
get (
"Maximum Restarts", maxRestarts_default_);
620 params_->set (
"Maximum Restarts", maxRestarts_);
625 maxIters_ = params->
get (
"Maximum Iterations", maxIters_default_);
628 params_->set (
"Maximum Iterations", maxIters_);
629 if (! maxIterTest_.is_null ()) {
630 maxIterTest_->setMaxIters (maxIters_);
636 blockSize_ = params->
get (
"Block Size", blockSize_default_);
638 blockSize_ <= 0, std::invalid_argument,
639 "Belos::PseudoBlockGmresSolMgr::setParameters: "
640 "The \"Block Size\" parameter must be strictly positive, "
641 "but you specified a value of " << blockSize_ <<
".");
644 params_->set (
"Block Size", blockSize_);
649 numBlocks_ = params->
get (
"Num Blocks", numBlocks_default_);
651 numBlocks_ <= 0, std::invalid_argument,
652 "Belos::PseudoBlockGmresSolMgr::setParameters: "
653 "The \"Num Blocks\" parameter must be strictly positive, "
654 "but you specified a value of " << numBlocks_ <<
".");
657 params_->set (
"Num Blocks", numBlocks_);
662 const std::string tempLabel = params->
get (
"Timer Label", label_default_);
665 if (tempLabel != label_) {
667 params_->set (
"Timer Label", label_);
668 const std::string solveLabel =
669 label_ +
": PseudoBlockGmresSolMgr total solve time";
670 #ifdef BELOS_TEUCHOS_TIME_MONITOR
672 #endif // BELOS_TEUCHOS_TIME_MONITOR
673 if (ortho_ != Teuchos::null) {
674 ortho_->setLabel( label_ );
682 if (Teuchos::isParameterType<int> (*params,
"Verbosity")) {
683 verbosity_ = params->
get (
"Verbosity", verbosity_default_);
685 verbosity_ = (int) Teuchos::getParameter<Belos::MsgType> (*params,
"Verbosity");
689 params_->set (
"Verbosity", verbosity_);
690 if (! printer_.is_null ()) {
691 printer_->setVerbosity (verbosity_);
697 if (Teuchos::isParameterType<int> (*params,
"Output Style")) {
698 outputStyle_ = params->
get (
"Output Style", outputStyle_default_);
700 outputStyle_ = (int) Teuchos::getParameter<Belos::OutputType> (*params,
"Output Style");
704 params_->set (
"Output Style", outputStyle_);
705 if (! outputTest_.is_null ()) {
712 if (params->
isSublist (
"User Status Tests")) {
714 if ( userStatusTestsList.
numParams() > 0 ) {
715 std::string userCombo_string = params->
get<std::string>(
"User Status Tests Combo Type",
"SEQ");
717 setUserConvStatusTest( testFactory->buildStatusTests(userStatusTestsList), testFactory->stringToComboType(userCombo_string) );
718 taggedTests_ = testFactory->getTaggedTests();
725 outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> > (*params,
"Output Stream");
728 params_->set(
"Output Stream", outputStream_);
729 if (! printer_.is_null ()) {
730 printer_->setOStream (outputStream_);
737 outputFreq_ = params->
get (
"Output Frequency", outputFreq_default_);
741 params_->set (
"Output Frequency", outputFreq_);
742 if (! outputTest_.is_null ()) {
743 outputTest_->setOutputFrequency (outputFreq_);
748 if (printer_.is_null ()) {
753 bool changedOrthoType =
false;
755 std::string tempOrthoType = params->
get (
"Orthogonalization", orthoType_default_);
756 if (tempOrthoType != orthoType_) {
757 orthoType_ = tempOrthoType;
758 changedOrthoType =
true;
761 params_->set(
"Orthogonalization", orthoType_);
764 if (params->
isParameter (
"Orthogonalization Constant")) {
765 if (params->
isType<MagnitudeType> (
"Orthogonalization Constant")) {
766 orthoKappa_ = params->
get (
"Orthogonalization Constant",
770 orthoKappa_ = params->
get (
"Orthogonalization Constant",
775 params_->set (
"Orthogonalization Constant", orthoKappa_);
776 if (orthoType_ ==
"DGKS") {
777 if (orthoKappa_ > 0 && ! ortho_.is_null() && !changedOrthoType) {
779 rcp_dynamic_cast<ortho_type> (ortho_)->setDepTol (orthoKappa_);
785 if (ortho_.is_null() || changedOrthoType) {
788 if (orthoType_==
"DGKS" && orthoKappa_ > 0) {
790 paramsOrtho->
set (
"depTol", orthoKappa_ );
793 ortho_ = factory.
makeMatOrthoManager (orthoType_, Teuchos::null, printer_, label_, paramsOrtho);
799 if (params->
isParameter (
"Convergence Tolerance")) {
800 if (params->
isType<MagnitudeType> (
"Convergence Tolerance")) {
801 convtol_ = params->
get (
"Convergence Tolerance",
809 params_->set (
"Convergence Tolerance", convtol_);
810 if (! impConvTest_.is_null ()) {
811 impConvTest_->setTolerance (convtol_);
813 if (! expConvTest_.is_null ()) {
814 expConvTest_->setTolerance (convtol_);
819 bool userDefinedResidualScalingUpdated =
false;
820 if (params->
isParameter (
"User Defined Residual Scaling")) {
822 if (params->
isType<MagnitudeType> (
"User Defined Residual Scaling")) {
823 tempResScaleFactor = params->
get (
"User Defined Residual Scaling",
827 tempResScaleFactor = params->
get (
"User Defined Residual Scaling",
832 if (resScaleFactor_ != tempResScaleFactor) {
833 resScaleFactor_ = tempResScaleFactor;
834 userDefinedResidualScalingUpdated =
true;
837 if(userDefinedResidualScalingUpdated)
839 if (! params->
isParameter (
"Implicit Residual Scaling") && ! impConvTest_.is_null ()) {
841 if(impResScale_ ==
"User Provided")
844 catch (std::exception& e) {
849 if (! params->
isParameter (
"Explicit Residual Scaling") && ! expConvTest_.is_null ()) {
851 if(expResScale_ ==
"User Provided")
854 catch (std::exception& e) {
863 if (params->
isParameter (
"Implicit Residual Scaling")) {
864 const std::string tempImpResScale =
865 Teuchos::getParameter<std::string> (*params,
"Implicit Residual Scaling");
868 if (impResScale_ != tempImpResScale) {
870 impResScale_ = tempImpResScale;
873 params_->set (
"Implicit Residual Scaling", impResScale_);
874 if (! impConvTest_.is_null ()) {
876 if(impResScale_ ==
"User Provided")
877 impConvTest_->defineScaleForm (impResScaleType,
Belos::TwoNorm, resScaleFactor_);
881 catch (std::exception& e) {
887 else if (userDefinedResidualScalingUpdated) {
890 if (! impConvTest_.is_null ()) {
892 if(impResScale_ ==
"User Provided")
893 impConvTest_->defineScaleForm (impResScaleType,
Belos::TwoNorm, resScaleFactor_);
895 catch (std::exception& e) {
903 if (params->
isParameter (
"Explicit Residual Scaling")) {
904 const std::string tempExpResScale =
905 Teuchos::getParameter<std::string> (*params,
"Explicit Residual Scaling");
908 if (expResScale_ != tempExpResScale) {
910 expResScale_ = tempExpResScale;
913 params_->set (
"Explicit Residual Scaling", expResScale_);
914 if (! expConvTest_.is_null ()) {
916 if(expResScale_ ==
"User Provided")
917 expConvTest_->defineScaleForm (expResScaleType,
Belos::TwoNorm, resScaleFactor_);
921 catch (std::exception& e) {
927 else if (userDefinedResidualScalingUpdated) {
930 if (! expConvTest_.is_null ()) {
932 if(expResScale_ ==
"User Provided")
933 expConvTest_->defineScaleForm (expResScaleType,
Belos::TwoNorm, resScaleFactor_);
935 catch (std::exception& e) {
944 if (params->
isParameter (
"Show Maximum Residual Norm Only")) {
945 showMaxResNormOnly_ =
946 Teuchos::getParameter<bool> (*params,
"Show Maximum Residual Norm Only");
949 params_->set (
"Show Maximum Residual Norm Only", showMaxResNormOnly_);
950 if (! impConvTest_.is_null ()) {
951 impConvTest_->setShowMaxResNormOnly (showMaxResNormOnly_);
953 if (! expConvTest_.is_null ()) {
954 expConvTest_->setShowMaxResNormOnly (showMaxResNormOnly_);
962 defQuorum_ = params->
get(
"Deflation Quorum", defQuorum_);
964 defQuorum_ > blockSize_, std::invalid_argument,
965 "Belos::PseudoBlockGmresSolMgr::setParameters: "
966 "The \"Deflation Quorum\" parameter (= " << defQuorum_ <<
") must not be "
967 "larger than \"Block Size\" (= " << blockSize_ <<
").");
968 params_->set (
"Deflation Quorum", defQuorum_);
969 if (! impConvTest_.is_null ()) {
970 impConvTest_->setQuorum (defQuorum_);
972 if (! expConvTest_.is_null ()) {
973 expConvTest_->setQuorum (defQuorum_);
978 if (timerSolve_ == Teuchos::null) {
979 std::string solveLabel = label_ +
": PseudoBlockGmresSolMgr total solve time";
980 #ifdef BELOS_TEUCHOS_TIME_MONITOR
990 template<
class ScalarType,
class MV,
class OP>
997 userConvStatusTest_ = userConvStatusTest;
998 comboType_ = comboType;
1001 template<
class ScalarType,
class MV,
class OP>
1007 debugStatusTest_ = debugStatusTest;
1012 template<
class ScalarType,
class MV,
class OP>
1025 "The relative residual tolerance that needs to be achieved by the\n"
1026 "iterative solver in order for the linear system to be declared converged.");
1027 pl->
set(
"Maximum Restarts", static_cast<int>(maxRestarts_default_),
1028 "The maximum number of restarts allowed for each\n"
1029 "set of RHS solved.");
1030 pl->
set(
"Maximum Iterations", static_cast<int>(maxIters_default_),
1031 "The maximum number of block iterations allowed for each\n"
1032 "set of RHS solved.");
1033 pl->
set(
"Num Blocks", static_cast<int>(numBlocks_default_),
1034 "The maximum number of vectors allowed in the Krylov subspace\n"
1035 "for each set of RHS solved.");
1036 pl->
set(
"Block Size", static_cast<int>(blockSize_default_),
1037 "The number of RHS solved simultaneously.");
1038 pl->
set(
"Verbosity", static_cast<int>(verbosity_default_),
1039 "What type(s) of solver information should be outputted\n"
1040 "to the output stream.");
1041 pl->
set(
"Output Style", static_cast<int>(outputStyle_default_),
1042 "What style is used for the solver information outputted\n"
1043 "to the output stream.");
1044 pl->
set(
"Output Frequency", static_cast<int>(outputFreq_default_),
1045 "How often convergence information should be outputted\n"
1046 "to the output stream.");
1047 pl->
set(
"Deflation Quorum", static_cast<int>(defQuorum_default_),
1048 "The number of linear systems that need to converge before\n"
1049 "they are deflated. This number should be <= block size.");
1050 pl->
set(
"Output Stream", Teuchos::rcpFromRef(std::cout),
1051 "A reference-counted pointer to the output stream where all\n"
1052 "solver output is sent.");
1053 pl->
set(
"Show Maximum Residual Norm Only", static_cast<bool>(showMaxResNormOnly_default_),
1054 "When convergence information is printed, only show the maximum\n"
1055 "relative residual norm when the block size is greater than one.");
1056 pl->
set(
"Implicit Residual Scaling", static_cast<const char *>(impResScale_default_),
1057 "The type of scaling used in the implicit residual convergence test.");
1058 pl->
set(
"Explicit Residual Scaling", static_cast<const char *>(expResScale_default_),
1059 "The type of scaling used in the explicit residual convergence test.");
1060 pl->
set(
"Timer Label", static_cast<const char *>(label_default_),
1061 "The string to use as a prefix for the timer labels.");
1062 pl->
set(
"Orthogonalization", static_cast<const char *>(orthoType_default_),
1063 "The type of orthogonalization to use.");
1065 "The constant used by DGKS orthogonalization to determine\n"
1066 "whether another step of classical Gram-Schmidt is necessary.");
1067 pl->
sublist(
"User Status Tests");
1068 pl->
set(
"User Status Tests Combo Type",
"SEQ",
1069 "Type of logical combination operation of user-defined\n"
1070 "and/or solver-specific status tests.");
1077 template<
class ScalarType,
class MV,
class OP>
1097 Teuchos::rcp(
new StatusTestGenResNorm_t( convtol_, defQuorum_ ) );
1098 if(impResScale_ ==
"User Provided")
1103 tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
1104 impConvTest_ = tmpImpConvTest;
1108 Teuchos::rcp(
new StatusTestGenResNorm_t( convtol_, defQuorum_ ) );
1109 tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit,
Belos::TwoNorm );
1110 if(expResScale_ ==
"User Provided")
1114 tmpExpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
1115 expConvTest_ = tmpExpConvTest;
1118 convTest_ =
Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
1125 Teuchos::rcp(
new StatusTestImpResNorm_t( convtol_, defQuorum_ ) );
1126 if(impResScale_ ==
"User Provided")
1130 tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
1131 impConvTest_ = tmpImpConvTest;
1134 expConvTest_ = impConvTest_;
1135 convTest_ = impConvTest_;
1138 if (
nonnull(userConvStatusTest_) ) {
1141 if (tmpComboTest != Teuchos::null) {
1142 std::vector<Teuchos::RCP<StatusTest<ScalarType,MV,OP> > > tmpVec = tmpComboTest->getStatusTests();
1143 comboType_ = tmpComboTest->getComboType();
1144 const int numResTests =
static_cast<int>(tmpVec.size());
1145 auto newConvTest =
Teuchos::rcp(
new StatusTestCombo_t(comboType_));
1146 newConvTest->addStatusTest(convTest_);
1147 for (
int j = 0; j < numResTests; ++j) {
1148 newConvTest->addStatusTest(tmpVec[j]);
1150 convTest_ = newConvTest;
1155 new StatusTestCombo_t( comboType_, convTest_, userConvStatusTest_ ) );
1163 sTest_ =
Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
1166 if (
nonnull(debugStatusTest_) ) {
1168 Teuchos::rcp_dynamic_cast<StatusTestCombo_t>(sTest_)->addStatusTest( debugStatusTest_ );
1173 StatusTestOutputFactory<ScalarType,MV,OP> stoFactory( outputStyle_, taggedTests_ );
1177 std::string solverDesc =
" Pseudo Block Gmres ";
1178 outputTest_->setSolverDesc( solverDesc );
1189 template<
class ScalarType,
class MV,
class OP>
1195 if (!isSet_) { setParameters( params_ ); }
1198 "Belos::PseudoBlockGmresSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
1201 if (!isSTSet_ || (!expResTest_ && !
Teuchos::is_null(problem_->getLeftPrec())) ) {
1203 "Belos::BlockGmresSolMgr::solve(): Linear problem and requested status tests are incompatible.");
1208 int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
1209 int numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
1211 std::vector<int> currIdx( numCurrRHS );
1212 blockSize_ = numCurrRHS;
1213 for (
int i=0; i<numCurrRHS; ++i)
1214 { currIdx[i] = startPtr+i; }
1217 problem_->setLSIndex( currIdx );
1222 plist.
set(
"Num Blocks",numBlocks_);
1225 outputTest_->reset();
1226 loaDetected_ =
false;
1229 bool isConverged =
true;
1239 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1243 while ( numRHS2Solve > 0 ) {
1246 std::vector<int> convRHSIdx;
1247 std::vector<int> currRHSIdx( currIdx );
1248 currRHSIdx.resize(numCurrRHS);
1251 block_gmres_iter->setNumBlocks( numBlocks_ );
1254 block_gmres_iter->resetNumIters();
1257 outputTest_->resetNumCalls();
1263 std::vector<int> index(1);
1264 Teuchos::RCP<MV> tmpV, R_0 = MVT::CloneCopy( *(problem_->getInitPrecResVec()), currIdx );
1265 newState.
V.resize( blockSize_ );
1266 newState.
Z.resize( blockSize_ );
1267 for (
int i=0; i<blockSize_; ++i) {
1269 tmpV = MVT::CloneViewNonConst( *R_0, index );
1276 int rank = ortho_->normalize( *tmpV, tmpZ );
1278 "Belos::PseudoBlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors.");
1280 newState.
V[i] = tmpV;
1281 newState.
Z[i] = tmpZ;
1285 block_gmres_iter->initialize(newState);
1286 int numRestarts = 0;
1292 block_gmres_iter->iterate();
1299 if ( convTest_->getStatus() ==
Passed ) {
1301 if ( expConvTest_->getLOADetected() ) {
1311 loaDetected_ =
true;
1313 "Belos::PseudoBlockGmresSolMgr::solve(): Warning! Solver has experienced a loss of accuracy!" << std::endl;
1314 isConverged =
false;
1318 std::vector<int> convIdx = expConvTest_->convIndices();
1322 if (convIdx.size() == currRHSIdx.size())
1329 problem_->setCurrLS();
1333 int curDim = oldState.
curDim;
1338 std::vector<int> oldRHSIdx( currRHSIdx );
1339 std::vector<int> defRHSIdx;
1340 for (
unsigned int i=0; i<currRHSIdx.size(); ++i) {
1342 for (
unsigned int j=0; j<convIdx.size(); ++j) {
1343 if (currRHSIdx[i] == convIdx[j]) {
1349 defRHSIdx.push_back( i );
1352 defState.
V.push_back( Teuchos::rcp_const_cast<MV>( oldState.
V[i] ) );
1357 currRHSIdx[have] = currRHSIdx[i];
1361 defRHSIdx.resize(currRHSIdx.size()-have);
1362 currRHSIdx.resize(have);
1370 problem_->setLSIndex( convIdx );
1373 problem_->updateSolution( defUpdate,
true );
1377 problem_->setLSIndex( currRHSIdx );
1380 defState.
curDim = curDim;
1383 block_gmres_iter->initialize(defState);
1390 else if ( maxIterTest_->getStatus() ==
Passed ) {
1392 isConverged =
false;
1400 else if ( block_gmres_iter->getCurSubspaceDim() == block_gmres_iter->getMaxSubspaceDim() ) {
1402 if ( numRestarts >= maxRestarts_ ) {
1403 isConverged =
false;
1408 printer_->stream(
Debug) <<
" Performing restart number " << numRestarts <<
" of " << maxRestarts_ << std::endl << std::endl;
1412 problem_->updateSolution( update,
true );
1419 newstate.
V.resize(currRHSIdx.size());
1420 newstate.
Z.resize(currRHSIdx.size());
1424 R_0 = MVT::Clone( *(problem_->getInitPrecResVec()), currRHSIdx.size() );
1425 problem_->computeCurrPrecResVec( &*R_0 );
1426 for (
unsigned int i=0; i<currRHSIdx.size(); ++i) {
1429 tmpV = MVT::CloneViewNonConst( *R_0, index );
1436 int rank = ortho_->normalize( *tmpV, tmpZ );
1438 "Belos::PseudoBlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors after the restart.");
1440 newstate.
V[i] = tmpV;
1441 newstate.
Z[i] = tmpZ;
1446 block_gmres_iter->initialize(newstate);
1459 "Belos::PseudoBlockGmresSolMgr::solve(): Invalid return from PseudoBlockGmresIter::iterate().");
1465 block_gmres_iter->updateLSQR( block_gmres_iter->getCurSubspaceDim() );
1468 sTest_->checkStatus( &*block_gmres_iter );
1469 if (convTest_->getStatus() !=
Passed)
1470 isConverged =
false;
1475 achievedTol_ = MT::one();
1477 MVT::MvInit( *X, SCT::zero() );
1478 printer_->stream(
Warnings) <<
"Belos::PseudoBlockGmresSolMgr::solve(): Warning! NaN has been detected!"
1482 catch (
const std::exception &e) {
1483 printer_->stream(
Errors) <<
"Error! Caught std::exception in PseudoBlockGmresIter::iterate() at iteration "
1484 << block_gmres_iter->getNumIters() << std::endl
1485 << e.what() << std::endl;
1492 if (
nonnull(userConvStatusTest_)) {
1495 problem_->updateSolution( update,
true );
1497 else if (
nonnull(expConvTest_->getSolution())) {
1501 MVT::MvAddMv( 0.0, *newX, 1.0, *newX, *curX );
1506 problem_->updateSolution( update,
true );
1510 problem_->setCurrLS();
1513 startPtr += numCurrRHS;
1514 numRHS2Solve -= numCurrRHS;
1515 if ( numRHS2Solve > 0 ) {
1516 numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
1518 blockSize_ = numCurrRHS;
1519 currIdx.resize( numCurrRHS );
1520 for (
int i=0; i<numCurrRHS; ++i)
1521 { currIdx[i] = startPtr+i; }
1524 if (defQuorum_ > blockSize_) {
1525 if (impConvTest_ != Teuchos::null)
1526 impConvTest_->setQuorum( blockSize_ );
1527 if (expConvTest_ != Teuchos::null)
1528 expConvTest_->setQuorum( blockSize_ );
1532 problem_->setLSIndex( currIdx );
1535 currIdx.resize( numRHS2Solve );
1546 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1555 numIters_ = maxIterTest_->getNumIters();
1568 const std::vector<MagnitudeType>* pTestValues = NULL;
1570 pTestValues = expConvTest_->getTestValue();
1571 if (pTestValues == NULL || pTestValues->size() < 1) {
1572 pTestValues = impConvTest_->getTestValue();
1577 pTestValues = impConvTest_->getTestValue();
1580 "Belos::PseudoBlockGmresSolMgr::solve(): The implicit convergence test's "
1581 "getTestValue() method returned NULL. Please report this bug to the "
1582 "Belos developers.");
1584 "Belos::PseudoBlockGmresSolMgr::solve(): The implicit convergence test's "
1585 "getTestValue() method returned a vector of length zero. Please report "
1586 "this bug to the Belos developers.");
1591 achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
1594 if (!isConverged || loaDetected_) {
1601 template<
class ScalarType,
class MV,
class OP>
1604 std::ostringstream out;
1606 out <<
"\"Belos::PseudoBlockGmresSolMgr\": {";
1607 if (this->getObjectLabel () !=
"") {
1608 out <<
"Label: " << this->getObjectLabel () <<
", ";
1610 out <<
"Num Blocks: " << numBlocks_
1611 <<
", Maximum Iterations: " << maxIters_
1612 <<
", Maximum Restarts: " << maxRestarts_
1613 <<
", Convergence Tolerance: " << convtol_
1619 template<
class ScalarType,
class MV,
class OP>
1641 out <<
"\"Belos::PseudoBlockGmresSolMgr\":" << endl;
1643 out <<
"Template parameters:" << endl;
1646 out <<
"ScalarType: " << TypeNameTraits<ScalarType>::name () << endl
1647 <<
"MV: " << TypeNameTraits<MV>::name () << endl
1648 <<
"OP: " << TypeNameTraits<OP>::name () << endl;
1650 if (this->getObjectLabel () !=
"") {
1651 out <<
"Label: " << this->getObjectLabel () << endl;
1653 out <<
"Num Blocks: " << numBlocks_ << endl
1654 <<
"Maximum Iterations: " << maxIters_ << endl
1655 <<
"Maximum Restarts: " << maxRestarts_ << endl
1656 <<
"Convergence Tolerance: " << convtol_ << endl;
static const double orthoKappa
DGKS orthogonalization constant.
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...
void setParameters(const Teuchos::RCP< Teuchos::ParameterList > ¶ms) override
Set the parameters the solver manager should use to solve the linear problem.
Class which manages the output and verbosity of the Belos solvers.
bool is_null(const boost::shared_ptr< T > &p)
PseudoBlockGmresSolMgr()
Empty constructor.
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...
T & get(const std::string &name, T def_value)
PseudoBlockGmresSolMgrOrthoFailure(const std::string &what_arg)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
bool is_null(const std::shared_ptr< T > &p)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void setDebugStatusTest(const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &debugStatusTest) override
Set a debug status test that will be checked at the same time as the top-level status test...
An abstract class of StatusTest for stopping criteria using residual norms.
std::vector< Teuchos::RCP< const MV > > V
The current Krylov basis.
An implementation of StatusTestResNorm using a family of residual norms.
Ordinal numParams() const
static const double convTol
Default convergence tolerance.
A pure virtual class for defining the status tests for the Belos iterative solvers.
Structure to contain pointers to PseudoBlockGmresIter state variables.
Teuchos::Array< Teuchos::RCP< Teuchos::Time > > getTimers() const
Return the timers for this object.
A factory class for generating StatusTestOutput objects.
Traits class which defines basic operations on multivectors.
bool isParameter(const std::string &name) const
A Belos::StatusTest class for specifying a maximum number of iterations.
ResetType
How to reset the solver.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Pure virtual base class which describes the basic interface for a solver manager. ...
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const override
Print the object with the given verbosity level to a FancyOStream.
static void summarize(Ptr< const Comm< int > > comm, std::ostream &out=std::cout, const bool alwaysWriteLocal=false, const bool writeGlobalStats=true, const bool writeZeroTimers=true, const ECounterSetOp setOp=Intersection, const std::string &filter="", const bool ignoreZeroTimers=false)
bool isSublist(const std::string &name) const
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.
void reset(const ResetType type) override
Performs a reset of the solver manager specified by the ResetType. This informs the solver manager th...
Interface to standard and "pseudoblock" GMRES.
int getNumIters() const override
Iteration count for the most recent call to solve().
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 ...
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
A list of valid default parameters for this solver.
void validateParameters(ParameterList const &validParamList, int const depth=1000, EValidateUsed const validateUsed=VALIDATE_USED_ENABLED, EValidateDefaults const validateDefaults=VALIDATE_DEFAULTS_ENABLED) const
static const EVerbosityLevel verbLevel_default
bool isLOADetected() const override
Whether a "loss of accuracy" was detected during the last solve().
bool nonnull(const boost::shared_ptr< T > &p)
std::vector< Teuchos::RCP< const Teuchos::SerialDenseVector< int, ScalarType > > > Z
The current right-hand side of the least squares system RY = Z.
void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem) override
Set the linear problem to solve.
Belos concrete class for performing the pseudo-block GMRES iteration.
bool isType(const std::string &name) const
Teuchos::RCP< SolverManager< ScalarType, MV, OP > > clone() const override
clone for Inverted Injection (DII)
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.
std::string description() const override
Return a one-line description of this object.
const StatusTestResNorm< ScalarType, MV, OP > * getResidualStatusTest() const
Return the residual status test.
Class which defines basic traits for the operator type.
ParameterList & sublist(const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
int curDim
The current dimension of the reduction.
Parent class to all Belos exceptions.
ReturnType solve() override
This method performs possibly repeated calls to the underlying linear solver's iterate() routine unti...
Default parameters common to most Belos solvers.
std::vector< Teuchos::RCP< const Teuchos::SerialDenseVector< int, ScalarType > > > sn
The current Given's rotation coefficients.
static const double resScaleFactor
User-defined residual scaling factor.
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.
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const override
The current parameters for this solver.
PseudoBlockGmresIterOrthoFailure is thrown when the orthogonalization manager is unable to generate o...
const LinearProblem< ScalarType, MV, OP > & getProblem() const override
Return a reference to the linear problem being solved by this solver manager.
PseudoBlockGmresSolMgrOrthoFailure is thrown when the orthogonalization manager is unable to generate...
virtual void setUserConvStatusTest(const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &userConvStatusTest, const typename StatusTestCombo< ScalarType, MV, OP >::ComboType &comboType=StatusTestCombo< ScalarType, MV, OP >::SEQ) override
Set a custom status test.
MagnitudeType achievedTol() const override
Tolerance achieved by the last solve() invocation.
Teuchos::RCP< Belos::MatOrthoManager< Scalar, MV, OP > > makeMatOrthoManager(const std::string &ortho, const Teuchos::RCP< const OP > &M, const Teuchos::RCP< OutputManager< Scalar > > &, const std::string &label, const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Return an instance of the specified MatOrthoManager subclass.