Belos Package Browser (Single Doxygen Collection)  Development
BelosPseudoBlockTFQMRSolMgr.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_PSEUDO_BLOCK_TFQMR_SOLMGR_HPP
43 #define BELOS_PSEUDO_BLOCK_TFQMR_SOLMGR_HPP
44 
49 #include "BelosConfigDefs.hpp"
50 #include "BelosTypes.hpp"
51 
52 #include "BelosLinearProblem.hpp"
53 #include "BelosSolverManager.hpp"
54 
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  PseudoBlockTFQMRSolMgrLinearProblemFailure(const std::string& what_arg) : BelosError(what_arg)
91  {}};
92 
100  PseudoBlockTFQMRSolMgrOrthoFailure(const std::string& what_arg) : BelosError(what_arg)
101  {}};
102 
103  template<class ScalarType, class MV, class OP>
104  class PseudoBlockTFQMRSolMgr : public SolverManager<ScalarType,MV,OP> {
105 
106  private:
112 
113  public:
114 
116 
117 
124 
143 
147 
149 
150 
152  return *problem_;
153  }
154 
158 
162 
169  return Teuchos::tuple(timerSolve_);
170  }
171 
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 
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.
256 
257  // Output manager.
260 
261  // Status test.
267 
268  // Current parameter list.
270 
271  // Default solver values.
274  static const int maxIters_default_;
275  static const bool expResTest_default_;
276  static const int verbosity_default_;
277  static const int outputStyle_default_;
278  static const int outputFreq_default_;
279  static const int defQuorum_default_;
280  static const std::string impResScale_default_;
281  static const std::string expResScale_default_;
282  static const std::string label_default_;
284 
285  // Current solver values.
291 
292  // Timers.
293  std::string label_;
295 
296  // Internal state variables.
298  };
299 
300 
301 // Default solver values.
302 template<class ScalarType, class MV, class OP>
304 
305 template<class ScalarType, class MV, class OP>
307 
308 template<class ScalarType, class MV, class OP>
310 
311 template<class ScalarType, class MV, class OP>
313 
314 template<class ScalarType, class MV, class OP>
316 
317 template<class ScalarType, class MV, class OP>
319 
320 template<class ScalarType, class MV, class OP>
322 
323 template<class ScalarType, class MV, class OP>
325 
326 template<class ScalarType, class MV, class OP>
327 const std::string PseudoBlockTFQMRSolMgr<ScalarType,MV,OP>::impResScale_default_ = "Norm of Preconditioned Initial Residual";
328 
329 template<class ScalarType, class MV, class OP>
330 const std::string PseudoBlockTFQMRSolMgr<ScalarType,MV,OP>::expResScale_default_ = "Norm of Initial Residual";
331 
332 template<class ScalarType, class MV, class OP>
334 
335 template<class ScalarType, class MV, class OP>
337 
338 
339 // Empty Constructor
340 template<class ScalarType, class MV, class OP>
342  outputStream_(outputStream_default_),
343  convtol_(convtol_default_),
344  impTolScale_(impTolScale_default_),
345  achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
346  maxIters_(maxIters_default_),
347  numIters_(0),
348  verbosity_(verbosity_default_),
349  outputStyle_(outputStyle_default_),
350  outputFreq_(outputFreq_default_),
351  defQuorum_(defQuorum_default_),
352  expResTest_(expResTest_default_),
353  impResScale_(impResScale_default_),
354  expResScale_(expResScale_default_),
355  label_(label_default_),
356  isSet_(false),
357  isSTSet_(false)
358 {}
359 
360 
361 // Basic Constructor
362 template<class ScalarType, class MV, class OP>
366  problem_(problem),
367  outputStream_(outputStream_default_),
368  convtol_(convtol_default_),
369  impTolScale_(impTolScale_default_),
370  achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
371  maxIters_(maxIters_default_),
372  numIters_(0),
373  verbosity_(verbosity_default_),
374  outputStyle_(outputStyle_default_),
375  outputFreq_(outputFreq_default_),
376  defQuorum_(defQuorum_default_),
377  expResTest_(expResTest_default_),
378  impResScale_(impResScale_default_),
379  expResScale_(expResScale_default_),
380  label_(label_default_),
381  isSet_(false),
382  isSTSet_(false)
383 {
384  TEUCHOS_TEST_FOR_EXCEPTION(problem_ == Teuchos::null, std::invalid_argument, "Problem not given to solver manager.");
385 
386  // If the parameter list pointer is null, then set the current parameters to the default parameter list.
387  if ( !is_null(pl) ) {
388  setParameters( pl );
389  }
390 }
391 
392 template<class ScalarType, class MV, class OP>
394 {
395  // Create the internal parameter list if ones doesn't already exist.
396  if (params_ == Teuchos::null) {
397  params_ = Teuchos::rcp( new Teuchos::ParameterList(*getValidParameters()) );
398  }
399  else {
400  params->validateParameters(*getValidParameters());
401  }
402 
403  // Check for maximum number of iterations
404  if (params->isParameter("Maximum Iterations")) {
405  maxIters_ = params->get("Maximum Iterations",maxIters_default_);
406 
407  // Update parameter in our list and in status test.
408  params_->set("Maximum Iterations", maxIters_);
409  if (maxIterTest_!=Teuchos::null)
410  maxIterTest_->setMaxIters( maxIters_ );
411  }
412 
413  // Check to see if the timer label changed.
414  if (params->isParameter("Timer Label")) {
415  std::string tempLabel = params->get("Timer Label", label_default_);
416 
417  // Update parameter in our list and solver timer
418  if (tempLabel != label_) {
419  label_ = tempLabel;
420  params_->set("Timer Label", label_);
421  std::string solveLabel = label_ + ": PseudoBlockTFQMRSolMgr total solve time";
422 #ifdef BELOS_TEUCHOS_TIME_MONITOR
423  timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
424 #endif
425  }
426  }
427 
428  // Check for a change in verbosity level
429  if (params->isParameter("Verbosity")) {
430  if (Teuchos::isParameterType<int>(*params,"Verbosity")) {
431  verbosity_ = params->get("Verbosity", verbosity_default_);
432  } else {
433  verbosity_ = (int)Teuchos::getParameter<Belos::MsgType>(*params,"Verbosity");
434  }
435 
436  // Update parameter in our list.
437  params_->set("Verbosity", verbosity_);
438  if (printer_ != Teuchos::null)
439  printer_->setVerbosity(verbosity_);
440  }
441 
442  // Check for a change in output style
443  if (params->isParameter("Output Style")) {
444  if (Teuchos::isParameterType<int>(*params,"Output Style")) {
445  outputStyle_ = params->get("Output Style", outputStyle_default_);
446  } else {
447  outputStyle_ = (int)Teuchos::getParameter<Belos::OutputType>(*params,"Output Style");
448  }
449 
450  // Reconstruct the convergence test if the explicit residual test is not being used.
451  params_->set("Output Style", outputStyle_);
452  isSTSet_ = false;
453  }
454 
455  // output stream
456  if (params->isParameter("Output Stream")) {
457  outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> >(*params,"Output Stream");
458 
459  // Update parameter in our list.
460  params_->set("Output Stream", outputStream_);
461  if (printer_ != Teuchos::null)
462  printer_->setOStream( outputStream_ );
463  }
464 
465  // frequency level
466  if (verbosity_ & Belos::StatusTestDetails) {
467  if (params->isParameter("Output Frequency")) {
468  outputFreq_ = params->get("Output Frequency", outputFreq_default_);
469  }
470 
471  // Update parameter in out list and output status test.
472  params_->set("Output Frequency", outputFreq_);
473  if (outputTest_ != Teuchos::null)
474  outputTest_->setOutputFrequency( outputFreq_ );
475  }
476 
477  // Create output manager if we need to.
478  if (printer_ == Teuchos::null) {
479  printer_ = Teuchos::rcp( new OutputManager<ScalarType>(verbosity_, outputStream_) );
480  }
481 
482  // Check for convergence tolerance
483  if (params->isParameter("Convergence Tolerance")) {
484  convtol_ = params->get("Convergence Tolerance",convtol_default_);
485 
486  // Update parameter in our list and residual tests.
487  params_->set("Convergence Tolerance", convtol_);
488  isSTSet_ = false;
489  }
490 
491  if (params->isParameter("Implicit Tolerance Scale Factor")) {
492  impTolScale_ = params->get("Implicit Tolerance Scale Factor",impTolScale_default_);
493 
494  // Update parameter in our list.
495  params_->set("Implicit Tolerance Scale Factor", impTolScale_);
496  isSTSet_ = false;
497  }
498 
499  if (params->isParameter("Implicit Residual Scaling")) {
500  std::string tempImpResScale = Teuchos::getParameter<std::string>( *params, "Implicit Residual Scaling" );
501 
502  // Only update the scaling if it's different.
503  if (impResScale_ != tempImpResScale) {
504  impResScale_ = tempImpResScale;
505 
506  // Update parameter in our list.
507  params_->set("Implicit Residual Scaling", impResScale_);
508  isSTSet_ = false;
509  }
510  }
511 
512  if (params->isParameter("Explicit Residual Scaling")) {
513  std::string tempExpResScale = Teuchos::getParameter<std::string>( *params, "Explicit Residual Scaling" );
514 
515  // Only update the scaling if it's different.
516  if (expResScale_ != tempExpResScale) {
517  expResScale_ = tempExpResScale;
518 
519  // Update parameter in our list.
520  params_->set("Explicit Residual Scaling", expResScale_);
521  isSTSet_ = false;
522  }
523  }
524 
525  if (params->isParameter("Explicit Residual Test")) {
526  expResTest_ = Teuchos::getParameter<bool>( *params,"Explicit Residual Test" );
527 
528  // Reconstruct the convergence test if the explicit residual test is not being used.
529  params_->set("Explicit Residual Test", expResTest_);
530  if (expConvTest_ == Teuchos::null) {
531  isSTSet_ = false;
532  }
533  }
534 
535  // Get the deflation quorum, or number of converged systems before deflation is allowed
536  if (params->isParameter("Deflation Quorum")) {
537  defQuorum_ = params->get("Deflation Quorum", defQuorum_);
538  params_->set ("Deflation Quorum", defQuorum_);
539  if (! impConvTest_.is_null ()) {
540  impConvTest_->setQuorum (defQuorum_);
541  }
542  if (! expConvTest_.is_null ()) {
543  expConvTest_->setQuorum (defQuorum_);
544  }
545  }
546 
547  // Create the timer if we need to.
548  if (timerSolve_ == Teuchos::null) {
549  std::string solveLabel = label_ + ": PseudoBlockTFQMRSolMgr total solve time";
550 #ifdef BELOS_TEUCHOS_TIME_MONITOR
551  timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
552 #endif
553  }
554 
555  // Inform the solver manager that the current parameters were set.
556  isSet_ = true;
557 }
558 
559 
560 // Check the status test versus the defined linear problem
561 template<class ScalarType, class MV, class OP>
563 
564  typedef Belos::StatusTestCombo<ScalarType,MV,OP> StatusTestCombo_t;
565  typedef Belos::StatusTestGenResNorm<ScalarType,MV,OP> StatusTestGenResNorm_t;
566 
567  // Basic test checks maximum iterations and native residual.
568  maxIterTest_ = Teuchos::rcp( new StatusTestMaxIters<ScalarType,MV,OP>( maxIters_ ) );
569 
570  if (expResTest_) {
571 
572  // Implicit residual test, using the native residual to determine if convergence was achieved.
573  Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
574  Teuchos::rcp( new StatusTestGenResNorm_t( impTolScale_*convtol_, defQuorum_ ) );
575  tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
576  impConvTest_ = tmpImpConvTest;
577 
578  // Explicit residual test once the native residual is below the tolerance
579  Teuchos::RCP<StatusTestGenResNorm_t> tmpExpConvTest =
580  Teuchos::rcp( new StatusTestGenResNorm_t( convtol_, defQuorum_ ) );
581  tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit, Belos::TwoNorm );
582  tmpExpConvTest->defineScaleForm( convertStringToScaleType(expResScale_), Belos::TwoNorm );
583  expConvTest_ = tmpExpConvTest;
584 
585  // The convergence test is a combination of the "cheap" implicit test and explicit test.
586  convTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
587  }
588  else {
589 
590  // Implicit residual test, using the native residual to determine if convergence was achieved.
591  Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
592  Teuchos::rcp( new StatusTestGenResNorm_t( convtol_, defQuorum_ ) );
593  tmpImpConvTest->defineScaleForm( convertStringToScaleType(impResScale_), Belos::TwoNorm );
594  impConvTest_ = tmpImpConvTest;
595 
596  // Set the explicit and total convergence test to this implicit test that checks for accuracy loss.
597  expConvTest_ = impConvTest_;
598  convTest_ = impConvTest_;
599  }
600  sTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
601 
602  // Create the status test output class.
603  // This class manages and formats the output from the status test.
604  StatusTestOutputFactory<ScalarType,MV,OP> stoFactory( outputStyle_ );
605  outputTest_ = stoFactory.create( printer_, sTest_, outputFreq_, Passed+Failed+Undefined );
606 
607  // Set the solver string for the output test
608  std::string solverDesc = " Pseudo Block TFQMR ";
609  outputTest_->setSolverDesc( solverDesc );
610 
611 
612  // The status test is now set.
613  isSTSet_ = true;
614 
615  return false;
616 }
617 
618 
619 template<class ScalarType, class MV, class OP>
622 {
624 
625  // Set all the valid parameters and their default values.
626  if(is_null(validPL)) {
627  Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
628  pl->set("Convergence Tolerance", convtol_default_,
629  "The relative residual tolerance that needs to be achieved by the\n"
630  "iterative solver in order for the linear system to be declared converged.");
631  pl->set("Implicit Tolerance Scale Factor", impTolScale_default_,
632  "The scale factor used by the implicit residual test when explicit residual\n"
633  "testing is used. May enable faster convergence when TFQMR bound is too loose.");
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("Deflation Quorum", defQuorum_default_,
647  "The number of linear systems that need to converge before they are deflated.");
648  pl->set("Output Stream", outputStream_default_,
649  "A reference-counted pointer to the output stream where all\n"
650  "solver output is sent.");
651  pl->set("Explicit Residual Test", expResTest_default_,
652  "Whether the explicitly computed residual should be used in the convergence test.");
653  pl->set("Implicit Residual Scaling", impResScale_default_,
654  "The type of scaling used in the implicit residual convergence test.");
655  pl->set("Explicit Residual Scaling", expResScale_default_,
656  "The type of scaling used in the explicit residual convergence test.");
657  pl->set("Timer Label", label_default_,
658  "The string to use as a prefix for the timer labels.");
659  // pl->set("Restart Timers", restartTimers_);
660  validPL = pl;
661  }
662  return validPL;
663 }
664 
665 
666 // solve()
667 template<class ScalarType, class MV, class OP>
669 
670  // Set the current parameters if they were not set before.
671  // NOTE: This may occur if the user generated the solver manager with the default constructor and
672  // then didn't set any parameters using setParameters().
673  if (!isSet_) {
674  setParameters(Teuchos::parameterList(*getValidParameters()));
675  }
676 
678  "Belos::PseudoBlockTFQMRSolMgr::solve(): Linear problem is not a valid object.");
679 
681  "Belos::PseudoBlockTFQMRSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
682 
683  if (!isSTSet_) {
685  "Belos::PseudoBlockTFQMRSolMgr::solve(): Linear problem and requested status tests are incompatible.");
686  }
687 
688  // Create indices for the linear systems to be solved.
689  int startPtr = 0;
690  int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
691  int numCurrRHS = numRHS2Solve;
692 
693  std::vector<int> currIdx( numRHS2Solve );
694  for (int i=0; i<numRHS2Solve; ++i) {
695  currIdx[i] = startPtr+i;
696  }
697 
698  // Inform the linear problem of the current linear system to solve.
699  problem_->setLSIndex( currIdx );
700 
702  // Parameter list
704 
705  // Reset the status test.
706  outputTest_->reset();
707 
708  // Assume convergence is achieved, then let any failed convergence set this to false.
709  bool isConverged = true;
710 
712  // TFQMR solver
713 
715  Teuchos::rcp( new PseudoBlockTFQMRIter<ScalarType,MV,OP>(problem_,printer_,outputTest_,plist) );
716 
717  // Enter solve() iterations
718  {
719 #ifdef BELOS_TEUCHOS_TIME_MONITOR
720  Teuchos::TimeMonitor slvtimer(*timerSolve_);
721 #endif
722 
723  while ( numRHS2Solve > 0 ) {
724  //
725  // Reset the active / converged vectors from this block
726  std::vector<int> convRHSIdx;
727  std::vector<int> currRHSIdx( currIdx );
728  currRHSIdx.resize(numCurrRHS);
729 
730  // Reset the number of iterations.
731  block_tfqmr_iter->resetNumIters();
732 
733  // Reset the number of calls that the status test output knows about.
734  outputTest_->resetNumCalls();
735 
736  // Get the current residual for this block of linear systems.
737  Teuchos::RCP<MV> R_0 = MVT::CloneViewNonConst( *(Teuchos::rcp_const_cast<MV>(problem_->getInitPrecResVec())), currIdx );
738 
739  // Set the new state and initialize the solver.
741  newstate.Rtilde = R_0;
742  block_tfqmr_iter->initializeTFQMR(newstate);
743 
744  while(1) {
745 
746  // tell block_tfqmr_iter to iterate
747  try {
748  block_tfqmr_iter->iterate();
749 
751  //
752  // check convergence first
753  //
755  if ( convTest_->getStatus() == Passed ) {
756 
757  // Figure out which linear systems converged.
758  std::vector<int> convIdx = expConvTest_->convIndices();
759 
760  // If the number of converged linear systems is equal to the
761  // number of current linear systems, then we are done with this block.
762  if (convIdx.size() == currRHSIdx.size())
763  break; // break from while(1){block_tfqmr_iter->iterate()}
764 
765  // Update the current solution with the update computed by the iteration object.
766  problem_->updateSolution( block_tfqmr_iter->getCurrentUpdate(), true );
767 
768  // Inform the linear problem that we are finished with this current linear system.
769  problem_->setCurrLS();
770 
771  // Reset currRHSIdx to have the right-hand sides that are left to converge for this block.
772  int have = 0;
773  std::vector<int> unconvIdx( currRHSIdx.size() );
774  for (unsigned int i=0; i<currRHSIdx.size(); ++i) {
775  bool found = false;
776  for (unsigned int j=0; j<convIdx.size(); ++j) {
777  if (currRHSIdx[i] == convIdx[j]) {
778  found = true;
779  break;
780  }
781  }
782  if (!found) {
783  unconvIdx[have] = i;
784  currRHSIdx[have++] = currRHSIdx[i];
785  }
786  }
787  unconvIdx.resize(have);
788  currRHSIdx.resize(have);
789 
790  // Set the remaining indices after deflation.
791  problem_->setLSIndex( currRHSIdx );
792 
793  // Get the current residual vector.
794  // Set the new state and initialize the solver.
795  PseudoBlockTFQMRIterState<ScalarType,MV> currentState = block_tfqmr_iter->getState();
796 
797  // Set the new state and initialize the solver.
799 
800  // Copy over the vectors.
801  defstate.Rtilde = MVT::CloneView( *currentState.Rtilde, unconvIdx);
802  defstate.U = MVT::CloneView( *currentState.U, unconvIdx );
803  defstate.AU = MVT::CloneView( *currentState.AU, unconvIdx );
804  defstate.V = MVT::CloneView( *currentState.V, unconvIdx );
805  defstate.W = MVT::CloneView( *currentState.W, unconvIdx );
806  defstate.D = MVT::CloneView( *currentState.D, unconvIdx );
807 
808  // Copy over the scalars.
809  for (std::vector<int>::iterator uIter = unconvIdx.begin(); uIter != unconvIdx.end(); uIter++)
810  {
811  defstate.alpha.push_back( currentState.alpha[ *uIter ] );
812  defstate.eta.push_back( currentState.eta[ *uIter ] );
813  defstate.rho.push_back( currentState.rho[ *uIter ] );
814  defstate.tau.push_back( currentState.tau[ *uIter ] );
815  defstate.theta.push_back( currentState.theta[ *uIter ] );
816  }
817 
818  block_tfqmr_iter->initializeTFQMR(defstate);
819  }
821  //
822  // check for maximum iterations
823  //
825  else if ( maxIterTest_->getStatus() == Passed ) {
826  // we don't have convergence
827  isConverged = false;
828  break; // break from while(1){block_tfqmr_iter->iterate()}
829  }
830 
832  //
833  // we returned from iterate(), but none of our status tests Passed.
834  // something is wrong, and it is probably our fault.
835  //
837 
838  else {
839  TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,
840  "Belos::PseudoBlockTFQMRSolMgr::solve(): Invalid return from PseudoBlockTFQMRIter::iterate().");
841  }
842  }
843  catch (const std::exception &e) {
844  printer_->stream(Errors) << "Error! Caught std::exception in PseudoBlockTFQMRIter::iterate() at iteration "
845  << block_tfqmr_iter->getNumIters() << std::endl
846  << e.what() << std::endl;
847  throw;
848  }
849  }
850 
851  // Update the current solution with the update computed by the iteration object.
852  problem_->updateSolution( block_tfqmr_iter->getCurrentUpdate(), true );
853 
854  // Inform the linear problem that we are finished with this block linear system.
855  problem_->setCurrLS();
856 
857  // Update indices for the linear systems to be solved.
858  startPtr += numCurrRHS;
859  numRHS2Solve -= numCurrRHS;
860  if ( numRHS2Solve > 0 ) {
861  numCurrRHS = numRHS2Solve;
862  currIdx.resize( numCurrRHS );
863  for (int i=0; i<numCurrRHS; ++i)
864  { currIdx[i] = startPtr+i; }
865 
866  // Adapt the status test quorum if we need to.
867  if (defQuorum_ > numCurrRHS) {
868  if (impConvTest_ != Teuchos::null)
869  impConvTest_->setQuorum( numCurrRHS );
870  if (expConvTest_ != Teuchos::null)
871  expConvTest_->setQuorum( numCurrRHS );
872  }
873 
874  // Set the next indices.
875  problem_->setLSIndex( currIdx );
876  }
877  else {
878  currIdx.resize( numRHS2Solve );
879  }
880 
881  }// while ( numRHS2Solve > 0 )
882 
883  }
884 
885  // print final summary
886  sTest_->print( printer_->stream(FinalSummary) );
887 
888  // print timing information
889 #ifdef BELOS_TEUCHOS_TIME_MONITOR
890  // Calling summarize() can be expensive, so don't call unless the
891  // user wants to print out timing details. summarize() will do all
892  // the work even if it's passed a "black hole" output stream.
893  if (verbosity_ & TimingDetails)
894  Teuchos::TimeMonitor::summarize( printer_->stream(TimingDetails) );
895 #endif
896 
897  // get iteration information for this solve
898  numIters_ = maxIterTest_->getNumIters();
899 
900  // Save the convergence test value ("achieved tolerance") for this
901  // solve. For this solver, convTest_ may either be a single
902  // (implicit) residual norm test, or a combination of two residual
903  // norm tests. In the latter case, the master convergence test
904  // convTest_ is a SEQ combo of the implicit resp. explicit tests.
905  // If the implicit test never passes, then the explicit test won't
906  // ever be executed. This manifests as
907  // expConvTest_->getTestValue()->size() < 1. We deal with this case
908  // by using the values returned by impConvTest_->getTestValue().
909  {
910  // We'll fetch the vector of residual norms one way or the other.
911  const std::vector<MagnitudeType>* pTestValues = NULL;
912  if (expResTest_) {
913  pTestValues = expConvTest_->getTestValue();
914  if (pTestValues == NULL || pTestValues->size() < 1) {
915  pTestValues = impConvTest_->getTestValue();
916  }
917  }
918  else {
919  // Only the implicit residual norm test is being used.
920  pTestValues = impConvTest_->getTestValue();
921  }
922  TEUCHOS_TEST_FOR_EXCEPTION(pTestValues == NULL, std::logic_error,
923  "Belos::PseudoBlockTFQMRSolMgr::solve(): The implicit convergence test's "
924  "getTestValue() method returned NULL. Please report this bug to the "
925  "Belos developers.");
926  TEUCHOS_TEST_FOR_EXCEPTION(pTestValues->size() < 1, std::logic_error,
927  "Belos::TMQMRSolMgr::solve(): The implicit convergence test's "
928  "getTestValue() method returned a vector of length zero. Please report "
929  "this bug to the Belos developers.");
930 
931  // FIXME (mfh 12 Dec 2011) Does pTestValues really contain the
932  // achieved tolerances for all vectors in the current solve(), or
933  // just for the vectors from the last deflation?
934  achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
935  }
936 
937  if (!isConverged) {
938  return Unconverged; // return from PseudoBlockTFQMRSolMgr::solve()
939  }
940  return Converged; // return from PseudoBlockTFQMRSolMgr::solve()
941 }
942 
943 // This method requires the solver manager to return a std::string that describes itself.
944 template<class ScalarType, class MV, class OP>
946 {
947  std::ostringstream oss;
948  oss << "Belos::PseudoBlockTFQMRSolMgr<...,"<<Teuchos::ScalarTraits<ScalarType>::name()<<">";
949  oss << "{}";
950  return oss.str();
951 }
952 
953 } // end Belos namespace
954 
955 #endif /* BELOS_PSEUDO_BLOCK_TFQMR_SOLMGR_HPP */
ScaleType convertStringToScaleType(const std::string &scaleType)
Convert the given string to its ScaleType enum value.
Definition: BelosTypes.cpp:87
int getNumIters() const
Get the iteration count for the most recent call to solve().
Collection of types and exceptions used within the Belos solvers.
Teuchos::Array< Teuchos::RCP< Teuchos::Time > > getTimers() const
Return the timers for this object.
PseudoBlockTFQMRSolMgrOrthoFailure is thrown when the orthogonalization manager is unable to generate...
Belos&#39;s basic output manager for sending information of select verbosity levels to the appropriate ou...
void reset(const ResetType type)
Performs a reset of the solver manager specified by the ResetType. This informs the solver manager th...
Class which manages the output and verbosity of the Belos solvers.
Teuchos::RCP< std::ostream > outputStream_
bool is_null(const boost::shared_ptr< T > &p)
This class implements the preconditioned transpose-free QMR algorithm for solving non-Hermitian linea...
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Get a parameter list containing the valid parameters for this object.
void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem)
Set the linear problem that needs to be solved.
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > convTest_
T & get(ParameterList &l, const std::string &name)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
Structure to contain pointers to PseudoBlockTFQMRIter state variables.
static RCP< Time > getNewCounter(const std::string &name)
bool is_null(const std::shared_ptr< T > &p)
Teuchos::RCP< Teuchos::Time > timerSolve_
const LinearProblem< ScalarType, MV, OP > & getProblem() const
Return a reference to the linear problem being solved by this solver manager.
Teuchos::RCP< StatusTestGenResNorm< ScalarType, MV, OP > > impConvTest_
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
A factory class for generating StatusTestOutput objects.
MultiVecTraits< ScalarType, MV > MVT
Belos concrete class for generating iterations with the preconditioned tranpose-free QMR (TFQMR) meth...
An implementation of StatusTestResNorm using a family of residual norms.
Belos::StatusTest class for specifying a maximum number of iterations.
static std::string name()
void setParameters(const Teuchos::RCP< Teuchos::ParameterList > &params)
Set the parameters the solver manager should use to solve the linear problem.
std::string description() const
Method to return description of the pseudo-block TFQMR solver manager.
A factory class for generating StatusTestOutput objects.
Traits class which defines basic operations on multivectors.
Belos::StatusTest for logically combining several status tests.
void validateParameters(ParameterList const &validParamList, int const depth=1000, EValidateUsed const validateUsed=VALIDATE_USED_ENABLED, EValidateDefaults const validateDefaults=VALIDATE_DEFAULTS_ENABLED) const
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > sTest_
A Belos::StatusTest class for specifying a maximum number of iterations.
PseudoBlockTFQMRSolMgr()
Empty constructor for PseudoBlockTFQMRSolMgr. This constructor takes no arguments and sets the defaul...
PseudoBlockTFQMRSolMgrLinearProblemFailure is thrown when the linear problem is not setup (i...
ResetType
How to reset the solver.
Definition: BelosTypes.hpp:205
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. ...
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)
OperatorTraits< ScalarType, MV, OP > OPT
PseudoBlockTFQMRSolMgrOrthoFailure(const std::string &what_arg)
static const MagnitudeType convtol_default_
ReturnType solve()
This method performs possibly repeated calls to the underlying linear solver&#39;s iterate() routine unti...
Teuchos::ScalarTraits< ScalarType > SCT
A linear system to solve, and its associated information.
Class which describes the linear problem to be solved by the iterative solver.
The Belos::PseudoBlockTFQMRSolMgr provides a powerful and fully-featured solver manager over the pseu...
ReturnType
Whether the Belos solve converged for all linear systems.
Definition: BelosTypes.hpp:154
Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > problem_
The Belos::SolverManager is a templated virtual base class that defines the basic interface that any ...
MagnitudeType achievedTol() const
Tolerance achieved by the last solve() invocation.
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const
Get a parameter list containing the current parameters for this object.
Teuchos::RCP< StatusTestGenResNorm< ScalarType, MV, OP > > expConvTest_
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< StatusTestOutput< ScalarType, MV, OP > > outputTest_
Belos::StatusTestResNorm for specifying general residual norm stopping criteria.
Teuchos::RCP< const MV > W
The current residual basis.
bool isLOADetected() const
Whether loss of accuracy was detected during the last solve() invocation.
A class for extending the status testing capabilities of Belos via logical combinations.
bool isParameter(const std::string &name) const
Class which defines basic traits for the operator type.
Parent class to all Belos exceptions.
Definition: BelosTypes.hpp:60
Teuchos::RCP< OutputManager< ScalarType > > printer_
Teuchos::ScalarTraits< MagnitudeType > MT
Teuchos::ScalarTraits< ScalarType >::magnitudeType MagnitudeType
Belos header file which uses auto-configuration information to include necessary C++ headers...
static const MagnitudeType impTolScale_default_
static const Teuchos::RCP< std::ostream > outputStream_default_
Teuchos::RCP< StatusTestMaxIters< ScalarType, MV, OP > > maxIterTest_