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>
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;
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
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")) {
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) {
789 paramsOrtho->
set (
"depTol", orthoKappa_ );
798 if (params->
isParameter (
"Convergence Tolerance")) {
800 convtol_ = params->
get (
"Convergence Tolerance",
808 params_->set (
"Convergence Tolerance", convtol_);
809 if (! impConvTest_.is_null ()) {
810 impConvTest_->setTolerance (convtol_);
812 if (! expConvTest_.is_null ()) {
813 expConvTest_->setTolerance (convtol_);
818 bool userDefinedResidualScalingUpdated =
false;
819 if (params->
isParameter (
"User Defined Residual Scaling")) {
822 tempResScaleFactor = params->
get (
"User Defined Residual Scaling",
826 tempResScaleFactor = params->
get (
"User Defined Residual Scaling",
831 if (resScaleFactor_ != tempResScaleFactor) {
832 resScaleFactor_ = tempResScaleFactor;
833 userDefinedResidualScalingUpdated =
true;
836 if(userDefinedResidualScalingUpdated)
838 if (! params->
isParameter (
"Implicit Residual Scaling") && ! impConvTest_.is_null ()) {
840 if(impResScale_ ==
"User Provided")
843 catch (std::exception& e) {
848 if (! params->
isParameter (
"Explicit Residual Scaling") && ! expConvTest_.is_null ()) {
850 if(expResScale_ ==
"User Provided")
853 catch (std::exception& e) {
862 if (params->
isParameter (
"Implicit Residual Scaling")) {
863 const std::string tempImpResScale =
864 Teuchos::getParameter<std::string> (*params,
"Implicit Residual Scaling");
867 if (impResScale_ != tempImpResScale) {
869 impResScale_ = tempImpResScale;
872 params_->set (
"Implicit Residual Scaling", impResScale_);
873 if (! impConvTest_.is_null ()) {
875 if(impResScale_ ==
"User Provided")
876 impConvTest_->defineScaleForm (impResScaleType,
Belos::TwoNorm, resScaleFactor_);
880 catch (std::exception& e) {
886 else if (userDefinedResidualScalingUpdated) {
889 if (! impConvTest_.is_null ()) {
891 if(impResScale_ ==
"User Provided")
892 impConvTest_->defineScaleForm (impResScaleType,
Belos::TwoNorm, resScaleFactor_);
894 catch (std::exception& e) {
902 if (params->
isParameter (
"Explicit Residual Scaling")) {
903 const std::string tempExpResScale =
904 Teuchos::getParameter<std::string> (*params,
"Explicit Residual Scaling");
907 if (expResScale_ != tempExpResScale) {
909 expResScale_ = tempExpResScale;
912 params_->set (
"Explicit Residual Scaling", expResScale_);
913 if (! expConvTest_.is_null ()) {
915 if(expResScale_ ==
"User Provided")
916 expConvTest_->defineScaleForm (expResScaleType,
Belos::TwoNorm, resScaleFactor_);
920 catch (std::exception& e) {
926 else if (userDefinedResidualScalingUpdated) {
929 if (! expConvTest_.is_null ()) {
931 if(expResScale_ ==
"User Provided")
932 expConvTest_->defineScaleForm (expResScaleType,
Belos::TwoNorm, resScaleFactor_);
934 catch (std::exception& e) {
943 if (params->
isParameter (
"Show Maximum Residual Norm Only")) {
944 showMaxResNormOnly_ =
945 Teuchos::getParameter<bool> (*params,
"Show Maximum Residual Norm Only");
948 params_->set (
"Show Maximum Residual Norm Only", showMaxResNormOnly_);
949 if (! impConvTest_.is_null ()) {
950 impConvTest_->setShowMaxResNormOnly (showMaxResNormOnly_);
952 if (! expConvTest_.is_null ()) {
953 expConvTest_->setShowMaxResNormOnly (showMaxResNormOnly_);
961 defQuorum_ = params->
get(
"Deflation Quorum", defQuorum_);
963 defQuorum_ > blockSize_, std::invalid_argument,
964 "Belos::PseudoBlockGmresSolMgr::setParameters: "
965 "The \"Deflation Quorum\" parameter (= " << defQuorum_ <<
") must not be "
966 "larger than \"Block Size\" (= " << blockSize_ <<
").");
967 params_->set (
"Deflation Quorum", defQuorum_);
968 if (! impConvTest_.is_null ()) {
969 impConvTest_->setQuorum (defQuorum_);
971 if (! expConvTest_.is_null ()) {
972 expConvTest_->setQuorum (defQuorum_);
978 std::string solveLabel = label_ +
": PseudoBlockGmresSolMgr total solve time";
979 #ifdef BELOS_TEUCHOS_TIME_MONITOR
989 template<
class ScalarType,
class MV,
class OP>
996 userConvStatusTest_ = userConvStatusTest;
997 comboType_ = comboType;
1000 template<
class ScalarType,
class MV,
class OP>
1006 debugStatusTest_ = debugStatusTest;
1011 template<
class ScalarType,
class MV,
class OP>
1024 "The relative residual tolerance that needs to be achieved by the\n"
1025 "iterative solver in order for the linear system to be declared converged.");
1026 pl->
set(
"Maximum Restarts", static_cast<int>(maxRestarts_default_),
1027 "The maximum number of restarts allowed for each\n"
1028 "set of RHS solved.");
1029 pl->
set(
"Maximum Iterations", static_cast<int>(maxIters_default_),
1030 "The maximum number of block iterations allowed for each\n"
1031 "set of RHS solved.");
1032 pl->
set(
"Num Blocks", static_cast<int>(numBlocks_default_),
1033 "The maximum number of vectors allowed in the Krylov subspace\n"
1034 "for each set of RHS solved.");
1035 pl->
set(
"Block Size", static_cast<int>(blockSize_default_),
1036 "The number of RHS solved simultaneously.");
1037 pl->
set(
"Verbosity", static_cast<int>(verbosity_default_),
1038 "What type(s) of solver information should be outputted\n"
1039 "to the output stream.");
1040 pl->
set(
"Output Style", static_cast<int>(outputStyle_default_),
1041 "What style is used for the solver information outputted\n"
1042 "to the output stream.");
1043 pl->
set(
"Output Frequency", static_cast<int>(outputFreq_default_),
1044 "How often convergence information should be outputted\n"
1045 "to the output stream.");
1046 pl->
set(
"Deflation Quorum", static_cast<int>(defQuorum_default_),
1047 "The number of linear systems that need to converge before\n"
1048 "they are deflated. This number should be <= block size.");
1049 pl->
set(
"Output Stream", Teuchos::rcpFromRef(std::cout),
1050 "A reference-counted pointer to the output stream where all\n"
1051 "solver output is sent.");
1052 pl->
set(
"Show Maximum Residual Norm Only", static_cast<bool>(showMaxResNormOnly_default_),
1053 "When convergence information is printed, only show the maximum\n"
1054 "relative residual norm when the block size is greater than one.");
1055 pl->
set(
"Implicit Residual Scaling", static_cast<const char *>(impResScale_default_),
1056 "The type of scaling used in the implicit residual convergence test.");
1057 pl->
set(
"Explicit Residual Scaling", static_cast<const char *>(expResScale_default_),
1058 "The type of scaling used in the explicit residual convergence test.");
1059 pl->
set(
"Timer Label", static_cast<const char *>(label_default_),
1060 "The string to use as a prefix for the timer labels.");
1061 pl->
set(
"Orthogonalization", static_cast<const char *>(orthoType_default_),
1062 "The type of orthogonalization to use.");
1064 "The constant used by DGKS orthogonalization to determine\n"
1065 "whether another step of classical Gram-Schmidt is necessary.");
1066 pl->
sublist(
"User Status Tests");
1067 pl->
set(
"User Status Tests Combo Type",
"SEQ",
1068 "Type of logical combination operation of user-defined\n"
1069 "and/or solver-specific status tests.");
1076 template<
class ScalarType,
class MV,
class OP>
1096 Teuchos::rcp(
new StatusTestGenResNorm_t( convtol_, defQuorum_ ) );
1097 if(impResScale_ ==
"User Provided")
1102 tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
1103 impConvTest_ = tmpImpConvTest;
1107 Teuchos::rcp(
new StatusTestGenResNorm_t( convtol_, defQuorum_ ) );
1108 tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit,
Belos::TwoNorm );
1109 if(expResScale_ ==
"User Provided")
1113 tmpExpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
1114 expConvTest_ = tmpExpConvTest;
1117 convTest_ =
Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
1124 Teuchos::rcp(
new StatusTestImpResNorm_t( convtol_, defQuorum_ ) );
1125 if(impResScale_ ==
"User Provided")
1129 tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
1130 impConvTest_ = tmpImpConvTest;
1133 expConvTest_ = impConvTest_;
1134 convTest_ = impConvTest_;
1137 if (
nonnull(userConvStatusTest_) ) {
1141 std::vector<Teuchos::RCP<StatusTest<ScalarType,MV,OP> > > tmpVec = tmpComboTest->getStatusTests();
1142 comboType_ = tmpComboTest->getComboType();
1143 const int numResTests =
static_cast<int>(tmpVec.size());
1144 auto newConvTest =
Teuchos::rcp(
new StatusTestCombo_t(comboType_));
1145 newConvTest->addStatusTest(convTest_);
1146 for (
int j = 0; j < numResTests; ++j) {
1147 newConvTest->addStatusTest(tmpVec[j]);
1149 convTest_ = newConvTest;
1154 new StatusTestCombo_t( comboType_, convTest_, userConvStatusTest_ ) );
1162 sTest_ =
Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
1165 if (
nonnull(debugStatusTest_) ) {
1167 Teuchos::rcp_dynamic_cast<StatusTestCombo_t>(sTest_)->addStatusTest( debugStatusTest_ );
1176 std::string solverDesc =
" Pseudo Block Gmres ";
1177 outputTest_->setSolverDesc( solverDesc );
1188 template<
class ScalarType,
class MV,
class OP>
1194 if (!isSet_) { setParameters( params_ ); }
1197 "Belos::PseudoBlockGmresSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
1200 if (!isSTSet_ || (!expResTest_ && !
Teuchos::is_null(problem_->getLeftPrec())) ) {
1202 "Belos::BlockGmresSolMgr::solve(): Linear problem and requested status tests are incompatible.");
1207 int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
1208 int numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
1210 std::vector<int> currIdx( numCurrRHS );
1211 blockSize_ = numCurrRHS;
1212 for (
int i=0; i<numCurrRHS; ++i)
1213 { currIdx[i] = startPtr+i; }
1216 problem_->setLSIndex( currIdx );
1221 plist.
set(
"Num Blocks",numBlocks_);
1224 outputTest_->reset();
1225 loaDetected_ =
false;
1228 bool isConverged =
true;
1238 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1242 while ( numRHS2Solve > 0 ) {
1245 std::vector<int> convRHSIdx;
1246 std::vector<int> currRHSIdx( currIdx );
1247 currRHSIdx.resize(numCurrRHS);
1250 block_gmres_iter->setNumBlocks( numBlocks_ );
1253 block_gmres_iter->resetNumIters();
1256 outputTest_->resetNumCalls();
1262 std::vector<int> index(1);
1263 Teuchos::RCP<MV> tmpV, R_0 = MVT::CloneCopy( *(problem_->getInitPrecResVec()), currIdx );
1264 newState.
V.resize( blockSize_ );
1265 newState.
Z.resize( blockSize_ );
1266 for (
int i=0; i<blockSize_; ++i) {
1268 tmpV = MVT::CloneViewNonConst( *R_0, index );
1275 int rank = ortho_->normalize( *tmpV, tmpZ );
1277 "Belos::PseudoBlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors.");
1279 newState.
V[i] = tmpV;
1280 newState.
Z[i] = tmpZ;
1284 block_gmres_iter->initialize(newState);
1285 int numRestarts = 0;
1291 block_gmres_iter->iterate();
1298 if ( convTest_->getStatus() ==
Passed ) {
1300 if ( expConvTest_->getLOADetected() ) {
1310 loaDetected_ =
true;
1312 "Belos::PseudoBlockGmresSolMgr::solve(): Warning! Solver has experienced a loss of accuracy!" << std::endl;
1313 isConverged =
false;
1317 std::vector<int> convIdx = expConvTest_->convIndices();
1321 if (convIdx.size() == currRHSIdx.size())
1328 problem_->setCurrLS();
1332 int curDim = oldState.
curDim;
1337 std::vector<int> oldRHSIdx( currRHSIdx );
1338 std::vector<int> defRHSIdx;
1339 for (
unsigned int i=0; i<currRHSIdx.size(); ++i) {
1341 for (
unsigned int j=0; j<convIdx.size(); ++j) {
1342 if (currRHSIdx[i] == convIdx[j]) {
1348 defRHSIdx.push_back( i );
1351 defState.
V.push_back( Teuchos::rcp_const_cast<MV>( oldState.
V[i] ) );
1356 currRHSIdx[have] = currRHSIdx[i];
1360 defRHSIdx.resize(currRHSIdx.size()-have);
1361 currRHSIdx.resize(have);
1369 problem_->setLSIndex( convIdx );
1372 problem_->updateSolution( defUpdate,
true );
1376 problem_->setLSIndex( currRHSIdx );
1379 defState.
curDim = curDim;
1382 block_gmres_iter->initialize(defState);
1389 else if ( maxIterTest_->getStatus() ==
Passed ) {
1391 isConverged =
false;
1399 else if ( block_gmres_iter->getCurSubspaceDim() == block_gmres_iter->getMaxSubspaceDim() ) {
1401 if ( numRestarts >= maxRestarts_ ) {
1402 isConverged =
false;
1407 printer_->stream(
Debug) <<
" Performing restart number " << numRestarts <<
" of " << maxRestarts_ << std::endl << std::endl;
1411 problem_->updateSolution( update,
true );
1418 newstate.
V.resize(currRHSIdx.size());
1419 newstate.
Z.resize(currRHSIdx.size());
1423 R_0 = MVT::Clone( *(problem_->getInitPrecResVec()), currRHSIdx.size() );
1424 problem_->computeCurrPrecResVec( &*R_0 );
1425 for (
unsigned int i=0; i<currRHSIdx.size(); ++i) {
1428 tmpV = MVT::CloneViewNonConst( *R_0, index );
1435 int rank = ortho_->normalize( *tmpV, tmpZ );
1437 "Belos::PseudoBlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors after the restart.");
1439 newstate.
V[i] = tmpV;
1440 newstate.
Z[i] = tmpZ;
1445 block_gmres_iter->initialize(newstate);
1458 "Belos::PseudoBlockGmresSolMgr::solve(): Invalid return from PseudoBlockGmresIter::iterate().");
1464 block_gmres_iter->updateLSQR( block_gmres_iter->getCurSubspaceDim() );
1467 sTest_->checkStatus( &*block_gmres_iter );
1468 if (convTest_->getStatus() !=
Passed)
1469 isConverged =
false;
1474 achievedTol_ = MT::one();
1476 MVT::MvInit( *X, SCT::zero() );
1477 printer_->stream(
Warnings) <<
"Belos::PseudoBlockGmresSolMgr::solve(): Warning! NaN has been detected!"
1481 catch (
const std::exception &e) {
1482 printer_->stream(
Errors) <<
"Error! Caught std::exception in PseudoBlockGmresIter::iterate() at iteration "
1483 << block_gmres_iter->getNumIters() << std::endl
1484 << e.what() << std::endl;
1491 if (
nonnull(userConvStatusTest_)) {
1494 problem_->updateSolution( update,
true );
1496 else if (
nonnull(expConvTest_->getSolution())) {
1500 MVT::MvAddMv( 0.0, *newX, 1.0, *newX, *curX );
1505 problem_->updateSolution( update,
true );
1509 problem_->setCurrLS();
1512 startPtr += numCurrRHS;
1513 numRHS2Solve -= numCurrRHS;
1514 if ( numRHS2Solve > 0 ) {
1515 numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
1517 blockSize_ = numCurrRHS;
1518 currIdx.resize( numCurrRHS );
1519 for (
int i=0; i<numCurrRHS; ++i)
1520 { currIdx[i] = startPtr+i; }
1523 if (defQuorum_ > blockSize_) {
1525 impConvTest_->setQuorum( blockSize_ );
1527 expConvTest_->setQuorum( blockSize_ );
1531 problem_->setLSIndex( currIdx );
1534 currIdx.resize( numRHS2Solve );
1545 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1554 numIters_ = maxIterTest_->getNumIters();
1567 const std::vector<MagnitudeType>* pTestValues = NULL;
1569 pTestValues = expConvTest_->getTestValue();
1570 if (pTestValues == NULL || pTestValues->size() < 1) {
1571 pTestValues = impConvTest_->getTestValue();
1576 pTestValues = impConvTest_->getTestValue();
1579 "Belos::PseudoBlockGmresSolMgr::solve(): The implicit convergence test's "
1580 "getTestValue() method returned NULL. Please report this bug to the "
1581 "Belos developers.");
1583 "Belos::PseudoBlockGmresSolMgr::solve(): The implicit convergence test's "
1584 "getTestValue() method returned a vector of length zero. Please report "
1585 "this bug to the Belos developers.");
1590 achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
1593 if (!isConverged || loaDetected_) {
1600 template<
class ScalarType,
class MV,
class OP>
1603 std::ostringstream out;
1605 out <<
"\"Belos::PseudoBlockGmresSolMgr\": {";
1606 if (this->getObjectLabel () !=
"") {
1607 out <<
"Label: " << this->getObjectLabel () <<
", ";
1609 out <<
"Num Blocks: " << numBlocks_
1610 <<
", Maximum Iterations: " << maxIters_
1611 <<
", Maximum Restarts: " << maxRestarts_
1612 <<
", Convergence Tolerance: " << convtol_
1618 template<
class ScalarType,
class MV,
class OP>
1640 out <<
"\"Belos::PseudoBlockGmresSolMgr\":" << endl;
1642 out <<
"Template parameters:" << endl;
1645 out <<
"ScalarType: " << TypeNameTraits<ScalarType>::name () << endl
1646 <<
"MV: " << TypeNameTraits<MV>::name () << endl
1647 <<
"OP: " << TypeNameTraits<OP>::name () << endl;
1649 if (this->getObjectLabel () !=
"") {
1650 out <<
"Label: " << this->getObjectLabel () << endl;
1652 out <<
"Num Blocks: " << numBlocks_ << endl
1653 <<
"Maximum Iterations: " << maxIters_ << endl
1654 <<
"Maximum Restarts: " << maxRestarts_ << endl
1655 <<
"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.
Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > problem_
The current linear problem to solve.
bool is_null(const boost::shared_ptr< T > &p)
Teuchos::RCP< Teuchos::ParameterList > params_
PseudoBlockGmresSolMgr()
Empty constructor.
PseudoBlockGmresSolMgrLinearProblemFailure(const std::string &what_arg)
MagnitudeType resScaleFactor_
ScaleType
The type of scaling to use on the residual norm value.
Teuchos::RCP< StatusTestMaxIters< ScalarType, MV, OP > > maxIterTest_
Convergence test using the implicit residual norm(s), with an explicit residual norm(s) check for los...
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > debugStatusTest_
T & get(ParameterList &l, const std::string &name)
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)
A factory class for generating StatusTestOutput objects.
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.
RCP< ParameterList > sublist(const RCP< ParameterList > ¶mList, const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
static constexpr int outputFreq_default_
Ordinal numParams() const
Teuchos::RCP< StatusTestOutput< ScalarType, MV, OP > > outputTest_
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.
static constexpr const char * label_default_
A factory class for generating StatusTestOutput objects.
Teuchos::RCP< std::ostream > outputStream_
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > convTest_
Teuchos::RCP< MatOrthoManager< ScalarType, MV, OP > > ortho_
Traits class which defines basic operations on multivectors.
MagnitudeType achievedTol_
bool isParameter(const std::string &name) const
A Belos::StatusTest class for specifying a maximum number of iterations.
MultiVecTraits< ScalarType, MV > MVT
static constexpr int blockSize_default_
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 constexpr int numBlocks_default_
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
static constexpr int defQuorum_default_
static constexpr bool showMaxResNormOnly_default_
static constexpr int verbosity_default_
Teuchos::ScalarTraits< ScalarType >::magnitudeType MagnitudeType
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.
Teuchos::RCP< std::map< std::string, Teuchos::RCP< StatusTest< ScalarType, MV, OP > > > > taggedTests_
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().
Teuchos::RCP< Teuchos::Time > timerSolve_
ReturnType
Whether the Belos solve converged for all linear systems.
bool checkStatusTest()
Check current status tests against current linear problem.
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.
static constexpr const char * orthoType_default_
static constexpr int maxRestarts_default_
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
static constexpr const char * impResScale_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.
Teuchos::RCP< StatusTestOutput< ScalarType, MV, OP > > create(const Teuchos::RCP< OutputManager< ScalarType > > &printer, Teuchos::RCP< StatusTest< ScalarType, MV, OP > > test, int mod, int printStates)
Create the StatusTestOutput object specified by the outputStyle.
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.
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > sTest_
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.
int curDim
The current dimension of the reduction.
Teuchos::ScalarTraits< MagnitudeType > MT
static constexpr const char * expResScale_default_
static constexpr int maxIters_default_
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.
static constexpr int outputStyle_default_
Belos header file which uses auto-configuration information to include necessary C++ headers...
Teuchos::RCP< StatusTestResNorm< ScalarType, MV, OP > > expConvTest_
OperatorTraits< ScalarType, MV, OP > OPT
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.
Teuchos::RCP< OutputManager< ScalarType > > printer_
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.
Teuchos::ScalarTraits< ScalarType > SCT
StatusTestCombo< ScalarType, MV, OP >::ComboType comboType_
Teuchos::RCP< StatusTestResNorm< ScalarType, MV, OP > > impConvTest_
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.
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > userConvStatusTest_
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.
MagnitudeType orthoKappa_