Belos Version of the Day
BelosMinresSolMgr.hpp
Go to the documentation of this file.
00001 //@HEADER
00002 // ************************************************************************
00003 //
00004 //                 Belos: Block Linear Solvers Package
00005 //                  Copyright 2004 Sandia Corporation
00006 //
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00008 // the U.S. Government retains certain rights in this software.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00038 //
00039 // ************************************************************************
00040 //@HEADER
00041 
00042 #ifndef BELOS_MINRES_SOLMGR_HPP
00043 #define BELOS_MINRES_SOLMGR_HPP
00044 
00047 
00048 #include "BelosConfigDefs.hpp"
00049 #include "BelosTypes.hpp"
00050 
00051 #include "BelosLinearProblem.hpp"
00052 #include "BelosSolverManager.hpp"
00053 
00054 #include "BelosMinresIter.hpp"
00055 #include "BelosStatusTestMaxIters.hpp"
00056 #include "BelosStatusTestGenResNorm.hpp"
00057 #include "BelosStatusTestCombo.hpp"
00058 #include "BelosStatusTestOutputFactory.hpp"
00059 #include "BelosOutputManager.hpp"
00060 #include "Teuchos_BLAS.hpp"
00061 #include "Teuchos_LAPACK.hpp"
00062 #ifdef BELOS_TEUCHOS_TIME_MONITOR
00063 #include "Teuchos_TimeMonitor.hpp"
00064 #endif
00065 
00066 #include "Teuchos_StandardParameterEntryValidators.hpp"
00067 // Teuchos::ScalarTraits<int> doesn't define rmax(), alas, so we get
00068 // INT_MAX from here.
00069 #include <climits>
00070 
00071 namespace Belos {
00072 
00074 
00075 
00085   //
00089   class MinresSolMgrLinearProblemFailure : public BelosError {
00090   public:
00091     MinresSolMgrLinearProblemFailure (const std::string& what_arg) :
00092       BelosError(what_arg)
00093     {}
00094   };
00095 
00114   template<class ScalarType, class MV, class OP>
00115   class MinresSolMgr : public SolverManager<ScalarType,MV,OP> {
00116 
00117   private:
00118     typedef MultiVecTraits<ScalarType,MV> MVT;
00119     typedef OperatorTraits<ScalarType,MV,OP> OPT;
00120     typedef Teuchos::ScalarTraits<ScalarType> SCT;
00121     typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
00122     typedef Teuchos::ScalarTraits< MagnitudeType > MT;
00123 
00124   public:
00125 
00137     static Teuchos::RCP<const Teuchos::ParameterList> defaultParameters();
00138 
00140 
00141 
00150     MinresSolMgr();
00151 
00183     MinresSolMgr (const Teuchos::RCP<LinearProblem< ScalarType, MV, OP> > &problem,
00184                   const Teuchos::RCP<Teuchos::ParameterList> &params);
00185 
00187     virtual ~MinresSolMgr() {};
00189 
00191 
00192 
00194     const LinearProblem<ScalarType,MV,OP>& getProblem() const {
00195       return *problem_;
00196     }
00197 
00199     Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const {
00200       if (defaultParams_.is_null()) {
00201         defaultParams_ = defaultParameters ();
00202       }
00203       return defaultParams_;
00204     }
00205 
00207     Teuchos::RCP<const Teuchos::ParameterList> getCurrentParameters() const {
00208       return params_;
00209     }
00210 
00220     Teuchos::Array<Teuchos::RCP<Teuchos::Time> > getTimers() const {
00221       return Teuchos::tuple (timerSolve_);
00222     }
00223 
00229     MagnitudeType achievedTol() const {
00230       return achievedTol_;
00231     }
00232 
00234     int getNumIters() const {
00235       return numIters_;
00236     }
00237 
00243     bool isLOADetected() const { return false; }
00244 
00246 
00248 
00249 
00250     void
00251     setProblem (const Teuchos::RCP<LinearProblem<ScalarType, MV, OP> > &problem)
00252     {
00253       problem_ = problem;
00254     }
00255 
00256     void
00257     setParameters (const Teuchos::RCP<Teuchos::ParameterList>& params);
00258 
00260 
00262 
00263 
00264     void
00265     reset (const ResetType type)
00266     {
00267       if ((type & Belos::Problem) && ! problem_.is_null()) {
00268         problem_->setProblem ();
00269       }
00270     }
00272 
00274 
00275 
00293     ReturnType solve();
00294 
00296 
00299 
00300     std::string description() const;
00301 
00303 
00304   private:
00306     Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > problem_;
00307 
00309     Teuchos::RCP<OutputManager<ScalarType> > printer_;
00310     Teuchos::RCP<std::ostream> outputStream_;
00311 
00318     Teuchos::RCP<StatusTest<ScalarType,MV,OP> > sTest_;
00319 
00323     Teuchos::RCP<StatusTestMaxIters<ScalarType,MV,OP> > maxIterTest_;
00324 
00328     Teuchos::RCP<StatusTest<ScalarType,MV,OP> > convTest_;
00329 
00333     Teuchos::RCP<StatusTestGenResNorm<ScalarType,MV,OP> > impConvTest_;
00334 
00338     Teuchos::RCP<StatusTestGenResNorm<ScalarType,MV,OP> > expConvTest_;
00339 
00344     Teuchos::RCP<StatusTestOutput<ScalarType,MV,OP> > outputTest_;
00345 
00349     mutable Teuchos::RCP<const Teuchos::ParameterList> defaultParams_;
00350 
00352     Teuchos::RCP<Teuchos::ParameterList> params_;
00353 
00355     MagnitudeType convtol_;
00356 
00358     MagnitudeType achievedTol_;
00359 
00361     int maxIters_;
00362 
00364     int numIters_;
00365 
00367     int blockSize_;
00368 
00370     int verbosity_;
00371 
00373     int outputStyle_;
00374 
00376     int outputFreq_;
00377 
00379     std::string label_;
00380 
00382     Teuchos::RCP<Teuchos::Time> timerSolve_;
00383 
00385     bool parametersSet_;
00386 
00392     static void
00393     validateProblem (const Teuchos::RCP<LinearProblem<ScalarType, MV, OP> >& problem);
00394   };
00395 
00396 
00397   template<class ScalarType, class MV, class OP>
00398   Teuchos::RCP<const Teuchos::ParameterList>
00399   MinresSolMgr<ScalarType, MV, OP>::defaultParameters()
00400   {
00401     using Teuchos::ParameterList;
00402     using Teuchos::parameterList;
00403     using Teuchos::RCP;
00404     using Teuchos::rcp;
00405     using Teuchos::rcpFromRef;
00406     using Teuchos::EnhancedNumberValidator;
00407     typedef MagnitudeType MT;
00408     typedef Teuchos::ScalarTraits<MT> MST;
00409 
00410     // List of parameters accepted by MINRES, and their default values.
00411     RCP<ParameterList> pl = parameterList ("MINRES");
00412 
00413     pl->set ("Convergence Tolerance", MST::squareroot (MST::eps()),
00414              "Relative residual tolerance that needs to be achieved by "
00415              "the iterative solver, in order for the linear system to be "
00416              "declared converged.",
00417              rcp (new EnhancedNumberValidator<MT> (MST::zero(), MST::rmax())));
00418     pl->set ("Maximum Iterations", static_cast<int>(1000),
00419              "Maximum number of iterations allowed for each right-hand "
00420              "side solved.",
00421              rcp (new EnhancedNumberValidator<int> (0, INT_MAX)));
00422     pl->set ("Num Blocks", static_cast<int> (-1),
00423              "Ignored, but permitted, for compatibility with other Belos "
00424              "solvers.");
00425     pl->set ("Block Size", static_cast<int> (1),
00426              "Number of vectors in each block.  WARNING: The current "
00427              "implementation of MINRES only accepts a block size of 1, "
00428              "since it can only solve for 1 right-hand side at a time.",
00429              rcp (new EnhancedNumberValidator<int> (1, 1)));
00430     pl->set ("Verbosity", (int) Belos::Errors,
00431              "The type(s) of solver information that should "
00432              "be written to the output stream.");
00433     pl->set ("Output Style", (int) Belos::General,
00434              "What style is used for the solver information written "
00435              "to the output stream.");
00436     pl->set ("Output Frequency", static_cast<int>(-1),
00437              "How often (in terms of number of iterations) intermediate "
00438              "convergence information should be written to the output stream."
00439              "  -1 means never.");
00440     pl->set ("Output Stream", rcpFromRef(std::cout),
00441              "A reference-counted pointer to the output stream where all "
00442              "solver output is sent.  The output stream defaults to stdout.");
00443     pl->set ("Timer Label", std::string("Belos"),
00444              "The string to use as a prefix for the timer labels.");
00445     return pl;
00446   }
00447 
00448   //
00449   // Empty Constructor
00450   //
00451   template<class ScalarType, class MV, class OP>
00452   MinresSolMgr<ScalarType,MV,OP>::MinresSolMgr () :
00453     numIters_ (0),
00454     parametersSet_ (false)
00455   {}
00456 
00457   //
00458   // Primary constructor (use this one)
00459   //
00460   template<class ScalarType, class MV, class OP>
00461   MinresSolMgr<ScalarType, MV, OP>::
00462   MinresSolMgr (const Teuchos::RCP<LinearProblem<ScalarType, MV, OP> > &problem,
00463                 const Teuchos::RCP<Teuchos::ParameterList>& params) :
00464     problem_ (problem),
00465     numIters_ (0),
00466     parametersSet_ (false)
00467   {
00468     TEUCHOS_TEST_FOR_EXCEPTION(problem_.is_null(), std::invalid_argument,
00469                                "MinresSolMgr: The version of the constructor "
00470                                "that takes a LinearProblem to solve was given a "
00471                                "null LinearProblem.");
00472     setParameters (params);
00473   }
00474 
00475   template<class ScalarType, class MV, class OP>
00476   void
00477   MinresSolMgr<ScalarType, MV, OP>::
00478   validateProblem (const Teuchos::RCP<LinearProblem<ScalarType, MV, OP> >& problem)
00479   {
00480     TEUCHOS_TEST_FOR_EXCEPTION(problem.is_null(),
00481       MinresSolMgrLinearProblemFailure,
00482       "MINRES requires that you have provided a nonnull LinearProblem to the "
00483       "solver manager, before you call the solve() method.");
00484     TEUCHOS_TEST_FOR_EXCEPTION(problem->getOperator().is_null(),
00485       MinresSolMgrLinearProblemFailure,
00486       "MINRES requires a LinearProblem object with a non-null operator (the "
00487       "matrix A).");
00488     TEUCHOS_TEST_FOR_EXCEPTION(problem->getRHS().is_null(),
00489       MinresSolMgrLinearProblemFailure,
00490       "MINRES requires a LinearProblem object with a non-null right-hand side.");
00491     TEUCHOS_TEST_FOR_EXCEPTION( ! problem->isProblemSet(),
00492       MinresSolMgrLinearProblemFailure,
00493       "MINRES requires that before you give it a LinearProblem to solve, you "
00494       "must first call the linear problem's setProblem() method.");
00495   }
00496 
00497   template<class ScalarType, class MV, class OP>
00498   void
00499   MinresSolMgr<ScalarType, MV, OP>::
00500   setParameters (const Teuchos::RCP<Teuchos::ParameterList>& params)
00501   {
00502     using Teuchos::ParameterList;
00503     using Teuchos::parameterList;
00504     using Teuchos::RCP;
00505     using Teuchos::rcp;
00506     using Teuchos::rcpFromRef;
00507     using Teuchos::null;
00508     using Teuchos::is_null;
00509     using std::string;
00510     using std::ostream;
00511     using std::endl;
00512 
00513     if (params_.is_null()) {
00514       params_ = parameterList (*getValidParameters());
00515     }
00516     RCP<ParameterList> pl = params;
00517     pl->validateParametersAndSetDefaults (*params_);
00518 
00519     //
00520     // Read parameters from the parameter list.  We have already
00521     // populated it with defaults.
00522     //
00523     blockSize_ = pl->get<int> ("Block Size");
00524     verbosity_ = pl->get<int> ("Verbosity");
00525     outputStyle_ = pl->get<int> ("Output Style");
00526     outputFreq_ = pl->get<int>("Output Frequency");
00527     outputStream_ = pl->get<RCP<std::ostream> > ("Output Stream");
00528     convtol_ = pl->get<MagnitudeType> ("Convergence Tolerance");
00529     maxIters_ = pl->get<int> ("Maximum Iterations");
00530     //
00531     // All done reading parameters from the parameter list.
00532     // Now we know it's valid and we can store it.
00533     //
00534     params_ = pl;
00535 
00536     // Change the timer label, and create the timer if necessary.
00537     const string newLabel = pl->get<string> ("Timer Label");
00538     {
00539       if (newLabel != label_ || timerSolve_.is_null()) {
00540         label_ = newLabel;
00541 #ifdef BELOS_TEUCHOS_TIME_MONITOR
00542         const string solveLabel = label_ + ": MinresSolMgr total solve time";
00543         // Unregister the old timer before creating a new one.
00544         if (! timerSolve_.is_null()) {
00545           Teuchos::TimeMonitor::clearCounter (label_);
00546           timerSolve_ = Teuchos::null;
00547         }
00548         timerSolve_ = Teuchos::TimeMonitor::getNewCounter (solveLabel);
00549 #endif // BELOS_TEUCHOS_TIME_MONITOR
00550       }
00551     }
00552 
00553     // Create output manager, if necessary; otherwise, set its parameters.
00554     bool recreatedPrinter = false;
00555     if (printer_.is_null()) {
00556       printer_ = rcp (new OutputManager<ScalarType> (verbosity_, outputStream_));
00557       recreatedPrinter = true;
00558     } else {
00559       // Set the output stream's verbosity level.
00560       printer_->setVerbosity (verbosity_);
00561       // Tell the output manager about the new output stream.
00562       printer_->setOStream (outputStream_);
00563     }
00564 
00565     //
00566     // Set up the convergence tests
00567     //
00568     typedef StatusTestGenResNorm<ScalarType, MV, OP> res_norm_type;
00569     typedef StatusTestCombo<ScalarType, MV, OP> combo_type;
00570 
00571     // Do we need to allocate at least one of the implicit or explicit
00572     // residual norm convergence tests?
00573     const bool allocatedConvergenceTests =
00574       impConvTest_.is_null() || expConvTest_.is_null();
00575 
00576     // Allocate or set the tolerance of the implicit residual norm
00577     // convergence test.
00578     if (impConvTest_.is_null()) {
00579       impConvTest_ = rcp (new res_norm_type (convtol_));
00580       impConvTest_->defineResForm (res_norm_type::Implicit, TwoNorm);
00581       // TODO (mfh 03 Nov 2011) Allow users to define the type of
00582       // scaling (or a custom scaling factor).
00583       impConvTest_->defineScaleForm (NormOfInitRes, TwoNorm);
00584     } else {
00585       impConvTest_->setTolerance (convtol_);
00586     }
00587 
00588     // Allocate or set the tolerance of the explicit residual norm
00589     // convergence test.
00590     if (expConvTest_.is_null()) {
00591       expConvTest_ = rcp (new res_norm_type (convtol_));
00592       expConvTest_->defineResForm (res_norm_type::Explicit, TwoNorm);
00593       // TODO (mfh 03 Nov 2011) Allow users to define the type of
00594       // scaling (or a custom scaling factor).
00595       expConvTest_->defineScaleForm (NormOfInitRes, TwoNorm);
00596     } else {
00597       expConvTest_->setTolerance (convtol_);
00598     }
00599 
00600     // Whether we need to recreate the full status test.  We only need
00601     // to do that if at least one of convTest_ or maxIterTest_ had to
00602     // be reallocated.
00603     bool needToRecreateFullStatusTest = sTest_.is_null();
00604 
00605     // Residual status test is a combo of the implicit and explicit
00606     // convergence tests.
00607     if (convTest_.is_null() || allocatedConvergenceTests) {
00608       convTest_ = rcp (new combo_type (combo_type::SEQ, impConvTest_, expConvTest_));
00609       needToRecreateFullStatusTest = true;
00610     }
00611 
00612     // Maximum number of iterations status test.  It tells the solver to
00613     // stop iteration, if the maximum number of iterations has been
00614     // exceeded.  Initialize it if we haven't yet done so, otherwise
00615     // tell it the new maximum number of iterations.
00616     if (maxIterTest_.is_null()) {
00617       maxIterTest_ = rcp (new StatusTestMaxIters<ScalarType,MV,OP> (maxIters_));
00618       needToRecreateFullStatusTest = true;
00619     } else {
00620       maxIterTest_->setMaxIters (maxIters_);
00621     }
00622 
00623     // Create the full status test if we need to.
00624     //
00625     // The full status test: the maximum number of iterations have
00626     // been reached, OR the residual has converged.
00627     //
00628     // "If we need to" means either that the status test was never
00629     // created before, or that its two component tests had to be
00630     // reallocated.
00631     if (needToRecreateFullStatusTest) {
00632       sTest_ = rcp (new combo_type (combo_type::OR, maxIterTest_, convTest_));
00633     }
00634 
00635     // If necessary, create the status test output class.  This class
00636     // manages and formats the output from the status test.  We have
00637     // to recreate the output test if we had to (re)allocate either
00638     // printer_ or sTest_.
00639     if (outputTest_.is_null() || needToRecreateFullStatusTest || recreatedPrinter) {
00640       StatusTestOutputFactory<ScalarType,MV,OP> stoFactory (outputStyle_);
00641       outputTest_ = stoFactory.create (printer_, sTest_, outputFreq_,
00642                                        Passed+Failed+Undefined);
00643     } else {
00644       outputTest_->setOutputFrequency (outputFreq_);
00645     }
00646     // Set the solver string for the output test.
00647     // StatusTestOutputFactory has no constructor argument for this.
00648     outputTest_->setSolverDesc (std::string (" MINRES "));
00649 
00650     // Inform the solver manager that the current parameters were set.
00651     parametersSet_ = true;
00652 
00653     if (verbosity_ & Debug) {
00654       using std::endl;
00655 
00656       std::ostream& dbg = printer_->stream (Debug);
00657       dbg << "MINRES parameters:" << endl << params_ << endl;
00658     }
00659   }
00660 
00661 
00662   template<class ScalarType, class MV, class OP>
00663   ReturnType MinresSolMgr<ScalarType,MV,OP>::solve()
00664   {
00665     using Teuchos::RCP;
00666     using Teuchos::rcp;
00667     using Teuchos::rcp_const_cast;
00668     using std::endl;
00669 
00670     if (! parametersSet_) {
00671       setParameters (params_);
00672     }
00673     std::ostream& dbg = printer_->stream (Debug);
00674 
00675 #ifdef BELOS_TEUCHOS_TIME_MONITOR
00676     Teuchos::TimeMonitor solveTimerMonitor (*timerSolve_);
00677 #endif // BELOS_TEUCHOS_TIME_MONITOR
00678 
00679     // We need a problem to solve, else we can't solve it.
00680     validateProblem (problem_);
00681 
00682     // Reset the status test for this solve.
00683     outputTest_->reset();
00684 
00685     // The linear problem has this many right-hand sides to solve.
00686     // MINRES can solve only one at a time, so we solve for each
00687     // right-hand side in succession.
00688     const int numRHS2Solve = MVT::GetNumberVecs (*(problem_->getRHS()));
00689 
00690     // Create MINRES iteration object.  Pass along the solver
00691     // manager's parameters, which have already been validated.
00692     typedef MinresIter<ScalarType, MV, OP> iter_type;
00693     RCP<iter_type> minres_iter =
00694       rcp (new iter_type (problem_, printer_, outputTest_, *params_));
00695 
00696     // The index/indices of the right-hand sides for which MINRES did
00697     // _not_ converge.  Hopefully this is empty after the for loop
00698     // below!  If it is not empty, at least one right-hand side did
00699     // not converge.
00700     std::vector<int> notConverged;
00701     std::vector<int> currentIndices(1);
00702 
00703     numIters_ = 0;
00704 
00705     // Solve for each right-hand side in turn.
00706     for (int currentRHS = 0; currentRHS < numRHS2Solve; ++currentRHS) {
00707       // Inform the linear problem of the right-hand side(s) currently
00708       // being solved.  MINRES only knows how to solve linear problems
00709       // with one right-hand side, so we only include one index, which
00710       // is the index of the current right-hand side.
00711       currentIndices[0] = currentRHS;
00712       problem_->setLSIndex (currentIndices);
00713 
00714       dbg << "-- Current right-hand side index being solved: "
00715           << currentRHS << endl;
00716 
00717       // Reset the number of iterations.
00718       minres_iter->resetNumIters();
00719       // Reset the number of calls that the status test output knows about.
00720       outputTest_->resetNumCalls();
00721       // Set the new state and initialize the solver.
00722       MinresIterationState<ScalarType, MV> newstate;
00723 
00724       // Get the residual vector for the current linear system
00725       // (that is, for the current right-hand side).
00726       newstate.Y = MVT::CloneViewNonConst (*(rcp_const_cast<MV> (problem_->getInitResVec())), currentIndices);
00727       minres_iter->initializeMinres (newstate);
00728 
00729       // Attempt to solve for the solution corresponding to the
00730       // current right-hand side.
00731       while (true) {
00732         try {
00733           minres_iter->iterate();
00734 
00735           // First check for convergence
00736           if (convTest_->getStatus() == Passed) {
00737             dbg << "---- Converged after " << maxIterTest_->getNumIters()
00738                 << " iterations" << endl;
00739             break;
00740           }
00741           // Now check for max # of iterations
00742           else if (maxIterTest_->getStatus() == Passed) {
00743             dbg << "---- Did not converge after " << maxIterTest_->getNumIters()
00744                 << " iterations" << endl;
00745             // This right-hand side didn't converge!
00746             notConverged.push_back (currentRHS);
00747             break;
00748           } else {
00749             // If we get here, we returned from iterate(), but none of
00750             // our status tests Passed.  Something is wrong, and it is
00751             // probably our fault.
00752             TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
00753               "Belos::MinresSolMgr::solve(): iterations neither converged, "
00754               "nor reached the maximum number of iterations " << maxIters_
00755               << ".  That means something went wrong.");
00756           }
00757         } catch (const std::exception &e) {
00758           printer_->stream (Errors)
00759             << "Error! Caught std::exception in MinresIter::iterate() at "
00760             << "iteration " << minres_iter->getNumIters() << endl
00761             << e.what() << endl;
00762           throw e;
00763         }
00764       }
00765 
00766       // Inform the linear problem that we are finished with the
00767       // current right-hand side.  It may or may not have converged,
00768       // but we don't try again if the first time didn't work.
00769       problem_->setCurrLS();
00770 
00771       // Get iteration information for this solve: total number of
00772       // iterations for all right-hand sides.
00773       numIters_ += maxIterTest_->getNumIters();
00774     }
00775 
00776     // Print final summary of the solution process
00777     sTest_->print (printer_->stream (FinalSummary));
00778 
00779     // Print timing information, if the corresponding compile-time and
00780     // run-time options are enabled.
00781 #ifdef BELOS_TEUCHOS_TIME_MONITOR
00782     // Calling summarize() can be expensive, so don't call unless the
00783     // user wants to print out timing details.  summarize() will do all
00784     // the work even if it's passed a "black hole" output stream.
00785     if (verbosity_ & TimingDetails) {
00786       Teuchos::TimeMonitor::summarize (printer_->stream (TimingDetails));
00787     }
00788 #endif // BELOS_TEUCHOS_TIME_MONITOR
00789 
00790     // Save the convergence test value ("achieved tolerance") for this
00791     // solve.  This solver always has two residual norm status tests:
00792     // an explicit and an implicit test.  The master convergence test
00793     // convTest_ is a SEQ combo of the implicit resp. explicit tests.
00794     // If the implicit test never passes, then the explicit test won't
00795     // ever be executed.  This manifests as
00796     // expConvTest_->getTestValue()->size() < 1.  We deal with this
00797     // case by using the values returned by
00798     // impConvTest_->getTestValue().
00799     {
00800       const std::vector<MagnitudeType>* pTestValues = expConvTest_->getTestValue();
00801       if (pTestValues == NULL || pTestValues->size() < 1) {
00802         pTestValues = impConvTest_->getTestValue();
00803       }
00804       TEUCHOS_TEST_FOR_EXCEPTION(pTestValues == NULL, std::logic_error,
00805         "Belos::MinresSolMgr::solve(): The implicit convergence test's getTestValue() "
00806         "method returned NULL.  Please report this bug to the Belos developers.");
00807       TEUCHOS_TEST_FOR_EXCEPTION(pTestValues->size() < 1, std::logic_error,
00808         "Belos::MinresSolMgr::solve(): The implicit convergence test's getTestValue() "
00809         "method returned a vector of length zero.  Please report this bug to the "
00810         "Belos developers.");
00811 
00812       // FIXME (mfh 12 Dec 2011) Does pTestValues really contain the
00813       // achieved tolerances for all vectors in the current solve(), or
00814       // just for the vectors from the last deflation?
00815       achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
00816     }
00817 
00818     if (notConverged.size() > 0) {
00819       return Unconverged;
00820     } else {
00821       return Converged;
00822     }
00823   }
00824 
00825   //  This method requires the solver manager to return a std::string that describes itself.
00826   template<class ScalarType, class MV, class OP>
00827   std::string MinresSolMgr<ScalarType,MV,OP>::description() const
00828   {
00829     std::ostringstream oss;
00830     oss << "Belos::MinresSolMgr< "
00831         << Teuchos::ScalarTraits<ScalarType>::name()
00832         <<", MV, OP >";
00833     // oss << "{";
00834     // oss << "Block Size=" << blockSize_;
00835     // oss << "}";
00836     return oss.str();
00837   }
00838 
00839 } // end Belos namespace
00840 
00841 #endif /* BELOS_MINRES_SOLMGR_HPP */
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines