Belos  Version of the Day
BelosPseudoBlockStochasticCGSolMgr.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_STOCHASTIC_CG_SOLMGR_HPP
43 #define BELOS_PSEUDO_BLOCK_STOCHASTIC_CG_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 #include "Teuchos_BLAS.hpp"
62 #ifdef BELOS_TEUCHOS_TIME_MONITOR
63 #include "Teuchos_TimeMonitor.hpp"
64 #endif
65 
75 namespace Belos {
76 
78 
79 
87  PseudoBlockStochasticCGSolMgrLinearProblemFailure(const std::string& what_arg) : BelosError(what_arg)
88  {}};
89 
97  PseudoBlockStochasticCGSolMgrOrthoFailure(const std::string& what_arg) : BelosError(what_arg)
98  {}};
99 
100  template<class ScalarType, class MV, class OP>
101  class PseudoBlockStochasticCGSolMgr : public SolverManager<ScalarType,MV,OP> {
102 
103  private:
107  typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
109 
110  public:
111 
113 
114 
121 
133 
137 
139 
140 
142  return *problem_;
143  }
144 
148 
152 
159  return Teuchos::tuple(timerSolve_);
160  }
161 
163  int getNumIters() const {
164  return numIters_;
165  }
166 
170  bool isLOADetected() const { return false; }
171 
173 
175 
176 
178  void setProblem( const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem ) { problem_ = problem; }
179 
182 
184 
186 
187 
191  void reset( const ResetType type ) { if ((type & Belos::Problem) && !Teuchos::is_null(problem_)) problem_->setProblem(); }
193 
195 
196 
214  ReturnType solve();
215 
217 
220 
223 
225  std::string description() const;
226 
228 
229  private:
230 
231  // Linear problem.
233 
234  // Output manager.
236  Teuchos::RCP<std::ostream> outputStream_;
237 
238  // Status test.
243 
244  // Orthogonalization manager.
246 
247  // Current parameter list.
249 
255  mutable Teuchos::RCP<const Teuchos::ParameterList> validParams_;
256 
257  // Default solver values.
258  static const MagnitudeType convtol_default_;
259  static const int maxIters_default_;
260  static const bool assertPositiveDefiniteness_default_;
261  static const bool showMaxResNormOnly_default_;
262  static const int verbosity_default_;
263  static const int outputStyle_default_;
264  static const int outputFreq_default_;
265  static const int defQuorum_default_;
266  static const std::string resScale_default_;
267  static const std::string label_default_;
268  static const Teuchos::RCP<std::ostream> outputStream_default_;
269 
270  // Current solver values.
271  MagnitudeType convtol_;
272  int maxIters_, numIters_;
273  int verbosity_, outputStyle_, outputFreq_, defQuorum_;
274  bool assertPositiveDefiniteness_, showMaxResNormOnly_;
275  std::string resScale_;
276 
277  // Timers.
278  std::string label_;
279  Teuchos::RCP<Teuchos::Time> timerSolve_;
280 
281  // Internal state variables.
282  bool isSet_;
283 
284  // Stashed copy of the stochastic vector
285  Teuchos::RCP<MV> Y_;
286 
287  };
288 
289 
290 // Default solver values.
291 template<class ScalarType, class MV, class OP>
292 const typename PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::MagnitudeType PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::convtol_default_ = 1e-8;
293 
294 template<class ScalarType, class MV, class OP>
295 const int PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::maxIters_default_ = 1000;
296 
297 template<class ScalarType, class MV, class OP>
298 const bool PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::assertPositiveDefiniteness_default_ = true;
299 
300 template<class ScalarType, class MV, class OP>
301 const bool PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::showMaxResNormOnly_default_ = false;
302 
303 template<class ScalarType, class MV, class OP>
304 const int PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::verbosity_default_ = Belos::Errors;
305 
306 template<class ScalarType, class MV, class OP>
307 const int PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::outputStyle_default_ = Belos::General;
308 
309 template<class ScalarType, class MV, class OP>
310 const int PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::outputFreq_default_ = -1;
311 
312 template<class ScalarType, class MV, class OP>
313 const int PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::defQuorum_default_ = 1;
314 
315 template<class ScalarType, class MV, class OP>
316 const std::string PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::resScale_default_ = "Norm of Initial Residual";
317 
318 template<class ScalarType, class MV, class OP>
319 const std::string PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::label_default_ = "Belos";
320 
321 template<class ScalarType, class MV, class OP>
322 const Teuchos::RCP<std::ostream> PseudoBlockStochasticCGSolMgr<ScalarType,MV,OP>::outputStream_default_ = Teuchos::rcp(&std::cout,false);
323 
324 
325 // Empty Constructor
326 template<class ScalarType, class MV, class OP>
328  outputStream_(outputStream_default_),
329  convtol_(convtol_default_),
330  maxIters_(maxIters_default_),
331  numIters_(0),
332  verbosity_(verbosity_default_),
333  outputStyle_(outputStyle_default_),
334  outputFreq_(outputFreq_default_),
335  defQuorum_(defQuorum_default_),
336  assertPositiveDefiniteness_(assertPositiveDefiniteness_default_),
337  showMaxResNormOnly_(showMaxResNormOnly_default_),
338  resScale_(resScale_default_),
339  label_(label_default_),
340  isSet_(false)
341 {}
342 
343 // Basic Constructor
344 template<class ScalarType, class MV, class OP>
348  problem_(problem),
349  outputStream_(outputStream_default_),
350  convtol_(convtol_default_),
351  maxIters_(maxIters_default_),
352  numIters_(0),
353  verbosity_(verbosity_default_),
354  outputStyle_(outputStyle_default_),
355  outputFreq_(outputFreq_default_),
356  defQuorum_(defQuorum_default_),
357  assertPositiveDefiniteness_(assertPositiveDefiniteness_default_),
358  showMaxResNormOnly_(showMaxResNormOnly_default_),
359  resScale_(resScale_default_),
360  label_(label_default_),
361  isSet_(false)
362 {
364  problem_.is_null (), std::invalid_argument,
365  "Belos::PseudoBlockStochasticCGSolMgr two-argument constructor: "
366  "'problem' is null. You must supply a non-null Belos::LinearProblem "
367  "instance when calling this constructor.");
368 
369  if (! pl.is_null ()) {
370  // Set the parameters using the list that was passed in.
371  setParameters (pl);
372  }
373 }
374 
375 template<class ScalarType, class MV, class OP>
377 {
379  using Teuchos::parameterList;
380  using Teuchos::RCP;
381 
382  RCP<const ParameterList> defaultParams = getValidParameters();
383 
384  // Create the internal parameter list if one doesn't already exist.
385  if (params_.is_null()) {
386  params_ = parameterList (*defaultParams);
387  } else {
388  params->validateParameters (*defaultParams);
389  }
390 
391  // Check for maximum number of iterations
392  if (params->isParameter("Maximum Iterations")) {
393  maxIters_ = params->get("Maximum Iterations",maxIters_default_);
394 
395  // Update parameter in our list and in status test.
396  params_->set("Maximum Iterations", maxIters_);
397  if (maxIterTest_!=Teuchos::null)
398  maxIterTest_->setMaxIters( maxIters_ );
399  }
400 
401  // Check if positive definiteness assertions are to be performed
402  if (params->isParameter("Assert Positive Definiteness")) {
403  assertPositiveDefiniteness_ = params->get("Assert Positive Definiteness",assertPositiveDefiniteness_default_);
404 
405  // Update parameter in our list.
406  params_->set("Assert Positive Definiteness", assertPositiveDefiniteness_);
407  }
408 
409  // Check to see if the timer label changed.
410  if (params->isParameter("Timer Label")) {
411  std::string tempLabel = params->get("Timer Label", label_default_);
412 
413  // Update parameter in our list and solver timer
414  if (tempLabel != label_) {
415  label_ = tempLabel;
416  params_->set("Timer Label", label_);
417  std::string solveLabel = label_ + ": PseudoBlockStochasticCGSolMgr total solve time";
418 #ifdef BELOS_TEUCHOS_TIME_MONITOR
419  timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
420 #endif
421  if (ortho_ != Teuchos::null) {
422  ortho_->setLabel( label_ );
423  }
424  }
425  }
426 
427  // Check for a change in verbosity level
428  if (params->isParameter("Verbosity")) {
429  if (Teuchos::isParameterType<int>(*params,"Verbosity")) {
430  verbosity_ = params->get("Verbosity", verbosity_default_);
431  } else {
432  verbosity_ = (int)Teuchos::getParameter<Belos::MsgType>(*params,"Verbosity");
433  }
434 
435  // Update parameter in our list.
436  params_->set("Verbosity", verbosity_);
437  if (printer_ != Teuchos::null)
438  printer_->setVerbosity(verbosity_);
439  }
440 
441  // Check for a change in output style
442  if (params->isParameter("Output Style")) {
443  if (Teuchos::isParameterType<int>(*params,"Output Style")) {
444  outputStyle_ = params->get("Output Style", outputStyle_default_);
445  } else {
446  outputStyle_ = (int)Teuchos::getParameter<Belos::OutputType>(*params,"Output Style");
447  }
448 
449  // Reconstruct the convergence test if the explicit residual test is not being used.
450  params_->set("Output Style", outputStyle_);
451  outputTest_ = Teuchos::null;
452  }
453 
454  // output stream
455  if (params->isParameter("Output Stream")) {
456  outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> >(*params,"Output Stream");
457 
458  // Update parameter in our list.
459  params_->set("Output Stream", outputStream_);
460  if (printer_ != Teuchos::null)
461  printer_->setOStream( outputStream_ );
462  }
463 
464  // frequency level
465  if (verbosity_ & Belos::StatusTestDetails) {
466  if (params->isParameter("Output Frequency")) {
467  outputFreq_ = params->get("Output Frequency", outputFreq_default_);
468  }
469 
470  // Update parameter in out list and output status test.
471  params_->set("Output Frequency", outputFreq_);
472  if (outputTest_ != Teuchos::null)
473  outputTest_->setOutputFrequency( outputFreq_ );
474  }
475 
476  // Create output manager if we need to.
477  if (printer_ == Teuchos::null) {
478  printer_ = Teuchos::rcp( new OutputManager<ScalarType>(verbosity_, outputStream_) );
479  }
480 
481  // Convergence
482  typedef Belos::StatusTestCombo<ScalarType,MV,OP> StatusTestCombo_t;
483  typedef Belos::StatusTestGenResNorm<ScalarType,MV,OP> StatusTestResNorm_t;
484 
485  // Check for convergence tolerance
486  if (params->isParameter("Convergence Tolerance")) {
487  convtol_ = params->get("Convergence Tolerance",convtol_default_);
488 
489  // Update parameter in our list and residual tests.
490  params_->set("Convergence Tolerance", convtol_);
491  if (convTest_ != Teuchos::null)
492  convTest_->setTolerance( convtol_ );
493  }
494 
495  if (params->isParameter("Show Maximum Residual Norm Only")) {
496  showMaxResNormOnly_ = Teuchos::getParameter<bool>(*params,"Show Maximum Residual Norm Only");
497 
498  // Update parameter in our list and residual tests
499  params_->set("Show Maximum Residual Norm Only", showMaxResNormOnly_);
500  if (convTest_ != Teuchos::null)
501  convTest_->setShowMaxResNormOnly( showMaxResNormOnly_ );
502  }
503 
504  // Check for a change in scaling, if so we need to build new residual tests.
505  bool newResTest = false;
506  {
507  // "Residual Scaling" is the old parameter name; "Implicit
508  // Residual Scaling" is the new name. We support both options for
509  // backwards compatibility.
510  std::string tempResScale = resScale_;
511  bool implicitResidualScalingName = false;
512  if (params->isParameter ("Residual Scaling")) {
513  tempResScale = params->get<std::string> ("Residual Scaling");
514  }
515  else if (params->isParameter ("Implicit Residual Scaling")) {
516  tempResScale = params->get<std::string> ("Implicit Residual Scaling");
517  implicitResidualScalingName = true;
518  }
519 
520  // Only update the scaling if it's different.
521  if (resScale_ != tempResScale) {
522  Belos::ScaleType resScaleType = convertStringToScaleType( tempResScale );
523  resScale_ = tempResScale;
524 
525  // Update parameter in our list and residual tests, using the
526  // given parameter name.
527  if (implicitResidualScalingName) {
528  params_->set ("Implicit Residual Scaling", resScale_);
529  }
530  else {
531  params_->set ("Residual Scaling", resScale_);
532  }
533 
534  if (! convTest_.is_null()) {
535  try {
536  convTest_->defineScaleForm( resScaleType, Belos::TwoNorm );
537  }
538  catch (std::exception& e) {
539  // Make sure the convergence test gets constructed again.
540  newResTest = true;
541  }
542  }
543  }
544  }
545 
546  // Get the deflation quorum, or number of converged systems before deflation is allowed
547  if (params->isParameter("Deflation Quorum")) {
548  defQuorum_ = params->get("Deflation Quorum", defQuorum_);
549  params_->set("Deflation Quorum", defQuorum_);
550  if (convTest_ != Teuchos::null)
551  convTest_->setQuorum( defQuorum_ );
552  }
553 
554  // Create status tests if we need to.
555 
556  // Basic test checks maximum iterations and native residual.
557  if (maxIterTest_ == Teuchos::null)
558  maxIterTest_ = Teuchos::rcp( new StatusTestMaxIters<ScalarType,MV,OP>( maxIters_ ) );
559 
560  // Implicit residual test, using the native residual to determine if convergence was achieved.
561  if (convTest_ == Teuchos::null || newResTest) {
562  convTest_ = Teuchos::rcp( new StatusTestResNorm_t( convtol_, defQuorum_, showMaxResNormOnly_ ) );
563  convTest_->defineScaleForm( convertStringToScaleType( resScale_ ), Belos::TwoNorm );
564  }
565 
566  if (sTest_ == Teuchos::null || newResTest)
567  sTest_ = Teuchos::rcp( new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
568 
569  if (outputTest_ == Teuchos::null || newResTest) {
570 
571  // Create the status test output class.
572  // This class manages and formats the output from the status test.
573  StatusTestOutputFactory<ScalarType,MV,OP> stoFactory( outputStyle_ );
574  outputTest_ = stoFactory.create( printer_, sTest_, outputFreq_, Passed+Failed+Undefined );
575 
576  // Set the solver string for the output test
577  std::string solverDesc = " Pseudo Block CG ";
578  outputTest_->setSolverDesc( solverDesc );
579 
580  }
581 
582  // Create the timer if we need to.
583  if (timerSolve_ == Teuchos::null) {
584  std::string solveLabel = label_ + ": PseudoBlockStochasticCGSolMgr total solve time";
585 #ifdef BELOS_TEUCHOS_TIME_MONITOR
586  timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
587 #endif
588  }
589 
590  // Inform the solver manager that the current parameters were set.
591  isSet_ = true;
592 }
593 
594 
595 template<class ScalarType, class MV, class OP>
598 {
600  using Teuchos::parameterList;
601  using Teuchos::RCP;
602 
603  if (validParams_.is_null()) {
604  // Set all the valid parameters and their default values.
605  RCP<ParameterList> pl = parameterList ();
606  pl->set("Convergence Tolerance", convtol_default_,
607  "The relative residual tolerance that needs to be achieved by the\n"
608  "iterative solver in order for the linera system to be declared converged.");
609  pl->set("Maximum Iterations", maxIters_default_,
610  "The maximum number of block iterations allowed for each\n"
611  "set of RHS solved.");
612  pl->set("Assert Positive Definiteness", assertPositiveDefiniteness_default_,
613  "Whether or not to assert that the linear operator\n"
614  "and the preconditioner are indeed positive definite.");
615  pl->set("Verbosity", verbosity_default_,
616  "What type(s) of solver information should be outputted\n"
617  "to the output stream.");
618  pl->set("Output Style", outputStyle_default_,
619  "What style is used for the solver information outputted\n"
620  "to the output stream.");
621  pl->set("Output Frequency", outputFreq_default_,
622  "How often convergence information should be outputted\n"
623  "to the output stream.");
624  pl->set("Deflation Quorum", defQuorum_default_,
625  "The number of linear systems that need to converge before\n"
626  "they are deflated. This number should be <= block size.");
627  pl->set("Output Stream", outputStream_default_,
628  "A reference-counted pointer to the output stream where all\n"
629  "solver output is sent.");
630  pl->set("Show Maximum Residual Norm Only", showMaxResNormOnly_default_,
631  "When convergence information is printed, only show the maximum\n"
632  "relative residual norm when the block size is greater than one.");
633  pl->set("Implicit Residual Scaling", resScale_default_,
634  "The type of scaling used in the residual convergence test.");
635  // We leave the old name as a valid parameter for backwards
636  // compatibility (so that validateParametersAndSetDefaults()
637  // doesn't raise an exception if it encounters "Residual
638  // Scaling"). The new name was added for compatibility with other
639  // solvers, none of which use "Residual Scaling".
640  pl->set("Residual Scaling", resScale_default_,
641  "The type of scaling used in the residual convergence test. This "
642  "name is deprecated; the new name is \"Implicit Residual Scaling\".");
643  pl->set("Timer Label", label_default_,
644  "The string to use as a prefix for the timer labels.");
645  // defaultParams_->set("Restart Timers", restartTimers_);
646  validParams_ = pl;
647  }
648  return validParams_;
649 }
650 
651 
652 // solve()
653 template<class ScalarType, class MV, class OP>
655 
656  // Set the current parameters if they were not set before.
657  // NOTE: This may occur if the user generated the solver manager with the default constructor and
658  // then didn't set any parameters using setParameters().
659  if (!isSet_) { setParameters( params_ ); }
660 
662 
664  "Belos::PseudoBlockStochasticCGSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
665 
666  // Create indices for the linear systems to be solved.
667  int startPtr = 0;
668  int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
669  int numCurrRHS = numRHS2Solve;
670 
671  std::vector<int> currIdx( numRHS2Solve ), currIdx2( numRHS2Solve );
672  for (int i=0; i<numRHS2Solve; ++i) {
673  currIdx[i] = startPtr+i;
674  currIdx2[i]=i;
675  }
676 
677  // Inform the linear problem of the current linear system to solve.
678  problem_->setLSIndex( currIdx );
679 
681  // Parameter list
683 
684  plist.set("Assert Positive Definiteness",assertPositiveDefiniteness_);
685 
686  // Reset the status test.
687  outputTest_->reset();
688 
689  // Assume convergence is achieved, then let any failed convergence set this to false.
690  bool isConverged = true;
691 
693  // Pseudo-Block CG solver
694 
696  = Teuchos::rcp( new PseudoBlockStochasticCGIter<ScalarType,MV,OP>(problem_,printer_,outputTest_,plist) );
697 
698  // Enter solve() iterations
699  {
700 #ifdef BELOS_TEUCHOS_TIME_MONITOR
701  Teuchos::TimeMonitor slvtimer(*timerSolve_);
702 #endif
703 
704  while ( numRHS2Solve > 0 ) {
705 
706  // Reset the active / converged vectors from this block
707  std::vector<int> convRHSIdx;
708  std::vector<int> currRHSIdx( currIdx );
709  currRHSIdx.resize(numCurrRHS);
710 
711  // Reset the number of iterations.
712  block_cg_iter->resetNumIters();
713 
714  // Reset the number of calls that the status test output knows about.
715  outputTest_->resetNumCalls();
716 
717  // Get the current residual for this block of linear systems.
718  Teuchos::RCP<MV> R_0 = MVT::CloneViewNonConst( *(Teuchos::rcp_const_cast<MV>(problem_->getInitResVec())), currIdx );
719 
720  // Get a new state struct and initialize the solver.
722  newState.R = R_0;
723  block_cg_iter->initializeCG(newState);
724 
725  while(1) {
726 
727  // tell block_gmres_iter to iterate
728  try {
729  block_cg_iter->iterate();
730 
732  //
733  // check convergence first
734  //
736  if ( convTest_->getStatus() == Passed ) {
737 
738  // Figure out which linear systems converged.
739  std::vector<int> convIdx = Teuchos::rcp_dynamic_cast<StatusTestGenResNorm<ScalarType,MV,OP> >(convTest_)->convIndices();
740 
741  // If the number of converged linear systems is equal to the
742  // number of current linear systems, then we are done with this block.
743  if (convIdx.size() == currRHSIdx.size())
744  break; // break from while(1){block_cg_iter->iterate()}
745 
746  // Inform the linear problem that we are finished with this current linear system.
747  problem_->setCurrLS();
748 
749  // Reset currRHSIdx to have the right-hand sides that are left to converge for this block.
750  int have = 0;
751  std::vector<int> unconvIdx(currRHSIdx.size());
752  for (unsigned int i=0; i<currRHSIdx.size(); ++i) {
753  bool found = false;
754  for (unsigned int j=0; j<convIdx.size(); ++j) {
755  if (currRHSIdx[i] == convIdx[j]) {
756  found = true;
757  break;
758  }
759  }
760  if (!found) {
761  currIdx2[have] = currIdx2[i];
762  currRHSIdx[have++] = currRHSIdx[i];
763  }
764  }
765  currRHSIdx.resize(have);
766  currIdx2.resize(have);
767 
768  // Set the remaining indices after deflation.
769  problem_->setLSIndex( currRHSIdx );
770 
771  // Get the current residual vector.
772  std::vector<MagnitudeType> norms;
773  R_0 = MVT::CloneCopy( *(block_cg_iter->getNativeResiduals(&norms)),currIdx2 );
774  for (int i=0; i<have; ++i) { currIdx2[i] = i; }
775 
776  // Set the new state and initialize the solver.
778  defstate.R = R_0;
779  block_cg_iter->initializeCG(defstate);
780  }
781 
783  //
784  // check for maximum iterations
785  //
787  else if ( maxIterTest_->getStatus() == Passed ) {
788  // we don't have convergence
789  isConverged = false;
790  break; // break from while(1){block_cg_iter->iterate()}
791  }
792 
794  //
795  // we returned from iterate(), but none of our status tests Passed.
796  // something is wrong, and it is probably our fault.
797  //
799 
800  else {
801  TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,
802  "Belos::PseudoBlockStochasticCGSolMgr::solve(): Invalid return from PseudoBlockStochasticCGIter::iterate().");
803  }
804  }
805  catch (const std::exception &e) {
806  printer_->stream(Errors) << "Error! Caught std::exception in PseudoBlockStochasticCGIter::iterate() at iteration "
807  << block_cg_iter->getNumIters() << std::endl
808  << e.what() << std::endl;
809  throw;
810  }
811  }
812 
813  // Inform the linear problem that we are finished with this block linear system.
814  problem_->setCurrLS();
815 
816  // Update indices for the linear systems to be solved.
817  startPtr += numCurrRHS;
818  numRHS2Solve -= numCurrRHS;
819 
820  if ( numRHS2Solve > 0 ) {
821 
822  numCurrRHS = numRHS2Solve;
823  currIdx.resize( numCurrRHS );
824  currIdx2.resize( numCurrRHS );
825  for (int i=0; i<numCurrRHS; ++i)
826  { currIdx[i] = startPtr+i; currIdx2[i] = i; }
827 
828  // Set the next indices.
829  problem_->setLSIndex( currIdx );
830  }
831  else {
832  currIdx.resize( numRHS2Solve );
833  }
834 
835  }// while ( numRHS2Solve > 0 )
836 
837  }
838 
839  // get the final stochastic vector
840  Y_=block_cg_iter->getStochasticVector();
841 
842 
843  // print final summary
844  sTest_->print( printer_->stream(FinalSummary) );
845 
846  // print timing information
847 #ifdef BELOS_TEUCHOS_TIME_MONITOR
848  // Calling summarize() can be expensive, so don't call unless the
849  // user wants to print out timing details. summarize() will do all
850  // the work even if it's passed a "black hole" output stream.
851  if (verbosity_ & TimingDetails)
852  Teuchos::TimeMonitor::summarize( printer_->stream(TimingDetails) );
853 #endif
854 
855  // get iteration information for this solve
856  numIters_ = maxIterTest_->getNumIters();
857 
858  if (!isConverged ) {
859  return Unconverged; // return from PseudoBlockStochasticCGSolMgr::solve()
860  }
861  return Converged; // return from PseudoBlockStochasticCGSolMgr::solve()
862 }
863 
864 // This method requires the solver manager to return a std::string that describes itself.
865 template<class ScalarType, class MV, class OP>
867 {
868  std::ostringstream oss;
869  oss << "Belos::PseudoBlockStochasticCGSolMgr<...,"<<Teuchos::ScalarTraits<ScalarType>::name()<<">";
870  oss << "{";
871  oss << "}";
872  return oss.str();
873 }
874 
875 } // end Belos namespace
876 
877 #endif /* BELOS_PSEUDO_BLOCK_CG_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.
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.
std::string description() const
Method to return description of the block CG solver manager.
ScaleType
The type of scaling to use on the residual norm value.
Definition: BelosTypes.hpp:119
T & get(const std::string &name, T def_value)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
static RCP< Time > getNewCounter(const std::string &name)
bool is_null(const std::shared_ptr< T > &p)
PseudoBlockStochasticCGSolMgrOrthoFailure is thrown when the orthogonalization manager is unable to g...
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
A factory class for generating StatusTestOutput objects.
An implementation of StatusTestResNorm using a family of residual norms.
bool isLOADetected() const
Return whether a loss of accuracy was detected by this solver during the most current solve...
Belos::StatusTest class for specifying a maximum number of iterations.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Get a parameter list containing the valid parameters for this object.
static std::string name()
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const
Get a parameter list containing the current parameters for this object.
A factory class for generating StatusTestOutput objects.
bool is_null() const
const LinearProblem< ScalarType, MV, OP > & getProblem() const
Return a reference to the linear problem being solved by this solver manager.
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
void reset(const ResetType type)
Performs a reset of the solver manager specified by the ResetType. This informs the solver manager th...
A Belos::StatusTest class for specifying a maximum number of iterations.
Belos concrete class for performing the stochastic pseudo-block CG iteration.
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)
Teuchos::Array< Teuchos::RCP< Teuchos::Time > > getTimers() const
Return the timers for this object.
void setParameters(const Teuchos::RCP< Teuchos::ParameterList > &params)
Set the parameters the solver manager should use to solve the linear problem.
int getNumIters() const
Get the iteration count for the most recent call to solve().
A linear system to solve, and its associated information.
Class which describes the linear problem to be solved by the iterative solver.
Teuchos::RCP< const MV > R
The current residual.
ReturnType
Whether the Belos solve converged for all linear systems.
Definition: BelosTypes.hpp:154
The Belos::SolverManager is a templated virtual base class that defines the basic interface that any ...
PseudoBlockStochasticCGSolMgr()
Empty constructor for BlockStochasticCGSolMgr. This constructor takes no arguments and sets the defau...
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.
Belos::StatusTestResNorm for specifying general residual norm stopping criteria.
This class implements the stochastic pseudo-block CG iteration, where the basic stochastic CG algorit...
The Belos::PseudoBlockStochasticCGSolMgr provides a powerful and fully-featured solver manager over t...
Teuchos::RCP< MV > getStochasticVector()
Get a copy of the final stochastic vector.
void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem)
Set the linear problem that needs to be solved.
A class for extending the status testing capabilities of Belos via logical combinations.
bool isParameter(const std::string &name) const
ReturnType solve()
This method performs possibly repeated calls to the underlying linear solver&#39;s iterate() routine unti...
Class which defines basic traits for the operator type.
Parent class to all Belos exceptions.
Definition: BelosTypes.hpp:60
Belos header file which uses auto-configuration information to include necessary C++ headers...
PseudoBlockStochasticCGSolMgrLinearProblemFailure is thrown when the linear problem is not setup (i...
Structure to contain pointers to CGIteration state variables.

Generated for Belos by doxygen 1.8.14