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

Generated on Tue Oct 20 12:48:33 2009 for Belos by doxygen 1.4.7