BelosBlockGmresSolMgr.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 terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 //
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00025 //
00026 // ***********************************************************************
00027 // @HEADER
00028 
00029 #ifndef BELOS_BLOCK_GMRES_SOLMGR_HPP
00030 #define BELOS_BLOCK_GMRES_SOLMGR_HPP
00031 
00036 #include "BelosConfigDefs.hpp"
00037 #include "BelosTypes.hpp"
00038 
00039 #include "BelosLinearProblem.hpp"
00040 #include "BelosSolverManager.hpp"
00041 
00042 #include "BelosGmresIteration.hpp"
00043 #include "BelosBlockGmresIter.hpp"
00044 #include "BelosBlockFGmresIter.hpp"
00045 #include "BelosDGKSOrthoManager.hpp"
00046 #include "BelosICGSOrthoManager.hpp"
00047 #include "BelosIMGSOrthoManager.hpp"
00048 #include "BelosStatusTestMaxIters.hpp"
00049 #include "BelosStatusTestGenResNorm.hpp"
00050 #include "BelosStatusTestImpResNorm.hpp"
00051 #include "BelosStatusTestCombo.hpp"
00052 #include "BelosStatusTestOutput.hpp"
00053 #include "BelosOutputManager.hpp"
00054 #include "Teuchos_BLAS.hpp"
00055 #include "Teuchos_LAPACK.hpp"
00056 #include "Teuchos_TimeMonitor.hpp"
00057 
00075 namespace Belos {
00076   
00078 
00079   
00086 class BlockGmresSolMgrLinearProblemFailure : public BelosError {public:
00087   BlockGmresSolMgrLinearProblemFailure(const std::string& what_arg) : BelosError(what_arg)
00088     {}};
00089   
00096 class BlockGmresSolMgrOrthoFailure : public BelosError {public:
00097   BlockGmresSolMgrOrthoFailure(const std::string& what_arg) : BelosError(what_arg)
00098     {}};
00099   
00100 template<class ScalarType, class MV, class OP>
00101 class BlockGmresSolMgr : public SolverManager<ScalarType,MV,OP> {
00102     
00103 private:
00104   typedef MultiVecTraits<ScalarType,MV> MVT;
00105   typedef OperatorTraits<ScalarType,MV,OP> OPT;
00106   typedef Teuchos::ScalarTraits<ScalarType> SCT;
00107   typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
00108   typedef Teuchos::ScalarTraits<MagnitudeType> MT;
00109     
00110 public:
00111     
00113 
00114    
00120   BlockGmresSolMgr();
00121  
00135   BlockGmresSolMgr( const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem,
00136     const Teuchos::RCP<Teuchos::ParameterList> &pl );
00137     
00139   virtual ~BlockGmresSolMgr() {};
00141     
00143 
00144     
00147   const LinearProblem<ScalarType,MV,OP>& getProblem() const {
00148     return *problem_;
00149   }
00150 
00153   Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const;
00154 
00157   Teuchos::RCP<const Teuchos::ParameterList> getCurrentParameters() const { return params_; }
00158  
00164   Teuchos::Array<Teuchos::RCP<Teuchos::Time> > getTimers() const {
00165     return tuple(timerSolve_);
00166   }
00167   
00169   int getNumIters() const {
00170     return numIters_;
00171   } 
00172  
00176   bool isLOADetected() const { return loaDetected_; }
00177  
00179     
00181 
00182     
00183   void setProblem( const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem ) { problem_ = problem; isSTSet_ = false; }
00184     
00185   void setParameters( const Teuchos::RCP<Teuchos::ParameterList> &params );
00186     
00188     
00190 
00191     
00209   ReturnType solve();
00210     
00212     
00215     
00217   std::string description() const;
00218     
00220     
00221 private:
00222 
00223   // Method to convert std::string to enumerated type for residual.
00224   Belos::ScaleType convertStringToScaleType( std::string& scaleType ) {
00225     if (scaleType == "Norm of Initial Residual") {
00226       return Belos::NormOfInitRes;
00227     } else if (scaleType == "Norm of Preconditioned Initial Residual") {
00228       return Belos::NormOfPrecInitRes;
00229     } else if (scaleType == "Norm of RHS") {
00230       return Belos::NormOfRHS;
00231     } else if (scaleType == "None") {
00232       return Belos::None;
00233     } else 
00234       TEST_FOR_EXCEPTION( true ,std::logic_error,
00235         "Belos::BlockGmresSolMgr(): Invalid residual scaling type.");
00236   }
00237   
00238   // Method for checking current status test against defined linear problem.
00239   bool checkStatusTest();
00240 
00241   // Linear problem.
00242   Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > problem_;
00243     
00244   // Output manager.
00245   Teuchos::RCP<OutputManager<ScalarType> > printer_;
00246   Teuchos::RCP<std::ostream> outputStream_;
00247 
00248   // Status test.
00249   Teuchos::RCP<StatusTest<ScalarType,MV,OP> > sTest_;
00250   Teuchos::RCP<StatusTestMaxIters<ScalarType,MV,OP> > maxIterTest_;
00251   Teuchos::RCP<StatusTest<ScalarType,MV,OP> > convTest_;
00252   Teuchos::RCP<StatusTestResNorm<ScalarType,MV,OP> > expConvTest_, impConvTest_;
00253   Teuchos::RCP<StatusTestOutput<ScalarType,MV,OP> > outputTest_;
00254 
00255   // Orthogonalization manager.
00256   Teuchos::RCP<MatOrthoManager<ScalarType,MV,OP> > ortho_; 
00257     
00258   // Current parameter list.
00259   Teuchos::RCP<ParameterList> params_;
00260 
00261   // Default solver values.
00262   static const MagnitudeType convtol_default_;
00263   static const MagnitudeType orthoKappa_default_;
00264   static const int maxRestarts_default_;
00265   static const int maxIters_default_;
00266   static const bool adaptiveBlockSize_default_;
00267   static const bool showMaxResNormOnly_default_;
00268   static const bool flexibleGmres_default_;
00269   static const bool expResTest_default_;
00270   static const int blockSize_default_;
00271   static const int numBlocks_default_;
00272   static const int verbosity_default_;
00273   static const int outputFreq_default_;
00274   static const std::string impResScale_default_; 
00275   static const std::string expResScale_default_; 
00276   static const std::string label_default_;
00277   static const std::string orthoType_default_;
00278   static const Teuchos::RCP<std::ostream> outputStream_default_;
00279 
00280   // Current solver values.
00281   MagnitudeType convtol_, orthoKappa_;
00282   int maxRestarts_, maxIters_, numIters_;
00283   int blockSize_, numBlocks_, verbosity_, outputFreq_;
00284   bool adaptiveBlockSize_, showMaxResNormOnly_, isFlexible_, expResTest_;
00285   std::string orthoType_; 
00286   std::string impResScale_, expResScale_;
00287     
00288   // Timers.
00289   std::string label_;
00290   Teuchos::RCP<Teuchos::Time> timerSolve_;
00291 
00292   // Internal state variables.
00293   bool isSet_, isSTSet_;
00294   bool loaDetected_;
00295 };
00296 
00297 
00298 // Default solver values.
00299 template<class ScalarType, class MV, class OP>
00300 const typename Teuchos::ScalarTraits<ScalarType>::magnitudeType BlockGmresSolMgr<ScalarType,MV,OP>::convtol_default_ = 1e-8;
00301 
00302 template<class ScalarType, class MV, class OP>
00303 const typename Teuchos::ScalarTraits<ScalarType>::magnitudeType BlockGmresSolMgr<ScalarType,MV,OP>::orthoKappa_default_ = -1.0;
00304 
00305 template<class ScalarType, class MV, class OP>
00306 const int BlockGmresSolMgr<ScalarType,MV,OP>::maxRestarts_default_ = 20;
00307 
00308 template<class ScalarType, class MV, class OP>
00309 const int BlockGmresSolMgr<ScalarType,MV,OP>::maxIters_default_ = 1000;
00310 
00311 template<class ScalarType, class MV, class OP>
00312 const bool BlockGmresSolMgr<ScalarType,MV,OP>::adaptiveBlockSize_default_ = true;
00313 
00314 template<class ScalarType, class MV, class OP>
00315 const bool BlockGmresSolMgr<ScalarType,MV,OP>::showMaxResNormOnly_default_ = false;
00316 
00317 template<class ScalarType, class MV, class OP>
00318 const bool BlockGmresSolMgr<ScalarType,MV,OP>::flexibleGmres_default_ = false;
00319 
00320 template<class ScalarType, class MV, class OP>
00321 const bool BlockGmresSolMgr<ScalarType,MV,OP>::expResTest_default_ = false;
00322 
00323 template<class ScalarType, class MV, class OP>
00324 const int BlockGmresSolMgr<ScalarType,MV,OP>::blockSize_default_ = 1;
00325 
00326 template<class ScalarType, class MV, class OP>
00327 const int BlockGmresSolMgr<ScalarType,MV,OP>::numBlocks_default_ = 300;
00328 
00329 template<class ScalarType, class MV, class OP>
00330 const int BlockGmresSolMgr<ScalarType,MV,OP>::verbosity_default_ = Belos::Errors;
00331 
00332 template<class ScalarType, class MV, class OP>
00333 const int BlockGmresSolMgr<ScalarType,MV,OP>::outputFreq_default_ = -1;
00334 
00335 template<class ScalarType, class MV, class OP>
00336 const std::string BlockGmresSolMgr<ScalarType,MV,OP>::impResScale_default_ = "Norm of Preconditioned Initial Residual";
00337 
00338 template<class ScalarType, class MV, class OP>
00339 const std::string BlockGmresSolMgr<ScalarType,MV,OP>::expResScale_default_ = "Norm of Initial Residual";
00340 
00341 template<class ScalarType, class MV, class OP>
00342 const std::string BlockGmresSolMgr<ScalarType,MV,OP>::label_default_ = "Belos";
00343 
00344 template<class ScalarType, class MV, class OP>
00345 const std::string BlockGmresSolMgr<ScalarType,MV,OP>::orthoType_default_ = "DGKS";
00346 
00347 template<class ScalarType, class MV, class OP>
00348 const Teuchos::RCP<std::ostream> BlockGmresSolMgr<ScalarType,MV,OP>::outputStream_default_ = Teuchos::rcp(&std::cout,false);
00349 
00350 
00351 // Empty Constructor
00352 template<class ScalarType, class MV, class OP>
00353 BlockGmresSolMgr<ScalarType,MV,OP>::BlockGmresSolMgr() :
00354   outputStream_(outputStream_default_),
00355   convtol_(convtol_default_),
00356   orthoKappa_(orthoKappa_default_),
00357   maxRestarts_(maxRestarts_default_),
00358   maxIters_(maxIters_default_),
00359   blockSize_(blockSize_default_),
00360   numBlocks_(numBlocks_default_),
00361   verbosity_(verbosity_default_),
00362   outputFreq_(outputFreq_default_),
00363   adaptiveBlockSize_(adaptiveBlockSize_default_),
00364   showMaxResNormOnly_(showMaxResNormOnly_default_),
00365   isFlexible_(flexibleGmres_default_),
00366   expResTest_(expResTest_default_),
00367   orthoType_(orthoType_default_),
00368   impResScale_(impResScale_default_),
00369   expResScale_(expResScale_default_),
00370   label_(label_default_),
00371   isSet_(false),
00372   isSTSet_(false),
00373   loaDetected_(false)
00374 {}
00375 
00376 
00377 // Basic Constructor
00378 template<class ScalarType, class MV, class OP>
00379 BlockGmresSolMgr<ScalarType,MV,OP>::BlockGmresSolMgr( 
00380   const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem,
00381   const Teuchos::RCP<Teuchos::ParameterList> &pl ) : 
00382   problem_(problem),
00383   outputStream_(outputStream_default_),
00384   convtol_(convtol_default_),
00385   orthoKappa_(orthoKappa_default_),
00386   maxRestarts_(maxRestarts_default_),
00387   maxIters_(maxIters_default_),
00388   blockSize_(blockSize_default_),
00389   numBlocks_(numBlocks_default_),
00390   verbosity_(verbosity_default_),
00391   outputFreq_(outputFreq_default_),
00392   adaptiveBlockSize_(adaptiveBlockSize_default_),
00393   showMaxResNormOnly_(showMaxResNormOnly_default_),
00394   isFlexible_(flexibleGmres_default_),
00395   expResTest_(expResTest_default_),
00396   orthoType_(orthoType_default_),
00397   impResScale_(impResScale_default_),
00398   expResScale_(expResScale_default_),
00399   label_(label_default_),
00400   isSet_(false),
00401   isSTSet_(false),
00402   loaDetected_(false)
00403 {
00404 
00405   TEST_FOR_EXCEPTION(problem_ == Teuchos::null, std::invalid_argument, "Problem not given to solver manager.");
00406 
00407   // If the parameter list pointer is null, then set the current parameters to the default parameter list.
00408   if ( !is_null(pl) ) {
00409     setParameters( pl );  
00410   }
00411 
00412 }
00413 
00414 
00415 template<class ScalarType, class MV, class OP>
00416 Teuchos::RCP<const Teuchos::ParameterList>
00417 BlockGmresSolMgr<ScalarType,MV,OP>::getValidParameters() const
00418 {
00419   static Teuchos::RCP<const Teuchos::ParameterList> validPL;
00420   if (is_null(validPL)) {
00421     Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
00422     pl->set("Convergence Tolerance", convtol_default_,
00423       "The relative residual tolerance that needs to be achieved by the\n"
00424       "iterative solver in order for the linear system to be declared converged." );
00425     pl->set("Maximum Restarts", maxRestarts_default_,
00426       "The maximum number of restarts allowed for each\n"
00427       "set of RHS solved.");
00428     pl->set(
00429       "Maximum Iterations", maxIters_default_,
00430       "The maximum number of block iterations allowed for each\n"
00431       "set of RHS solved.");
00432     pl->set("Num Blocks", numBlocks_default_,
00433       "The maximum number of blocks allowed in the Krylov subspace\n"
00434       "for each set of RHS solved.");
00435     pl->set("Block Size", blockSize_default_,
00436       "The number of vectors in each block.  This number times the\n"
00437       "number of blocks is the total Krylov subspace dimension.");
00438     pl->set("Adaptive Block Size", adaptiveBlockSize_default_,
00439       "Whether the solver manager should adapt the block size\n"
00440       "based on the number of RHS to solve.");
00441     pl->set("Verbosity", verbosity_default_,
00442       "What type(s) of solver information should be outputted\n"
00443       "to the output stream.");
00444     pl->set("Output Frequency", outputFreq_default_,
00445       "How often convergence information should be outputted\n"
00446       "to the output stream.");  
00447     pl->set("Output Stream", outputStream_default_,
00448       "A reference-counted pointer to the output stream where all\n"
00449       "solver output is sent.");
00450     pl->set("Show Maximum Residual Norm Only", showMaxResNormOnly_default_,
00451       "When convergence information is printed, only show the maximum\n"
00452       "relative residual norm when the block size is greater than one.");
00453     pl->set("Flexible Gmres", flexibleGmres_default_,
00454       "Whether the solver manager should use the flexible variant\n"
00455       "of GMRES.");
00456     pl->set("Explicit Residual Test", expResTest_default_,
00457       "Whether the explicitly computed residual should be used in the convergence test.");
00458     pl->set("Implicit Residual Scaling", impResScale_default_,
00459       "The type of scaling used in the implicit residual convergence test.");
00460     pl->set("Explicit Residual Scaling", expResScale_default_,
00461       "The type of scaling used in the explicit residual convergence test.");
00462     pl->set("Timer Label", label_default_,
00463       "The string to use as a prefix for the timer labels.");
00464     //  pl->set("Restart Timers", restartTimers_);
00465     pl->set("Orthogonalization", orthoType_default_,
00466       "The type of orthogonalization to use: DGKS, ICGS, or IMGS.");
00467     pl->set("Orthogonalization Constant",orthoKappa_default_,
00468       "The constant used by DGKS orthogonalization to determine\n"
00469       "whether another step of classical Gram-Schmidt is necessary.");
00470     validPL = pl;
00471   }
00472   return validPL;
00473 }
00474 
00475 
00476 template<class ScalarType, class MV, class OP>
00477 void BlockGmresSolMgr<ScalarType,MV,OP>::setParameters( const Teuchos::RCP<Teuchos::ParameterList> &params )
00478 {
00479 
00480   // Create the internal parameter list if ones doesn't already exist.
00481   if (params_ == Teuchos::null) {
00482     params_ = Teuchos::rcp( new Teuchos::ParameterList(*getValidParameters()) );
00483   }
00484   else {
00485     params->validateParameters(*getValidParameters());
00486   }
00487 
00488   // Check for maximum number of restarts
00489   if (params->isParameter("Maximum Restarts")) {
00490     maxRestarts_ = params->get("Maximum Restarts",maxRestarts_default_);
00491 
00492     // Update parameter in our list.
00493     params_->set("Maximum Restarts", maxRestarts_);
00494   }
00495 
00496   // Check for maximum number of iterations
00497   if (params->isParameter("Maximum Iterations")) {
00498     maxIters_ = params->get("Maximum Iterations",maxIters_default_);
00499 
00500     // Update parameter in our list and in status test.
00501     params_->set("Maximum Iterations", maxIters_);
00502     if (maxIterTest_!=Teuchos::null)
00503       maxIterTest_->setMaxIters( maxIters_ );
00504   }
00505 
00506   // Check for blocksize
00507   if (params->isParameter("Block Size")) {
00508     blockSize_ = params->get("Block Size",blockSize_default_);    
00509     TEST_FOR_EXCEPTION(blockSize_ <= 0, std::invalid_argument,
00510       "Belos::BlockGmresSolMgr: \"Block Size\" must be strictly positive.");
00511 
00512     // Update parameter in our list.
00513     params_->set("Block Size", blockSize_);
00514   }
00515 
00516   // Check if the blocksize should be adaptive
00517   if (params->isParameter("Adapative Block Size")) {
00518     adaptiveBlockSize_ = params->get("Adaptive Block Size",adaptiveBlockSize_default_);
00519     
00520     // Update parameter in our list.
00521     params_->set("Adaptive Block Size", adaptiveBlockSize_);
00522   }
00523 
00524   // Check for the maximum number of blocks.
00525   if (params->isParameter("Num Blocks")) {
00526     numBlocks_ = params->get("Num Blocks",numBlocks_default_);
00527     TEST_FOR_EXCEPTION(numBlocks_ <= 0, std::invalid_argument,
00528       "Belos::BlockGmresSolMgr: \"Num Blocks\" must be strictly positive.");
00529 
00530     // Update parameter in our list.
00531     params_->set("Num Blocks", numBlocks_);
00532   }
00533 
00534   // Check to see if the timer label changed.
00535   if (params->isParameter("Timer Label")) {
00536     std::string tempLabel = params->get("Timer Label", label_default_);
00537 
00538     // Update parameter in our list and solver timer
00539     if (tempLabel != label_) {
00540       label_ = tempLabel;
00541       params_->set("Timer Label", label_);
00542       std::string solveLabel = label_ + ": BlockGmresSolMgr total solve time";
00543       timerSolve_ = Teuchos::TimeMonitor::getNewTimer(solveLabel);
00544     }
00545   }
00546 
00547   // Determine whether this solver should be "flexible".
00548   if (params->isParameter("Flexible Gmres")) {
00549     isFlexible_ = Teuchos::getParameter<bool>(*params,"Flexible Gmres");
00550     params_->set("Flexible Gmres", isFlexible_);
00551     if (isFlexible_ && expResTest_) {
00552       // Use an implicit convergence test if the Gmres solver is flexible
00553       isSTSet_ = false;
00554     }
00555   }
00556 
00557 
00558   // Check if the orthogonalization changed.
00559   if (params->isParameter("Orthogonalization")) {
00560     std::string tempOrthoType = params->get("Orthogonalization",orthoType_default_);
00561     TEST_FOR_EXCEPTION( tempOrthoType != "DGKS" && tempOrthoType != "ICGS" && tempOrthoType != "IMGS", 
00562                         std::invalid_argument,
00563       "Belos::BlockGmresSolMgr: \"Orthogonalization\" must be either \"DGKS\", \"ICGS\", or \"IMGS\".");
00564     if (tempOrthoType != orthoType_) {
00565       orthoType_ = tempOrthoType;
00566       // Create orthogonalization manager
00567       if (orthoType_=="DGKS") {
00568         if (orthoKappa_ <= 0) {
00569           ortho_ = Teuchos::rcp( new DGKSOrthoManager<ScalarType,MV,OP>( label_ ) );
00570         }
00571         else {
00572           ortho_ = Teuchos::rcp( new DGKSOrthoManager<ScalarType,MV,OP>( label_ ) );
00573           Teuchos::rcp_dynamic_cast<DGKSOrthoManager<ScalarType,MV,OP> >(ortho_)->setDepTol( orthoKappa_ );
00574         }
00575       }
00576       else if (orthoType_=="ICGS") {
00577         ortho_ = Teuchos::rcp( new ICGSOrthoManager<ScalarType,MV,OP>( label_ ) );
00578       } 
00579       else if (orthoType_=="IMGS") {
00580         ortho_ = Teuchos::rcp( new IMGSOrthoManager<ScalarType,MV,OP>( label_ ) );
00581       } 
00582     }  
00583   }
00584 
00585   // Check which orthogonalization constant to use.
00586   if (params->isParameter("Orthogonalization Constant")) {
00587     orthoKappa_ = params->get("Orthogonalization Constant",orthoKappa_default_);
00588 
00589     // Update parameter in our list.
00590     params_->set("Orthogonalization Constant",orthoKappa_);
00591     if (orthoType_=="DGKS") {
00592       if (orthoKappa_ > 0 && ortho_ != Teuchos::null) {
00593         Teuchos::rcp_dynamic_cast<DGKSOrthoManager<ScalarType,MV,OP> >(ortho_)->setDepTol( orthoKappa_ );
00594       }
00595     } 
00596   }
00597 
00598   // Check for a change in verbosity level
00599   if (params->isParameter("Verbosity")) {
00600     if (Teuchos::isParameterType<int>(*params,"Verbosity")) {
00601       verbosity_ = params->get("Verbosity", verbosity_default_);
00602     } else {
00603       verbosity_ = (int)Teuchos::getParameter<Belos::MsgType>(*params,"Verbosity");
00604     }
00605 
00606     // Update parameter in our list.
00607     params_->set("Verbosity", verbosity_);
00608     if (printer_ != Teuchos::null)
00609       printer_->setVerbosity(verbosity_);
00610   }
00611 
00612   // output stream
00613   if (params->isParameter("Output Stream")) {
00614     outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> >(*params,"Output Stream");
00615 
00616     // Update parameter in our list.
00617     params_->set("Output Stream", outputStream_);
00618     if (printer_ != Teuchos::null)
00619       printer_->setOStream( outputStream_ );
00620   }
00621 
00622   // frequency level
00623   if (verbosity_ & Belos::StatusTestDetails) {
00624     if (params->isParameter("Output Frequency")) {
00625       outputFreq_ = params->get("Output Frequency", outputFreq_default_);
00626     }
00627 
00628     // Update parameter in out list and output status test.
00629     params_->set("Output Frequency", outputFreq_);
00630     if (outputTest_ != Teuchos::null)
00631       outputTest_->setOutputFrequency( outputFreq_ );
00632   }
00633 
00634   // Create output manager if we need to.
00635   if (printer_ == Teuchos::null) {
00636     printer_ = Teuchos::rcp( new OutputManager<ScalarType>(verbosity_, outputStream_) );
00637   }  
00638   
00639   // Convergence
00640   typedef Belos::StatusTestCombo<ScalarType,MV,OP>  StatusTestCombo_t;
00641   typedef Belos::StatusTestGenResNorm<ScalarType,MV,OP>  StatusTestResNorm_t;
00642 
00643   // Check for convergence tolerance
00644   if (params->isParameter("Convergence Tolerance")) {
00645     convtol_ = params->get("Convergence Tolerance",convtol_default_);
00646 
00647     // Update parameter in our list and residual tests.
00648     params_->set("Convergence Tolerance", convtol_);
00649     if (impConvTest_ != Teuchos::null)
00650       impConvTest_->setTolerance( convtol_ );
00651     if (expConvTest_ != Teuchos::null)
00652       expConvTest_->setTolerance( convtol_ );
00653   }
00654  
00655   // Check for a change in scaling, if so we need to build new residual tests.
00656   if (params->isParameter("Implicit Residual Scaling")) {
00657     std::string tempImpResScale = Teuchos::getParameter<std::string>( *params, "Implicit Residual Scaling" );
00658 
00659     // Only update the scaling if it's different.
00660     if (impResScale_ != tempImpResScale) {
00661       Belos::ScaleType impResScaleType = convertStringToScaleType( tempImpResScale );
00662       impResScale_ = tempImpResScale;
00663 
00664       // Update parameter in our list and residual tests
00665       params_->set("Implicit Residual Scaling", impResScale_);
00666       if (impConvTest_ != Teuchos::null) {
00667         try { 
00668           impConvTest_->defineScaleForm( impResScaleType, Belos::TwoNorm );
00669         }
00670         catch (std::exception& e) { 
00671           // Make sure the convergence test gets constructed again.
00672           isSTSet_ = false;
00673         }
00674       }
00675     }      
00676   }
00677   
00678   if (params->isParameter("Explicit Residual Scaling")) {
00679     std::string tempExpResScale = Teuchos::getParameter<std::string>( *params, "Explicit Residual Scaling" );
00680 
00681     // Only update the scaling if it's different.
00682     if (expResScale_ != tempExpResScale) {
00683       Belos::ScaleType expResScaleType = convertStringToScaleType( tempExpResScale );
00684       expResScale_ = tempExpResScale;
00685 
00686       // Update parameter in our list and residual tests
00687       params_->set("Explicit Residual Scaling", expResScale_);
00688       if (expConvTest_ != Teuchos::null) {
00689         try { 
00690           expConvTest_->defineScaleForm( expResScaleType, Belos::TwoNorm );
00691         }
00692         catch (std::exception& e) {
00693           // Make sure the convergence test gets constructed again.
00694           isSTSet_ = false;
00695         }
00696       }
00697     }      
00698   }
00699 
00700   if (params->isParameter("Explicit Residual Test")) {
00701     expResTest_ = Teuchos::getParameter<bool>( *params,"Explicit Residual Test" );
00702 
00703     // Reconstruct the convergence test if the explicit residual test is not being used.
00704     params_->set("Explicit Residual Test", expResTest_);
00705     if (expConvTest_ == Teuchos::null) {
00706       isSTSet_ = false;
00707     }
00708   }
00709 
00710 
00711   if (params->isParameter("Show Maximum Residual Norm Only")) {
00712     showMaxResNormOnly_ = Teuchos::getParameter<bool>(*params,"Show Maximum Residual Norm Only");
00713 
00714     // Update parameter in our list and residual tests
00715     params_->set("Show Maximum Residual Norm Only", showMaxResNormOnly_);
00716     if (impConvTest_ != Teuchos::null)
00717       impConvTest_->setShowMaxResNormOnly( showMaxResNormOnly_ );
00718     if (expConvTest_ != Teuchos::null)
00719       expConvTest_->setShowMaxResNormOnly( showMaxResNormOnly_ );
00720   }
00721 
00722   // Create orthogonalization manager if we need to.
00723   if (ortho_ == Teuchos::null) {
00724     if (orthoType_=="DGKS") {
00725       if (orthoKappa_ <= 0) {
00726         ortho_ = Teuchos::rcp( new DGKSOrthoManager<ScalarType,MV,OP>( label_ ) );
00727       }
00728       else {
00729         ortho_ = Teuchos::rcp( new DGKSOrthoManager<ScalarType,MV,OP>( label_ ) );
00730         Teuchos::rcp_dynamic_cast<DGKSOrthoManager<ScalarType,MV,OP> >(ortho_)->setDepTol( orthoKappa_ );
00731       }
00732     }
00733     else if (orthoType_=="ICGS") {
00734       ortho_ = Teuchos::rcp( new ICGSOrthoManager<ScalarType,MV,OP>( label_ ) );
00735     } 
00736     else if (orthoType_=="IMGS") {
00737       ortho_ = Teuchos::rcp( new IMGSOrthoManager<ScalarType,MV,OP>( label_ ) );
00738     } 
00739     else {
00740       TEST_FOR_EXCEPTION(orthoType_!="ICGS"&&orthoType_!="DGKS"&&orthoType_!="IMGS",std::logic_error,
00741         "Belos::BlockGmresSolMgr(): Invalid orthogonalization type.");
00742     }  
00743   }
00744 
00745   // Create the timer if we need to.
00746   if (timerSolve_ == Teuchos::null) {
00747     std::string solveLabel = label_ + ": BlockGmresSolMgr total solve time";
00748     timerSolve_ = Teuchos::TimeMonitor::getNewTimer(solveLabel);
00749   }
00750 
00751   // Inform the solver manager that the current parameters were set.
00752   isSet_ = true;
00753 }
00754 
00755 // Check the status test versus the defined linear problem
00756 template<class ScalarType, class MV, class OP>
00757 bool BlockGmresSolMgr<ScalarType,MV,OP>::checkStatusTest() {
00758 
00759   typedef Belos::StatusTestCombo<ScalarType,MV,OP>  StatusTestCombo_t;
00760   typedef Belos::StatusTestGenResNorm<ScalarType,MV,OP>  StatusTestGenResNorm_t;
00761   typedef Belos::StatusTestImpResNorm<ScalarType,MV,OP>  StatusTestImpResNorm_t;
00762 
00763   // Basic test checks maximum iterations and native residual.
00764   maxIterTest_ = Teuchos::rcp( new StatusTestMaxIters<ScalarType,MV,OP>( maxIters_ ) );
00765 
00766   // If there is a left preconditioner, we create a combined status test that checks the implicit
00767   // and then explicit residual norm to see if we have convergence.
00768   if (!Teuchos::is_null(problem_->getLeftPrec()) && !isFlexible_) {
00769     expResTest_ = true;
00770   }
00771 
00772   if (expResTest_) {
00773    
00774     // Implicit residual test, using the native residual to determine if convergence was achieved.
00775     Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
00776       Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
00777     tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
00778     tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
00779     impConvTest_ = tmpImpConvTest;
00780 
00781     // Explicit residual test once the native residual is below the tolerance
00782     Teuchos::RCP<StatusTestGenResNorm_t> tmpExpConvTest =
00783       Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
00784     tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit, Belos::TwoNorm );
00785     tmpExpConvTest->defineScaleForm( convertStringToScaleType(expResScale_), Belos::TwoNorm );
00786     tmpExpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
00787     expConvTest_ = tmpExpConvTest;
00788 
00789     // The convergence test is a combination of the "cheap" implicit test and explicit test.
00790     convTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
00791   }
00792   else {
00793 
00794     if (isFlexible_) {
00795       // Implicit residual test, using the native residual to determine if convergence was achieved.
00796       Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
00797         Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
00798       tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
00799       tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
00800       impConvTest_ = tmpImpConvTest;
00801     }
00802     else {
00803       // Implicit residual test, using the native residual to determine if convergence was achieved.
00804       // Use test that checks for loss of accuracy.
00805       Teuchos::RCP<StatusTestImpResNorm_t> tmpImpConvTest =
00806         Teuchos::rcp( new StatusTestImpResNorm_t( convtol_ ) );
00807       tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
00808       tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
00809       impConvTest_ = tmpImpConvTest;
00810     }
00811 
00812     // Set the explicit and total convergence test to this implicit test that checks for accuracy loss.
00813     expConvTest_ = impConvTest_;
00814     convTest_ = impConvTest_;
00815   }
00816 
00817   sTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
00818 
00819   if (outputFreq_ > 0) {
00820     outputTest_ = Teuchos::rcp( new StatusTestOutput<ScalarType,MV,OP>( printer_,
00821         sTest_,
00822         outputFreq_,
00823         Passed+Failed+Undefined ) );
00824   }
00825   else {
00826     outputTest_ = Teuchos::rcp( new StatusTestOutput<ScalarType,MV,OP>( printer_,
00827         sTest_, 1 ) );
00828   }
00829 
00830   // The status test is now set.
00831   isSTSet_ = true;
00832 
00833   return false;
00834 }
00835   
00836 // solve()
00837 template<class ScalarType, class MV, class OP>
00838 ReturnType BlockGmresSolMgr<ScalarType,MV,OP>::solve() {
00839 
00840   // Set the current parameters if they were not set before.
00841   // NOTE:  This may occur if the user generated the solver manager with the default constructor and 
00842   // then didn't set any parameters using setParameters().
00843   if (!isSet_) {
00844     setParameters(Teuchos::parameterList(*getValidParameters()));
00845   }
00846 
00847   Teuchos::BLAS<int,ScalarType> blas;
00848   Teuchos::LAPACK<int,ScalarType> lapack;
00849   
00850   TEST_FOR_EXCEPTION(problem_ == Teuchos::null,BlockGmresSolMgrLinearProblemFailure,
00851     "Belos::BlockGmresSolMgr::solve(): Linear problem is not a valid object.");
00852 
00853   TEST_FOR_EXCEPTION(!problem_->isProblemSet(),BlockGmresSolMgrLinearProblemFailure,
00854     "Belos::BlockGmresSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
00855 
00856   if (isFlexible_) {
00857     TEST_FOR_EXCEPTION(problem_->getRightPrec()==Teuchos::null,BlockGmresSolMgrLinearProblemFailure,
00858       "Belos::BlockGmresSolMgr::solve(): Linear problem does not have a preconditioner required for flexible GMRES, call setRightPrec().");
00859   }
00860  
00861   if (!isSTSet_ || (!expResTest_ && !Teuchos::is_null(problem_->getLeftPrec())) ) {
00862     TEST_FOR_EXCEPTION( checkStatusTest(),BlockGmresSolMgrLinearProblemFailure,
00863       "Belos::BlockGmresSolMgr::solve(): Linear problem and requested status tests are incompatible.");
00864   }
00865 
00866   // Create indices for the linear systems to be solved.
00867   int startPtr = 0;
00868   int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
00869   int numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
00870 
00871   std::vector<int> currIdx;
00872   //  If an adaptive block size is allowed then only the linear systems that need to be solved are solved.
00873   //  Otherwise, the index set is generated that informs the linear problem that some linear systems are augmented.
00874   if ( adaptiveBlockSize_ ) {
00875     blockSize_ = numCurrRHS;
00876     currIdx.resize( numCurrRHS  );
00877     for (int i=0; i<numCurrRHS; ++i) 
00878     { currIdx[i] = startPtr+i; }
00879     
00880   }
00881   else {
00882     currIdx.resize( blockSize_ );
00883     for (int i=0; i<numCurrRHS; ++i) 
00884     { currIdx[i] = startPtr+i; }
00885     for (int i=numCurrRHS; i<blockSize_; ++i) 
00886     { currIdx[i] = -1; }
00887   }
00888 
00889   // Inform the linear problem of the current linear system to solve.
00890   problem_->setLSIndex( currIdx );
00891 
00893   // Parameter list
00894   Teuchos::ParameterList plist;
00895   plist.set("Block Size",blockSize_);
00896 
00897   int dim = MVT::GetVecLength( *(problem_->getRHS()) );  
00898   if (blockSize_*numBlocks_ > dim) {
00899     int tmpNumBlocks = 0;
00900     if (blockSize_ == 1)
00901       tmpNumBlocks = dim / blockSize_;  // Allow for a good breakdown.
00902     else
00903       tmpNumBlocks = ( dim - blockSize_) / blockSize_;  // Allow for restarting.
00904     printer_->stream(Warnings) << 
00905       "Warning! Requested Krylov subspace dimension is larger that operator dimension!" << std::endl <<
00906       " The maximum number of blocks allowed for the Krylov subspace will be adjusted to " << tmpNumBlocks << std::endl;
00907     plist.set("Num Blocks",tmpNumBlocks);
00908   } 
00909   else 
00910     plist.set("Num Blocks",numBlocks_);
00911   
00912   // Reset the status test.  
00913   outputTest_->reset();
00914   loaDetected_ = false;
00915 
00916   // Assume convergence is achieved, then let any failed convergence set this to false.
00917   bool isConverged = true;  
00918 
00920   // BlockGmres solver
00921 
00922   Teuchos::RCP<GmresIteration<ScalarType,MV,OP> > block_gmres_iter;
00923 
00924   if (isFlexible_)
00925     block_gmres_iter = Teuchos::rcp( new BlockFGmresIter<ScalarType,MV,OP>(problem_,printer_,outputTest_,ortho_,plist) );
00926   else
00927     block_gmres_iter = Teuchos::rcp( new BlockGmresIter<ScalarType,MV,OP>(problem_,printer_,outputTest_,ortho_,plist) );
00928   
00929   // Enter solve() iterations
00930   {
00931     Teuchos::TimeMonitor slvtimer(*timerSolve_);
00932 
00933     while ( numRHS2Solve > 0 ) {
00934 
00935       // Set the current number of blocks and blocksize with the Gmres iteration.
00936       if (blockSize_*numBlocks_ > dim) {
00937         int tmpNumBlocks = 0;
00938         if (blockSize_ == 1)
00939           tmpNumBlocks = dim / blockSize_;  // Allow for a good breakdown.
00940         else
00941           tmpNumBlocks = ( dim - blockSize_) / blockSize_;  // Allow for restarting.
00942         block_gmres_iter->setSize( blockSize_, tmpNumBlocks );
00943       }
00944       else
00945         block_gmres_iter->setSize( blockSize_, numBlocks_ );
00946 
00947       // Reset the number of iterations.
00948       block_gmres_iter->resetNumIters();
00949 
00950       // Reset the number of calls that the status test output knows about.
00951       outputTest_->resetNumCalls();
00952 
00953       // Create the first block in the current Krylov basis.
00954       Teuchos::RCP<MV> V_0;
00955       if (isFlexible_)
00956         V_0 = MVT::CloneCopy( *(problem_->getInitResVec()), currIdx );
00957       else 
00958         V_0 = MVT::CloneCopy( *(problem_->getInitPrecResVec()), currIdx );
00959 
00960 
00961       // Get a matrix to hold the orthonormalization coefficients.
00962       Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > z_0 = 
00963         rcp( new Teuchos::SerialDenseMatrix<int,ScalarType>( blockSize_, blockSize_ ) );
00964  
00965       // Orthonormalize the new V_0
00966       int rank = ortho_->normalize( *V_0, z_0 );
00967       TEST_FOR_EXCEPTION(rank != blockSize_,BlockGmresSolMgrOrthoFailure,
00968         "Belos::BlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors.");
00969       
00970       // Set the new state and initialize the solver.
00971       GmresIterationState<ScalarType,MV> newstate;
00972       newstate.V = V_0;
00973       newstate.z = z_0;
00974       newstate.curDim = 0;
00975       block_gmres_iter->initializeGmres(newstate);
00976       int numRestarts = 0;
00977 
00978       while(1) {
00979         // tell block_gmres_iter to iterate
00980         try {
00981           block_gmres_iter->iterate();
00982     
00984           //
00985           // check convergence first
00986           //
00988           if ( convTest_->getStatus() == Passed ) {
00989             if ( expConvTest_->getLOADetected() ) {
00990               // we don't have convergence
00991               loaDetected_ = true;
00992               isConverged = false;
00993             }
00994             break;  // break from while(1){block_gmres_iter->iterate()}
00995           }
00997           //
00998           // check for maximum iterations
00999           //
01001           else if ( maxIterTest_->getStatus() == Passed ) {
01002             // we don't have convergence
01003             isConverged = false;
01004             break;  // break from while(1){block_gmres_iter->iterate()}
01005           }
01007           //
01008           // check for restarting, i.e. the subspace is full
01009           //
01011           else if ( block_gmres_iter->getCurSubspaceDim() == block_gmres_iter->getMaxSubspaceDim() ) {
01012       
01013             if ( numRestarts >= maxRestarts_ ) {
01014               isConverged = false;
01015               break; // break from while(1){block_gmres_iter->iterate()}
01016             }
01017             numRestarts++;
01018       
01019             printer_->stream(Debug) << " Performing restart number " << numRestarts << " of " << maxRestarts_ << std::endl << std::endl;
01020       
01021             // Update the linear problem.
01022             Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
01023             if (isFlexible_) {
01024               // Update the solution manually, since the preconditioning doesn't need to be undone.
01025               Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
01026               MVT::MvAddMv( 1.0, *curX, 1.0, *update, *curX );
01027             }
01028             else 
01029               problem_->updateSolution( update, true );
01030 
01031             // Get the state.
01032             GmresIterationState<ScalarType,MV> oldState = block_gmres_iter->getState();
01033       
01034             // Compute the restart std::vector.
01035             // Get a view of the current Krylov basis.
01036             V_0  = MVT::Clone( *(oldState.V), blockSize_ );
01037             if (isFlexible_)
01038               problem_->computeCurrResVec( &*V_0 );
01039             else
01040               problem_->computeCurrPrecResVec( &*V_0 );
01041 
01042             // Get a view of the first block of the Krylov basis.
01043             z_0 = rcp( new Teuchos::SerialDenseMatrix<int,ScalarType>( blockSize_, blockSize_ ) );
01044       
01045             // Orthonormalize the new V_0
01046             rank = ortho_->normalize( *V_0, z_0 );
01047             TEST_FOR_EXCEPTION(rank != blockSize_,BlockGmresSolMgrOrthoFailure,
01048               "Belos::BlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors after restart.");
01049 
01050             // Set the new state and initialize the solver.
01051             newstate.V = V_0;
01052             newstate.z = z_0;
01053             newstate.curDim = 0;
01054             block_gmres_iter->initializeGmres(newstate);
01055 
01056           } // end of restarting
01057 
01059           //
01060           // we returned from iterate(), but none of our status tests Passed.
01061           // something is wrong, and it is probably our fault.
01062           //
01064 
01065           else {
01066             TEST_FOR_EXCEPTION(true,std::logic_error,
01067               "Belos::BlockGmresSolMgr::solve(): Invalid return from BlockGmresIter::iterate().");
01068           }
01069         }
01070         catch (GmresIterationOrthoFailure e) {
01071           // If the block size is not one, it's not considered a lucky breakdown.
01072           if (blockSize_ != 1) {
01073             printer_->stream(Errors) << "Error! Caught std::exception in BlockGmresIter::iterate() at iteration " 
01074                                      << block_gmres_iter->getNumIters() << std::endl 
01075                                      << e.what() << std::endl;
01076             if (convTest_->getStatus() != Passed)
01077               isConverged = false;
01078             break;
01079           } 
01080           else {
01081             // If the block size is one, try to recover the most recent least-squares solution
01082             block_gmres_iter->updateLSQR( block_gmres_iter->getCurSubspaceDim() );
01083 
01084             // Check to see if the most recent least-squares solution yielded convergence.
01085             sTest_->checkStatus( &*block_gmres_iter );
01086             if (convTest_->getStatus() != Passed)
01087               isConverged = false;
01088             break;
01089           }
01090         }
01091         catch (std::exception e) {
01092           printer_->stream(Errors) << "Error! Caught std::exception in BlockGmresIter::iterate() at iteration " 
01093                                    << block_gmres_iter->getNumIters() << std::endl 
01094                                    << e.what() << std::endl;
01095           throw;
01096         }
01097       }
01098       
01099       // Compute the current solution.
01100       // Update the linear problem.
01101       if (isFlexible_) {
01102         // Update the solution manually, since the preconditioning doesn't need to be undone.
01103         Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
01104         Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
01105         MVT::MvAddMv( 1.0, *curX, 1.0, *update, *curX );
01106       }
01107       else {
01108         // Attempt to get the current solution from the residual status test, if it has one.
01109         if ( !Teuchos::is_null(expConvTest_->getSolution()) ) {
01110           Teuchos::RCP<MV> newX = expConvTest_->getSolution();
01111           Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
01112           MVT::MvAddMv( 0.0, *newX, 1.0, *newX, *curX );
01113         }
01114         else {
01115           Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
01116           problem_->updateSolution( update, true );
01117         }  
01118       }
01119 
01120       // Inform the linear problem that we are finished with this block linear system.
01121       problem_->setCurrLS();
01122       
01123       // Update indices for the linear systems to be solved.
01124       startPtr += numCurrRHS;
01125       numRHS2Solve -= numCurrRHS;
01126       if ( numRHS2Solve > 0 ) {
01127         numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
01128 
01129         if ( adaptiveBlockSize_ ) {
01130           blockSize_ = numCurrRHS;
01131           currIdx.resize( numCurrRHS  );
01132           for (int i=0; i<numCurrRHS; ++i) 
01133           { currIdx[i] = startPtr+i; }    
01134         }
01135         else {
01136           currIdx.resize( blockSize_ );
01137           for (int i=0; i<numCurrRHS; ++i) 
01138           { currIdx[i] = startPtr+i; }
01139           for (int i=numCurrRHS; i<blockSize_; ++i) 
01140           { currIdx[i] = -1; }
01141         }
01142         // Set the next indices.
01143         problem_->setLSIndex( currIdx );
01144       }
01145       else {
01146         currIdx.resize( numRHS2Solve );
01147       }
01148       
01149     }// while ( numRHS2Solve > 0 )
01150     
01151   }
01152 
01153   // print final summary
01154   sTest_->print( printer_->stream(FinalSummary) );
01155  
01156   // print timing information
01157   Teuchos::TimeMonitor::summarize( printer_->stream(TimingDetails) );
01158 
01159   // get iteration information for this solve
01160   numIters_ = maxIterTest_->getNumIters();
01161   
01162   if (!isConverged || loaDetected_) {
01163     return Unconverged; // return from BlockGmresSolMgr::solve() 
01164   }
01165   return Converged; // return from BlockGmresSolMgr::solve() 
01166 }
01167 
01168 //  This method requires the solver manager to return a std::string that describes itself.
01169 template<class ScalarType, class MV, class OP>
01170 std::string BlockGmresSolMgr<ScalarType,MV,OP>::description() const
01171 {
01172   std::ostringstream oss;
01173   oss << "Belos::BlockGmresSolMgr<...,"<<Teuchos::ScalarTraits<ScalarType>::name()<<">";
01174   oss << "{";
01175   if (isFlexible_) {
01176     oss << "Variant=\'Flexible\', ";
01177   }
01178   oss << "Ortho Type='"<<orthoType_<<"\', Block Size=" <<blockSize_;
01179   oss << ", Num Blocks=" <<numBlocks_<< ", Max Restarts=" << maxRestarts_;
01180   oss << "}";
01181   return oss.str();
01182 }
01183   
01184 } // end Belos namespace
01185 
01186 #endif /* BELOS_BLOCK_GMRES_SOLMGR_HPP */

Generated on Wed May 12 21:45:50 2010 for Belos by  doxygen 1.4.7