Belos Package Browser (Single Doxygen Collection) Development
BelosTFQMRSolMgr.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_TFQMR_SOLMGR_HPP
00043 #define BELOS_TFQMR_SOLMGR_HPP
00044 
00049 #include "BelosConfigDefs.hpp"
00050 #include "BelosTypes.hpp"
00051 
00052 #include "BelosLinearProblem.hpp"
00053 #include "BelosSolverManager.hpp"
00054 
00055 #include "BelosTFQMRIter.hpp"
00056 #include "BelosStatusTestMaxIters.hpp"
00057 #include "BelosStatusTestGenResNorm.hpp"
00058 #include "BelosStatusTestCombo.hpp"
00059 #include "BelosStatusTestOutputFactory.hpp"
00060 #include "BelosOutputManager.hpp"
00061 #include "Teuchos_BLAS.hpp"
00062 #include "Teuchos_LAPACK.hpp"
00063 #include "Teuchos_TimeMonitor.hpp"
00064 
00078 namespace Belos {
00079   
00081 
00082   
00089   class TFQMRSolMgrLinearProblemFailure : public BelosError {public:
00090     TFQMRSolMgrLinearProblemFailure(const std::string& what_arg) : BelosError(what_arg)
00091     {}};
00092   
00099   class TFQMRSolMgrOrthoFailure : public BelosError {public:
00100     TFQMRSolMgrOrthoFailure(const std::string& what_arg) : BelosError(what_arg)
00101     {}};
00102   
00103   template<class ScalarType, class MV, class OP>
00104   class TFQMRSolMgr : public SolverManager<ScalarType,MV,OP> {
00105     
00106   private:
00107     typedef MultiVecTraits<ScalarType,MV> MVT;
00108     typedef OperatorTraits<ScalarType,MV,OP> OPT;
00109     typedef Teuchos::ScalarTraits<ScalarType> SCT;
00110     typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
00111     typedef Teuchos::ScalarTraits<MagnitudeType> MT;
00112     
00113   public:
00114     
00116 
00117 
00123      TFQMRSolMgr();
00124 
00141     TFQMRSolMgr( const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem,
00142      const Teuchos::RCP<Teuchos::ParameterList> &pl );
00143     
00145     virtual ~TFQMRSolMgr() {};
00147     
00149 
00150     
00151     const LinearProblem<ScalarType,MV,OP>& getProblem() const {
00152       return *problem_;
00153     }
00154 
00157     Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const;
00158     
00161     Teuchos::RCP<const Teuchos::ParameterList> getCurrentParameters() const { return params_; }
00162     
00168     Teuchos::Array<Teuchos::RCP<Teuchos::Time> > getTimers() const {
00169       return tuple(timerSolve_);
00170     }
00171     
00173     int getNumIters() const {
00174       return numIters_;
00175     }
00176     
00179     bool isLOADetected() const { return false; }
00180     
00182     
00184 
00185     
00187     void setProblem( const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem ) { problem_ = problem; }
00188     
00190     void setParameters( const Teuchos::RCP<Teuchos::ParameterList> &params );
00191     
00193     
00195 
00196 
00200     void reset( const ResetType type ) { if ((type & Belos::Problem) && !Teuchos::is_null(problem_)) problem_->setProblem(); }
00202     
00204 
00205     
00223     ReturnType solve();
00224     
00226     
00229     
00231     std::string description() const;
00232     
00234     
00235   private:
00236 
00237     // Method to convert std::string to enumerated type for residual.
00238     Belos::ScaleType convertStringToScaleType( std::string& scaleType ) {
00239       if (scaleType == "Norm of Initial Residual") {
00240   return Belos::NormOfInitRes;
00241       } else if (scaleType == "Norm of Preconditioned Initial Residual") {
00242   return Belos::NormOfPrecInitRes;
00243       } else if (scaleType == "Norm of RHS") {
00244   return Belos::NormOfRHS;
00245       } else if (scaleType == "None") {
00246   return Belos::None;
00247       } else 
00248   TEST_FOR_EXCEPTION( true ,std::logic_error,
00249           "Belos::TFQMRSolMgr(): Invalid residual scaling type.");
00250     }
00251 
00252     // Method for checking current status test against defined linear problem.
00253     bool checkStatusTest();
00254     
00255     // Linear problem.
00256     Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > problem_;
00257     
00258     // Output manager.
00259     Teuchos::RCP<OutputManager<ScalarType> > printer_;
00260     Teuchos::RCP<std::ostream> outputStream_;
00261     
00262     // Status test.
00263     Teuchos::RCP<StatusTest<ScalarType,MV,OP> > sTest_;
00264     Teuchos::RCP<StatusTestMaxIters<ScalarType,MV,OP> > maxIterTest_;
00265     Teuchos::RCP<StatusTest<ScalarType,MV,OP> > convTest_;
00266     Teuchos::RCP<StatusTestGenResNorm<ScalarType,MV,OP> > expConvTest_, impConvTest_;
00267     Teuchos::RCP<StatusTestOutput<ScalarType,MV,OP> > outputTest_;
00268     
00269     // Current parameter list.
00270     Teuchos::RCP<ParameterList> params_;
00271     
00272     // Default solver values.
00273     static const MagnitudeType convtol_default_;
00274     static const int maxIters_default_;
00275     static const bool expResTest_default_;
00276     static const int verbosity_default_;
00277     static const int outputStyle_default_;
00278     static const int outputFreq_default_;
00279     static const std::string impResScale_default_; 
00280     static const std::string expResScale_default_; 
00281     static const std::string label_default_;
00282     static const Teuchos::RCP<std::ostream> outputStream_default_;
00283 
00284     // Current solver values.
00285     MagnitudeType convtol_;
00286     int maxIters_, numIters_;
00287     int verbosity_, outputStyle_, outputFreq_;
00288     int blockSize_;
00289     bool expResTest_;
00290     std::string impResScale_, expResScale_;
00291     
00292     // Timers.
00293     std::string label_;
00294     Teuchos::RCP<Teuchos::Time> timerSolve_;
00295 
00296     // Internal state variables.
00297     bool isSet_, isSTSet_;
00298   };
00299 
00300 
00301 // Default solver values.
00302 template<class ScalarType, class MV, class OP>
00303 const typename TFQMRSolMgr<ScalarType,MV,OP>::MagnitudeType TFQMRSolMgr<ScalarType,MV,OP>::convtol_default_ = 1e-8;
00304 
00305 template<class ScalarType, class MV, class OP>
00306 const int TFQMRSolMgr<ScalarType,MV,OP>::maxIters_default_ = 1000;
00307 
00308 template<class ScalarType, class MV, class OP>
00309 const bool TFQMRSolMgr<ScalarType,MV,OP>::expResTest_default_ = false;
00310 
00311 template<class ScalarType, class MV, class OP>
00312 const int TFQMRSolMgr<ScalarType,MV,OP>::verbosity_default_ = Belos::Errors;
00313 
00314 template<class ScalarType, class MV, class OP>
00315 const int TFQMRSolMgr<ScalarType,MV,OP>::outputStyle_default_ = Belos::General;
00316 
00317 template<class ScalarType, class MV, class OP>
00318 const int TFQMRSolMgr<ScalarType,MV,OP>::outputFreq_default_ = -1;
00319 
00320 template<class ScalarType, class MV, class OP>
00321 const std::string TFQMRSolMgr<ScalarType,MV,OP>::impResScale_default_ = "Norm of Preconditioned Initial Residual";
00322 
00323 template<class ScalarType, class MV, class OP>
00324 const std::string TFQMRSolMgr<ScalarType,MV,OP>::expResScale_default_ = "Norm of Initial Residual";
00325 
00326 template<class ScalarType, class MV, class OP>
00327 const std::string TFQMRSolMgr<ScalarType,MV,OP>::label_default_ = "Belos";
00328 
00329 template<class ScalarType, class MV, class OP>
00330 const Teuchos::RCP<std::ostream> TFQMRSolMgr<ScalarType,MV,OP>::outputStream_default_ = Teuchos::rcp(&std::cout,false);
00331 
00332 
00333 // Empty Constructor
00334 template<class ScalarType, class MV, class OP>
00335 TFQMRSolMgr<ScalarType,MV,OP>::TFQMRSolMgr() :
00336   outputStream_(outputStream_default_),
00337   convtol_(convtol_default_),
00338   maxIters_(maxIters_default_),
00339   verbosity_(verbosity_default_),
00340   outputStyle_(outputStyle_default_),
00341   outputFreq_(outputFreq_default_),
00342   blockSize_(1),
00343   expResTest_(expResTest_default_),
00344   impResScale_(impResScale_default_),
00345   expResScale_(expResScale_default_),
00346   label_(label_default_),
00347   isSet_(false),
00348   isSTSet_(false)
00349 {}
00350 
00351 
00352 // Basic Constructor
00353 template<class ScalarType, class MV, class OP>
00354 TFQMRSolMgr<ScalarType,MV,OP>::TFQMRSolMgr( 
00355                const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem,
00356                const Teuchos::RCP<Teuchos::ParameterList> &pl ) : 
00357   problem_(problem),
00358   outputStream_(outputStream_default_),
00359   convtol_(convtol_default_),
00360   maxIters_(maxIters_default_),
00361   verbosity_(verbosity_default_),
00362   outputStyle_(outputStyle_default_),
00363   outputFreq_(outputFreq_default_),
00364   blockSize_(1),
00365   expResTest_(expResTest_default_),
00366   impResScale_(impResScale_default_),
00367   expResScale_(expResScale_default_),
00368   label_(label_default_),
00369   isSet_(false),
00370   isSTSet_(false)
00371 {
00372   TEST_FOR_EXCEPTION(problem_ == Teuchos::null, std::invalid_argument, "Problem not given to solver manager.");
00373   
00374   // If the parameter list pointer is null, then set the current parameters to the default parameter list.
00375   if ( !is_null(pl) ) {
00376     setParameters( pl );  
00377   }
00378 }
00379   
00380 template<class ScalarType, class MV, class OP>
00381 void TFQMRSolMgr<ScalarType,MV,OP>::setParameters( const Teuchos::RCP<Teuchos::ParameterList> &params )
00382 {
00383   // Create the internal parameter list if ones doesn't already exist.
00384   if (params_ == Teuchos::null) {
00385     params_ = Teuchos::rcp( new Teuchos::ParameterList(*getValidParameters()) );
00386   }
00387   else {
00388     params->validateParameters(*getValidParameters());
00389   }
00390 
00391   // Check for maximum number of iterations
00392   if (params->isParameter("Maximum Iterations")) {
00393     maxIters_ = params->get("Maximum Iterations",maxIters_default_);
00394 
00395     // Update parameter in our list and in status test.
00396     params_->set("Maximum Iterations", maxIters_);
00397     if (maxIterTest_!=Teuchos::null)
00398       maxIterTest_->setMaxIters( maxIters_ );
00399   }
00400 
00401   // Check for blocksize
00402   if (params->isParameter("Block Size")) {
00403     blockSize_ = params->get("Block Size",1);    
00404     TEST_FOR_EXCEPTION(blockSize_ != 1, std::invalid_argument,
00405            "Belos::TFQMRSolMgr: \"Block Size\" must be 1.");
00406 
00407     // Update parameter in our list.
00408     params_->set("Block Size", blockSize_);
00409   }
00410 
00411   // Check to see if the timer label changed.
00412   if (params->isParameter("Timer Label")) {
00413     std::string tempLabel = params->get("Timer Label", label_default_);
00414 
00415     // Update parameter in our list and solver timer
00416     if (tempLabel != label_) {
00417       label_ = tempLabel;
00418       params_->set("Timer Label", label_);
00419       std::string solveLabel = label_ + ": TFQMRSolMgr total solve time";
00420       timerSolve_ = Teuchos::TimeMonitor::getNewTimer(solveLabel);
00421     }
00422   }
00423 
00424   // Check for a change in verbosity level
00425   if (params->isParameter("Verbosity")) {
00426     if (Teuchos::isParameterType<int>(*params,"Verbosity")) {
00427       verbosity_ = params->get("Verbosity", verbosity_default_);
00428     } else {
00429       verbosity_ = (int)Teuchos::getParameter<Belos::MsgType>(*params,"Verbosity");
00430     }
00431 
00432     // Update parameter in our list.
00433     params_->set("Verbosity", verbosity_);
00434     if (printer_ != Teuchos::null)
00435       printer_->setVerbosity(verbosity_);
00436   }
00437 
00438   // Check for a change in output style
00439   if (params->isParameter("Output Style")) {
00440     if (Teuchos::isParameterType<int>(*params,"Output Style")) {
00441       outputStyle_ = params->get("Output Style", outputStyle_default_);
00442     } else {
00443       outputStyle_ = (int)Teuchos::getParameter<Belos::OutputType>(*params,"Output Style");
00444     }
00445 
00446     // Reconstruct the convergence test if the explicit residual test is not being used.
00447     params_->set("Output Style", outputStyle_);
00448     isSTSet_ = false;
00449   }
00450 
00451   // output stream
00452   if (params->isParameter("Output Stream")) {
00453     outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> >(*params,"Output Stream");
00454 
00455     // Update parameter in our list.
00456     params_->set("Output Stream", outputStream_);
00457     if (printer_ != Teuchos::null)
00458       printer_->setOStream( outputStream_ );
00459   }
00460 
00461   // frequency level
00462   if (verbosity_ & Belos::StatusTestDetails) {
00463     if (params->isParameter("Output Frequency")) {
00464       outputFreq_ = params->get("Output Frequency", outputFreq_default_);
00465     }
00466 
00467     // Update parameter in out list and output status test.
00468     params_->set("Output Frequency", outputFreq_);
00469     if (outputTest_ != Teuchos::null)
00470       outputTest_->setOutputFrequency( outputFreq_ );
00471   }
00472 
00473   // Create output manager if we need to.
00474   if (printer_ == Teuchos::null) {
00475     printer_ = Teuchos::rcp( new OutputManager<ScalarType>(verbosity_, outputStream_) );
00476   }  
00477   
00478   // Check for convergence tolerance
00479   if (params->isParameter("Convergence Tolerance")) {
00480     convtol_ = params->get("Convergence Tolerance",convtol_default_);
00481 
00482     // Update parameter in our list and residual tests.
00483     params_->set("Convergence Tolerance", convtol_);
00484     if (impConvTest_ != Teuchos::null)
00485       impConvTest_->setTolerance( convtol_ );
00486     if (expConvTest_ != Teuchos::null)
00487       expConvTest_->setTolerance( convtol_ );
00488   }
00489   
00490   // Check for a change in scaling, if so we need to build new residual tests.
00491   if (params->isParameter("Implicit Residual Scaling")) {
00492     std::string tempImpResScale = Teuchos::getParameter<std::string>( *params, "Implicit Residual Scaling" );
00493 
00494     // Only update the scaling if it's different.
00495     if (impResScale_ != tempImpResScale) {
00496       Belos::ScaleType impResScaleType = convertStringToScaleType( tempImpResScale );
00497       impResScale_ = tempImpResScale;
00498 
00499       // Update parameter in our list and residual tests
00500       params_->set("Implicit Residual Scaling", impResScale_);
00501       if (impConvTest_ != Teuchos::null) {
00502         try { 
00503           impConvTest_->defineScaleForm( impResScaleType, Belos::TwoNorm );
00504         }
00505         catch (std::exception& e) { 
00506           // Make sure the convergence test gets constructed again.
00507           isSTSet_ = false;
00508         }
00509       }
00510     }      
00511   }
00512   
00513   if (params->isParameter("Explicit Residual Scaling")) {
00514     std::string tempExpResScale = Teuchos::getParameter<std::string>( *params, "Explicit Residual Scaling" );
00515 
00516     // Only update the scaling if it's different.
00517     if (expResScale_ != tempExpResScale) {
00518       Belos::ScaleType expResScaleType = convertStringToScaleType( tempExpResScale );
00519       expResScale_ = tempExpResScale;
00520 
00521       // Update parameter in our list and residual tests
00522       params_->set("Explicit Residual Scaling", expResScale_);
00523       if (expConvTest_ != Teuchos::null) {
00524         try { 
00525           expConvTest_->defineScaleForm( expResScaleType, Belos::TwoNorm );
00526         }
00527         catch (std::exception& e) {
00528           // Make sure the convergence test gets constructed again.
00529           isSTSet_ = false;
00530         }
00531       }
00532     }      
00533   }
00534 
00535   if (params->isParameter("Explicit Residual Test")) {
00536     expResTest_ = Teuchos::getParameter<bool>( *params,"Explicit Residual Test" );
00537 
00538     // Reconstruct the convergence test if the explicit residual test is not being used.
00539     params_->set("Explicit Residual Test", expResTest_);
00540     if (expConvTest_ == Teuchos::null) {
00541       isSTSet_ = false;
00542     }
00543   }
00544 
00545   // Create the timer if we need to.
00546   if (timerSolve_ == Teuchos::null) {
00547     std::string solveLabel = label_ + ": TFQMRSolMgr total solve time";
00548     timerSolve_ = Teuchos::TimeMonitor::getNewTimer(solveLabel);
00549   }
00550 
00551   // Inform the solver manager that the current parameters were set.
00552   isSet_ = true;
00553 }
00554 
00555 
00556 // Check the status test versus the defined linear problem
00557 template<class ScalarType, class MV, class OP>
00558 bool TFQMRSolMgr<ScalarType,MV,OP>::checkStatusTest() {
00559 
00560   typedef Belos::StatusTestCombo<ScalarType,MV,OP>  StatusTestCombo_t;
00561   typedef Belos::StatusTestGenResNorm<ScalarType,MV,OP>  StatusTestGenResNorm_t;
00562 
00563   // Basic test checks maximum iterations and native residual.
00564   maxIterTest_ = Teuchos::rcp( new StatusTestMaxIters<ScalarType,MV,OP>( maxIters_ ) );
00565 
00566   if (expResTest_) {
00567    
00568     // Implicit residual test, using the native residual to determine if convergence was achieved.
00569     Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
00570       Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
00571     tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
00572     impConvTest_ = tmpImpConvTest;
00573 
00574     // Explicit residual test once the native residual is below the tolerance
00575     Teuchos::RCP<StatusTestGenResNorm_t> tmpExpConvTest =
00576       Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
00577     tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit, Belos::TwoNorm );
00578     tmpExpConvTest->defineScaleForm( convertStringToScaleType(expResScale_), Belos::TwoNorm );
00579     expConvTest_ = tmpExpConvTest;
00580 
00581     // The convergence test is a combination of the "cheap" implicit test and explicit test.
00582     convTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
00583   }
00584   else {
00585 
00586     // Implicit residual test, using the native residual to determine if convergence was achieved.
00587     Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
00588       Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
00589     tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
00590     impConvTest_ = tmpImpConvTest;
00591 
00592     // Set the explicit and total convergence test to this implicit test that checks for accuracy loss.
00593     expConvTest_ = impConvTest_;
00594     convTest_ = impConvTest_;
00595   }
00596   sTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
00597 
00598   // Create the status test output class.
00599   // This class manages and formats the output from the status test.
00600   StatusTestOutputFactory<ScalarType,MV,OP> stoFactory( outputStyle_ );
00601   outputTest_ = stoFactory.create( printer_, sTest_, outputFreq_, Passed+Failed+Undefined );
00602 
00603   // Set the solver string for the output test
00604   std::string solverDesc = " TFQMR ";
00605   outputTest_->setSolverDesc( solverDesc );
00606 
00607   
00608   // The status test is now set.
00609   isSTSet_ = true;
00610 
00611   return false;
00612 }
00613 
00614     
00615 template<class ScalarType, class MV, class OP>
00616 Teuchos::RCP<const Teuchos::ParameterList> 
00617 TFQMRSolMgr<ScalarType,MV,OP>::getValidParameters() const
00618 {
00619   static Teuchos::RCP<const Teuchos::ParameterList> validPL;
00620   
00621   // Set all the valid parameters and their default values.
00622   if(is_null(validPL)) {
00623     Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
00624     pl->set("Convergence Tolerance", convtol_default_,
00625       "The relative residual tolerance that needs to be achieved by the\n"
00626       "iterative solver in order for the linear system to be declared converged.");
00627     pl->set("Maximum Iterations", maxIters_default_,
00628       "The maximum number of block iterations allowed for each\n"
00629       "set of RHS solved.");
00630     pl->set("Verbosity", verbosity_default_,
00631       "What type(s) of solver information should be outputted\n"
00632       "to the output stream.");
00633     pl->set("Output Style", outputStyle_default_,
00634       "What style is used for the solver information outputted\n"
00635       "to the output stream.");
00636     pl->set("Output Frequency", outputFreq_default_,
00637       "How often convergence information should be outputted\n"
00638       "to the output stream.");  
00639     pl->set("Output Stream", outputStream_default_,
00640       "A reference-counted pointer to the output stream where all\n"
00641       "solver output is sent.");
00642     pl->set("Explicit Residual Test", expResTest_default_,
00643       "Whether the explicitly computed residual should be used in the convergence test.");
00644     pl->set("Implicit Residual Scaling", impResScale_default_,
00645       "The type of scaling used in the implicit residual convergence test.");
00646     pl->set("Explicit Residual Scaling", expResScale_default_,
00647       "The type of scaling used in the explicit residual convergence test.");
00648     pl->set("Timer Label", label_default_,
00649       "The string to use as a prefix for the timer labels.");
00650     //  pl->set("Restart Timers", restartTimers_);
00651     validPL = pl;
00652   }
00653   return validPL;
00654 }
00655 
00656   
00657 // solve()
00658 template<class ScalarType, class MV, class OP>
00659 ReturnType TFQMRSolMgr<ScalarType,MV,OP>::solve() {
00660 
00661   // Set the current parameters if they were not set before.
00662   // NOTE:  This may occur if the user generated the solver manager with the default constructor and 
00663   // then didn't set any parameters using setParameters().
00664   if (!isSet_) {
00665     setParameters(Teuchos::parameterList(*getValidParameters()));
00666   }
00667 
00668   Teuchos::BLAS<int,ScalarType> blas;
00669   Teuchos::LAPACK<int,ScalarType> lapack;
00670 
00671   TEST_FOR_EXCEPTION(problem_ == Teuchos::null,TFQMRSolMgrLinearProblemFailure,
00672          "Belos::TFQMRSolMgr::solve(): Linear problem is not a valid object.");
00673 
00674   TEST_FOR_EXCEPTION(!problem_->isProblemSet(),TFQMRSolMgrLinearProblemFailure,
00675                      "Belos::TFQMRSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
00676 
00677   if (!isSTSet_) {
00678     TEST_FOR_EXCEPTION( checkStatusTest(),TFQMRSolMgrLinearProblemFailure,
00679       "Belos::TFQMRSolMgr::solve(): Linear problem and requested status tests are incompatible.");
00680   }
00681 
00682   // Create indices for the linear systems to be solved.
00683   int startPtr = 0;
00684   int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
00685   int numCurrRHS = blockSize_;
00686 
00687   std::vector<int> currIdx, currIdx2;
00688 
00689   //  The index set is generated that informs the linear problem that some linear systems are augmented.
00690   currIdx.resize( blockSize_ );
00691   currIdx2.resize( blockSize_ );
00692   for (int i=0; i<numCurrRHS; ++i) 
00693     { currIdx[i] = startPtr+i; currIdx2[i]=i; }
00694 
00695   // Inform the linear problem of the current linear system to solve.
00696   problem_->setLSIndex( currIdx );
00697 
00699   // Parameter list
00700   Teuchos::ParameterList plist;
00701   plist.set("Block Size",blockSize_);
00702   
00703   // Reset the status test.  
00704   outputTest_->reset();
00705 
00706   // Assume convergence is achieved, then let any failed convergence set this to false.
00707   bool isConverged = true;  
00708 
00710   // TFQMR solver
00711 
00712   Teuchos::RCP<TFQMRIter<ScalarType,MV,OP> > tfqmr_iter = 
00713     Teuchos::rcp( new TFQMRIter<ScalarType,MV,OP>(problem_,printer_,outputTest_,plist) );
00714 
00715   // Enter solve() iterations
00716   {
00717     Teuchos::TimeMonitor slvtimer(*timerSolve_);
00718 
00719     while ( numRHS2Solve > 0 ) {
00720       //
00721       // Reset the active / converged vectors from this block
00722       std::vector<int> convRHSIdx;
00723       std::vector<int> currRHSIdx( currIdx );
00724       currRHSIdx.resize(numCurrRHS);
00725 
00726       // Reset the number of iterations.
00727       tfqmr_iter->resetNumIters();
00728 
00729       // Reset the number of calls that the status test output knows about.
00730       outputTest_->resetNumCalls();
00731 
00732       // Get the current residual for this block of linear systems.
00733       Teuchos::RCP<MV> R_0 = MVT::CloneViewNonConst( *(Teuchos::rcp_const_cast<MV>(problem_->getInitPrecResVec())), currIdx );
00734 
00735       // Set the new state and initialize the solver.
00736       TFQMRIterState<ScalarType,MV> newstate;
00737       newstate.R = R_0;
00738       tfqmr_iter->initializeTFQMR(newstate);
00739 
00740       while(1) {
00741   
00742   // tell tfqmr_iter to iterate
00743   try {
00744     tfqmr_iter->iterate();
00745     
00747     //
00748     // check convergence first
00749     //
00751     if ( convTest_->getStatus() == Passed ) {
00752       // We have convergence of the linear system.
00753       break;  // break from while(1){tfqmr_iter->iterate()}
00754     }
00756     //
00757     // check for maximum iterations
00758     //
00760     else if ( maxIterTest_->getStatus() == Passed ) {
00761       // we don't have convergence
00762       isConverged = false;
00763       break;  // break from while(1){tfqmr_iter->iterate()}
00764     }
00765 
00767     //
00768     // we returned from iterate(), but none of our status tests Passed.
00769     // something is wrong, and it is probably our fault.
00770     //
00772 
00773     else {
00774       TEST_FOR_EXCEPTION(true,std::logic_error,
00775              "Belos::TFQMRSolMgr::solve(): Invalid return from TFQMRIter::iterate().");
00776     }
00777   }
00778   catch (const std::exception &e) {
00779     printer_->stream(Errors) << "Error! Caught std::exception in TFQMRIter::iterate() at iteration " 
00780            << tfqmr_iter->getNumIters() << std::endl 
00781            << e.what() << std::endl;
00782     throw;
00783   }
00784       }
00785       
00786       // Inform the linear problem that we are finished with this block linear system.
00787       problem_->setCurrLS();
00788       
00789       // Update indices for the linear systems to be solved.
00790       startPtr += numCurrRHS;
00791       numRHS2Solve -= numCurrRHS;
00792       if ( numRHS2Solve > 0 ) {
00793   numCurrRHS = blockSize_;
00794 
00795   currIdx.resize( blockSize_ );
00796   currIdx2.resize( blockSize_ );
00797   for (int i=0; i<numCurrRHS; ++i) 
00798     { currIdx[i] = startPtr+i; currIdx2[i] = i; }
00799   // Set the next indices.
00800   problem_->setLSIndex( currIdx );
00801 
00802   // Set the new blocksize for the solver.
00803   tfqmr_iter->setBlockSize( blockSize_ ); 
00804       }
00805       else {
00806         currIdx.resize( numRHS2Solve );
00807       }
00808       
00809     }// while ( numRHS2Solve > 0 )
00810     
00811   }
00812 
00813   // print final summary
00814   sTest_->print( printer_->stream(FinalSummary) );
00815  
00816   // print timing information
00817   Teuchos::TimeMonitor::summarize( printer_->stream(TimingDetails) );
00818  
00819   // get iteration information for this solve
00820   numIters_ = maxIterTest_->getNumIters();
00821  
00822   if (!isConverged) {
00823     return Unconverged; // return from TFQMRSolMgr::solve() 
00824   }
00825   return Converged; // return from TFQMRSolMgr::solve() 
00826 }
00827 
00828 //  This method requires the solver manager to return a std::string that describes itself.
00829 template<class ScalarType, class MV, class OP>
00830 std::string TFQMRSolMgr<ScalarType,MV,OP>::description() const
00831 {
00832   std::ostringstream oss;
00833   oss << "Belos::TFQMRSolMgr<...,"<<Teuchos::ScalarTraits<ScalarType>::name()<<">";
00834   oss << "{}";
00835   return oss.str();
00836 }
00837   
00838 } // end Belos namespace
00839 
00840 #endif /* BELOS_TFQMR_SOLMGR_HPP */
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines