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