42 #ifndef BELOS_BLOCK_GMRES_SOLMGR_HPP 43 #define BELOS_BLOCK_GMRES_SOLMGR_HPP 68 #include "Teuchos_BLAS.hpp" 69 #ifdef BELOS_TEUCHOS_TIME_MONITOR 70 #include "Teuchos_TimeMonitor.hpp" 124 template<
class ScalarType,
class MV,
class OP>
130 typedef Teuchos::ScalarTraits<ScalarType>
SCT;
131 typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType
MagnitudeType;
132 typedef Teuchos::ScalarTraits<MagnitudeType>
MT;
167 const Teuchos::RCP<Teuchos::ParameterList> &pl );
195 Teuchos::Array<Teuchos::RCP<Teuchos::Time> >
getTimers()
const {
232 void setParameters(
const Teuchos::RCP<Teuchos::ParameterList> ¶ms );
279 describe (Teuchos::FancyOStream& out,
280 const Teuchos::EVerbosityLevel verbLevel =
281 Teuchos::Describable::verbLevel_default)
const;
294 Teuchos::RCP<LinearProblem<ScalarType,MV,OP> >
problem_;
301 Teuchos::RCP<StatusTest<ScalarType,MV,OP> >
sTest_;
308 Teuchos::RCP<MatOrthoManager<ScalarType,MV,OP> >
ortho_;
352 template<
class ScalarType,
class MV,
class OP>
355 template<
class ScalarType,
class MV,
class OP>
358 template<
class ScalarType,
class MV,
class OP>
361 template<
class ScalarType,
class MV,
class OP>
364 template<
class ScalarType,
class MV,
class OP>
367 template<
class ScalarType,
class MV,
class OP>
370 template<
class ScalarType,
class MV,
class OP>
373 template<
class ScalarType,
class MV,
class OP>
376 template<
class ScalarType,
class MV,
class OP>
379 template<
class ScalarType,
class MV,
class OP>
382 template<
class ScalarType,
class MV,
class OP>
385 template<
class ScalarType,
class MV,
class OP>
388 template<
class ScalarType,
class MV,
class OP>
391 template<
class ScalarType,
class MV,
class OP>
394 template<
class ScalarType,
class MV,
class OP>
397 template<
class ScalarType,
class MV,
class OP>
400 template<
class ScalarType,
class MV,
class OP>
403 template<
class ScalarType,
class MV,
class OP>
408 template<
class ScalarType,
class MV,
class OP>
410 outputStream_(outputStream_default_),
411 convtol_(convtol_default_),
412 orthoKappa_(orthoKappa_default_),
413 achievedTol_(
Teuchos::ScalarTraits<typename
Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
414 maxRestarts_(maxRestarts_default_),
415 maxIters_(maxIters_default_),
417 blockSize_(blockSize_default_),
418 numBlocks_(numBlocks_default_),
419 verbosity_(verbosity_default_),
420 outputStyle_(outputStyle_default_),
421 outputFreq_(outputFreq_default_),
422 adaptiveBlockSize_(adaptiveBlockSize_default_),
423 showMaxResNormOnly_(showMaxResNormOnly_default_),
424 isFlexible_(flexibleGmres_default_),
425 expResTest_(expResTest_default_),
426 orthoType_(orthoType_default_),
427 impResScale_(impResScale_default_),
428 expResScale_(expResScale_default_),
429 label_(label_default_),
437 template<
class ScalarType,
class MV,
class OP>
440 const Teuchos::RCP<Teuchos::ParameterList> &pl) :
442 outputStream_(outputStream_default_),
443 convtol_(convtol_default_),
444 orthoKappa_(orthoKappa_default_),
445 achievedTol_(
Teuchos::ScalarTraits<typename
Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
446 maxRestarts_(maxRestarts_default_),
447 maxIters_(maxIters_default_),
449 blockSize_(blockSize_default_),
450 numBlocks_(numBlocks_default_),
451 verbosity_(verbosity_default_),
452 outputStyle_(outputStyle_default_),
453 outputFreq_(outputFreq_default_),
454 adaptiveBlockSize_(adaptiveBlockSize_default_),
455 showMaxResNormOnly_(showMaxResNormOnly_default_),
456 isFlexible_(flexibleGmres_default_),
457 expResTest_(expResTest_default_),
458 orthoType_(orthoType_default_),
459 impResScale_(impResScale_default_),
460 expResScale_(expResScale_default_),
461 label_(label_default_),
467 TEUCHOS_TEST_FOR_EXCEPTION(
problem_ == Teuchos::null, std::invalid_argument,
"Problem not given to solver manager.");
470 if ( !is_null(pl) ) {
477 template<
class ScalarType,
class MV,
class OP>
478 Teuchos::RCP<const Teuchos::ParameterList>
481 static Teuchos::RCP<const Teuchos::ParameterList> validPL;
482 if (is_null(validPL)) {
483 Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
484 pl->set(
"Convergence Tolerance", convtol_default_,
485 "The relative residual tolerance that needs to be achieved by the\n" 486 "iterative solver in order for the linear system to be declared converged." );
487 pl->set(
"Maximum Restarts", maxRestarts_default_,
488 "The maximum number of restarts allowed for each\n" 489 "set of RHS solved.");
491 "Maximum Iterations", maxIters_default_,
492 "The maximum number of block iterations allowed for each\n" 493 "set of RHS solved.");
494 pl->set(
"Num Blocks", numBlocks_default_,
495 "The maximum number of blocks allowed in the Krylov subspace\n" 496 "for each set of RHS solved.");
497 pl->set(
"Block Size", blockSize_default_,
498 "The number of vectors in each block. This number times the\n" 499 "number of blocks is the total Krylov subspace dimension.");
500 pl->set(
"Adaptive Block Size", adaptiveBlockSize_default_,
501 "Whether the solver manager should adapt the block size\n" 502 "based on the number of RHS to solve.");
503 pl->set(
"Verbosity", verbosity_default_,
504 "What type(s) of solver information should be outputted\n" 505 "to the output stream.");
506 pl->set(
"Output Style", outputStyle_default_,
507 "What style is used for the solver information outputted\n" 508 "to the output stream.");
509 pl->set(
"Output Frequency", outputFreq_default_,
510 "How often convergence information should be outputted\n" 511 "to the output stream.");
512 pl->set(
"Output Stream", outputStream_default_,
513 "A reference-counted pointer to the output stream where all\n" 514 "solver output is sent.");
515 pl->set(
"Show Maximum Residual Norm Only", showMaxResNormOnly_default_,
516 "When convergence information is printed, only show the maximum\n" 517 "relative residual norm when the block size is greater than one.");
518 pl->set(
"Flexible Gmres", flexibleGmres_default_,
519 "Whether the solver manager should use the flexible variant\n" 521 pl->set(
"Explicit Residual Test", expResTest_default_,
522 "Whether the explicitly computed residual should be used in the convergence test.");
523 pl->set(
"Implicit Residual Scaling", impResScale_default_,
524 "The type of scaling used in the implicit residual convergence test.");
525 pl->set(
"Explicit Residual Scaling", expResScale_default_,
526 "The type of scaling used in the explicit residual convergence test.");
527 pl->set(
"Timer Label", label_default_,
528 "The string to use as a prefix for the timer labels.");
530 pl->set(
"Orthogonalization", orthoType_default_,
531 "The type of orthogonalization to use: DGKS, ICGS, or IMGS.");
532 pl->set(
"Orthogonalization Constant",orthoKappa_default_,
533 "The constant used by DGKS orthogonalization to determine\n" 534 "whether another step of classical Gram-Schmidt is necessary.");
541 template<
class ScalarType,
class MV,
class OP>
546 if (params_ == Teuchos::null) {
547 params_ = Teuchos::rcp(
new Teuchos::ParameterList(*getValidParameters()) );
550 params->validateParameters(*getValidParameters());
554 if (params->isParameter(
"Maximum Restarts")) {
555 maxRestarts_ = params->get(
"Maximum Restarts",maxRestarts_default_);
558 params_->set(
"Maximum Restarts", maxRestarts_);
562 if (params->isParameter(
"Maximum Iterations")) {
563 maxIters_ = params->get(
"Maximum Iterations",maxIters_default_);
566 params_->set(
"Maximum Iterations", maxIters_);
567 if (maxIterTest_!=Teuchos::null)
568 maxIterTest_->setMaxIters( maxIters_ );
572 if (params->isParameter(
"Block Size")) {
573 blockSize_ = params->get(
"Block Size",blockSize_default_);
574 TEUCHOS_TEST_FOR_EXCEPTION(blockSize_ <= 0, std::invalid_argument,
575 "Belos::BlockGmresSolMgr: \"Block Size\" must be strictly positive.");
578 params_->set(
"Block Size", blockSize_);
582 if (params->isParameter(
"Adaptive Block Size")) {
583 adaptiveBlockSize_ = params->get(
"Adaptive Block Size",adaptiveBlockSize_default_);
586 params_->set(
"Adaptive Block Size", adaptiveBlockSize_);
590 if (params->isParameter(
"Num Blocks")) {
591 numBlocks_ = params->get(
"Num Blocks",numBlocks_default_);
592 TEUCHOS_TEST_FOR_EXCEPTION(numBlocks_ <= 0, std::invalid_argument,
593 "Belos::BlockGmresSolMgr: \"Num Blocks\" must be strictly positive.");
596 params_->set(
"Num Blocks", numBlocks_);
600 if (params->isParameter(
"Timer Label")) {
601 std::string tempLabel = params->get(
"Timer Label", label_default_);
604 if (tempLabel != label_) {
606 params_->set(
"Timer Label", label_);
607 std::string solveLabel = label_ +
": BlockGmresSolMgr total solve time";
608 #ifdef BELOS_TEUCHOS_TIME_MONITOR 609 timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
611 if (ortho_ != Teuchos::null) {
612 ortho_->setLabel( label_ );
618 if (params->isParameter(
"Flexible Gmres")) {
619 isFlexible_ = Teuchos::getParameter<bool>(*params,
"Flexible Gmres");
620 params_->set(
"Flexible Gmres", isFlexible_);
621 if (isFlexible_ && expResTest_) {
629 if (params->isParameter(
"Orthogonalization")) {
630 std::string tempOrthoType = params->get(
"Orthogonalization",orthoType_default_);
631 TEUCHOS_TEST_FOR_EXCEPTION( tempOrthoType !=
"DGKS" && tempOrthoType !=
"ICGS" && tempOrthoType !=
"IMGS",
632 std::invalid_argument,
633 "Belos::BlockGmresSolMgr: \"Orthogonalization\" must be either \"DGKS\", \"ICGS\", or \"IMGS\".");
634 if (tempOrthoType != orthoType_) {
635 orthoType_ = tempOrthoType;
637 if (orthoType_==
"DGKS") {
638 if (orthoKappa_ <= 0) {
646 else if (orthoType_==
"ICGS") {
649 else if (orthoType_==
"IMGS") {
656 if (params->isParameter(
"Orthogonalization Constant")) {
657 orthoKappa_ = params->get(
"Orthogonalization Constant",orthoKappa_default_);
660 params_->set(
"Orthogonalization Constant",orthoKappa_);
661 if (orthoType_==
"DGKS") {
662 if (orthoKappa_ > 0 && ortho_ != Teuchos::null) {
669 if (params->isParameter(
"Verbosity")) {
670 if (Teuchos::isParameterType<int>(*params,
"Verbosity")) {
671 verbosity_ = params->get(
"Verbosity", verbosity_default_);
673 verbosity_ = (int)Teuchos::getParameter<Belos::MsgType>(*params,
"Verbosity");
677 params_->set(
"Verbosity", verbosity_);
678 if (printer_ != Teuchos::null)
679 printer_->setVerbosity(verbosity_);
683 if (params->isParameter(
"Output Style")) {
684 if (Teuchos::isParameterType<int>(*params,
"Output Style")) {
685 outputStyle_ = params->get(
"Output Style", outputStyle_default_);
687 outputStyle_ = (int)Teuchos::getParameter<Belos::OutputType>(*params,
"Output Style");
691 params_->set(
"Output Style", outputStyle_);
692 if (outputTest_ != Teuchos::null) {
698 if (params->isParameter(
"Output Stream")) {
699 outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> >(*params,
"Output Stream");
702 params_->set(
"Output Stream", outputStream_);
703 if (printer_ != Teuchos::null)
704 printer_->setOStream( outputStream_ );
709 if (params->isParameter(
"Output Frequency")) {
710 outputFreq_ = params->get(
"Output Frequency", outputFreq_default_);
714 params_->set(
"Output Frequency", outputFreq_);
715 if (outputTest_ != Teuchos::null)
716 outputTest_->setOutputFrequency( outputFreq_ );
720 if (printer_ == Teuchos::null) {
725 if (params->isParameter(
"Convergence Tolerance")) {
726 convtol_ = params->get(
"Convergence Tolerance",convtol_default_);
729 params_->set(
"Convergence Tolerance", convtol_);
730 if (impConvTest_ != Teuchos::null)
731 impConvTest_->setTolerance( convtol_ );
732 if (expConvTest_ != Teuchos::null)
733 expConvTest_->setTolerance( convtol_ );
737 if (params->isParameter(
"Implicit Residual Scaling")) {
738 std::string tempImpResScale = Teuchos::getParameter<std::string>( *params,
"Implicit Residual Scaling" );
741 if (impResScale_ != tempImpResScale) {
743 impResScale_ = tempImpResScale;
746 params_->set(
"Implicit Residual Scaling", impResScale_);
747 if (impConvTest_ != Teuchos::null) {
751 catch (std::exception& e) {
759 if (params->isParameter(
"Explicit Residual Scaling")) {
760 std::string tempExpResScale = Teuchos::getParameter<std::string>( *params,
"Explicit Residual Scaling" );
763 if (expResScale_ != tempExpResScale) {
765 expResScale_ = tempExpResScale;
768 params_->set(
"Explicit Residual Scaling", expResScale_);
769 if (expConvTest_ != Teuchos::null) {
773 catch (std::exception& e) {
781 if (params->isParameter(
"Explicit Residual Test")) {
782 expResTest_ = Teuchos::getParameter<bool>( *params,
"Explicit Residual Test" );
785 params_->set(
"Explicit Residual Test", expResTest_);
786 if (expConvTest_ == Teuchos::null) {
792 if (params->isParameter(
"Show Maximum Residual Norm Only")) {
793 showMaxResNormOnly_ = Teuchos::getParameter<bool>(*params,
"Show Maximum Residual Norm Only");
796 params_->set(
"Show Maximum Residual Norm Only", showMaxResNormOnly_);
797 if (impConvTest_ != Teuchos::null)
798 impConvTest_->setShowMaxResNormOnly( showMaxResNormOnly_ );
799 if (expConvTest_ != Teuchos::null)
800 expConvTest_->setShowMaxResNormOnly( showMaxResNormOnly_ );
804 if (ortho_ == Teuchos::null) {
805 if (orthoType_==
"DGKS") {
806 if (orthoKappa_ <= 0) {
814 else if (orthoType_==
"ICGS") {
817 else if (orthoType_==
"IMGS") {
821 TEUCHOS_TEST_FOR_EXCEPTION(orthoType_!=
"ICGS"&&orthoType_!=
"DGKS"&&orthoType_!=
"IMGS",std::logic_error,
822 "Belos::BlockGmresSolMgr(): Invalid orthogonalization type.");
827 if (timerSolve_ == Teuchos::null) {
828 std::string solveLabel = label_ +
": BlockGmresSolMgr total solve time";
829 #ifdef BELOS_TEUCHOS_TIME_MONITOR 830 timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
839 template<
class ScalarType,
class MV,
class OP>
851 if (!Teuchos::is_null(problem_->getLeftPrec()) && !isFlexible_) {
858 Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
859 Teuchos::rcp(
new StatusTestGenResNorm_t( convtol_ ) );
861 tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
862 impConvTest_ = tmpImpConvTest;
865 Teuchos::RCP<StatusTestGenResNorm_t> tmpExpConvTest =
866 Teuchos::rcp(
new StatusTestGenResNorm_t( convtol_ ) );
867 tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit,
Belos::TwoNorm );
869 tmpExpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
870 expConvTest_ = tmpExpConvTest;
873 convTest_ = Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
879 Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
880 Teuchos::rcp(
new StatusTestGenResNorm_t( convtol_ ) );
882 tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
883 impConvTest_ = tmpImpConvTest;
888 Teuchos::RCP<StatusTestImpResNorm_t> tmpImpConvTest =
889 Teuchos::rcp(
new StatusTestImpResNorm_t( convtol_ ) );
891 tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
892 impConvTest_ = tmpImpConvTest;
896 expConvTest_ = impConvTest_;
897 convTest_ = impConvTest_;
901 sTest_ = Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
909 std::string solverDesc =
" Block Gmres ";
911 solverDesc =
"Flexible" + solverDesc;
912 outputTest_->setSolverDesc( solverDesc );
921 template<
class ScalarType,
class MV,
class OP>
928 setParameters(Teuchos::parameterList(*getValidParameters()));
931 Teuchos::BLAS<int,ScalarType> blas;
934 "Belos::BlockGmresSolMgr::solve(): Linear problem is not a valid object.");
937 "Belos::BlockGmresSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
941 "Belos::BlockGmresSolMgr::solve(): Linear problem does not have a preconditioner required for flexible GMRES, call setRightPrec().");
944 if (!isSTSet_ || (!expResTest_ && !Teuchos::is_null(problem_->getLeftPrec())) ) {
946 "Belos::BlockGmresSolMgr::solve(): Linear problem and requested status tests are incompatible.");
951 int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
952 int numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
954 std::vector<int> currIdx;
957 if ( adaptiveBlockSize_ ) {
958 blockSize_ = numCurrRHS;
959 currIdx.resize( numCurrRHS );
960 for (
int i=0; i<numCurrRHS; ++i)
961 { currIdx[i] = startPtr+i; }
965 currIdx.resize( blockSize_ );
966 for (
int i=0; i<numCurrRHS; ++i)
967 { currIdx[i] = startPtr+i; }
968 for (
int i=numCurrRHS; i<blockSize_; ++i)
973 problem_->setLSIndex( currIdx );
977 Teuchos::ParameterList plist;
978 plist.set(
"Block Size",blockSize_);
980 ptrdiff_t dim = MVT::GetGlobalLength( *(problem_->getRHS()) );
981 if (blockSize_*static_cast<ptrdiff_t>(numBlocks_) > dim) {
982 int tmpNumBlocks = 0;
984 tmpNumBlocks = dim / blockSize_;
986 tmpNumBlocks = ( dim - blockSize_) / blockSize_;
988 "Belos::BlockGmresSolMgr::solve(): Warning! Requested Krylov subspace dimension is larger than operator dimension!" 989 << std::endl <<
" The maximum number of blocks allowed for the Krylov subspace will be adjusted to " << tmpNumBlocks << std::endl;
990 plist.set(
"Num Blocks",tmpNumBlocks);
993 plist.set(
"Num Blocks",numBlocks_);
996 outputTest_->reset();
997 loaDetected_ =
false;
1000 bool isConverged =
true;
1005 Teuchos::RCP<GmresIteration<ScalarType,MV,OP> > block_gmres_iter;
1014 #ifdef BELOS_TEUCHOS_TIME_MONITOR 1015 Teuchos::TimeMonitor slvtimer(*timerSolve_);
1018 while ( numRHS2Solve > 0 ) {
1021 if (blockSize_*numBlocks_ > dim) {
1022 int tmpNumBlocks = 0;
1023 if (blockSize_ == 1)
1024 tmpNumBlocks = dim / blockSize_;
1026 tmpNumBlocks = ( dim - blockSize_) / blockSize_;
1027 block_gmres_iter->setSize( blockSize_, tmpNumBlocks );
1030 block_gmres_iter->setSize( blockSize_, numBlocks_ );
1033 block_gmres_iter->resetNumIters();
1036 outputTest_->resetNumCalls();
1039 Teuchos::RCP<MV> V_0;
1042 if (currIdx[blockSize_-1] == -1) {
1043 V_0 = MVT::Clone( *(problem_->getInitResVec()), blockSize_ );
1044 problem_->computeCurrResVec( &*V_0 );
1047 V_0 = MVT::CloneCopy( *(problem_->getInitResVec()), currIdx );
1052 if (currIdx[blockSize_-1] == -1) {
1053 V_0 = MVT::Clone( *(problem_->getInitPrecResVec()), blockSize_ );
1054 problem_->computeCurrPrecResVec( &*V_0 );
1057 V_0 = MVT::CloneCopy( *(problem_->getInitPrecResVec()), currIdx );
1062 Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > z_0 =
1063 Teuchos::rcp(
new Teuchos::SerialDenseMatrix<int,ScalarType>( blockSize_, blockSize_ ) );
1066 int rank = ortho_->normalize( *V_0, z_0 );
1068 "Belos::BlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors.");
1075 block_gmres_iter->initializeGmres(newstate);
1076 int numRestarts = 0;
1081 block_gmres_iter->iterate();
1088 if ( convTest_->getStatus() ==
Passed ) {
1089 if ( expConvTest_->getLOADetected() ) {
1091 loaDetected_ =
true;
1093 "Belos::BlockGmresSolMgr::solve(): Warning! Solver has experienced a loss of accuracy!" << std::endl;
1094 isConverged =
false;
1103 else if ( maxIterTest_->getStatus() ==
Passed ) {
1105 isConverged =
false;
1113 else if ( block_gmres_iter->getCurSubspaceDim() == block_gmres_iter->getMaxSubspaceDim() ) {
1115 if ( numRestarts >= maxRestarts_ ) {
1116 isConverged =
false;
1121 printer_->stream(
Debug) <<
" Performing restart number " << numRestarts <<
" of " << maxRestarts_ << std::endl << std::endl;
1124 Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1127 Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
1128 MVT::MvAddMv( 1.0, *curX, 1.0, *update, *curX );
1131 problem_->updateSolution( update,
true );
1138 V_0 = MVT::Clone( *(oldState.
V), blockSize_ );
1140 problem_->computeCurrResVec( &*V_0 );
1142 problem_->computeCurrPrecResVec( &*V_0 );
1145 z_0 = Teuchos::rcp(
new Teuchos::SerialDenseMatrix<int,ScalarType>( blockSize_, blockSize_ ) );
1148 rank = ortho_->normalize( *V_0, z_0 );
1150 "Belos::BlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors after restart.");
1156 block_gmres_iter->initializeGmres(newstate);
1168 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
1169 "Belos::BlockGmresSolMgr::solve(): Invalid return from BlockGmresIter::iterate().");
1174 if (blockSize_ != 1) {
1175 printer_->stream(
Errors) <<
"Error! Caught std::exception in BlockGmresIter::iterate() at iteration " 1176 << block_gmres_iter->getNumIters() << std::endl
1177 << e.what() << std::endl;
1178 if (convTest_->getStatus() !=
Passed)
1179 isConverged =
false;
1184 block_gmres_iter->updateLSQR( block_gmres_iter->getCurSubspaceDim() );
1187 sTest_->checkStatus( &*block_gmres_iter );
1188 if (convTest_->getStatus() !=
Passed)
1189 isConverged =
false;
1193 catch (
const std::exception &e) {
1194 printer_->stream(
Errors) <<
"Error! Caught std::exception in BlockGmresIter::iterate() at iteration " 1195 << block_gmres_iter->getNumIters() << std::endl
1196 << e.what() << std::endl;
1205 Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1206 Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
1208 if (update != Teuchos::null)
1209 MVT::MvAddMv( 1.0, *curX, 1.0, *update, *curX );
1213 if ( !Teuchos::is_null(expConvTest_->getSolution()) ) {
1214 Teuchos::RCP<MV> newX = expConvTest_->getSolution();
1215 Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
1216 MVT::MvAddMv( 0.0, *newX, 1.0, *newX, *curX );
1219 Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1220 problem_->updateSolution( update,
true );
1225 problem_->setCurrLS();
1228 startPtr += numCurrRHS;
1229 numRHS2Solve -= numCurrRHS;
1230 if ( numRHS2Solve > 0 ) {
1231 numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
1233 if ( adaptiveBlockSize_ ) {
1234 blockSize_ = numCurrRHS;
1235 currIdx.resize( numCurrRHS );
1236 for (
int i=0; i<numCurrRHS; ++i)
1237 { currIdx[i] = startPtr+i; }
1240 currIdx.resize( blockSize_ );
1241 for (
int i=0; i<numCurrRHS; ++i)
1242 { currIdx[i] = startPtr+i; }
1243 for (
int i=numCurrRHS; i<blockSize_; ++i)
1244 { currIdx[i] = -1; }
1247 problem_->setLSIndex( currIdx );
1250 currIdx.resize( numRHS2Solve );
1261 #ifdef BELOS_TEUCHOS_TIME_MONITOR 1266 Teuchos::TimeMonitor::summarize( printer_->stream(
TimingDetails) );
1270 numIters_ = maxIterTest_->getNumIters();
1284 const std::vector<MagnitudeType>* pTestValues = NULL;
1286 pTestValues = expConvTest_->getTestValue();
1287 if (pTestValues == NULL || pTestValues->size() < 1) {
1288 pTestValues = impConvTest_->getTestValue();
1293 pTestValues = impConvTest_->getTestValue();
1295 TEUCHOS_TEST_FOR_EXCEPTION(pTestValues == NULL, std::logic_error,
1296 "Belos::BlockGmresSolMgr::solve(): The implicit convergence test's " 1297 "getTestValue() method returned NULL. Please report this bug to the " 1298 "Belos developers.");
1299 TEUCHOS_TEST_FOR_EXCEPTION(pTestValues->size() < 1, std::logic_error,
1300 "Belos::BlockGmresSolMgr::solve(): The implicit convergence test's " 1301 "getTestValue() method returned a vector of length zero. Please report " 1302 "this bug to the Belos developers.");
1307 achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
1310 if (!isConverged || loaDetected_) {
1317 template<
class ScalarType,
class MV,
class OP>
1320 std::ostringstream out;
1321 out <<
"\"Belos::BlockGmresSolMgr\": {";
1322 if (this->getObjectLabel () !=
"") {
1323 out <<
"Label: " << this->getObjectLabel () <<
", ";
1325 out <<
"Flexible: " << (isFlexible_ ?
"true" :
"false")
1326 <<
", Num Blocks: " << numBlocks_
1327 <<
", Maximum Iterations: " << maxIters_
1328 <<
", Maximum Restarts: " << maxRestarts_
1329 <<
", Convergence Tolerance: " << convtol_
1335 template<
class ScalarType,
class MV,
class OP>
1339 const Teuchos::EVerbosityLevel verbLevel)
const 1341 using Teuchos::TypeNameTraits;
1342 using Teuchos::VERB_DEFAULT;
1343 using Teuchos::VERB_NONE;
1344 using Teuchos::VERB_LOW;
1351 const Teuchos::EVerbosityLevel vl =
1352 (verbLevel == VERB_DEFAULT) ? VERB_LOW : verbLevel;
1354 if (vl != VERB_NONE) {
1355 Teuchos::OSTab tab0 (out);
1357 out <<
"\"Belos::BlockGmresSolMgr\":" << endl;
1358 Teuchos::OSTab tab1 (out);
1359 out <<
"Template parameters:" << endl;
1361 Teuchos::OSTab tab2 (out);
1362 out <<
"ScalarType: " << TypeNameTraits<ScalarType>::name () << endl
1363 <<
"MV: " << TypeNameTraits<MV>::name () << endl
1364 <<
"OP: " << TypeNameTraits<OP>::name () << endl;
1366 if (this->getObjectLabel () !=
"") {
1367 out <<
"Label: " << this->getObjectLabel () << endl;
1369 out <<
"Flexible: " << (isFlexible_ ?
"true" :
"false") << endl
1370 <<
"Num Blocks: " << numBlocks_ << endl
1371 <<
"Maximum Iterations: " << maxIters_ << endl
1372 <<
"Maximum Restarts: " << maxRestarts_ << endl
1373 <<
"Convergence Tolerance: " << convtol_ << endl;
ScaleType convertStringToScaleType(const std::string &scaleType)
Convert the given string to its ScaleType enum value.
Teuchos::RCP< OutputManager< ScalarType > > printer_
Collection of types and exceptions used within the Belos solvers.
Teuchos::ScalarTraits< MagnitudeType > MT
static const std::string impResScale_default_
Teuchos::ScalarTraits< ScalarType >::magnitudeType MagnitudeType
Belos's basic output manager for sending information of select verbosity levels to the appropriate ou...
Class which manages the output and verbosity of the Belos solvers.
static const bool flexibleGmres_default_
Teuchos::Array< Teuchos::RCP< Teuchos::Time > > getTimers() const
Return the timers for this object.
void reset(const ResetType type)
Performs a reset of the solver manager specified by the ResetType. This informs the solver manager th...
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > convTest_
Teuchos::RCP< StatusTestOutput< ScalarType, MV, OP > > outputTest_
static const std::string expResScale_default_
static const int maxRestarts_default_
ScaleType
The type of scaling to use on the residual norm value.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Get a parameter list containing the valid parameters for this object.
Teuchos::RCP< Teuchos::Time > timerSolve_
This class implements the block GMRES iteration, where a block Krylov subspace is constructed...
Convergence test using the implicit residual norm(s), with an explicit residual norm(s) check for los...
Teuchos::RCP< const MV > V
The current Krylov basis.
Virtual base class for StatusTest that printing status tests.
Teuchos::RCP< MatOrthoManager< ScalarType, MV, OP > > ortho_
static const int blockSize_default_
static const bool expResTest_default_
A factory class for generating StatusTestOutput objects.
static const bool adaptiveBlockSize_default_
Teuchos::RCP< const Teuchos::SerialDenseMatrix< int, ScalarType > > z
The current right-hand side of the least squares system RY = Z.
An implementation of StatusTestResNorm using a family of residual norms.
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const
Get a parameter list containing the current parameters for this object.
BlockGmresSolMgrLinearProblemFailure is thrown when the linear problem is not setup (i...
Structure to contain pointers to GmresIteration state variables.
Belos::StatusTest class for specifying a maximum number of iterations.
Teuchos::ScalarTraits< ScalarType > SCT
Interface to Block GMRES and Flexible GMRES.
virtual ~BlockGmresSolMgr()
Destructor.
A factory class for generating StatusTestOutput objects.
Iterated Modified Gram-Schmidt (IMGS) implementation of the Belos::OrthoManager class.
static const std::string orthoType_default_
Traits class which defines basic operations on multivectors.
Belos::StatusTest for logically combining several status tests.
static const int verbosity_default_
Classical Gram-Schmidt (with DGKS correction) implementation of the Belos::OrthoManager class...
A Belos::StatusTest class for specifying a maximum number of iterations.
Teuchos::RCP< StatusTestMaxIters< ScalarType, MV, OP > > maxIterTest_
ResetType
How to reset the solver.
static const std::string label_default_
Pure virtual base class which describes the basic interface for a solver manager. ...
int curDim
The current dimension of the reduction.
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.
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > sTest_
BlockGmresSolMgrLinearProblemFailure(const std::string &what_arg)
static const int outputStyle_default_
Teuchos::RCP< std::ostream > outputStream_
A linear system to solve, and its associated information.
Class which describes the linear problem to be solved by the iterative solver.
This class implements the block flexible GMRES iteration, where a block Krylov subspace is constructe...
std::string description() const
Return a one-line description of this object.
static const int outputFreq_default_
bool isLOADetected() const
Return whether a loss of accuracy was detected by this solver during the most current solve...
static const int numBlocks_default_
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 ...
ReturnType solve()
This method performs possibly repeated calls to the underlying linear solver's iterate() routine unti...
Belos concrete class for performing the block GMRES iteration.
void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem)
Set the linear problem that needs to be solved.
Iterated Classical Gram-Schmidt (ICGS) implementation of the Belos::OrthoManager class.
An implementation of the Belos::MatOrthoManager that performs orthogonalization using multiple steps ...
static const MagnitudeType orthoKappa_default_
An implementation of the Belos::MatOrthoManager that performs orthogonalization using multiple steps ...
static const Teuchos::RCP< std::ostream > outputStream_default_
MagnitudeType achievedTol() const
Tolerance achieved by the last solve() invocation.
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.
static const MagnitudeType convtol_default_
Belos::StatusTestResNorm for specifying general residual norm stopping criteria.
Belos::StatusTest for specifying an implicit residual norm stopping criteria that checks for loss of ...
BlockGmresSolMgr()
Empty constructor for BlockGmresSolMgr. This constructor takes no arguments and sets the default valu...
int getNumIters() const
Get the iteration count for the most recent call to solve().
An implementation of the Belos::MatOrthoManager that performs orthogonalization using (potentially) m...
BlockGmresSolMgrOrthoFailure(const std::string &what_arg)
A class for extending the status testing capabilities of Belos via logical combinations.
static const bool showMaxResNormOnly_default_
Teuchos::RCP< StatusTestResNorm< ScalarType, MV, OP > > expConvTest_
Class which defines basic traits for the operator type.
Parent class to all Belos exceptions.
BlockGmresSolMgrOrthoFailure is thrown when the orthogonalization manager is unable to generate ortho...
Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > problem_
static const int maxIters_default_
MagnitudeType achievedTol_
Pure virtual base class which augments the basic interface for a Gmres linear solver iteration...
GmresIterationOrthoFailure is thrown when the GmresIteration object is unable to compute independent ...
OperatorTraits< ScalarType, MV, OP > OPT
Belos header file which uses auto-configuration information to include necessary C++ headers...
Teuchos::RCP< Teuchos::ParameterList > params_
Teuchos::RCP< StatusTestResNorm< ScalarType, MV, OP > > impConvTest_
MagnitudeType orthoKappa_
const LinearProblem< ScalarType, MV, OP > & getProblem() const
Get current linear problem being solved for in this object.
Belos concrete class for performing the block, flexible GMRES iteration.
MultiVecTraits< ScalarType, MV > MVT
void setParameters(const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Set the parameters the solver manager should use to solve the linear problem.