Belos  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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:
131  typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
133 
134 public:
135 
137 
138 
145 
168 
170  virtual ~BlockGmresSolMgr() {};
171 
175  }
177 
179 
180 
183  const LinearProblem<ScalarType,MV,OP>& getProblem() const override {
184  return *problem_;
185  }
186 
190 
194 
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
288  const Teuchos::EVerbosityLevel verbLevel =
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.
303 
304  // Output manager.
306  Teuchos::RCP<std::ostream> outputStream_;
307 
308  // Status test.
313  Teuchos::RCP<StatusTestResNorm<ScalarType,MV,OP> > expConvTest_, impConvTest_;
315 
316  // Orthogonalization manager.
318 
319  // Current parameter list.
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>
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>
431 {
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>
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  // Perform sanity checking for flexible Gmres here.
817  // NOTE: If the user requests that the solver manager use flexible GMRES, but there is no right preconditioner, don't use flexible GMRES.
818  // Throw an error is the user provided a left preconditioner, as that is inconsistent with flexible GMRES.
819  if (isFlexible_ && Teuchos::is_null(problem_->getRightPrec())) {
820  isFlexible_ = false;
821  params_->set("Flexible Gmres", isFlexible_);
822 
823  // If the user specified the preconditioner as a left preconditioner, throw an error.
825  "Belos::BlockGmresSolMgr::solve(): Linear problem has a left preconditioner, not a right preconditioner, which is required for flexible GMRES.");
826  }
827 
828  // If there is a left preconditioner, we create a combined status test that checks the implicit
829  // and then explicit residual norm to see if we have convergence.
830  if (!Teuchos::is_null(problem_->getLeftPrec()) && !isFlexible_) {
831  expResTest_ = true;
832  }
833 
834  if (expResTest_) {
835 
836  // Implicit residual test, using the native residual to determine if convergence was achieved.
837  Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
838  Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
839  tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
840  tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
841  impConvTest_ = tmpImpConvTest;
842 
843  // Explicit residual test once the native residual is below the tolerance
844  Teuchos::RCP<StatusTestGenResNorm_t> tmpExpConvTest =
845  Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
846  tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit, Belos::TwoNorm );
847  tmpExpConvTest->defineScaleForm( convertStringToScaleType(expResScale_), Belos::TwoNorm );
848  tmpExpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
849  expConvTest_ = tmpExpConvTest;
850 
851  // The convergence test is a combination of the "cheap" implicit test and explicit test.
852  convTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
853  }
854  else {
855 
856  if (isFlexible_) {
857  // Implicit residual test, using the native residual to determine if convergence was achieved.
858  Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
859  Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
860  tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
861  tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
862  impConvTest_ = tmpImpConvTest;
863  }
864  else {
865  // Implicit residual test, using the native residual to determine if convergence was achieved.
866  // Use test that checks for loss of accuracy.
867  Teuchos::RCP<StatusTestImpResNorm_t> tmpImpConvTest =
868  Teuchos::rcp( new StatusTestImpResNorm_t( convtol_ ) );
869  tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
870  tmpImpConvTest->setShowMaxResNormOnly( showMaxResNormOnly_ );
871  impConvTest_ = tmpImpConvTest;
872  }
873 
874  // Set the explicit and total convergence test to this implicit test that checks for accuracy loss.
875  expConvTest_ = impConvTest_;
876  convTest_ = impConvTest_;
877  }
878 
879  // Create the status test.
880  sTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
881 
882  // Add debug status test, if one is provided by the user
883  if (nonnull(debugStatusTest_) ) {
884  // Add debug convergence test
885  Teuchos::rcp_dynamic_cast<StatusTestCombo_t>(sTest_)->addStatusTest( debugStatusTest_ );
886  }
887 
888  // Create the status test output class.
889  // This class manages and formats the output from the status test.
890  StatusTestOutputFactory<ScalarType,MV,OP> stoFactory( outputStyle_ );
891  outputTest_ = stoFactory.create( printer_, sTest_, outputFreq_, Passed+Failed+Undefined );
892 
893  // Set the solver string for the output test
894  std::string solverDesc = " Block Gmres ";
895  if (isFlexible_)
896  solverDesc = "Flexible" + solverDesc;
897  outputTest_->setSolverDesc( solverDesc );
898 
899  // The status test is now set.
900  isSTSet_ = true;
901 
902  return false;
903 }
904 
905 template<class ScalarType, class MV, class OP>
907  const Teuchos::RCP<StatusTest<ScalarType,MV,OP> > &debugStatusTest
908  )
909 {
910  debugStatusTest_ = debugStatusTest;
911 }
912 
913 
914 // solve()
915 template<class ScalarType, class MV, class OP>
917 
918  // Set the current parameters if they were not set before.
919  // NOTE: This may occur if the user generated the solver manager with the default constructor and
920  // then didn't set any parameters using setParameters().
921  if (!isSet_) {
922  setParameters(Teuchos::parameterList(*getValidParameters()));
923  }
924 
926 
928  "Belos::BlockGmresSolMgr::solve(): Linear problem is not a valid object.");
929 
931  "Belos::BlockGmresSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
932 
933  if (!isSTSet_ || (!expResTest_ && !Teuchos::is_null(problem_->getLeftPrec())) ) {
935  "Belos::BlockGmresSolMgr::solve(): Linear problem and requested status tests are incompatible.");
936  }
937 
938  // Create indices for the linear systems to be solved.
939  int startPtr = 0;
940  int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
941  int numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
942 
943  std::vector<int> currIdx;
944  // If an adaptive block size is allowed then only the linear systems that need to be solved are solved.
945  // Otherwise, the index set is generated that informs the linear problem that some linear systems are augmented.
946  if ( adaptiveBlockSize_ ) {
947  blockSize_ = numCurrRHS;
948  currIdx.resize( numCurrRHS );
949  for (int i=0; i<numCurrRHS; ++i)
950  { currIdx[i] = startPtr+i; }
951 
952  }
953  else {
954  currIdx.resize( blockSize_ );
955  for (int i=0; i<numCurrRHS; ++i)
956  { currIdx[i] = startPtr+i; }
957  for (int i=numCurrRHS; i<blockSize_; ++i)
958  { currIdx[i] = -1; }
959  }
960 
961  // Inform the linear problem of the current linear system to solve.
962  problem_->setLSIndex( currIdx );
963 
965  // Parameter list
967  plist.set("Block Size",blockSize_);
968 
969  ptrdiff_t dim = MVT::GetGlobalLength( *(problem_->getRHS()) );
970  if (blockSize_*static_cast<ptrdiff_t>(numBlocks_) > dim) {
971  int tmpNumBlocks = 0;
972  if (blockSize_ == 1)
973  tmpNumBlocks = dim / blockSize_; // Allow for a good breakdown.
974  else
975  tmpNumBlocks = ( dim - blockSize_) / blockSize_; // Allow for restarting.
976  printer_->stream(Warnings) <<
977  "Belos::BlockGmresSolMgr::solve(): Warning! Requested Krylov subspace dimension is larger than operator dimension!"
978  << std::endl << " The maximum number of blocks allowed for the Krylov subspace will be adjusted to " << tmpNumBlocks << std::endl;
979  plist.set("Num Blocks",tmpNumBlocks);
980  }
981  else
982  plist.set("Num Blocks",numBlocks_);
983 
984  // Reset the status test.
985  outputTest_->reset();
986  loaDetected_ = false;
987 
988  // Assume convergence is achieved, then let any failed convergence set this to false.
989  bool isConverged = true;
990 
992  // BlockGmres solver
993 
995 
996  if (isFlexible_)
997  block_gmres_iter = Teuchos::rcp( new BlockFGmresIter<ScalarType,MV,OP>(problem_,printer_,outputTest_,ortho_,plist) );
998  else
999  block_gmres_iter = Teuchos::rcp( new BlockGmresIter<ScalarType,MV,OP>(problem_,printer_,outputTest_,ortho_,plist) );
1000 
1001  // Enter solve() iterations
1002  {
1003 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1004  Teuchos::TimeMonitor slvtimer(*timerSolve_);
1005 #endif
1006 
1007  while ( numRHS2Solve > 0 ) {
1008 
1009  // Set the current number of blocks and blocksize with the Gmres iteration.
1010  if (blockSize_*numBlocks_ > dim) {
1011  int tmpNumBlocks = 0;
1012  if (blockSize_ == 1)
1013  tmpNumBlocks = dim / blockSize_; // Allow for a good breakdown.
1014  else
1015  tmpNumBlocks = ( dim - blockSize_) / blockSize_; // Allow for restarting.
1016  block_gmres_iter->setSize( blockSize_, tmpNumBlocks );
1017  }
1018  else
1019  block_gmres_iter->setSize( blockSize_, numBlocks_ );
1020 
1021  // Reset the number of iterations.
1022  block_gmres_iter->resetNumIters();
1023 
1024  // Reset the number of calls that the status test output knows about.
1025  outputTest_->resetNumCalls();
1026 
1027  // Create the first block in the current Krylov basis.
1028  Teuchos::RCP<MV> V_0;
1029  if (isFlexible_) {
1030  // Load the correct residual if the system is augmented
1031  if (currIdx[blockSize_-1] == -1) {
1032  V_0 = MVT::Clone( *(problem_->getInitResVec()), blockSize_ );
1033  problem_->computeCurrResVec( &*V_0 );
1034  }
1035  else {
1036  V_0 = MVT::CloneCopy( *(problem_->getInitResVec()), currIdx );
1037  }
1038  }
1039  else {
1040  // Load the correct residual if the system is augmented
1041  if (currIdx[blockSize_-1] == -1) {
1042  V_0 = MVT::Clone( *(problem_->getInitPrecResVec()), blockSize_ );
1043  problem_->computeCurrPrecResVec( &*V_0 );
1044  }
1045  else {
1046  V_0 = MVT::CloneCopy( *(problem_->getInitPrecResVec()), currIdx );
1047  }
1048  }
1049 
1050  // Get a matrix to hold the orthonormalization coefficients.
1052  Teuchos::rcp( new Teuchos::SerialDenseMatrix<int,ScalarType>( blockSize_, blockSize_ ) );
1053 
1054  // Orthonormalize the new V_0
1055  int rank = ortho_->normalize( *V_0, z_0 );
1057  "Belos::BlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors.");
1058 
1059  // Set the new state and initialize the solver.
1061  newstate.V = V_0;
1062  newstate.z = z_0;
1063  newstate.curDim = 0;
1064  block_gmres_iter->initializeGmres(newstate);
1065  int numRestarts = 0;
1066 
1067  while(1) {
1068  // tell block_gmres_iter to iterate
1069  try {
1070  block_gmres_iter->iterate();
1071 
1073  //
1074  // check convergence first
1075  //
1077  if ( convTest_->getStatus() == Passed ) {
1078  if ( expConvTest_->getLOADetected() ) {
1079  // we don't have convergence
1080  loaDetected_ = true;
1081  printer_->stream(Warnings) <<
1082  "Belos::BlockGmresSolMgr::solve(): Warning! Solver has experienced a loss of accuracy!" << std::endl;
1083  isConverged = false;
1084  }
1085  break; // break from while(1){block_gmres_iter->iterate()}
1086  }
1088  //
1089  // check for maximum iterations
1090  //
1092  else if ( maxIterTest_->getStatus() == Passed ) {
1093  // we don't have convergence
1094  isConverged = false;
1095  break; // break from while(1){block_gmres_iter->iterate()}
1096  }
1098  //
1099  // check for restarting, i.e. the subspace is full
1100  //
1102  else if ( block_gmres_iter->getCurSubspaceDim() == block_gmres_iter->getMaxSubspaceDim() ) {
1103 
1104  if ( numRestarts >= maxRestarts_ ) {
1105  isConverged = false;
1106  break; // break from while(1){block_gmres_iter->iterate()}
1107  }
1108  numRestarts++;
1109 
1110  printer_->stream(Debug) << " Performing restart number " << numRestarts << " of " << maxRestarts_ << std::endl << std::endl;
1111 
1112  // Update the linear problem.
1113  Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1114  if (isFlexible_) {
1115  // Update the solution manually, since the preconditioning doesn't need to be undone.
1116  Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
1117  MVT::MvAddMv( 1.0, *curX, 1.0, *update, *curX );
1118  }
1119  else
1120  problem_->updateSolution( update, true );
1121 
1122  // Get the state.
1123  GmresIterationState<ScalarType,MV> oldState = block_gmres_iter->getState();
1124 
1125  // Compute the restart std::vector.
1126  // Get a view of the current Krylov basis.
1127  V_0 = MVT::Clone( *(oldState.V), blockSize_ );
1128  if (isFlexible_)
1129  problem_->computeCurrResVec( &*V_0 );
1130  else
1131  problem_->computeCurrPrecResVec( &*V_0 );
1132 
1133  // Get a view of the first block of the Krylov basis.
1134  z_0 = Teuchos::rcp( new Teuchos::SerialDenseMatrix<int,ScalarType>( blockSize_, blockSize_ ) );
1135 
1136  // Orthonormalize the new V_0
1137  rank = ortho_->normalize( *V_0, z_0 );
1139  "Belos::BlockGmresSolMgr::solve(): Failed to compute initial block of orthonormal vectors after restart.");
1140 
1141  // Set the new state and initialize the solver.
1142  newstate.V = V_0;
1143  newstate.z = z_0;
1144  newstate.curDim = 0;
1145  block_gmres_iter->initializeGmres(newstate);
1146 
1147  } // end of restarting
1148 
1150  //
1151  // we returned from iterate(), but none of our status tests Passed.
1152  // something is wrong, and it is probably our fault.
1153  //
1155 
1156  else {
1157  TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,
1158  "Belos::BlockGmresSolMgr::solve(): Invalid return from BlockGmresIter::iterate().");
1159  }
1160  }
1161  catch (const GmresIterationOrthoFailure &e) {
1162  // If the block size is not one, it's not considered a lucky breakdown.
1163  if (blockSize_ != 1) {
1164  printer_->stream(Errors) << "Error! Caught std::exception in BlockGmresIter::iterate() at iteration "
1165  << block_gmres_iter->getNumIters() << std::endl
1166  << e.what() << std::endl;
1167  if (convTest_->getStatus() != Passed)
1168  isConverged = false;
1169  break;
1170  }
1171  else {
1172  // If the block size is one, try to recover the most recent least-squares solution
1173  block_gmres_iter->updateLSQR( block_gmres_iter->getCurSubspaceDim() );
1174 
1175  // Check to see if the most recent least-squares solution yielded convergence.
1176  sTest_->checkStatus( &*block_gmres_iter );
1177  if (convTest_->getStatus() != Passed)
1178  isConverged = false;
1179  break;
1180  }
1181  }
1182  catch (const std::exception &e) {
1183  printer_->stream(Errors) << "Error! Caught std::exception in BlockGmresIter::iterate() at iteration "
1184  << block_gmres_iter->getNumIters() << std::endl
1185  << e.what() << std::endl;
1186  throw;
1187  }
1188  }
1189 
1190  // Compute the current solution.
1191  // Update the linear problem.
1192  if (isFlexible_) {
1193  // Update the solution manually, since the preconditioning doesn't need to be undone.
1194  Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1195  Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
1196  // Update the solution only if there is a valid update from the iteration
1197  if (update != Teuchos::null)
1198  MVT::MvAddMv( 1.0, *curX, 1.0, *update, *curX );
1199  }
1200  else {
1201  // Attempt to get the current solution from the residual status test, if it has one.
1202  if ( !Teuchos::is_null(expConvTest_->getSolution()) ) {
1203  Teuchos::RCP<MV> newX = expConvTest_->getSolution();
1204  Teuchos::RCP<MV> curX = problem_->getCurrLHSVec();
1205  MVT::MvAddMv( 0.0, *newX, 1.0, *newX, *curX );
1206  }
1207  else {
1208  Teuchos::RCP<MV> update = block_gmres_iter->getCurrentUpdate();
1209  problem_->updateSolution( update, true );
1210  }
1211  }
1212 
1213  // Inform the linear problem that we are finished with this block linear system.
1214  problem_->setCurrLS();
1215 
1216  // Update indices for the linear systems to be solved.
1217  startPtr += numCurrRHS;
1218  numRHS2Solve -= numCurrRHS;
1219  if ( numRHS2Solve > 0 ) {
1220  numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
1221 
1222  if ( adaptiveBlockSize_ ) {
1223  blockSize_ = numCurrRHS;
1224  currIdx.resize( numCurrRHS );
1225  for (int i=0; i<numCurrRHS; ++i)
1226  { currIdx[i] = startPtr+i; }
1227  }
1228  else {
1229  currIdx.resize( blockSize_ );
1230  for (int i=0; i<numCurrRHS; ++i)
1231  { currIdx[i] = startPtr+i; }
1232  for (int i=numCurrRHS; i<blockSize_; ++i)
1233  { currIdx[i] = -1; }
1234  }
1235  // Set the next indices.
1236  problem_->setLSIndex( currIdx );
1237  }
1238  else {
1239  currIdx.resize( numRHS2Solve );
1240  }
1241 
1242  }// while ( numRHS2Solve > 0 )
1243 
1244  }
1245 
1246  // print final summary
1247  sTest_->print( printer_->stream(FinalSummary) );
1248 
1249  // print timing information
1250 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1251  // Calling summarize() can be expensive, so don't call unless the
1252  // user wants to print out timing details. summarize() will do all
1253  // the work even if it's passed a "black hole" output stream.
1254  if (verbosity_ & TimingDetails)
1255  Teuchos::TimeMonitor::summarize( printer_->stream(TimingDetails) );
1256 #endif
1257 
1258  // get iteration information for this solve
1259  numIters_ = maxIterTest_->getNumIters();
1260 
1261  // Save the convergence test value ("achieved tolerance") for this
1262  // solve. This requires a bit more work than for BlockCGSolMgr,
1263  // since for this solver, convTest_ may either be a single residual
1264  // norm test, or a combination of two residual norm tests. In the
1265  // latter case, the master convergence test convTest_ is a SEQ combo
1266  // of the implicit resp. explicit tests. If the implicit test never
1267  // passes, then the explicit test won't ever be executed. This
1268  // manifests as expConvTest_->getTestValue()->size() < 1. We deal
1269  // with this case by using the values returned by
1270  // impConvTest_->getTestValue().
1271  {
1272  // We'll fetch the vector of residual norms one way or the other.
1273  const std::vector<MagnitudeType>* pTestValues = NULL;
1274  if (expResTest_) {
1275  pTestValues = expConvTest_->getTestValue();
1276  if (pTestValues == NULL || pTestValues->size() < 1) {
1277  pTestValues = impConvTest_->getTestValue();
1278  }
1279  }
1280  else {
1281  // Only the implicit residual norm test is being used.
1282  pTestValues = impConvTest_->getTestValue();
1283  }
1284  TEUCHOS_TEST_FOR_EXCEPTION(pTestValues == NULL, std::logic_error,
1285  "Belos::BlockGmresSolMgr::solve(): The implicit convergence test's "
1286  "getTestValue() method returned NULL. Please report this bug to the "
1287  "Belos developers.");
1288  TEUCHOS_TEST_FOR_EXCEPTION(pTestValues->size() < 1, std::logic_error,
1289  "Belos::BlockGmresSolMgr::solve(): The implicit convergence test's "
1290  "getTestValue() method returned a vector of length zero. Please report "
1291  "this bug to the Belos developers.");
1292 
1293  // FIXME (mfh 12 Dec 2011) Does pTestValues really contain the
1294  // achieved tolerances for all vectors in the current solve(), or
1295  // just for the vectors from the last deflation?
1296  achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
1297  }
1298 
1299  if (!isConverged || loaDetected_) {
1300  return Unconverged; // return from BlockGmresSolMgr::solve()
1301  }
1302  return Converged; // return from BlockGmresSolMgr::solve()
1303 }
1304 
1305 
1306 template<class ScalarType, class MV, class OP>
1308 {
1309  std::ostringstream out;
1310  out << "\"Belos::BlockGmresSolMgr\": {";
1311  if (this->getObjectLabel () != "") {
1312  out << "Label: " << this->getObjectLabel () << ", ";
1313  }
1314  out << "Flexible: " << (isFlexible_ ? "true" : "false")
1315  << ", Num Blocks: " << numBlocks_
1316  << ", Maximum Iterations: " << maxIters_
1317  << ", Maximum Restarts: " << maxRestarts_
1318  << ", Convergence Tolerance: " << convtol_
1319  << "}";
1320  return out.str ();
1321 }
1322 
1323 
1324 template<class ScalarType, class MV, class OP>
1325 void
1328  const Teuchos::EVerbosityLevel verbLevel) const
1329 {
1331  using Teuchos::VERB_DEFAULT;
1332  using Teuchos::VERB_NONE;
1333  using Teuchos::VERB_LOW;
1334  // using Teuchos::VERB_MEDIUM;
1335  // using Teuchos::VERB_HIGH;
1336  // using Teuchos::VERB_EXTREME;
1337  using std::endl;
1338 
1339  // Set default verbosity if applicable.
1340  const Teuchos::EVerbosityLevel vl =
1341  (verbLevel == VERB_DEFAULT) ? VERB_LOW : verbLevel;
1342 
1343  if (vl != VERB_NONE) {
1344  Teuchos::OSTab tab0 (out);
1345 
1346  out << "\"Belos::BlockGmresSolMgr\":" << endl;
1347  Teuchos::OSTab tab1 (out);
1348  out << "Template parameters:" << endl;
1349  {
1350  Teuchos::OSTab tab2 (out);
1351  out << "ScalarType: " << TypeNameTraits<ScalarType>::name () << endl
1352  << "MV: " << TypeNameTraits<MV>::name () << endl
1353  << "OP: " << TypeNameTraits<OP>::name () << endl;
1354  }
1355  if (this->getObjectLabel () != "") {
1356  out << "Label: " << this->getObjectLabel () << endl;
1357  }
1358  out << "Flexible: " << (isFlexible_ ? "true" : "false") << endl
1359  << "Num Blocks: " << numBlocks_ << endl
1360  << "Maximum Iterations: " << maxIters_ << endl
1361  << "Maximum Restarts: " << maxRestarts_ << endl
1362  << "Convergence Tolerance: " << convtol_ << endl;
1363  }
1364 }
1365 
1366 
1367 } // end Belos namespace
1368 
1369 #endif /* BELOS_BLOCK_GMRES_SOLMGR_HPP */
static const double orthoKappa
DGKS orthogonalization constant.
Definition: BelosTypes.hpp:307
ScaleType convertStringToScaleType(const std::string &scaleType)
Convert the given string to its ScaleType enum value.
Definition: BelosTypes.cpp:106
Teuchos::Array< Teuchos::RCP< Teuchos::Time > > getTimers() const
Return the timers for this object.
Collection of types and exceptions used within the Belos solvers.
Belos&#39;s basic output manager for sending information of select verbosity levels to the appropriate ou...
Class which manages the output and verbosity of the Belos solvers.
bool is_null(const boost::shared_ptr< T > &p)
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...
ScaleType
The type of scaling to use on the residual norm value.
Definition: BelosTypes.hpp:120
This class implements the block GMRES iteration, where a block Krylov subspace is constructed...
Convergence test using the implicit residual norm(s), with an explicit residual norm(s) check for los...
Teuchos::RCP< const MV > V
The current Krylov basis.
T & get(const std::string &name, T def_value)
Virtual base class for StatusTest that printing status tests.
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
static RCP< Time > getNewCounter(const std::string &name)
bool is_null(const std::shared_ptr< T > &p)
MagnitudeType achievedTol() const override
Tolerance achieved by the last solve() invocation.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
Get a parameter list containing the valid parameters for this object.
Teuchos::RCP< const Teuchos::SerialDenseMatrix< int, ScalarType > > z
The current right-hand side of the least squares system RY = Z.
An implementation of StatusTestResNorm using a family of residual norms.
BlockGmresSolMgrLinearProblemFailure is thrown when the linear problem is not setup (i...
Structure to contain pointers to GmresIteration state variables.
static const double convTol
Default convergence tolerance.
Definition: BelosTypes.hpp:301
Belos::StatusTest class for specifying a maximum number of iterations.
Interface to Block GMRES and Flexible GMRES.
A pure virtual class for defining the status tests for the Belos iterative solvers.
virtual ~BlockGmresSolMgr()
Destructor.
void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem) override
Set the linear problem that needs to be solved.
A factory class for generating StatusTestOutput objects.
Iterated Modified Gram-Schmidt (IMGS) implementation of the Belos::OrthoManager class.
ReturnType solve() override
This method performs possibly repeated calls to the underlying linear solver&#39;s iterate() routine unti...
Traits class which defines basic operations on multivectors.
void setParameters(const Teuchos::RCP< Teuchos::ParameterList > &params) override
Set the parameters the solver manager should use to solve the linear problem.
Belos::StatusTest for logically combining several status tests.
bool isParameter(const std::string &name) const
Classical Gram-Schmidt (with DGKS correction) implementation of the Belos::OrthoManager class...
A Belos::StatusTest class for specifying a maximum number of iterations.
ResetType
How to reset the solver.
Definition: BelosTypes.hpp:206
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Pure virtual base class which describes the basic interface for a solver manager. ...
int curDim
The current dimension of the reduction.
static void summarize(Ptr< const Comm< int > > comm, std::ostream &out=std::cout, const bool alwaysWriteLocal=false, const bool writeGlobalStats=true, const bool writeZeroTimers=true, const ECounterSetOp setOp=Intersection, const std::string &filter="", const bool ignoreZeroTimers=false)
BlockGmresSolMgrLinearProblemFailure(const std::string &what_arg)
void reset(const ResetType type) override
Performs a reset of the solver manager specified by the ResetType. This informs the solver manager th...
A linear system to solve, and its associated information.
Class which describes the linear problem to be solved by the iterative solver.
This class implements the block flexible GMRES iteration, where a block Krylov subspace is constructe...
int getNumIters() const override
Get the iteration count for the most recent call to solve().
Teuchos::RCP< SolverManager< ScalarType, MV, OP > > clone() const override
clone for Inverted Injection (DII)
ReturnType
Whether the Belos solve converged for all linear systems.
Definition: BelosTypes.hpp:155
The Belos::SolverManager is a templated virtual base class that defines the basic interface that any ...
Belos concrete class for performing the block GMRES iteration.
Iterated Classical Gram-Schmidt (ICGS) implementation of the Belos::OrthoManager class.
void validateParameters(ParameterList const &validParamList, int const depth=1000, EValidateUsed const validateUsed=VALIDATE_USED_ENABLED, EValidateDefaults const validateDefaults=VALIDATE_DEFAULTS_ENABLED) const
An implementation of the Belos::MatOrthoManager that performs orthogonalization using multiple steps ...
static const EVerbosityLevel verbLevel_default
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.
An implementation of the Belos::MatOrthoManager that performs orthogonalization using multiple steps ...
bool nonnull(const boost::shared_ptr< T > &p)
bool isLOADetected() const override
Return whether a loss of accuracy was detected by this solver during the most current solve...
Belos::StatusTestResNorm for specifying general residual norm stopping criteria.
Belos::StatusTest for specifying an implicit residual norm stopping criteria that checks for loss of ...
BlockGmresSolMgr()
Empty constructor for BlockGmresSolMgr. This constructor takes no arguments and sets the default valu...
bool isType(const std::string &name) const
An implementation of the Belos::MatOrthoManager that performs orthogonalization using (potentially) m...
BlockGmresSolMgrOrthoFailure(const std::string &what_arg)
A class for extending the status testing capabilities of Belos via logical combinations.
std::string description() const override
Return a one-line description of this object.
Class which defines basic traits for the operator type.
Parent class to all Belos exceptions.
Definition: BelosTypes.hpp:60
BlockGmresSolMgrOrthoFailure is thrown when the orthogonalization manager is unable to generate ortho...
Default parameters common to most Belos solvers.
Definition: BelosTypes.hpp:291
Pure virtual base class which augments the basic interface for a Gmres linear solver iteration...
GmresIterationOrthoFailure is thrown when the GmresIteration object is unable to compute independent ...
Belos header file which uses auto-configuration information to include necessary C++ headers...
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const override
Get a parameter list containing the current parameters for this object.
const LinearProblem< ScalarType, MV, OP > & getProblem() const override
Get current linear problem being solved for in this object.
Belos concrete class for performing the block, flexible GMRES iteration.

Generated on Fri Jun 5 2020 10:20:38 for Belos by doxygen 1.8.5