Belos Version of the Day
BelosGmresInnerSolver.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_GmresInnerSolver_hpp
00043 #define __Belos_GmresInnerSolver_hpp
00044 
00045 #include <BelosInnerSolver.hpp>
00046 #include <BelosGmresSolMgr.hpp>
00047 #include <Teuchos_TypeNameTraits.hpp>
00048 
00049 namespace Belos {
00050 
00054   template<class Scalar, class MV, class OP>
00055   class GmresInnerSolver : public InnerSolver<Scalar, MV, OP> {
00056   public:
00057     typedef InnerSolver<Scalar, MV, OP> base_type;
00058     typedef typename base_type::scalar_type scalar_type;
00059     typedef typename base_type::magnitude_type magnitude_type;
00060     typedef typename base_type::multivector_type multivector_type;
00061     typedef typename base_type::operator_type operator_type;
00062 
00063   private:
00067     static Teuchos::RCP<LinearProblem<Scalar, MV, OP> >
00068     validProblem (const Teuchos::RCP<LinearProblem<Scalar, MV, OP> >& problem)
00069     {
00070       if (problem.is_null())
00071   throw std::invalid_argument("GmresInnerSolver constructor: The linear "
00072             "problem to solve is null.");
00073       else if (problem->getOperator().is_null())
00074   throw std::invalid_argument("GmresInnerSolver constructor: the matrix A"
00075             " in the linear system AX=B to solve is "
00076             "null.");
00077       else
00078   return problem;
00079     }
00080 
00081   public:
00102     GmresInnerSolver (const Teuchos::RCP<LinearProblem<Scalar,MV,OP> >& problem,
00103           const Teuchos::RCP<const Teuchos::ParameterList>& params,
00104           const bool debug = false)
00105       : solMgr_ (validProblem (problem), params, debug)
00106     {}
00107 
00114     Teuchos::RCP<const Teuchos::ParameterList> getCurrentParameters() const {
00115       return solMgr_.getCurrentParameters();
00116     }
00117 
00139     InnerSolveResult
00140     solve (const Teuchos::RCP<MV>& X,
00141      const Teuchos::RCP<const MV>& B,
00142      const magnitude_type convTol,
00143      const int maxItersPerRestart,
00144      const int maxNumRestarts)
00145     {
00146       solMgr_.setRHS (B);
00147       solMgr_.setLHS (X);
00148       // FIXME (mfh 07 Mar 2011) 
00149       //
00150       // I have a vague feeling that the stopping criteria should be
00151       // changed _after_ setting the RHS, but I'm not sure if this is
00152       // true.  (Does it have something to do with one of the stopping
00153       // criteria computing or storing its own residual vector(s)?)
00154       solMgr_.changeStoppingCriteria (convTol, maxItersPerRestart, 
00155               maxNumRestarts);
00156       return invokeSolver ();
00157     }
00158 
00176     InnerSolveResult
00177     solve (const Teuchos::RCP<MV>& X,
00178      const Teuchos::RCP<const MV>& B)
00179     {
00180       // The solver manager retains the last configuration of
00181       // convergence tolerance, max number of iterations per restart
00182       // cycle, and max number of restart cycles. Just set the
00183       // left-hand side X and right-hand side B, and call solve().
00184       solMgr_.setLHS (X);
00185       solMgr_.setRHS (B);
00186       return invokeSolver ();
00187     }
00188 
00189   private:
00195     GmresInnerSolver();
00196 
00202     GmresSolMgr<Scalar, MV, OP> solMgr_;
00203 
00211     InnerSolveResult
00212     invokeSolver ()
00213     {
00214       using std::pair;
00215       using std::vector;
00216 
00217       ReturnType result = solMgr_.solve ();
00218       TEST_FOR_EXCEPTION(result != Converged && result != Unconverged,
00219        std::logic_error,
00220        "The solver manager returned an invalid ResultType "
00221        "enum value " << result << ".  Valid values are "
00222        "Converged=" << Converged << " and Unconverged=" 
00223        << Unconverged << ".");
00224       //
00225       // Compute max (total number of iterations, number of restart cycles)
00226       // for all right-hand side(s).
00227       //
00228       const vector<pair<int, int> > iterInfo = solMgr_.totalNumIters();
00229       // If we could use lambdas, or if C++ had better iteration
00230       // primitives, we wouldn't need a for loop here.
00231       typedef vector<pair<int, int>> >::const_iterator iter_type;
00232       int totalNumIters = 0, numRestartCycles = 0;
00233       for (iter_type it = iterInfo.begin(); it != iterInfo.end(); ++it)
00234   {
00235     numRestartCycles = std::max (it->first, numRestartCycles);
00236     totalNumIters = std::max (it->second, totalNumIters);
00237   }
00238       return InnerSolveResult (result, numRestartCycles, totalNumIters);
00239     }
00240 
00241   };
00242 
00243 } // namespace Belos
00244 
00245 #endif // __Belos_GmresInnerSolver_hpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines