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

Generated on Wed May 12 21:30:08 2010 for Belos by  doxygen 1.4.7