Belos Package Browser (Single Doxygen Collection)  Development
BelosTFQMRSolMgr.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_TFQMR_SOLMGR_HPP
43 #define BELOS_TFQMR_SOLMGR_HPP
44 
49 #include "BelosConfigDefs.hpp"
50 #include "BelosTypes.hpp"
51 
52 #include "BelosLinearProblem.hpp"
53 #include "BelosSolverManager.hpp"
54 
55 #include "BelosTFQMRIter.hpp"
58 #include "BelosStatusTestCombo.hpp"
60 #include "BelosOutputManager.hpp"
61 #ifdef BELOS_TEUCHOS_TIME_MONITOR
62 #include "Teuchos_TimeMonitor.hpp"
63 #endif
64 
78 namespace Belos {
79 
81 
82 
90  TFQMRSolMgrLinearProblemFailure(const std::string& what_arg) : BelosError(what_arg)
91  {}};
92 
93  template<class ScalarType, class MV, class OP>
94  class TFQMRSolMgr : public SolverManager<ScalarType,MV,OP> {
95 
96  private:
99  typedef Teuchos::ScalarTraits<ScalarType> SCT;
100  typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
101  typedef Teuchos::ScalarTraits<MagnitudeType> MT;
102 
103  public:
104 
106 
107 
113  TFQMRSolMgr();
114 
131  TFQMRSolMgr( const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem,
132  const Teuchos::RCP<Teuchos::ParameterList> &pl );
133 
135  virtual ~TFQMRSolMgr() {};
136 
138  Teuchos::RCP<SolverManager<ScalarType, MV, OP> > clone () const override {
139  return Teuchos::rcp(new TFQMRSolMgr<ScalarType,MV,OP>);
140  }
142 
144 
145 
146  const LinearProblem<ScalarType,MV,OP>& getProblem() const override {
147  return *problem_;
148  }
149 
152  Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const override;
153 
156  Teuchos::RCP<const Teuchos::ParameterList> getCurrentParameters() const override { return params_; }
157 
163  Teuchos::Array<Teuchos::RCP<Teuchos::Time> > getTimers() const {
164  return Teuchos::tuple(timerSolve_);
165  }
166 
172  MagnitudeType achievedTol() const override {
173  return achievedTol_;
174  }
175 
177  int getNumIters() const override {
178  return numIters_;
179  }
180 
188  bool isLOADetected() const override { return false; }
190 
192 
193 
195  void setProblem( const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem ) override { problem_ = problem; }
196 
198  void setParameters( const Teuchos::RCP<Teuchos::ParameterList> &params ) override;
199 
201 
203 
208  void reset( const ResetType type ) override { if ((type & Belos::Problem) && !Teuchos::is_null(problem_)) problem_->setProblem(); }
210 
212 
213 
231  ReturnType solve() override;
232 
234 
236 
238  std::string description() const override;
240 
241  private:
242 
243  // Method for checking current status test against defined linear problem.
244  bool checkStatusTest();
245 
246  // Linear problem.
247  Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > problem_;
248 
249  // Output manager.
250  Teuchos::RCP<OutputManager<ScalarType> > printer_;
251  Teuchos::RCP<std::ostream> outputStream_;
252 
253  // Status test.
254  Teuchos::RCP<StatusTest<ScalarType,MV,OP> > sTest_;
255  Teuchos::RCP<StatusTestMaxIters<ScalarType,MV,OP> > maxIterTest_;
256  Teuchos::RCP<StatusTest<ScalarType,MV,OP> > convTest_;
257  Teuchos::RCP<StatusTestGenResNorm<ScalarType,MV,OP> > expConvTest_, impConvTest_;
258  Teuchos::RCP<StatusTestOutput<ScalarType,MV,OP> > outputTest_;
259 
260  // Current parameter list.
261  Teuchos::RCP<Teuchos::ParameterList> params_;
262 
263  // Default solver values.
264  static constexpr int maxIters_default_ = 1000;
265  static constexpr bool expResTest_default_ = false;
266  static constexpr int verbosity_default_ = Belos::Errors;
267  static constexpr int outputStyle_default_ = Belos::General;
268  static constexpr int outputFreq_default_ = -1;
269  static constexpr const char * impResScale_default_ = "Norm of Preconditioned Initial Residual";
270  static constexpr const char * expResScale_default_ = "Norm of Initial Residual";
271  static constexpr const char * label_default_ = "Belos";
272 // https://stackoverflow.com/questions/24398102/constexpr-and-initialization-of-a-static-const-void-pointer-with-reinterpret-cas
273 #if defined(_WIN32) && defined(__clang__)
274  static constexpr std::ostream * outputStream_default_ =
275  __builtin_constant_p(reinterpret_cast<const std::ostream*>(&std::cout));
276 #else
277  static constexpr std::ostream * outputStream_default_ = &std::cout;
278 #endif
279 
280  // Current solver values.
287 
288  // Timers.
289  std::string label_;
290  Teuchos::RCP<Teuchos::Time> timerSolve_;
291 
292  // Internal state variables.
294  };
295 
296 
297 // Empty Constructor
298 template<class ScalarType, class MV, class OP>
300  outputStream_(Teuchos::rcp(outputStream_default_,false)),
301  convtol_(DefaultSolverParameters::convTol),
302  impTolScale_(DefaultSolverParameters::impTolScale),
303  achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
304  maxIters_(maxIters_default_),
305  numIters_(0),
306  verbosity_(verbosity_default_),
307  outputStyle_(outputStyle_default_),
308  outputFreq_(outputFreq_default_),
309  blockSize_(1),
310  expResTest_(expResTest_default_),
311  impResScale_(impResScale_default_),
312  expResScale_(expResScale_default_),
313  label_(label_default_),
314  isSet_(false),
315  isSTSet_(false)
316 {}
317 
318 
319 // Basic Constructor
320 template<class ScalarType, class MV, class OP>
322  const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem,
323  const Teuchos::RCP<Teuchos::ParameterList> &pl ) :
324  problem_(problem),
325  outputStream_(Teuchos::rcp(outputStream_default_,false)),
326  convtol_(DefaultSolverParameters::convTol),
327  impTolScale_(DefaultSolverParameters::impTolScale),
328  achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
329  maxIters_(maxIters_default_),
330  numIters_(0),
331  verbosity_(verbosity_default_),
332  outputStyle_(outputStyle_default_),
333  outputFreq_(outputFreq_default_),
334  blockSize_(1),
335  expResTest_(expResTest_default_),
336  impResScale_(impResScale_default_),
337  expResScale_(expResScale_default_),
338  label_(label_default_),
339  isSet_(false),
340  isSTSet_(false)
341 {
342  TEUCHOS_TEST_FOR_EXCEPTION(problem_ == Teuchos::null, std::invalid_argument, "Problem not given to solver manager.");
343 
344  // If the parameter list pointer is null, then set the current parameters to the default parameter list.
345  if ( !is_null(pl) ) {
346  setParameters( pl );
347  }
348 }
349 
350 template<class ScalarType, class MV, class OP>
351 void TFQMRSolMgr<ScalarType,MV,OP>::setParameters( const Teuchos::RCP<Teuchos::ParameterList> &params )
352 {
353  // Create the internal parameter list if ones doesn't already exist.
354  if (params_ == Teuchos::null) {
355  params_ = Teuchos::rcp( new Teuchos::ParameterList(*getValidParameters()) );
356  }
357  else {
358  params->validateParameters(*getValidParameters());
359  }
360 
361  // Check for maximum number of iterations
362  if (params->isParameter("Maximum Iterations")) {
363  maxIters_ = params->get("Maximum Iterations",maxIters_default_);
364 
365  // Update parameter in our list and in status test.
366  params_->set("Maximum Iterations", maxIters_);
367  if (maxIterTest_!=Teuchos::null)
368  maxIterTest_->setMaxIters( maxIters_ );
369  }
370 
371  // Check for blocksize
372  if (params->isParameter("Block Size")) {
373  blockSize_ = params->get("Block Size",1);
374  TEUCHOS_TEST_FOR_EXCEPTION(blockSize_ != 1, std::invalid_argument,
375  "Belos::TFQMRSolMgr: \"Block Size\" must be 1.");
376 
377  // Update parameter in our list.
378  params_->set("Block Size", blockSize_);
379  }
380 
381  // Check to see if the timer label changed.
382  if (params->isParameter("Timer Label")) {
383  std::string tempLabel = params->get("Timer Label", label_default_);
384 
385  // Update parameter in our list and solver timer
386  if (tempLabel != label_) {
387  label_ = tempLabel;
388  params_->set("Timer Label", label_);
389  std::string solveLabel = label_ + ": TFQMRSolMgr total solve time";
390 #ifdef BELOS_TEUCHOS_TIME_MONITOR
391  timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
392 #endif
393  }
394  }
395 
396  // Check for a change in verbosity level
397  if (params->isParameter("Verbosity")) {
398  if (Teuchos::isParameterType<int>(*params,"Verbosity")) {
399  verbosity_ = params->get("Verbosity", verbosity_default_);
400  } else {
401  verbosity_ = (int)Teuchos::getParameter<Belos::MsgType>(*params,"Verbosity");
402  }
403 
404  // Update parameter in our list.
405  params_->set("Verbosity", verbosity_);
406  if (printer_ != Teuchos::null)
407  printer_->setVerbosity(verbosity_);
408  }
409 
410  // Check for a change in output style
411  if (params->isParameter("Output Style")) {
412  if (Teuchos::isParameterType<int>(*params,"Output Style")) {
413  outputStyle_ = params->get("Output Style", outputStyle_default_);
414  } else {
415  outputStyle_ = (int)Teuchos::getParameter<Belos::OutputType>(*params,"Output Style");
416  }
417 
418  // Reconstruct the convergence test if the explicit residual test is not being used.
419  params_->set("Output Style", outputStyle_);
420  isSTSet_ = false;
421  }
422 
423  // output stream
424  if (params->isParameter("Output Stream")) {
425  outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> >(*params,"Output Stream");
426 
427  // Update parameter in our list.
428  params_->set("Output Stream", outputStream_);
429  if (printer_ != Teuchos::null)
430  printer_->setOStream( outputStream_ );
431  }
432 
433  // frequency level
434  if (verbosity_ & Belos::StatusTestDetails) {
435  if (params->isParameter("Output Frequency")) {
436  outputFreq_ = params->get("Output Frequency", outputFreq_default_);
437  }
438 
439  // Update parameter in out list and output status test.
440  params_->set("Output Frequency", outputFreq_);
441  if (outputTest_ != Teuchos::null)
442  outputTest_->setOutputFrequency( outputFreq_ );
443  }
444 
445  // Create output manager if we need to.
446  if (printer_ == Teuchos::null) {
447  printer_ = Teuchos::rcp( new OutputManager<ScalarType>(verbosity_, outputStream_) );
448  }
449 
450  // Check for convergence tolerance
451  if (params->isParameter("Convergence Tolerance")) {
452  if (params->isType<MagnitudeType> ("Convergence Tolerance")) {
453  convtol_ = params->get ("Convergence Tolerance",
454  static_cast<MagnitudeType> (DefaultSolverParameters::convTol));
455  }
456  else {
457  convtol_ = params->get ("Convergence Tolerance", DefaultSolverParameters::convTol);
458  }
459 
460  // Update parameter in our list.
461  params_->set("Convergence Tolerance", convtol_);
462  isSTSet_ = false;
463  }
464 
465  // Check for implicit residual scaling
466  if (params->isParameter("Implicit Tolerance Scale Factor")) {
467  if (params->isType<MagnitudeType> ("Implicit Tolerance Scale Factor")) {
468  impTolScale_ = params->get ("Implicit Tolerance Scale Factor",
469  static_cast<MagnitudeType> (DefaultSolverParameters::impTolScale));
470 
471  }
472  else {
473  impTolScale_ = params->get ("Implicit Tolerance Scale Factor",
475  }
476 
477  // Update parameter in our list.
478  params_->set("Implicit Tolerance Scale Factor", impTolScale_);
479  isSTSet_ = false;
480  }
481 
482  // Check for a change in scaling, if so we need to build new residual tests.
483  if (params->isParameter("Implicit Residual Scaling")) {
484  std::string tempImpResScale = Teuchos::getParameter<std::string>( *params, "Implicit Residual Scaling" );
485 
486  // Only update the scaling if it's different.
487  if (impResScale_ != tempImpResScale) {
488  impResScale_ = tempImpResScale;
489 
490  // Update parameter in our list and residual tests
491  params_->set("Implicit Residual Scaling", impResScale_);
492 
493  // Make sure the convergence test gets constructed again.
494  isSTSet_ = false;
495  }
496  }
497 
498  if (params->isParameter("Explicit Residual Scaling")) {
499  std::string tempExpResScale = Teuchos::getParameter<std::string>( *params, "Explicit Residual Scaling" );
500 
501  // Only update the scaling if it's different.
502  if (expResScale_ != tempExpResScale) {
503  expResScale_ = tempExpResScale;
504 
505  // Update parameter in our list and residual tests
506  params_->set("Explicit Residual Scaling", expResScale_);
507 
508  // Make sure the convergence test gets constructed again.
509  isSTSet_ = false;
510  }
511  }
512 
513  if (params->isParameter("Explicit Residual Test")) {
514  expResTest_ = Teuchos::getParameter<bool>( *params,"Explicit Residual Test" );
515 
516  // Reconstruct the convergence test if the explicit residual test is not being used.
517  params_->set("Explicit Residual Test", expResTest_);
518  if (expConvTest_ == Teuchos::null) {
519  isSTSet_ = false;
520  }
521  }
522 
523  // Create the timer if we need to.
524  if (timerSolve_ == Teuchos::null) {
525  std::string solveLabel = label_ + ": TFQMRSolMgr total solve time";
526 #ifdef BELOS_TEUCHOS_TIME_MONITOR
527  timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
528 #endif
529  }
530 
531  // Inform the solver manager that the current parameters were set.
532  isSet_ = true;
533 }
534 
535 
536 // Check the status test versus the defined linear problem
537 template<class ScalarType, class MV, class OP>
539 
540  typedef Belos::StatusTestCombo<ScalarType,MV,OP> StatusTestCombo_t;
541  typedef Belos::StatusTestGenResNorm<ScalarType,MV,OP> StatusTestGenResNorm_t;
542 
543  // Basic test checks maximum iterations and native residual.
544  maxIterTest_ = Teuchos::rcp( new StatusTestMaxIters<ScalarType,MV,OP>( maxIters_ ) );
545 
546  if (expResTest_) {
547 
548  // Implicit residual test, using the native residual to determine if convergence was achieved.
549  Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
550  Teuchos::rcp( new StatusTestGenResNorm_t( impTolScale_*convtol_ ) );
551  tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
552  impConvTest_ = tmpImpConvTest;
553 
554  // Explicit residual test once the native residual is below the tolerance
555  Teuchos::RCP<StatusTestGenResNorm_t> tmpExpConvTest =
556  Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
557  tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit, Belos::TwoNorm );
558  tmpExpConvTest->defineScaleForm( convertStringToScaleType(expResScale_), Belos::TwoNorm );
559  expConvTest_ = tmpExpConvTest;
560 
561  // The convergence test is a combination of the "cheap" implicit test and explicit test.
562  convTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
563  }
564  else {
565 
566  // Implicit residual test, using the native residual to determine if convergence was achieved.
567  Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
568  Teuchos::rcp( new StatusTestGenResNorm_t( convtol_ ) );
569  tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
570  impConvTest_ = tmpImpConvTest;
571 
572  // Set the explicit and total convergence test to this implicit test that checks for accuracy loss.
573  expConvTest_ = impConvTest_;
574  convTest_ = impConvTest_;
575  }
576  sTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
577 
578  // Create the status test output class.
579  // This class manages and formats the output from the status test.
580  StatusTestOutputFactory<ScalarType,MV,OP> stoFactory( outputStyle_ );
581  outputTest_ = stoFactory.create( printer_, sTest_, outputFreq_, Passed+Failed+Undefined );
582 
583  // Set the solver string for the output test
584  std::string solverDesc = " TFQMR ";
585  outputTest_->setSolverDesc( solverDesc );
586 
587 
588  // The status test is now set.
589  isSTSet_ = true;
590 
591  return false;
592 }
593 
594 
595 template<class ScalarType, class MV, class OP>
596 Teuchos::RCP<const Teuchos::ParameterList>
598 {
599  static Teuchos::RCP<const Teuchos::ParameterList> validPL;
600 
601  // Set all the valid parameters and their default values.
602  if(is_null(validPL)) {
603  Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
604 
605  // The static_cast is to resolve an issue with older clang versions which
606  // would cause the constexpr to link fail. With c++17 the problem is resolved.
607  pl->set("Convergence Tolerance", static_cast<MagnitudeType>(DefaultSolverParameters::convTol),
608  "The relative residual tolerance that needs to be achieved by the\n"
609  "iterative solver in order for the linear system to be declared converged.");
610  pl->set("Implicit Tolerance Scale Factor", static_cast<MagnitudeType>(DefaultSolverParameters::impTolScale),
611  "The scale factor used by the implicit residual test when explicit residual\n"
612  "testing is used. May enable faster convergence when TFQMR bound is too loose.");
613  pl->set("Maximum Iterations", static_cast<int>(maxIters_default_),
614  "The maximum number of block iterations allowed for each\n"
615  "set of RHS solved.");
616  pl->set("Verbosity", static_cast<int>(verbosity_default_),
617  "What type(s) of solver information should be outputted\n"
618  "to the output stream.");
619  pl->set("Output Style", static_cast<int>(outputStyle_default_),
620  "What style is used for the solver information outputted\n"
621  "to the output stream.");
622  pl->set("Output Frequency", static_cast<int>(outputFreq_default_),
623  "How often convergence information should be outputted\n"
624  "to the output stream.");
625  pl->set("Output Stream", Teuchos::rcp(outputStream_default_,false),
626  "A reference-counted pointer to the output stream where all\n"
627  "solver output is sent.");
628  pl->set("Explicit Residual Test", static_cast<bool>(expResTest_default_),
629  "Whether the explicitly computed residual should be used in the convergence test.");
630  pl->set("Implicit Residual Scaling", static_cast<const char *>(impResScale_default_),
631  "The type of scaling used in the implicit residual convergence test.");
632  pl->set("Explicit Residual Scaling", static_cast<const char *>(expResScale_default_),
633  "The type of scaling used in the explicit residual convergence test.");
634  pl->set("Timer Label", static_cast<const char *>(label_default_),
635  "The string to use as a prefix for the timer labels.");
636  validPL = pl;
637  }
638  return validPL;
639 }
640 
641 
642 // solve()
643 template<class ScalarType, class MV, class OP>
645 
646  // Set the current parameters if they were not set before.
647  // NOTE: This may occur if the user generated the solver manager with the default constructor and
648  // then didn't set any parameters using setParameters().
649  if (!isSet_) {
650  setParameters(Teuchos::parameterList(*getValidParameters()));
651  }
652 
653  TEUCHOS_TEST_FOR_EXCEPTION(problem_ == Teuchos::null,TFQMRSolMgrLinearProblemFailure,
654  "Belos::TFQMRSolMgr::solve(): Linear problem is not a valid object.");
655 
656  TEUCHOS_TEST_FOR_EXCEPTION(!problem_->isProblemSet(),TFQMRSolMgrLinearProblemFailure,
657  "Belos::TFQMRSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
658 
659  if (!isSTSet_) {
660  TEUCHOS_TEST_FOR_EXCEPTION( checkStatusTest(),TFQMRSolMgrLinearProblemFailure,
661  "Belos::TFQMRSolMgr::solve(): Linear problem and requested status tests are incompatible.");
662  }
663 
664  // Create indices for the linear systems to be solved.
665  int startPtr = 0;
666  int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
667  int numCurrRHS = blockSize_;
668 
669  std::vector<int> currIdx, currIdx2;
670 
671  // The index set is generated that informs the linear problem that some linear systems are augmented.
672  currIdx.resize( blockSize_ );
673  currIdx2.resize( blockSize_ );
674  for (int i=0; i<numCurrRHS; ++i)
675  { currIdx[i] = startPtr+i; currIdx2[i]=i; }
676 
677  // Inform the linear problem of the current linear system to solve.
678  problem_->setLSIndex( currIdx );
679 
681  // Parameter list
682  Teuchos::ParameterList plist;
683  plist.set("Block Size",blockSize_);
684 
685  // Reset the status test.
686  outputTest_->reset();
687 
688  // Assume convergence is achieved, then let any failed convergence set this to false.
689  bool isConverged = true;
690 
692  // TFQMR solver
693 
694  Teuchos::RCP<TFQMRIter<ScalarType,MV,OP> > tfqmr_iter =
695  Teuchos::rcp( new TFQMRIter<ScalarType,MV,OP>(problem_,printer_,outputTest_,plist) );
696 
697  // Enter solve() iterations
698  {
699 #ifdef BELOS_TEUCHOS_TIME_MONITOR
700  Teuchos::TimeMonitor slvtimer(*timerSolve_);
701 #endif
702 
703  while ( numRHS2Solve > 0 ) {
704  //
705  // Reset the active / converged vectors from this block
706  std::vector<int> convRHSIdx;
707  std::vector<int> currRHSIdx( currIdx );
708  currRHSIdx.resize(numCurrRHS);
709 
710  // Reset the number of iterations.
711  tfqmr_iter->resetNumIters();
712 
713  // Reset the number of calls that the status test output knows about.
714  outputTest_->resetNumCalls();
715 
716  // Get the current residual for this block of linear systems.
717  Teuchos::RCP<MV> R_0 = MVT::CloneViewNonConst( *(Teuchos::rcp_const_cast<MV>(problem_->getInitPrecResVec())), currIdx );
718 
719  // Set the new state and initialize the solver.
721  newstate.R = R_0;
722  tfqmr_iter->initializeTFQMR(newstate);
723 
724  while(1) {
725 
726  // tell tfqmr_iter to iterate
727  try {
728  tfqmr_iter->iterate();
729 
731  //
732  // check convergence first
733  //
735  if ( convTest_->getStatus() == Passed ) {
736  // We have convergence of the linear system.
737  break; // break from while(1){tfqmr_iter->iterate()}
738  }
740  //
741  // check for maximum iterations
742  //
744  else if ( maxIterTest_->getStatus() == Passed ) {
745  // we don't have convergence
746  isConverged = false;
747  break; // break from while(1){tfqmr_iter->iterate()}
748  }
749 
751  //
752  // we returned from iterate(), but none of our status tests Passed.
753  // something is wrong, and it is probably our fault.
754  //
756 
757  else {
758  TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,
759  "Belos::TFQMRSolMgr::solve(): Invalid return from TFQMRIter::iterate().");
760  }
761  }
762  catch (const std::exception &e) {
763  printer_->stream(Errors) << "Error! Caught std::exception in TFQMRIter::iterate() at iteration "
764  << tfqmr_iter->getNumIters() << std::endl
765  << e.what() << std::endl;
766  throw;
767  }
768  }
769 
770  // Update the current solution with the update computed by the iteration object.
771  problem_->updateSolution( tfqmr_iter->getCurrentUpdate(), true );
772 
773  // Inform the linear problem that we are finished with this block linear system.
774  problem_->setCurrLS();
775 
776  // Update indices for the linear systems to be solved.
777  startPtr += numCurrRHS;
778  numRHS2Solve -= numCurrRHS;
779  if ( numRHS2Solve > 0 ) {
780  numCurrRHS = blockSize_;
781 
782  currIdx.resize( blockSize_ );
783  currIdx2.resize( blockSize_ );
784  for (int i=0; i<numCurrRHS; ++i)
785  { currIdx[i] = startPtr+i; currIdx2[i] = i; }
786  // Set the next indices.
787  problem_->setLSIndex( currIdx );
788 
789  // Set the new blocksize for the solver.
790  tfqmr_iter->setBlockSize( blockSize_ );
791  }
792  else {
793  currIdx.resize( numRHS2Solve );
794  }
795 
796  }// while ( numRHS2Solve > 0 )
797 
798  }
799 
800  // print final summary
801  sTest_->print( printer_->stream(FinalSummary) );
802 
803  // print timing information
804 #ifdef BELOS_TEUCHOS_TIME_MONITOR
805  // Calling summarize() can be expensive, so don't call unless the
806  // user wants to print out timing details. summarize() will do all
807  // the work even if it's passed a "black hole" output stream.
808  if (verbosity_ & TimingDetails)
809  Teuchos::TimeMonitor::summarize( printer_->stream(TimingDetails) );
810 #endif
811 
812  // get iteration information for this solve
813  numIters_ = maxIterTest_->getNumIters();
814 
815  // Save the convergence test value ("achieved tolerance") for this
816  // solve. For this solver, convTest_ may either be a single
817  // (implicit) residual norm test, or a combination of two residual
818  // norm tests. In the latter case, the master convergence test
819  // convTest_ is a SEQ combo of the implicit resp. explicit tests.
820  // If the implicit test never passes, then the explicit test won't
821  // ever be executed. This manifests as
822  // expConvTest_->getTestValue()->size() < 1. We deal with this case
823  // by using the values returned by impConvTest_->getTestValue().
824  {
825  // We'll fetch the vector of residual norms one way or the other.
826  const std::vector<MagnitudeType>* pTestValues = NULL;
827  if (expResTest_) {
828  pTestValues = expConvTest_->getTestValue();
829  if (pTestValues == NULL || pTestValues->size() < 1) {
830  pTestValues = impConvTest_->getTestValue();
831  }
832  }
833  else {
834  // Only the implicit residual norm test is being used.
835  pTestValues = impConvTest_->getTestValue();
836  }
837  TEUCHOS_TEST_FOR_EXCEPTION(pTestValues == NULL, std::logic_error,
838  "Belos::TFQMRSolMgr::solve(): The implicit convergence test's "
839  "getTestValue() method returned NULL. Please report this bug to the "
840  "Belos developers.");
841  TEUCHOS_TEST_FOR_EXCEPTION(pTestValues->size() < 1, std::logic_error,
842  "Belos::TMQMRSolMgr::solve(): The implicit convergence test's "
843  "getTestValue() method returned a vector of length zero. Please report "
844  "this bug to the Belos developers.");
845 
846  // FIXME (mfh 12 Dec 2011) Does pTestValues really contain the
847  // achieved tolerances for all vectors in the current solve(), or
848  // just for the vectors from the last deflation?
849  achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
850  }
851 
852  if (!isConverged) {
853  return Unconverged; // return from TFQMRSolMgr::solve()
854  }
855  return Converged; // return from TFQMRSolMgr::solve()
856 }
857 
858 // This method requires the solver manager to return a std::string that describes itself.
859 template<class ScalarType, class MV, class OP>
861 {
862  std::ostringstream oss;
863  oss << "Belos::TFQMRSolMgr<...,"<<Teuchos::ScalarTraits<ScalarType>::name()<<">";
864  oss << "{}";
865  return oss.str();
866 }
867 
868 } // end Belos namespace
869 
870 #endif /* BELOS_TFQMR_SOLMGR_HPP */
ScaleType convertStringToScaleType(const std::string &scaleType)
Convert the given string to its ScaleType enum value.
Definition: BelosTypes.cpp:106
Collection of types and exceptions used within the Belos solvers.
This class implements the preconditioned transpose-free QMR algorithm for solving non-Hermitian linea...
int getNumIters() const override
Get the iteration count for the most recent call to solve().
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.
void setParameters(const Teuchos::RCP< Teuchos::ParameterList > &params) override
Set the parameters the solver manager should use to solve the linear problem.
static constexpr int outputStyle_default_
Teuchos::RCP< StatusTestOutput< ScalarType, MV, OP > > outputTest_
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > convTest_
static constexpr bool expResTest_default_
bool isLOADetected() const override
Whether loss of accuracy was detected during the last solve() invocation.
static constexpr int outputFreq_default_
MultiVecTraits< ScalarType, MV > MVT
A factory class for generating StatusTestOutput objects.
const LinearProblem< ScalarType, MV, OP > & getProblem() const override
Return a reference to the linear problem being solved by this solver manager.
Teuchos::RCP< Teuchos::Time > timerSolve_
Teuchos::RCP< SolverManager< ScalarType, MV, OP > > clone() const override
clone for Inverted Injection (DII)
static constexpr const char * label_default_
An implementation of StatusTestResNorm using a family of residual norms.
static constexpr int maxIters_default_
Belos concrete class for generating iterations with the preconditioned tranpose-free QMR (TFQMR) meth...
Teuchos::RCP< const MV > R
The current residual basis.
static const double convTol
Default convergence tolerance.
Definition: BelosTypes.hpp:293
Belos::StatusTest class for specifying a maximum number of iterations.
Teuchos::RCP< StatusTestGenResNorm< ScalarType, MV, OP > > expConvTest_
Teuchos::RCP< OutputManager< ScalarType > > printer_
std::string description() const override
Method to return description of the TFQMR solver manager.
Teuchos::ScalarTraits< ScalarType > SCT
A factory class for generating StatusTestOutput objects.
Teuchos::RCP< StatusTestGenResNorm< ScalarType, MV, OP > > impConvTest_
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.
Belos::StatusTest for logically combining several status tests.
static constexpr int verbosity_default_
MagnitudeType achievedTol() const override
Tolerance achieved by the last solve() invocation.
A Belos::StatusTest class for specifying a maximum number of iterations.
ResetType
How to reset the solver.
Definition: BelosTypes.hpp:206
Pure virtual base class which describes the basic interface for a solver manager. ...
static const double impTolScale
"Implicit Tolerance Scale Factor"
Definition: BelosTypes.hpp:305
Teuchos::RCP< StatusTestMaxIters< ScalarType, MV, OP > > maxIterTest_
Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > problem_
Teuchos::ScalarTraits< ScalarType >::magnitudeType MagnitudeType
A linear system to solve, and its associated information.
Structure to contain pointers to TFQMRIter state variables.
Class which describes the linear problem to be solved by the iterative solver.
TFQMRSolMgr()
Empty constructor for TFQMRSolMgr. This constructor takes no arguments and sets the default values fo...
static constexpr std::ostream * outputStream_default_
Teuchos::RCP< std::ostream > outputStream_
virtual ~TFQMRSolMgr()
Destructor.
ReturnType
Whether the Belos solve converged for all linear systems.
Definition: BelosTypes.hpp:155
void reset(const ResetType type) override
Performs a reset of the solver manager specified by the ResetType. This informs the solver manager th...
The Belos::SolverManager is a templated virtual base class that defines the basic interface that any ...
The Belos::TFQMRSolMgr provides a powerful and fully-featured solver manager over the TFQMR linear so...
MagnitudeType impTolScale_
static constexpr const char * expResScale_default_
Teuchos::RCP< Teuchos::ParameterList > params_
Teuchos::RCP< StatusTestOutput< ScalarType, MV, OP > > create(const Teuchos::RCP< OutputManager< ScalarType > > &printer, Teuchos::RCP< StatusTest< ScalarType, MV, OP > > test, int mod, int printStates)
Create the StatusTestOutput object specified by the outputStyle.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
Get a parameter list containing the valid parameters for this object.
Belos::StatusTestResNorm for specifying general residual norm stopping criteria.
A class for extending the status testing capabilities of Belos via logical combinations.
MagnitudeType achievedTol_
TFQMRSolMgrLinearProblemFailure(const std::string &what_arg)
static constexpr const char * impResScale_default_
Class which defines basic traits for the operator type.
Parent class to all Belos exceptions.
Definition: BelosTypes.hpp:60
Teuchos::Array< Teuchos::RCP< Teuchos::Time > > getTimers() const
Return the timers for this object.
Default parameters common to most Belos solvers.
Definition: BelosTypes.hpp:283
TFQMRSolMgrLinearProblemFailure is thrown when the linear problem is not setup (i.e.
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > sTest_
Belos header file which uses auto-configuration information to include necessary C++ headers...
OperatorTraits< ScalarType, MV, OP > OPT
void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem) override
Set the linear problem that needs to be solved.
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const override
Get a parameter list containing the current parameters for this object.
Teuchos::ScalarTraits< MagnitudeType > MT