Belos  Version of the Day
BelosBlockGmresSolMgr.hpp
Go to the documentation of this file.
1 //@HEADER
2 // ************************************************************************
3 //
4 // Belos: Block Linear Solvers Package
5 // Copyright 2004 Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ************************************************************************
40 //@HEADER
41 
42 #ifndef BELOS_BLOCK_GMRES_SOLMGR_HPP
43 #define BELOS_BLOCK_GMRES_SOLMGR_HPP
44 
49 #include "BelosConfigDefs.hpp"
50 #include "BelosTypes.hpp"
51 
52 #include "BelosLinearProblem.hpp"
53 #include "BelosSolverManager.hpp"
54 
55 #include "BelosGmresIteration.hpp"
56 #include "BelosBlockGmresIter.hpp"
57 #include "BelosBlockFGmresIter.hpp"
64 #include "BelosStatusTestCombo.hpp"
67 #include "BelosOutputManager.hpp"
68 #include "Teuchos_BLAS.hpp"
69 #ifdef BELOS_TEUCHOS_TIME_MONITOR
70 #include "Teuchos_TimeMonitor.hpp"
71 #endif
72 
83 namespace Belos {
84 
86 
87 
95  BlockGmresSolMgrLinearProblemFailure(const std::string& what_arg) : BelosError(what_arg)
96  {}};
97 
105  BlockGmresSolMgrOrthoFailure(const std::string& what_arg) : BelosError(what_arg)
106  {}};
107 
124 template<class ScalarType, class MV, class OP>
125 class BlockGmresSolMgr : public SolverManager<ScalarType,MV,OP> {
126 
127 private:
130  typedef Teuchos::ScalarTraits<ScalarType> SCT;
131  typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
132  typedef Teuchos::ScalarTraits<MagnitudeType> MT;
133 
134 public:
135 
137 
138 
145 
166  BlockGmresSolMgr( const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem,
167  const Teuchos::RCP<Teuchos::ParameterList> &pl );
168 
170  virtual ~BlockGmresSolMgr() {};
171 
173  Teuchos::RCP<SolverManager<ScalarType, MV, OP> > clone () const override {
174  return Teuchos::rcp(new BlockGmresSolMgr<ScalarType,MV,OP>);
175  }
177 
179 
180 
183  const LinearProblem<ScalarType,MV,OP>& getProblem() const override {
184  return *problem_;
185  }
186 
189  Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const override;
190 
193  Teuchos::RCP<const Teuchos::ParameterList> getCurrentParameters() const override { return params_; }
194 
200  Teuchos::Array<Teuchos::RCP<Teuchos::Time> > getTimers() const {
201  return Teuchos::tuple(timerSolve_);
202  }
203 
214  MagnitudeType achievedTol() const override {
215  return achievedTol_;
216  }
217 
219  int getNumIters() const override {
220  return numIters_;
221  }
222 
226  bool isLOADetected() const override { return loaDetected_; }
227 
229 
231 
232 
234  void setProblem( const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem ) override { problem_ = problem; isSTSet_ = false; }
235 
237  void setParameters( const Teuchos::RCP<Teuchos::ParameterList> &params ) override;
238 
240  void setDebugStatusTest( const Teuchos::RCP<StatusTest<ScalarType,MV,OP> > &debugStatusTest ) override;
241 
243 
245 
246 
250  void reset( const ResetType type ) override { if ((type & Belos::Problem) && !Teuchos::is_null(problem_)) problem_->setProblem(); }
252 
254 
255 
273  ReturnType solve() override;
274 
276 
279 
286  void
287  describe (Teuchos::FancyOStream& out,
288  const Teuchos::EVerbosityLevel verbLevel =
289  Teuchos::Describable::verbLevel_default) const override;
290 
292  std::string description () const override;
293 
295 
296 private:
297 
298  // Method for checking current status test against defined linear problem.
299  bool checkStatusTest();
300 
301  // Linear problem.
302  Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > problem_;
303 
304  // Output manager.
305  Teuchos::RCP<OutputManager<ScalarType> > printer_;
306  Teuchos::RCP<std::ostream> outputStream_;
307 
308  // Status test.
309  Teuchos::RCP<StatusTest<ScalarType,MV,OP> > debugStatusTest_;
310  Teuchos::RCP<StatusTest<ScalarType,MV,OP> > sTest_;
311  Teuchos::RCP<StatusTestMaxIters<ScalarType,MV,OP> > maxIterTest_;
312  Teuchos::RCP<StatusTest<ScalarType,MV,OP> > convTest_;
313  Teuchos::RCP<StatusTestResNorm<ScalarType,MV,OP> > expConvTest_, impConvTest_;
314  Teuchos::RCP<StatusTestOutput<ScalarType,MV,OP> > outputTest_;
315 
316  // Orthogonalization manager.
317  Teuchos::RCP<MatOrthoManager<ScalarType,MV,OP> > ortho_;
318 
319  // Current parameter list.
320  Teuchos::RCP<Teuchos::ParameterList> params_;
321 
322  // Default solver values.
323  static constexpr int maxRestarts_default_ = 20;
324  static constexpr int maxIters_default_ = 1000;
325  static constexpr bool adaptiveBlockSize_default_ = true;
326  static constexpr bool showMaxResNormOnly_default_ = false;
327  static constexpr bool flexibleGmres_default_ = false;
328  static constexpr bool expResTest_default_ = false;
329  static constexpr int blockSize_default_ = 1;
330  static constexpr int numBlocks_default_ = 300;
331  static constexpr int verbosity_default_ = Belos::Errors;
332  static constexpr int outputStyle_default_ = Belos::General;
333  static constexpr int outputFreq_default_ = -1;
334  static constexpr const char * impResScale_default_ = "Norm of Preconditioned Initial Residual";
335  static constexpr const char * expResScale_default_ = "Norm of Initial Residual";
336  static constexpr const char * label_default_ = "Belos";
337  static constexpr const char * orthoType_default_ = "DGKS";
338  static constexpr std::ostream * outputStream_default_ = &std::cout;
339 
340  // Current solver values.
341  MagnitudeType convtol_, orthoKappa_, achievedTol_;
342  int maxRestarts_, maxIters_, numIters_;
343  int blockSize_, numBlocks_, verbosity_, outputStyle_, outputFreq_;
344  bool adaptiveBlockSize_, showMaxResNormOnly_, isFlexible_, expResTest_;
345  std::string orthoType_;
346  std::string impResScale_, expResScale_;
347 
348  // Timers.
349  std::string label_;
350  Teuchos::RCP<Teuchos::Time> timerSolve_;
351 
352  // Internal state variables.
353  bool isSet_, isSTSet_;
354  bool loaDetected_;
355 };
356 
357 
358 // Empty Constructor
359 template<class ScalarType, class MV, class OP>
361  outputStream_(Teuchos::rcp(outputStream_default_,false)),
362  convtol_(DefaultSolverParameters::convTol),
363  orthoKappa_(DefaultSolverParameters::orthoKappa),
364  achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
365  maxRestarts_(maxRestarts_default_),
366  maxIters_(maxIters_default_),
367  numIters_(0),
368  blockSize_(blockSize_default_),
369  numBlocks_(numBlocks_default_),
370  verbosity_(verbosity_default_),
371  outputStyle_(outputStyle_default_),
372  outputFreq_(outputFreq_default_),
373  adaptiveBlockSize_(adaptiveBlockSize_default_),
374  showMaxResNormOnly_(showMaxResNormOnly_default_),
375  isFlexible_(flexibleGmres_default_),
376  expResTest_(expResTest_default_),
377  orthoType_(orthoType_default_),
378  impResScale_(impResScale_default_),
379  expResScale_(expResScale_default_),
380  label_(label_default_),
381  isSet_(false),
382  isSTSet_(false),
383  loaDetected_(false)
384 {}
385 
386 
387 // Basic Constructor
388 template<class ScalarType, class MV, class OP>
390 BlockGmresSolMgr (const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem,
391  const Teuchos::RCP<Teuchos::ParameterList> &pl) :
392  problem_(problem),
393  outputStream_(Teuchos::rcp(outputStream_default_,false)),
394  convtol_(DefaultSolverParameters::convTol),
395  orthoKappa_(DefaultSolverParameters::orthoKappa),
396  achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
397  maxRestarts_(maxRestarts_default_),
398  maxIters_(maxIters_default_),
399  numIters_(0),
400  blockSize_(blockSize_default_),
401  numBlocks_(numBlocks_default_),
402  verbosity_(verbosity_default_),
403  outputStyle_(outputStyle_default_),
404  outputFreq_(outputFreq_default_),
405  adaptiveBlockSize_(adaptiveBlockSize_default_),
406  showMaxResNormOnly_(showMaxResNormOnly_default_),
407  isFlexible_(flexibleGmres_default_),
408  expResTest_(expResTest_default_),
409  orthoType_(orthoType_default_),
410  impResScale_(impResScale_default_),
411  expResScale_(expResScale_default_),
412  label_(label_default_),
413  isSet_(false),
414  isSTSet_(false),
415  loaDetected_(false)
416 {
417 
418  TEUCHOS_TEST_FOR_EXCEPTION(problem_ == Teuchos::null, std::invalid_argument, "Problem not given to solver manager.");
419 
420  // If the parameter list pointer is null, then set the current parameters to the default parameter list.
421  if ( !is_null(pl) ) {
422  setParameters( pl );
423  }
424 
425 }
426 
427 
428 template<class ScalarType, class MV, class OP>
429 Teuchos::RCP<const Teuchos::ParameterList>
431 {
432  static Teuchos::RCP<const Teuchos::ParameterList> validPL;
433  if (is_null(validPL)) {
434  Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
435 
436  // The static_cast is to resolve an issue with older clang versions which
437  // would cause the constexpr to link fail. With c++17 the problem is resolved.
438  pl->set("Convergence Tolerance", static_cast<MagnitudeType>(DefaultSolverParameters::convTol),
439  "The relative residual tolerance that needs to be achieved by the\n"
440  "iterative solver in order for the linear system to be declared converged." );
441  pl->set("Maximum Restarts", static_cast<int>(maxRestarts_default_),
442  "The maximum number of restarts allowed for each\n"
443  "set of RHS solved.");
444  pl->set("Maximum Iterations", static_cast<int>(maxIters_default_),
445  "The maximum number of block iterations allowed for each\n"
446  "set of RHS solved.");
447  pl->set("Num Blocks", static_cast<int>(numBlocks_default_),
448  "The maximum number of blocks allowed in the Krylov subspace\n"
449  "for each set of RHS solved.");
450  pl->set("Block Size", static_cast<int>(blockSize_default_),
451  "The number of vectors in each block. This number times the\n"
452  "number of blocks is the total Krylov subspace dimension.");
453  pl->set("Adaptive Block Size", static_cast<bool>(adaptiveBlockSize_default_),
454  "Whether the solver manager should adapt the block size\n"
455  "based on the number of RHS to solve.");
456  pl->set("Verbosity", static_cast<int>(verbosity_default_),
457  "What type(s) of solver information should be outputted\n"
458  "to the output stream.");
459  pl->set("Output Style", static_cast<int>(outputStyle_default_),
460  "What style is used for the solver information outputted\n"
461  "to the output stream.");
462  pl->set("Output Frequency", static_cast<int>(outputFreq_default_),
463  "How often convergence information should be outputted\n"
464  "to the output stream.");
465  pl->set("Output Stream", Teuchos::rcp(outputStream_default_,false),
466  "A reference-counted pointer to the output stream where all\n"
467  "solver output is sent.");
468  pl->set("Show Maximum Residual Norm Only", static_cast<bool>(showMaxResNormOnly_default_),
469  "When convergence information is printed, only show the maximum\n"
470  "relative residual norm when the block size is greater than one.");
471  pl->set("Flexible Gmres", static_cast<bool>(flexibleGmres_default_),
472  "Whether the solver manager should use the flexible variant\n"
473  "of GMRES.");
474  pl->set("Explicit Residual Test", static_cast<bool>(expResTest_default_),
475  "Whether the explicitly computed residual should be used in the convergence test.");
476  pl->set("Implicit Residual Scaling", static_cast<const char *>(impResScale_default_),
477  "The type of scaling used in the implicit residual convergence test.");
478  pl->set("Explicit Residual Scaling", static_cast<const char *>(expResScale_default_),
479  "The type of scaling used in the explicit residual convergence test.");
480  pl->set("Timer Label", static_cast<const char *>(label_default_),
481  "The string to use as a prefix for the timer labels.");
482  pl->set("Orthogonalization", static_cast<const char *>(orthoType_default_),
483  "The type of orthogonalization to use: DGKS, ICGS, or IMGS.");
484  pl->set("Orthogonalization Constant",static_cast<MagnitudeType>(DefaultSolverParameters::orthoKappa),
485  "The constant used by DGKS orthogonalization to determine\n"
486  "whether another step of classical Gram-Schmidt is necessary.");
487  validPL = pl;
488  }
489  return validPL;
490 }
491 
492 
493 template<class ScalarType, class MV, class OP>
494 void BlockGmresSolMgr<ScalarType,MV,OP>::setParameters( const Teuchos::RCP<Teuchos::ParameterList> &params )
495 {
496 
497  // Create the internal parameter list if ones doesn't already exist.
498  if (params_ == Teuchos::null) {
499  params_ = Teuchos::rcp( new Teuchos::ParameterList(*getValidParameters()) );
500  }
501  else {
502  params->validateParameters(*getValidParameters());
503  }
504 
505  // Check for maximum number of restarts
506  if (params->isParameter("Maximum Restarts")) {
507  maxRestarts_ = params->get("Maximum Restarts",maxRestarts_default_);
508 
509  // Update parameter in our list.
510  params_->set("Maximum Restarts", maxRestarts_);
511  }
512 
513  // Check for maximum number of iterations
514  if (params->isParameter("Maximum Iterations")) {
515  maxIters_ = params->get("Maximum Iterations",maxIters_default_);
516 
517  // Update parameter in our list and in status test.
518  params_->set("Maximum Iterations", maxIters_);
519  if (maxIterTest_!=Teuchos::null)
520  maxIterTest_->setMaxIters( maxIters_ );
521  }
522 
523  // Check for blocksize
524  if (params->isParameter("Block Size")) {
525  blockSize_ = params->get("Block Size",blockSize_default_);
526  TEUCHOS_TEST_FOR_EXCEPTION(blockSize_ <= 0, std::invalid_argument,
527  "Belos::BlockGmresSolMgr: \"Block Size\" must be strictly positive.");
528 
529  // Update parameter in our list.
530  params_->set("Block Size", blockSize_);
531  }
532 
533  // Check if the blocksize should be adaptive
534  if (params->isParameter("Adaptive Block Size")) {
535  adaptiveBlockSize_ = params->get("Adaptive Block Size",adaptiveBlockSize_default_);
536 
537  // Update parameter in our list.
538  params_->set("Adaptive Block Size", adaptiveBlockSize_);
539  }
540 
541  // Check for the maximum number of blocks.
542  if (params->isParameter("Num Blocks")) {
543  numBlocks_ = params->get("Num Blocks",numBlocks_default_);
544  TEUCHOS_TEST_FOR_EXCEPTION(numBlocks_ <= 0, std::invalid_argument,
545  "Belos::BlockGmresSolMgr: \"Num Blocks\" must be strictly positive.");
546 
547  // Update parameter in our list.
548  params_->set("Num Blocks", numBlocks_);
549  }
550 
551  // Check to see if the timer label changed.
552  if (params->isParameter("Timer Label")) {
553  std::string tempLabel = params->get("Timer Label", label_default_);
554 
555  // Update parameter in our list, solver timer, and orthogonalization label
556  if (tempLabel != label_) {
557  label_ = tempLabel;
558  params_->set("Timer Label", label_);
559  std::string solveLabel = label_ + ": BlockGmresSolMgr total solve time";
560 #ifdef BELOS_TEUCHOS_TIME_MONITOR
561  timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
562 #endif
563  if (ortho_ != Teuchos::null) {
564  ortho_->setLabel( label_ );
565  }
566  }
567  }
568 
569  // Determine whether this solver should be "flexible".
570  if (params->isParameter("Flexible Gmres")) {
571  isFlexible_ = Teuchos::getParameter<bool>(*params,"Flexible Gmres");
572  params_->set("Flexible Gmres", isFlexible_);
573  if (isFlexible_ && expResTest_) {
574  // Use an implicit convergence test if the Gmres solver is flexible
575  isSTSet_ = false;
576  }
577  }
578 
579 
580  // Check if the orthogonalization changed.
581  if (params->isParameter("Orthogonalization")) {
582  std::string tempOrthoType = params->get("Orthogonalization",orthoType_default_);
583  TEUCHOS_TEST_FOR_EXCEPTION( tempOrthoType != "DGKS" && tempOrthoType != "ICGS" && tempOrthoType != "IMGS",
584  std::invalid_argument,
585  "Belos::BlockGmresSolMgr: \"Orthogonalization\" must be either \"DGKS\", \"ICGS\", or \"IMGS\".");
586  if (tempOrthoType != orthoType_) {
587  orthoType_ = tempOrthoType;
588  params_->set("Orthogonalization", orthoType_);
589  // Create orthogonalization manager
590  if (orthoType_=="DGKS") {
591  if (orthoKappa_ <= 0) {
592  ortho_ = Teuchos::rcp( new DGKSOrthoManager<ScalarType,MV,OP>( label_ ) );
593  }
594  else {
595  ortho_ = Teuchos::rcp( new DGKSOrthoManager<ScalarType,MV,OP>( label_ ) );
596  Teuchos::rcp_dynamic_cast<DGKSOrthoManager<ScalarType,MV,OP> >(ortho_)->setDepTol( orthoKappa_ );
597  }
598  }
599  else if (orthoType_=="ICGS") {
600  ortho_ = Teuchos::rcp( new ICGSOrthoManager<ScalarType,MV,OP>( label_ ) );
601  }
602  else if (orthoType_=="IMGS") {
603  ortho_ = Teuchos::rcp( new IMGSOrthoManager<ScalarType,MV,OP>( label_ ) );
604  }
605  }
606  }
607 
608  // Check which orthogonalization constant to use.
609  if (params->isParameter("Orthogonalization Constant")) {
610  if (params->isType<MagnitudeType> ("Orthogonalization Constant")) {
611  orthoKappa_ = params->get ("Orthogonalization Constant",
612  static_cast<MagnitudeType> (DefaultSolverParameters::orthoKappa));
613  }
614  else {
615  orthoKappa_ = params->get ("Orthogonalization Constant",
617  }
618 
619  // Update parameter in our list.
620  params_->set("Orthogonalization Constant",orthoKappa_);
621  if (orthoType_=="DGKS") {
622  if (orthoKappa_ > 0 && ortho_ != Teuchos::null) {
623  Teuchos::rcp_dynamic_cast<DGKSOrthoManager<ScalarType,MV,OP> >(ortho_)->setDepTol( orthoKappa_ );
624  }
625  }
626  }
627 
628  // Check for a change in verbosity level
629  if (params->isParameter("Verbosity")) {
630  if (Teuchos::isParameterType<int>(*params,"Verbosity")) {
631  verbosity_ = params->get("Verbosity", verbosity_default_);
632  } else {
633  verbosity_ = (int)Teuchos::getParameter<Belos::MsgType>(*params,"Verbosity");
634  }
635 
636  // Update parameter in our list.
637  params_->set("Verbosity", verbosity_);
638  if (printer_ != Teuchos::null)
639  printer_->setVerbosity(verbosity_);
640  }
641 
642  // Check for a change in output style
643  if (params->isParameter("Output Style")) {
644  if (Teuchos::isParameterType<int>(*params,"Output Style")) {
645  outputStyle_ = params->get("Output Style", outputStyle_default_);
646  } else {
647  outputStyle_ = (int)Teuchos::getParameter<Belos::OutputType>(*params,"Output Style");
648  }
649 
650  // Reconstruct the convergence test if the explicit residual test is not being used.
651  params_->set("Output Style", outputStyle_);
652  if (outputTest_ != Teuchos::null) {
653  isSTSet_ = false;
654  }
655  }
656 
657  // output stream
658  if (params->isParameter("Output Stream")) {
659  outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> >(*params,"Output Stream");
660 
661  // Update parameter in our list.
662  params_->set("Output Stream", outputStream_);
663  if (printer_ != Teuchos::null)
664  printer_->setOStream( outputStream_ );
665  }
666 
667  // frequency level
668  if (verbosity_ & Belos::StatusTestDetails) {
669  if (params->isParameter("Output Frequency")) {
670  outputFreq_ = params->get("Output Frequency", outputFreq_default_);
671  }
672 
673  // Update parameter in out list and output status test.
674  params_->set("Output Frequency", outputFreq_);
675  if (outputTest_ != Teuchos::null)
676  outputTest_->setOutputFrequency( outputFreq_ );
677  }
678 
679  // Create output manager if we need to.
680  if (printer_ == Teuchos::null) {
681  printer_ = Teuchos::rcp( new OutputManager<ScalarType>(verbosity_, outputStream_) );
682  }
683 
684  // Check for convergence tolerance
685  if (params->isParameter("Convergence Tolerance")) {
686  if (params->isType<MagnitudeType> ("Convergence Tolerance")) {
687  convtol_ = params->get ("Convergence Tolerance",
688  static_cast<MagnitudeType> (DefaultSolverParameters::convTol));
689  }
690  else {
691  convtol_ = params->get ("Convergence Tolerance", DefaultSolverParameters::convTol);
692  }
693 
694  // Update parameter in our list and residual tests.
695  params_->set("Convergence Tolerance", convtol_);
696  if (impConvTest_ != Teuchos::null)
697  impConvTest_->setTolerance( convtol_ );
698  if (expConvTest_ != Teuchos::null)
699  expConvTest_->setTolerance( convtol_ );
700  }
701 
702  // Check for a change in scaling, if so we need to build new residual tests.
703  if (params->isParameter("Implicit Residual Scaling")) {
704  std::string tempImpResScale = Teuchos::getParameter<std::string>( *params, "Implicit Residual Scaling" );
705 
706  // Only update the scaling if it's different.
707  if (impResScale_ != tempImpResScale) {
708  Belos::ScaleType impResScaleType = convertStringToScaleType( tempImpResScale );
709  impResScale_ = tempImpResScale;
710 
711  // Update parameter in our list and residual tests
712  params_->set("Implicit Residual Scaling", impResScale_);
713  if (impConvTest_ != Teuchos::null) {
714  try {
715  impConvTest_->defineScaleForm( impResScaleType, Belos::TwoNorm );
716  }
717  catch (std::exception& e) {
718  // Make sure the convergence test gets constructed again.
719  isSTSet_ = false;
720  }
721  }
722  }
723  }
724 
725  if (params->isParameter("Explicit Residual Scaling")) {
726  std::string tempExpResScale = Teuchos::getParameter<std::string>( *params, "Explicit Residual Scaling" );
727 
728  // Only update the scaling if it's different.
729  if (expResScale_ != tempExpResScale) {
730  Belos::ScaleType expResScaleType = convertStringToScaleType( tempExpResScale );
731  expResScale_ = tempExpResScale;
732 
733  // Update parameter in our list and residual tests
734  params_->set("Explicit Residual Scaling", expResScale_);
735  if (expConvTest_ != Teuchos::null) {
736  try {
737  expConvTest_->defineScaleForm( expResScaleType, Belos::TwoNorm );
738  }
739  catch (std::exception& e) {
740  // Make sure the convergence test gets constructed again.
741  isSTSet_ = false;
742  }
743  }
744  }
745  }
746 
747  if (params->isParameter("Explicit Residual Test")) {
748  expResTest_ = Teuchos::getParameter<bool>( *params,"Explicit Residual Test" );
749 
750  // Reconstruct the convergence test if the explicit residual test is not being used.
751  params_->set("Explicit Residual Test", expResTest_);
752  if (expConvTest_ == Teuchos::null) {
753  isSTSet_ = false;
754  }
755  }
756 
757 
758  if (params->isParameter("Show Maximum Residual Norm Only")) {
759  showMaxResNormOnly_ = Teuchos::getParameter<bool>(*params,"Show Maximum Residual Norm Only");
760 
761  // Update parameter in our list and residual tests
762  params_->set("Show Maximum Residual Norm Only", showMaxResNormOnly_);
763  if (impConvTest_ != Teuchos::null)
764  impConvTest_->setShowMaxResNormOnly( showMaxResNormOnly_ );
765  if (expConvTest_ != Teuchos::null)
766  expConvTest_->setShowMaxResNormOnly( showMaxResNormOnly_ );
767  }
768 
769  // Create orthogonalization manager if we need to.
770  if (ortho_ == Teuchos::null) {
771  params_->set("Orthogonalization", orthoType_);
772  if (orthoType_=="DGKS") {
773  if (orthoKappa_ <= 0) {
774  ortho_ = Teuchos::rcp( new DGKSOrthoManager<ScalarType,MV,OP>( label_ ) );
775  }
776  else {
777  ortho_ = Teuchos::rcp( new DGKSOrthoManager<ScalarType,MV,OP>( label_ ) );
778  Teuchos::rcp_dynamic_cast<DGKSOrthoManager<ScalarType,MV,OP> >(ortho_)->setDepTol( orthoKappa_ );
779  }
780  }
781  else if (orthoType_=="ICGS") {
782  ortho_ = Teuchos::rcp( new ICGSOrthoManager<ScalarType,MV,OP>( label_ ) );
783  }
784  else if (orthoType_=="IMGS") {
785  ortho_ = Teuchos::rcp( new IMGSOrthoManager<ScalarType,MV,OP>( label_ ) );
786  }
787  else {
788  TEUCHOS_TEST_FOR_EXCEPTION(orthoType_!="ICGS"&&orthoType_!="DGKS"&&orthoType_!="IMGS",std::logic_error,
789  "Belos::BlockGmresSolMgr(): Invalid orthogonalization type.");
790  }
791  }
792 
793  // Create the timer if we need to.
794  if (timerSolve_ == Teuchos::null) {
795  std::string solveLabel = label_ + ": BlockGmresSolMgr total solve time";
796 #ifdef BELOS_TEUCHOS_TIME_MONITOR
797  timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
798 #endif
799  }
800 
801  // Inform the solver manager that the current parameters were set.
802  isSet_ = true;
803 }
804 
805 // Check the status test versus the defined linear problem
806 template<class ScalarType, class MV, class OP>
808 
809  typedef Belos::StatusTestCombo<ScalarType,MV,OP> StatusTestCombo_t;
810  typedef Belos::StatusTestGenResNorm<ScalarType,MV,OP> StatusTestGenResNorm_t;
811  typedef Belos::StatusTestImpResNorm<ScalarType,MV,OP> StatusTestImpResNorm_t;
812 
813  // Basic test checks maximum iterations and native residual.
814  maxIterTest_ = Teuchos::rcp( new StatusTestMaxIters<ScalarType,MV,OP>( maxIters_ ) );
815 
816  // If there is a left preconditioner, we create a combined status test that checks the implicit
817  // and then explicit residual norm to see if we have convergence.
818  if (!Teuchos::is_null(problem_->getLeftPrec()) && !isFlexible_) {
819  expResTest_ = true;
820  }
821 
822  if (expResTest_) {
823 
824  // Implicit residual test, using the native residual to determine if convergence was achieved.
825  Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
826  Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
827  tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
828  tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
829  impConvTest_ = tmpImpConvTest;
830 
831  // Explicit residual test once the native residual is below the tolerance
832  Teuchos::RCP<StatusTestGenResNorm_t> tmpExpConvTest =
833  Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
834  tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit, Belos::TwoNorm );
835  tmpExpConvTest->defineScaleForm( convertStringToScaleType(expResScale_), Belos::TwoNorm );
836  tmpExpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
837  expConvTest_ = tmpExpConvTest;
838 
839  // The convergence test is a combination of the "cheap" implicit test and explicit test.
840  convTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
841  }
842  else {
843 
844  if (isFlexible_) {
845  // Implicit residual test, using the native residual to determine if convergence was achieved.
846  Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
847  Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
848  tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
849  tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
850  impConvTest_ = tmpImpConvTest;
851  }
852  else {
853  // Implicit residual test, using the native residual to determine if convergence was achieved.
854  // Use test that checks for loss of accuracy.
855  Teuchos::RCP<StatusTestImpResNorm_t> tmpImpConvTest =
856  Teuchos::rcp( new StatusTestImpResNorm_t( convtol_ ) );
857  tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
858  tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
859  impConvTest_ = tmpImpConvTest;
860  }
861 
862  // Set the explicit and total convergence test to this implicit test that checks for accuracy loss.
863  expConvTest_ = impConvTest_;
864  convTest_ = impConvTest_;
865  }
866 
867  // Create the status test.
868  sTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
869 
870  // Add debug status test, if one is provided by the user
871  if (nonnull(debugStatusTest_) ) {
872  // Add debug convergence test
873  Teuchos::rcp_dynamic_cast<StatusTestCombo_t>(sTest_)->addStatusTest( debugStatusTest_ );
874  }
875 
876  // Create the status test output class.
877  // This class manages and formats the output from the status test.
878  StatusTestOutputFactory<ScalarType,MV,OP> stoFactory( outputStyle_ );
879  outputTest_ = stoFactory.create( printer_, sTest_, outputFreq_, Passed+Failed+Undefined );
880 
881  // Set the solver string for the output test
882  std::string solverDesc = " Block Gmres ";
883  if (isFlexible_)
884  solverDesc = "Flexible" + solverDesc;
885  outputTest_->setSolverDesc( solverDesc );
886 
887  // The status test is now set.
888  isSTSet_ = true;
889 
890  return false;
891 }
892 
893 template<class ScalarType, class MV, class OP>
895  const Teuchos::RCP<StatusTest<ScalarType,MV,OP> > &debugStatusTest
896  )
897 {
898  debugStatusTest_ = debugStatusTest;
899 }
900 
901 
902 // solve()
903 template<class ScalarType, class MV, class OP>
905 
906  // Set the current parameters if they were not set before.
907  // NOTE: This may occur if the user generated the solver manager with the default constructor and
908  // then didn't set any parameters using setParameters().
909  if (!isSet_) {
910  setParameters(Teuchos::parameterList(*getValidParameters()));
911  }
912 
913  Teuchos::BLAS<int,ScalarType> blas;
914 
915  TEUCHOS_TEST_FOR_EXCEPTION(problem_ == Teuchos::null,BlockGmresSolMgrLinearProblemFailure,
916  "Belos::BlockGmresSolMgr::solve(): Linear problem is not a valid object.");
917 
918  TEUCHOS_TEST_FOR_EXCEPTION(!problem_->isProblemSet(),BlockGmresSolMgrLinearProblemFailure,
919  "Belos::BlockGmresSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
920 
921  if (isFlexible_) {
922  TEUCHOS_TEST_FOR_EXCEPTION(problem_->getRightPrec()==Teuchos::null,BlockGmresSolMgrLinearProblemFailure,
923  "Belos::BlockGmresSolMgr::solve(): Linear problem does not have a preconditioner required for flexible GMRES, call setRightPrec().");
924  }
925 
926  if (!isSTSet_ || (!expResTest_ && !Teuchos::is_null(problem_->getLeftPrec())) ) {
927  TEUCHOS_TEST_FOR_EXCEPTION( checkStatusTest(),BlockGmresSolMgrLinearProblemFailure,
928  "Belos::BlockGmresSolMgr::solve(): Linear problem and requested status tests are incompatible.");
929  }
930 
931  // Create indices for the linear systems to be solved.
932  int startPtr = 0;
933  int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
934  int numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
935 
936  std::vector<int> currIdx;
937  // If an adaptive block size is allowed then only the linear systems that need to be solved are solved.
938  // Otherwise, the index set is generated that informs the linear problem that some linear systems are augmented.
939  if ( adaptiveBlockSize_ ) {
940  blockSize_ = numCurrRHS;
941  currIdx.resize( numCurrRHS );
942  for (int i=0; i<numCurrRHS; ++i)
943  { currIdx[i] = startPtr+i; }
944 
945  }
946  else {
947  currIdx.resize( blockSize_ );
948  for (int i=0; i<numCurrRHS; ++i)
949  { currIdx[i] = startPtr+i; }
950  for (int i=numCurrRHS; i<blockSize_; ++i)
951  { currIdx[i] = -1; }
952  }
953 
954  // Inform the linear problem of the current linear system to solve.
955  problem_->setLSIndex( currIdx );
956 
958  // Parameter list
959  Teuchos::ParameterList plist;
960  plist.set("Block Size",blockSize_);
961 
962  ptrdiff_t dim = MVT::GetGlobalLength( *(problem_->getRHS()) );
963  if (blockSize_*static_cast<ptrdiff_t>(numBlocks_) > dim) {
964  int tmpNumBlocks = 0;
965  if (blockSize_ == 1)
966  tmpNumBlocks = dim / blockSize_; // Allow for a good breakdown.
967  else
968  tmpNumBlocks = ( dim - blockSize_) / blockSize_; // Allow for restarting.
969  printer_->stream(Warnings) <<
970  "Belos::BlockGmresSolMgr::solve(): Warning! Requested Krylov subspace dimension is larger than operator dimension!"
971  << std::endl << " The maximum number of blocks allowed for the Krylov subspace will be adjusted to " << tmpNumBlocks << std::endl;
972  plist.set("Num Blocks",tmpNumBlocks);
973  }
974  else
975  plist.set("Num Blocks",numBlocks_);
976 
977  // Reset the status test.
978  outputTest_->reset();
979  loaDetected_ = false;
980 
981  // Assume convergence is achieved, then let any failed convergence set this to false.
982  bool isConverged = true;
983 
985  // BlockGmres solver
986 
987  Teuchos::RCP<GmresIteration<ScalarType,MV,OP> > block_gmres_iter;
988 
989  if (isFlexible_)
990  block_gmres_iter = Teuchos::rcp( new BlockFGmresIter<ScalarType,MV,OP>(problem_,printer_,outputTest_,ortho_,plist) );
991  else
992  block_gmres_iter = Teuchos::rcp( new BlockGmresIter<ScalarType,MV,OP>(problem_,printer_,outputTest_,ortho_,plist) );
993 
994  // Enter solve() iterations
995  {
996 #ifdef BELOS_TEUCHOS_TIME_MONITOR
997  Teuchos::TimeMonitor slvtimer(*timerSolve_);
998 #endif
999 
1000  while ( numRHS2Solve > 0 ) {
1001 
1002  // Set the current number of blocks and blocksize with the Gmres iteration.
1003  if (blockSize_*numBlocks_ > dim) {
1004  int tmpNumBlocks = 0;
1005  if (blockSize_ == 1)
1006  tmpNumBlocks = dim / blockSize_; // Allow for a good breakdown.
1007  else
1008  tmpNumBlocks = ( dim - blockSize_) / blockSize_; // Allow for restarting.
1009  block_gmres_iter->setSize( blockSize_, tmpNumBlocks );
1010  }
1011  else
1012  block_gmres_iter->setSize( blockSize_, numBlocks_ );
1013 
1014  // Reset the number of iterations.
1015  block_gmres_iter->resetNumIters();
1016 
1017  // Reset the number of calls that the status test output knows about.
1018  outputTest_->resetNumCalls();
1019 
1020  // Create the first block in the current Krylov basis.
1021  Teuchos::RCP<MV> V_0;
1022  if (isFlexible_) {
1023  // Load the correct residual if the system is augmented
1024  if (currIdx[blockSize_-1] == -1) {
1025  V_0 = MVT::Clone( *(problem_->getInitResVec()), blockSize_ );
1026  problem_->computeCurrResVec( &*V_0 );
1027  }
1028  else {
1029  V_0 = MVT::CloneCopy( *(problem_->getInitResVec()), currIdx );
1030  }
1031  }
1032  else {
1033  // Load the correct residual if the system is augmented
1034  if (currIdx[blockSize_-1] == -1) {
1035  V_0 = MVT::Clone( *(problem_->getInitPrecResVec()), blockSize_ );
1036  problem_->computeCurrPrecResVec( &*V_0 );
1037  }
1038  else {
1039  V_0 = MVT::CloneCopy( *(problem_->getInitPrecResVec()), currIdx );
1040  }
1041  }
1042 
1043  // Get a matrix to hold the orthonormalization coefficients.
1044  Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > z_0 =
1045  Teuchos::rcp( new Teuchos::SerialDenseMatrix<int,ScalarType>( blockSize_, blockSize_ ) );
1046 
1047  // Orthonormalize the new V_0
1048  int rank = ortho_->normalize( *V_0, z_0 );
1049  TEUCHOS_TEST_FOR_EXCEPTION(rank != blockSize_,BlockGmresSolMgrOrthoFailure,
1050  "Belos::BlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors.");
1051 
1052  // Set the new state and initialize the solver.
1054  newstate.V = V_0;
1055  newstate.z = z_0;
1056  newstate.curDim = 0;
1057  block_gmres_iter->initializeGmres(newstate);
1058  int numRestarts = 0;
1059 
1060  while(1) {
1061  // tell block_gmres_iter to iterate
1062  try {
1063  block_gmres_iter->iterate();
1064 
1066  //
1067  // check convergence first
1068  //
1070  if ( convTest_->getStatus() == Passed ) {
1071  if ( expConvTest_->getLOADetected() ) {
1072  // we don't have convergence
1073  loaDetected_ = true;
1074  printer_->stream(Warnings) <<
1075  "Belos::BlockGmresSolMgr::solve(): Warning! Solver has experienced a loss of accuracy!" << std::endl;
1076  isConverged = false;
1077  }
1078  break; // break from while(1){block_gmres_iter->iterate()}
1079  }
1081  //
1082  // check for maximum iterations
1083  //
1085  else if ( maxIterTest_->getStatus() == Passed ) {
1086  // we don't have convergence
1087  isConverged = false;
1088  break; // break from while(1){block_gmres_iter->iterate()}
1089  }
1091  //
1092  // check for restarting, i.e. the subspace is full
1093  //
1095  else if ( block_gmres_iter->getCurSubspaceDim() == block_gmres_iter->getMaxSubspaceDim() ) {
1096 
1097  if ( numRestarts >= maxRestarts_ ) {
1098  isConverged = false;
1099  break; // break from while(1){block_gmres_iter->iterate()}
1100  }
1101  numRestarts++;
1102 
1103  printer_->stream(Debug) << " Performing restart number " << numRestarts << " of " << maxRestarts_ << std::endl << std::endl;
1104 
1105  // Update the linear problem.
1106  Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1107  if (isFlexible_) {
1108  // Update the solution manually, since the preconditioning doesn't need to be undone.
1109  Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
1110  MVT::MvAddMv( 1.0, *curX, 1.0, *update, *curX );
1111  }
1112  else
1113  problem_->updateSolution( update, true );
1114 
1115  // Get the state.
1116  GmresIterationState<ScalarType,MV> oldState = block_gmres_iter->getState();
1117 
1118  // Compute the restart std::vector.
1119  // Get a view of the current Krylov basis.
1120  V_0 = MVT::Clone( *(oldState.V), blockSize_ );
1121  if (isFlexible_)
1122  problem_->computeCurrResVec( &*V_0 );
1123  else
1124  problem_->computeCurrPrecResVec( &*V_0 );
1125 
1126  // Get a view of the first block of the Krylov basis.
1127  z_0 = Teuchos::rcp( new Teuchos::SerialDenseMatrix<int,ScalarType>( blockSize_, blockSize_ ) );
1128 
1129  // Orthonormalize the new V_0
1130  rank = ortho_->normalize( *V_0, z_0 );
1131  TEUCHOS_TEST_FOR_EXCEPTION(rank != blockSize_,BlockGmresSolMgrOrthoFailure,
1132  "Belos::BlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors after restart.");
1133 
1134  // Set the new state and initialize the solver.
1135  newstate.V = V_0;
1136  newstate.z = z_0;
1137  newstate.curDim = 0;
1138  block_gmres_iter->initializeGmres(newstate);
1139 
1140  } // end of restarting
1141 
1143  //
1144  // we returned from iterate(), but none of our status tests Passed.
1145  // something is wrong, and it is probably our fault.
1146  //
1148 
1149  else {
1150  TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,
1151  "Belos::BlockGmresSolMgr::solve(): Invalid return from BlockGmresIter::iterate().");
1152  }
1153  }
1154  catch (const GmresIterationOrthoFailure &e) {
1155  // If the block size is not one, it's not considered a lucky breakdown.
1156  if (blockSize_ != 1) {
1157  printer_->stream(Errors) << "Error! Caught std::exception in BlockGmresIter::iterate() at iteration "
1158  << block_gmres_iter->getNumIters() << std::endl
1159  << e.what() << std::endl;
1160  if (convTest_->getStatus() != Passed)
1161  isConverged = false;
1162  break;
1163  }
1164  else {
1165  // If the block size is one, try to recover the most recent least-squares solution
1166  block_gmres_iter->updateLSQR( block_gmres_iter->getCurSubspaceDim() );
1167 
1168  // Check to see if the most recent least-squares solution yielded convergence.
1169  sTest_->checkStatus( &*block_gmres_iter );
1170  if (convTest_->getStatus() != Passed)
1171  isConverged = false;
1172  break;
1173  }
1174  }
1175  catch (const std::exception &e) {
1176  printer_->stream(Errors) << "Error! Caught std::exception in BlockGmresIter::iterate() at iteration "
1177  << block_gmres_iter->getNumIters() << std::endl
1178  << e.what() << std::endl;
1179  throw;
1180  }
1181  }
1182 
1183  // Compute the current solution.
1184  // Update the linear problem.
1185  if (isFlexible_) {
1186  // Update the solution manually, since the preconditioning doesn't need to be undone.
1187  Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1188  Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
1189  // Update the solution only if there is a valid update from the iteration
1190  if (update != Teuchos::null)
1191  MVT::MvAddMv( 1.0, *curX, 1.0, *update, *curX );
1192  }
1193  else {
1194  // Attempt to get the current solution from the residual status test, if it has one.
1195  if ( !Teuchos::is_null(expConvTest_->getSolution()) ) {
1196  Teuchos::RCP<MV> newX = expConvTest_->getSolution();
1197  Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
1198  MVT::MvAddMv( 0.0, *newX, 1.0, *newX, *curX );
1199  }
1200  else {
1201  Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1202  problem_->updateSolution( update, true );
1203  }
1204  }
1205 
1206  // Inform the linear problem that we are finished with this block linear system.
1207  problem_->setCurrLS();
1208 
1209  // Update indices for the linear systems to be solved.
1210  startPtr += numCurrRHS;
1211  numRHS2Solve -= numCurrRHS;
1212  if ( numRHS2Solve > 0 ) {
1213  numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
1214 
1215  if ( adaptiveBlockSize_ ) {
1216  blockSize_ = numCurrRHS;
1217  currIdx.resize( numCurrRHS );
1218  for (int i=0; i<numCurrRHS; ++i)
1219  { currIdx[i] = startPtr+i; }
1220  }
1221  else {
1222  currIdx.resize( blockSize_ );
1223  for (int i=0; i<numCurrRHS; ++i)
1224  { currIdx[i] = startPtr+i; }
1225  for (int i=numCurrRHS; i<blockSize_; ++i)
1226  { currIdx[i] = -1; }
1227  }
1228  // Set the next indices.
1229  problem_->setLSIndex( currIdx );
1230  }
1231  else {
1232  currIdx.resize( numRHS2Solve );
1233  }
1234 
1235  }// while ( numRHS2Solve > 0 )
1236 
1237  }
1238 
1239  // print final summary
1240  sTest_->print( printer_->stream(FinalSummary) );
1241 
1242  // print timing information
1243 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1244  // Calling summarize() can be expensive, so don't call unless the
1245  // user wants to print out timing details. summarize() will do all
1246  // the work even if it's passed a "black hole" output stream.
1247  if (verbosity_ & TimingDetails)
1248  Teuchos::TimeMonitor::summarize( printer_->stream(TimingDetails) );
1249 #endif
1250 
1251  // get iteration information for this solve
1252  numIters_ = maxIterTest_->getNumIters();
1253 
1254  // Save the convergence test value ("achieved tolerance") for this
1255  // solve. This requires a bit more work than for BlockCGSolMgr,
1256  // since for this solver, convTest_ may either be a single residual
1257  // norm test, or a combination of two residual norm tests. In the
1258  // latter case, the master convergence test convTest_ is a SEQ combo
1259  // of the implicit resp. explicit tests. If the implicit test never
1260  // passes, then the explicit test won't ever be executed. This
1261  // manifests as expConvTest_->getTestValue()->size() < 1. We deal
1262  // with this case by using the values returned by
1263  // impConvTest_->getTestValue().
1264  {
1265  // We'll fetch the vector of residual norms one way or the other.
1266  const std::vector<MagnitudeType>* pTestValues = NULL;
1267  if (expResTest_) {
1268  pTestValues = expConvTest_->getTestValue();
1269  if (pTestValues == NULL || pTestValues->size() < 1) {
1270  pTestValues = impConvTest_->getTestValue();
1271  }
1272  }
1273  else {
1274  // Only the implicit residual norm test is being used.
1275  pTestValues = impConvTest_->getTestValue();
1276  }
1277  TEUCHOS_TEST_FOR_EXCEPTION(pTestValues == NULL, std::logic_error,
1278  "Belos::BlockGmresSolMgr::solve(): The implicit convergence test's "
1279  "getTestValue() method returned NULL. Please report this bug to the "
1280  "Belos developers.");
1281  TEUCHOS_TEST_FOR_EXCEPTION(pTestValues->size() < 1, std::logic_error,
1282  "Belos::BlockGmresSolMgr::solve(): The implicit convergence test's "
1283  "getTestValue() method returned a vector of length zero. Please report "
1284  "this bug to the Belos developers.");
1285 
1286  // FIXME (mfh 12 Dec 2011) Does pTestValues really contain the
1287  // achieved tolerances for all vectors in the current solve(), or
1288  // just for the vectors from the last deflation?
1289  achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
1290  }
1291 
1292  if (!isConverged || loaDetected_) {
1293  return Unconverged; // return from BlockGmresSolMgr::solve()
1294  }
1295  return Converged; // return from BlockGmresSolMgr::solve()
1296 }
1297 
1298 
1299 template<class ScalarType, class MV, class OP>
1301 {
1302  std::ostringstream out;
1303  out << "\"Belos::BlockGmresSolMgr\": {";
1304  if (this->getObjectLabel () != "") {
1305  out << "Label: " << this->getObjectLabel () << ", ";
1306  }
1307  out << "Flexible: " << (isFlexible_ ? "true" : "false")
1308  << ", Num Blocks: " << numBlocks_
1309  << ", Maximum Iterations: " << maxIters_
1310  << ", Maximum Restarts: " << maxRestarts_
1311  << ", Convergence Tolerance: " << convtol_
1312  << "}";
1313  return out.str ();
1314 }
1315 
1316 
1317 template<class ScalarType, class MV, class OP>
1318 void
1320 describe (Teuchos::FancyOStream &out,
1321  const Teuchos::EVerbosityLevel verbLevel) const
1322 {
1323  using Teuchos::TypeNameTraits;
1324  using Teuchos::VERB_DEFAULT;
1325  using Teuchos::VERB_NONE;
1326  using Teuchos::VERB_LOW;
1327  // using Teuchos::VERB_MEDIUM;
1328  // using Teuchos::VERB_HIGH;
1329  // using Teuchos::VERB_EXTREME;
1330  using std::endl;
1331 
1332  // Set default verbosity if applicable.
1333  const Teuchos::EVerbosityLevel vl =
1334  (verbLevel == VERB_DEFAULT) ? VERB_LOW : verbLevel;
1335 
1336  if (vl != VERB_NONE) {
1337  Teuchos::OSTab tab0 (out);
1338 
1339  out << "\"Belos::BlockGmresSolMgr\":" << endl;
1340  Teuchos::OSTab tab1 (out);
1341  out << "Template parameters:" << endl;
1342  {
1343  Teuchos::OSTab tab2 (out);
1344  out << "ScalarType: " << TypeNameTraits<ScalarType>::name () << endl
1345  << "MV: " << TypeNameTraits<MV>::name () << endl
1346  << "OP: " << TypeNameTraits<OP>::name () << endl;
1347  }
1348  if (this->getObjectLabel () != "") {
1349  out << "Label: " << this->getObjectLabel () << endl;
1350  }
1351  out << "Flexible: " << (isFlexible_ ? "true" : "false") << endl
1352  << "Num Blocks: " << numBlocks_ << endl
1353  << "Maximum Iterations: " << maxIters_ << endl
1354  << "Maximum Restarts: " << maxRestarts_ << endl
1355  << "Convergence Tolerance: " << convtol_ << endl;
1356  }
1357 }
1358 
1359 
1360 } // end Belos namespace
1361 
1362 #endif /* BELOS_BLOCK_GMRES_SOLMGR_HPP */
Belos::BlockGmresSolMgr::getTimers
Teuchos::Array< Teuchos::RCP< Teuchos::Time > > getTimers() const
Return the timers for this object.
Definition: BelosBlockGmresSolMgr.hpp:200
BelosConfigDefs.hpp
Belos header file which uses auto-configuration information to include necessary C++ headers.
Belos::BlockGmresSolMgr
Interface to Block GMRES and Flexible GMRES.
Definition: BelosBlockGmresSolMgr.hpp:125
Belos::ScaleType
ScaleType
The type of scaling to use on the residual norm value.
Definition: BelosTypes.hpp:119
Belos::BlockGmresSolMgr::solve
ReturnType solve() override
This method performs possibly repeated calls to the underlying linear solver's iterate() routine unti...
Definition: BelosBlockGmresSolMgr.hpp:904
Belos::BlockGmresSolMgr::description
std::string description() const override
Return a one-line description of this object.
Definition: BelosBlockGmresSolMgr.hpp:1300
Belos::BlockGmresSolMgr::isLOADetected
bool isLOADetected() const override
Return whether a loss of accuracy was detected by this solver during the most current solve.
Definition: BelosBlockGmresSolMgr.hpp:226
Belos::GmresIterationState::z
Teuchos::RCP< const Teuchos::SerialDenseMatrix< int, ScalarType > > z
The current right-hand side of the least squares system RY = Z.
Definition: BelosGmresIteration.hpp:88
Belos::GmresIterationState::V
Teuchos::RCP< const MV > V
The current Krylov basis.
Definition: BelosGmresIteration.hpp:71
Belos::DGKSOrthoManager
An implementation of the Belos::MatOrthoManager that performs orthogonalization using (potentially) m...
Definition: BelosDGKSOrthoManager.hpp:81
Belos::OperatorTraits
Class which defines basic traits for the operator type.
Definition: BelosOperatorTraits.hpp:109
Belos::DefaultSolverParameters
Default parameters common to most Belos solvers.
Definition: BelosTypes.hpp:286
Belos::Warnings
Definition: BelosTypes.hpp:259
Belos::OutputManager
Belos's basic output manager for sending information of select verbosity levels to the appropriate ou...
Definition: BelosIteration.hpp:64
BelosLinearProblem.hpp
Class which describes the linear problem to be solved by the iterative solver.
Belos::BlockGmresSolMgr::achievedTol
MagnitudeType achievedTol() const override
Tolerance achieved by the last solve() invocation.
Definition: BelosBlockGmresSolMgr.hpp:214
BelosIMGSOrthoManager.hpp
Iterated Modified Gram-Schmidt (IMGS) implementation of the Belos::OrthoManager class.
Belos::General
Definition: BelosTypes.hpp:140
Belos::BlockGmresSolMgr::getNumIters
int getNumIters() const override
Get the iteration count for the most recent call to solve().
Definition: BelosBlockGmresSolMgr.hpp:219
Belos::TimingDetails
Definition: BelosTypes.hpp:263
BelosDGKSOrthoManager.hpp
Classical Gram-Schmidt (with DGKS correction) implementation of the Belos::OrthoManager class.
Belos::StatusTestMaxIters
A Belos::StatusTest class for specifying a maximum number of iterations.
Definition: BelosStatusTestMaxIters.hpp:63
BelosSolverManager.hpp
Pure virtual base class which describes the basic interface for a solver manager.
Belos::Passed
Definition: BelosTypes.hpp:188
Belos::StatusTestDetails
Definition: BelosTypes.hpp:264
Belos::TwoNorm
Definition: BelosTypes.hpp:98
Belos::FinalSummary
Definition: BelosTypes.hpp:262
BelosStatusTestOutput.hpp
Virtual base class for StatusTest that printing status tests.
Belos::LinearProblem
A linear system to solve, and its associated information.
Definition: BelosIteration.hpp:61
BelosStatusTestGenResNorm.hpp
Belos::StatusTestResNorm for specifying general residual norm stopping criteria.
BelosStatusTestImpResNorm.hpp
Belos::StatusTest for specifying an implicit residual norm stopping criteria that checks for loss of ...
BelosBlockFGmresIter.hpp
Belos concrete class for performing the block, flexible GMRES iteration.
Belos::BlockGmresSolMgrOrthoFailure
BlockGmresSolMgrOrthoFailure is thrown when the orthogonalization manager is unable to generate ortho...
Definition: BelosBlockGmresSolMgr.hpp:104
Belos::DefaultSolverParameters::orthoKappa
static constexpr double orthoKappa
DGKS orthogonalization constant.
Definition: BelosTypes.hpp:300
Belos::BlockGmresSolMgr::getCurrentParameters
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const override
Get a parameter list containing the current parameters for this object.
Definition: BelosBlockGmresSolMgr.hpp:193
Belos::BlockGmresSolMgr::setDebugStatusTest
void setDebugStatusTest(const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &debugStatusTest) override
Set a debug status test that will be checked at the same time as the top-level status test.
Definition: BelosBlockGmresSolMgr.hpp:894
Belos::BlockFGmresIter
This class implements the block flexible GMRES iteration, where a block Krylov subspace is constructe...
Definition: BelosBlockFGmresIter.hpp:83
Belos::GmresIterationOrthoFailure
GmresIterationOrthoFailure is thrown when the GmresIteration object is unable to compute independent ...
Definition: BelosGmresIteration.hpp:123
BelosOutputManager.hpp
Class which manages the output and verbosity of the Belos solvers.
Belos
Definition: Belos_Details_EBelosSolverType.cpp:45
BelosStatusTestMaxIters.hpp
Belos::StatusTest class for specifying a maximum number of iterations.
Belos::BlockGmresSolMgrLinearProblemFailure
BlockGmresSolMgrLinearProblemFailure is thrown when the linear problem is not setup (i....
Definition: BelosBlockGmresSolMgr.hpp:94
Belos::Undefined
Definition: BelosTypes.hpp:190
Belos::BlockGmresSolMgr::BlockGmresSolMgr
BlockGmresSolMgr()
Empty constructor for BlockGmresSolMgr. This constructor takes no arguments and sets the default valu...
Definition: BelosBlockGmresSolMgr.hpp:360
Belos::StatusTestImpResNorm
Convergence test using the implicit residual norm(s), with an explicit residual norm(s) check for los...
Definition: BelosStatusTestImpResNorm.hpp:104
Belos::BlockGmresSolMgr::reset
void reset(const ResetType type) override
Performs a reset of the solver manager specified by the ResetType. This informs the solver manager th...
Definition: BelosBlockGmresSolMgr.hpp:250
Belos::BlockGmresSolMgr::setProblem
void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem) override
Set the linear problem that needs to be solved.
Definition: BelosBlockGmresSolMgr.hpp:234
Belos::BlockGmresSolMgrOrthoFailure::BlockGmresSolMgrOrthoFailure
BlockGmresSolMgrOrthoFailure(const std::string &what_arg)
Definition: BelosBlockGmresSolMgr.hpp:105
BelosICGSOrthoManager.hpp
Iterated Classical Gram-Schmidt (ICGS) implementation of the Belos::OrthoManager class.
Belos::BlockGmresSolMgr::getProblem
const LinearProblem< ScalarType, MV, OP > & getProblem() const override
Get current linear problem being solved for in this object.
Definition: BelosBlockGmresSolMgr.hpp:183
Belos::StatusTestGenResNorm
An implementation of StatusTestResNorm using a family of residual norms.
Definition: BelosStatusTestGenResNorm.hpp:79
Belos::Problem
Definition: BelosTypes.hpp:205
Belos::BlockGmresSolMgrLinearProblemFailure::BlockGmresSolMgrLinearProblemFailure
BlockGmresSolMgrLinearProblemFailure(const std::string &what_arg)
Definition: BelosBlockGmresSolMgr.hpp:95
Belos::ReturnType
ReturnType
Whether the Belos solve converged for all linear systems.
Definition: BelosTypes.hpp:154
Belos::Unconverged
Definition: BelosTypes.hpp:156
Belos::IMGSOrthoManager
An implementation of the Belos::MatOrthoManager that performs orthogonalization using multiple steps ...
Definition: BelosIMGSOrthoManager.hpp:83
Belos::BlockGmresSolMgr::~BlockGmresSolMgr
virtual ~BlockGmresSolMgr()
Destructor.
Definition: BelosBlockGmresSolMgr.hpp:170
Belos::ICGSOrthoManager
An implementation of the Belos::MatOrthoManager that performs orthogonalization using multiple steps ...
Definition: BelosICGSOrthoManager.hpp:81
Belos::GmresIterationState::curDim
int curDim
The current dimension of the reduction.
Definition: BelosGmresIteration.hpp:68
Belos::BlockGmresSolMgr::getValidParameters
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
Get a parameter list containing the valid parameters for this object.
Definition: BelosBlockGmresSolMgr.hpp:430
Belos::BlockGmresSolMgr::setParameters
void setParameters(const Teuchos::RCP< Teuchos::ParameterList > &params) override
Set the parameters the solver manager should use to solve the linear problem.
Definition: BelosBlockGmresSolMgr.hpp:494
Belos::convertStringToScaleType
ScaleType convertStringToScaleType(const std::string &scaleType)
Convert the given string to its ScaleType enum value.
Definition: BelosTypes.cpp:88
Belos::ResetType
ResetType
How to reset the solver.
Definition: BelosTypes.hpp:205
Belos::StatusTest
A pure virtual class for defining the status tests for the Belos iterative solvers.
Definition: BelosIteration.hpp:67
Belos::BlockGmresIter
This class implements the block GMRES iteration, where a block Krylov subspace is constructed....
Definition: BelosBlockGmresIter.hpp:83
BelosBlockGmresIter.hpp
Belos concrete class for performing the block GMRES iteration.
Belos::BelosError
Parent class to all Belos exceptions.
Definition: BelosTypes.hpp:60
Belos::StatusTestCombo
A class for extending the status testing capabilities of Belos via logical combinations.
Definition: BelosStatusTestCombo.hpp:91
BelosTypes.hpp
Collection of types and exceptions used within the Belos solvers.
Belos::GmresIterationState
Structure to contain pointers to GmresIteration state variables.
Definition: BelosGmresIteration.hpp:63
BelosGmresIteration.hpp
Pure virtual base class which augments the basic interface for a Gmres linear solver iteration.
Belos::Debug
Definition: BelosTypes.hpp:265
Belos::BlockGmresSolMgr::describe
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const override
Print the object with the given verbosity level to a FancyOStream.
Definition: BelosBlockGmresSolMgr.hpp:1320
Belos::DefaultSolverParameters::convTol
static constexpr double convTol
Default convergence tolerance.
Definition: BelosTypes.hpp:294
Belos::BlockGmresSolMgr::clone
Teuchos::RCP< SolverManager< ScalarType, MV, OP > > clone() const override
clone for Inverted Injection (DII)
Definition: BelosBlockGmresSolMgr.hpp:173
Belos::Converged
Definition: BelosTypes.hpp:155
Belos::Errors
Definition: BelosTypes.hpp:258
Belos::SolverManager
The Belos::SolverManager is a templated virtual base class that defines the basic interface that any ...
Definition: BelosSolverManager.hpp:72
BelosStatusTestOutputFactory.hpp
A factory class for generating StatusTestOutput objects.
Belos::Failed
Definition: BelosTypes.hpp:189
Belos::MultiVecTraits
Traits class which defines basic operations on multivectors.
Definition: BelosMultiVecTraits.hpp:129
BelosStatusTestCombo.hpp
Belos::StatusTest for logically combining several status tests.

Generated on Thu Feb 27 2020 16:06:46 for Belos by doxygen 1.8.16