42 #ifndef BELOS_MINRES_SOLMGR_HPP
43 #define BELOS_MINRES_SOLMGR_HPP
62 #ifdef BELOS_TEUCHOS_TIME_MONITOR
66 #include "Teuchos_StandardParameterEntryValidators.hpp"
114 template<
class ScalarType,
class MV,
class OP>
205 if (defaultParams_.
is_null()) {
208 return defaultParams_;
226 return Teuchos::tuple (timerSolve_);
273 problem_->setProblem ();
360 MagnitudeType convtol_;
363 MagnitudeType achievedTol_;
402 template<
class ScalarType,
class MV,
class OP>
407 using Teuchos::parameterList;
410 using Teuchos::rcpFromRef;
412 typedef MagnitudeType
MT;
416 RCP<ParameterList> pl = parameterList (
"MINRES");
418 pl->set (
"Convergence Tolerance", MST::squareroot (MST::eps()),
419 "Relative residual tolerance that needs to be achieved by "
420 "the iterative solver, in order for the linear system to be "
421 "declared converged.",
422 rcp (
new EnhancedNumberValidator<MT> (MST::zero(), MST::rmax())));
423 pl->set (
"Maximum Iterations", static_cast<int>(1000),
424 "Maximum number of iterations allowed for each right-hand "
426 rcp (
new EnhancedNumberValidator<int> (0, INT_MAX)));
427 pl->set (
"Num Blocks", static_cast<int> (-1),
428 "Ignored, but permitted, for compatibility with other Belos "
430 pl->set (
"Block Size", static_cast<int> (1),
431 "Number of vectors in each block. WARNING: The current "
432 "implementation of MINRES only accepts a block size of 1, "
433 "since it can only solve for 1 right-hand side at a time.",
434 rcp (
new EnhancedNumberValidator<int> (1, 1)));
436 "The type(s) of solver information that should "
437 "be written to the output stream.");
439 "What style is used for the solver information written "
440 "to the output stream.");
441 pl->set (
"Output Frequency", static_cast<int>(-1),
442 "How often (in terms of number of iterations) intermediate "
443 "convergence information should be written to the output stream."
445 pl->set (
"Output Stream", rcpFromRef(std::cout),
446 "A reference-counted pointer to the output stream where all "
447 "solver output is sent. The output stream defaults to stdout.");
448 pl->set (
"Timer Label", std::string(
"Belos"),
449 "The string to use as a prefix for the timer labels.");
456 template<
class ScalarType,
class MV,
class OP>
466 parametersSet_ (false)
472 template<
class ScalarType,
class MV,
class OP>
478 parametersSet_ (false)
481 "MinresSolMgr: The version of the constructor "
482 "that takes a LinearProblem to solve was given a "
483 "null LinearProblem.");
487 template<
class ScalarType,
class MV,
class OP>
494 "MINRES requires that you have provided a nonnull LinearProblem to the "
495 "solver manager, before you call the solve() method.");
498 "MINRES requires a LinearProblem object with a non-null operator (the "
502 "MINRES requires a LinearProblem object with a non-null right-hand side.");
505 "MINRES requires that before you give it a LinearProblem to solve, you "
506 "must first call the linear problem's setProblem() method.");
509 template<
class ScalarType,
class MV,
class OP>
515 using Teuchos::parameterList;
518 using Teuchos::rcpFromRef;
525 if (params_.is_null()) {
526 params_ = parameterList (*getValidParameters());
528 RCP<ParameterList> pl = params;
535 blockSize_ = pl->get<
int> (
"Block Size");
536 verbosity_ = pl->get<
int> (
"Verbosity");
537 outputStyle_ = pl->get<
int> (
"Output Style");
538 outputFreq_ = pl->get<
int>(
"Output Frequency");
539 outputStream_ = pl->get<RCP<std::ostream> > (
"Output Stream");
540 convtol_ = pl->get<MagnitudeType> (
"Convergence Tolerance");
541 maxIters_ = pl->get<
int> (
"Maximum Iterations");
549 const string newLabel = pl->get<
string> (
"Timer Label");
551 if (newLabel != label_ || timerSolve_.is_null()) {
553 #ifdef BELOS_TEUCHOS_TIME_MONITOR
554 const string solveLabel = label_ +
": MinresSolMgr total solve time";
556 if (! timerSolve_.is_null()) {
558 timerSolve_ = Teuchos::null;
561 #endif // BELOS_TEUCHOS_TIME_MONITOR
566 bool recreatedPrinter =
false;
567 if (printer_.is_null()) {
569 recreatedPrinter =
true;
572 printer_->setVerbosity (verbosity_);
574 printer_->setOStream (outputStream_);
585 const bool allocatedConvergenceTests =
586 impConvTest_.is_null() || expConvTest_.is_null();
590 if (impConvTest_.is_null()) {
591 impConvTest_ =
rcp (
new res_norm_type (convtol_));
592 impConvTest_->defineResForm (res_norm_type::Implicit,
TwoNorm);
597 impConvTest_->setTolerance (convtol_);
602 if (expConvTest_.is_null()) {
603 expConvTest_ =
rcp (
new res_norm_type (convtol_));
604 expConvTest_->defineResForm (res_norm_type::Explicit,
TwoNorm);
609 expConvTest_->setTolerance (convtol_);
615 bool needToRecreateFullStatusTest = sTest_.is_null();
619 if (convTest_.is_null() || allocatedConvergenceTests) {
620 convTest_ =
rcp (
new combo_type (combo_type::SEQ, impConvTest_, expConvTest_));
621 needToRecreateFullStatusTest =
true;
628 if (maxIterTest_.is_null()) {
630 needToRecreateFullStatusTest =
true;
632 maxIterTest_->setMaxIters (maxIters_);
643 if (needToRecreateFullStatusTest) {
644 sTest_ =
rcp (
new combo_type (combo_type::OR, maxIterTest_, convTest_));
651 if (outputTest_.is_null() || needToRecreateFullStatusTest || recreatedPrinter) {
653 outputTest_ = stoFactory.
create (printer_, sTest_, outputFreq_,
656 outputTest_->setOutputFrequency (outputFreq_);
660 outputTest_->setSolverDesc (std::string (
" MINRES "));
663 parametersSet_ =
true;
665 if (verbosity_ &
Debug) {
668 std::ostream& dbg = printer_->stream (
Debug);
669 dbg <<
"MINRES parameters:" << endl << params_ << endl;
674 template<
class ScalarType,
class MV,
class OP>
679 using Teuchos::rcp_const_cast;
682 if (! parametersSet_) {
683 setParameters (params_);
685 std::ostream& dbg = printer_->stream (
Debug);
687 #ifdef BELOS_TEUCHOS_TIME_MONITOR
689 #endif // BELOS_TEUCHOS_TIME_MONITOR
692 validateProblem (problem_);
695 outputTest_->reset();
700 const int numRHS2Solve = MVT::GetNumberVecs (*(problem_->getRHS()));
705 RCP<iter_type> minres_iter =
706 rcp (
new iter_type (problem_, printer_, outputTest_, *params_));
712 std::vector<int> notConverged;
713 std::vector<int> currentIndices(1);
718 for (
int currentRHS = 0; currentRHS < numRHS2Solve; ++currentRHS) {
723 currentIndices[0] = currentRHS;
724 problem_->setLSIndex (currentIndices);
726 dbg <<
"-- Current right-hand side index being solved: "
727 << currentRHS << endl;
730 minres_iter->resetNumIters();
732 outputTest_->resetNumCalls();
738 newstate.
Y = MVT::CloneViewNonConst (*(rcp_const_cast<MV> (problem_->getInitResVec())), currentIndices);
739 minres_iter->initializeMinres (newstate);
745 minres_iter->iterate();
748 if (convTest_->getStatus() ==
Passed) {
749 dbg <<
"---- Converged after " << maxIterTest_->getNumIters()
750 <<
" iterations" << endl;
754 else if (maxIterTest_->getStatus() ==
Passed) {
755 dbg <<
"---- Did not converge after " << maxIterTest_->getNumIters()
756 <<
" iterations" << endl;
758 notConverged.push_back (currentRHS);
765 "Belos::MinresSolMgr::solve(): iterations neither converged, "
766 "nor reached the maximum number of iterations " << maxIters_
767 <<
". That means something went wrong.");
769 }
catch (
const std::exception &e) {
771 <<
"Error! Caught std::exception in MinresIter::iterate() at "
772 <<
"iteration " << minres_iter->getNumIters() << endl
781 problem_->setCurrLS();
785 numIters_ += maxIterTest_->getNumIters();
793 #ifdef BELOS_TEUCHOS_TIME_MONITOR
800 #endif // BELOS_TEUCHOS_TIME_MONITOR
812 const std::vector<MagnitudeType>* pTestValues = expConvTest_->getTestValue();
813 if (pTestValues == NULL || pTestValues->size() < 1) {
814 pTestValues = impConvTest_->getTestValue();
817 "Belos::MinresSolMgr::solve(): The implicit convergence test's getTestValue() "
818 "method returned NULL. Please report this bug to the Belos developers.");
820 "Belos::MinresSolMgr::solve(): The implicit convergence test's getTestValue() "
821 "method returned a vector of length zero. Please report this bug to the "
822 "Belos developers.");
827 achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
830 if (notConverged.size() > 0) {
838 template<
class ScalarType,
class MV,
class OP>
841 std::ostringstream oss;
842 oss <<
"Belos::MinresSolMgr< "