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

Generated on Thu Jul 21 2016 14:43:57 for Belos by doxygen 1.8.11