Belos Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
BelosStatusTestGenResNorm.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_STATUS_TEST_GEN_RESNORM_H
43 #define BELOS_STATUS_TEST_GEN_RESNORM_H
44 
51 #include "BelosLinearProblem.hpp"
52 #include "BelosMultiVecTraits.hpp"
53 
76 namespace Belos {
77 
78 template <class ScalarType, class MV, class OP>
79 class StatusTestGenResNorm: public StatusTestResNorm<ScalarType,MV,OP> {
80 
81  public:
82 
83  // Convenience typedefs
87 
89 
90 
94  enum ResType {Implicit,
97  };
98 
100 
102 
103 
117  StatusTestGenResNorm( MagnitudeType Tolerance, int quorum = -1, bool showMaxResNormOnly = false );
118 
120  virtual ~StatusTestGenResNorm();
122 
124 
125 
127 
134  int defineResForm( ResType TypeOfResidual, NormType TypeOfNorm);
135 
137 
157  int defineScaleForm( ScaleType TypeOfScaling, NormType TypeOfNorm, MagnitudeType ScaleValue = Teuchos::ScalarTraits<MagnitudeType>::one());
158 
160 
162 
165  int setTolerance(MagnitudeType tolerance) {tolerance_ = tolerance; return(0);}
166 
169  int setQuorum(int quorum) {quorum_ = quorum; return(0);}
170 
172  int setShowMaxResNormOnly(bool showMaxResNormOnly) {showMaxResNormOnly_ = showMaxResNormOnly; return(0);}
173 
175 
177 
178 
186 
188  StatusType getStatus() const {return(status_);};
190 
192 
193 
195  void reset();
196 
198 
200 
201 
203  void print(std::ostream& os, int indent = 0) const;
204 
206  void printStatus(std::ostream& os, StatusType type) const;
208 
210 
211 
215  Teuchos::RCP<MV> getSolution() { if (restype_==Implicit) { return Teuchos::null; } else { return curSoln_; } }
216 
219  int getQuorum() const { return quorum_; }
220 
223 
225  std::vector<int> convIndices() { return ind_; }
226 
229 
231  const std::vector<MagnitudeType>* getTestValue() const {return(&testvector_);};
232 
234  const std::vector<MagnitudeType>* getResNormValue() const {return(&resvector_);};
235 
237  const std::vector<MagnitudeType>* getScaledNormValue() const {return(&scalevector_);};
238 
241  bool getLOADetected() const { return false; }
242 
244 
245 
248 
256 
259 
261  std::string description() const
262  {
263  std::ostringstream oss;
264  oss << "Belos::StatusTestGenResNorm<>: " << resFormStr();
265  oss << ", tol = " << tolerance_;
266  return oss.str();
267  }
269 
270  protected:
271 
272  private:
273 
275 
276 
277  std::string resFormStr() const
278  {
279  std::ostringstream oss;
280  oss << "(";
281  oss << ((resnormtype_==OneNorm) ? "1-Norm" : (resnormtype_==TwoNorm) ? "2-Norm" : "Inf-Norm");
282  oss << ((restype_==Explicit) ? " Exp" : " Imp");
283  oss << " Res Vec) ";
284 
285  // If there is no residual scaling, return current string.
286  if (scaletype_!=None)
287  {
288  // Insert division sign.
289  oss << "/ ";
290 
291  // Determine output string for scaling, if there is any.
293  oss << " (User Scale)";
294  else {
295  oss << "(";
296  oss << ((scalenormtype_==OneNorm) ? "1-Norm" : (resnormtype_==TwoNorm) ? "2-Norm" : "Inf-Norm");
298  oss << " Res0";
300  oss << " Prec Res0";
301  else
302  oss << " RHS ";
303  oss << ")";
304  }
305  }
306 
307  return oss.str();
308  }
309 
311 
313 
314 
317 
319  int quorum_;
320 
323 
326 
329 
332 
335 
338 
340  std::vector<MagnitudeType> scalevector_;
341 
343  std::vector<MagnitudeType> resvector_;
344 
346  std::vector<MagnitudeType> testvector_;
347 
349  std::vector<int> ind_;
350 
353 
356 
359 
362 
364  std::vector<int> curLSIdx_;
365 
368 
370  int numrhs_;
371 
374 
377 
380 
382 
383 };
384 
385 template <class ScalarType, class MV, class OP>
387 StatusTestGenResNorm (MagnitudeType Tolerance, int quorum, bool showMaxResNormOnly)
388  : tolerance_(Tolerance),
389  quorum_(quorum),
390  showMaxResNormOnly_(showMaxResNormOnly),
391  restype_(Implicit),
392  resnormtype_(TwoNorm),
393  scaletype_(NormOfInitRes),
394  scalenormtype_(TwoNorm),
395  scalevalue_(Teuchos::ScalarTraits<MagnitudeType>::one ()),
396  status_(Undefined),
397  curBlksz_(0),
398  curNumRHS_(0),
399  curLSNum_(0),
400  numrhs_(0),
401  firstcallCheckStatus_(true),
402  firstcallDefineResForm_(true),
403  firstcallDefineScaleForm_(true)
404 {
405  // This constructor will compute the residual ||r_i||/||r0_i|| <= tolerance using the 2-norm of
406  // the implicit residual std::vector.
407 }
408 
409 template <class ScalarType, class MV, class OP>
411 {}
412 
413 template <class ScalarType, class MV, class OP>
415 {
416  status_ = Undefined;
417  curBlksz_ = 0;
418  curLSNum_ = 0;
419  curLSIdx_.resize(0);
420  numrhs_ = 0;
421  ind_.resize(0);
422  firstcallCheckStatus_ = true;
423  curSoln_ = Teuchos::null;
424 }
425 
426 template <class ScalarType, class MV, class OP>
428 {
429  TEUCHOS_TEST_FOR_EXCEPTION(firstcallDefineResForm_==false,StatusTestError,
430  "StatusTestGenResNorm::defineResForm(): The residual form has already been defined.");
431  firstcallDefineResForm_ = false;
432 
433  restype_ = TypeOfResidual;
434  resnormtype_ = TypeOfNorm;
435 
436  return(0);
437 }
438 
439 template <class ScalarType, class MV, class OP>
441  MagnitudeType ScaleValue )
442 {
443  TEUCHOS_TEST_FOR_EXCEPTION(firstcallDefineScaleForm_==false,StatusTestError,
444  "StatusTestGenResNorm::defineScaleForm(): The scaling type has already been defined.");
445  firstcallDefineScaleForm_ = false;
446 
447  scaletype_ = TypeOfScaling;
448  scalenormtype_ = TypeOfNorm;
449  scalevalue_ = ScaleValue;
450 
451  return(0);
452 }
453 
454 template <class ScalarType, class MV, class OP>
456 {
458  const LinearProblem<ScalarType,MV,OP>& lp = iSolver->getProblem();
459  // Compute scaling term (done once for each block that's being solved)
460  if (firstcallCheckStatus_) {
461  StatusType status = firstCallCheckStatusSetup(iSolver);
462  if(status==Failed) {
463  status_ = Failed;
464  return(status_);
465  }
466  }
467  //
468  // This section computes the norm of the residual std::vector
469  //
470  if ( curLSNum_ != lp.getLSNumber() ) {
471  //
472  // We have moved on to the next rhs block
473  //
474  curLSNum_ = lp.getLSNumber();
475  curLSIdx_ = lp.getLSIndex();
476  curBlksz_ = (int)curLSIdx_.size();
477  int validLS = 0;
478  for (int i=0; i<curBlksz_; ++i) {
479  if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
480  validLS++;
481  }
482  curNumRHS_ = validLS;
483  curSoln_ = Teuchos::null;
484  //
485  } else {
486  //
487  // We are in the same rhs block, return if we are converged
488  //
489  if (status_==Passed) { return status_; }
490  }
491  if (restype_==Implicit) {
492  //
493  // get the native residual norms from the solver for this block of right-hand sides.
494  // If the residual is returned in multivector form, use the resnormtype to compute the residual norms.
495  // Otherwise the native residual is assumed to be stored in the resvector_.
496  //
497  std::vector<MagnitudeType> tmp_resvector( curBlksz_ );
498  Teuchos::RCP<const MV> residMV = iSolver->getNativeResiduals( &tmp_resvector );
499  if ( residMV != Teuchos::null ) {
500  tmp_resvector.resize( MVT::GetNumberVecs( *residMV ) );
501  MVT::MvNorm( *residMV, tmp_resvector, resnormtype_ );
502  typename std::vector<int>::iterator p = curLSIdx_.begin();
503  for (int i=0; p<curLSIdx_.end(); ++p, ++i) {
504  // Check if this index is valid
505  if (*p != -1)
506  resvector_[*p] = tmp_resvector[i];
507  }
508  } else {
509  typename std::vector<int>::iterator p = curLSIdx_.begin();
510  for (int i=0; p<curLSIdx_.end(); ++p, ++i) {
511  // Check if this index is valid
512  if (*p != -1)
513  resvector_[*p] = tmp_resvector[i];
514  }
515  }
516  }
517  else if (restype_==Explicit) {
518  //
519  // Request the true residual for this block of right-hand sides.
520  //
521  Teuchos::RCP<MV> cur_update = iSolver->getCurrentUpdate();
522  curSoln_ = lp.updateSolution( cur_update );
523  Teuchos::RCP<MV> cur_res = MVT::Clone( *curSoln_, MVT::GetNumberVecs( *curSoln_ ) );
524  lp.computeCurrResVec( &*cur_res, &*curSoln_ );
525  std::vector<MagnitudeType> tmp_resvector( MVT::GetNumberVecs( *cur_res ) );
526  MVT::MvNorm( *cur_res, tmp_resvector, resnormtype_ );
527  typename std::vector<int>::iterator p = curLSIdx_.begin();
528  for (int i=0; p<curLSIdx_.end(); ++p, ++i) {
529  // Check if this index is valid
530  if (*p != -1)
531  resvector_[*p] = tmp_resvector[i];
532  }
533  }
534  //
535  // Compute the new linear system residuals for testing.
536  // (if any of them don't meet the tolerance or are NaN, then we exit with that status)
537  //
538  if ( scalevector_.size() > 0 ) {
539  typename std::vector<int>::iterator p = curLSIdx_.begin();
540  for (; p<curLSIdx_.end(); ++p) {
541  // Check if this index is valid
542  if (*p != -1) {
543  // Scale the std::vector accordingly
544  if ( scalevector_[ *p ] != zero ) {
545  // Don't intentionally divide by zero.
546  testvector_[ *p ] = resvector_[ *p ] / scalevector_[ *p ] / scalevalue_;
547  } else {
548  testvector_[ *p ] = resvector_[ *p ] / scalevalue_;
549  }
550  }
551  }
552  }
553  else {
554  typename std::vector<int>::iterator p = curLSIdx_.begin();
555  for (; p<curLSIdx_.end(); ++p) {
556  // Check if this index is valid
557  if (*p != -1)
558  testvector_[ *p ] = resvector_[ *p ] / scalevalue_;
559  }
560  }
561 
562  // Check status of new linear system residuals and see if we have the quorum.
563  int have = 0;
564  ind_.resize( curLSIdx_.size() );
565  typename std::vector<int>::iterator p = curLSIdx_.begin();
566  for (; p<curLSIdx_.end(); ++p) {
567  // Check if this index is valid
568  if (*p != -1) {
569  // Check if any of the residuals are larger than the tolerance.
570  if (testvector_[ *p ] > tolerance_) {
571  // do nothing.
572  } else if (testvector_[ *p ] <= tolerance_) {
573  ind_[have] = *p;
574  have++;
575  } else {
576  // Throw an std::exception if a NaN is found.
577  status_ = Failed;
578  TEUCHOS_TEST_FOR_EXCEPTION(true,StatusTestError,"StatusTestGenResNorm::checkStatus(): NaN has been detected.");
579  }
580  }
581  }
582  ind_.resize(have);
583  int need = (quorum_ == -1) ? curNumRHS_: quorum_;
584  status_ = (have >= need) ? Passed : Failed;
585 
586  // Return the current status
587  return status_;
588 }
589 
590 template <class ScalarType, class MV, class OP>
591 void StatusTestGenResNorm<ScalarType,MV,OP>::print(std::ostream& os, int indent) const
592 {
593  for (int j = 0; j < indent; j ++)
594  os << ' ';
595  printStatus(os, status_);
596  os << resFormStr();
597  if (status_==Undefined)
598  os << ", tol = " << tolerance_ << std::endl;
599  else {
600  os << std::endl;
601  if(showMaxResNormOnly_ && curBlksz_ > 1) {
602  const MagnitudeType maxRelRes = *std::max_element(
603  testvector_.begin()+curLSIdx_[0],testvector_.begin()+curLSIdx_[curBlksz_-1]
604  );
605  for (int j = 0; j < indent + 13; j ++)
606  os << ' ';
607  os << "max{residual["<<curLSIdx_[0]<<"..."<<curLSIdx_[curBlksz_-1]<<"]} = " << maxRelRes
608  << ( maxRelRes <= tolerance_ ? " <= " : " > " ) << tolerance_ << std::endl;
609  }
610  else {
611  for ( int i=0; i<numrhs_; i++ ) {
612  for (int j = 0; j < indent + 13; j ++)
613  os << ' ';
614  os << "residual [ " << i << " ] = " << testvector_[ i ];
615  os << ((testvector_[i]<tolerance_) ? " < " : (testvector_[i]==tolerance_) ? " == " : (testvector_[i]>tolerance_) ? " > " : " " ) << tolerance_ << std::endl;
616  }
617  }
618  }
619  os << std::endl;
620 }
621 
622 template <class ScalarType, class MV, class OP>
624 {
625  os << std::left << std::setw(13) << std::setfill('.');
626  switch (type) {
627  case Passed:
628  os << "Converged";
629  break;
630  case Failed:
631  os << "Unconverged";
632  break;
633  case Undefined:
634  default:
635  os << "**";
636  break;
637  }
638  os << std::left << std::setfill(' ');
639  return;
640 }
641 
642 template <class ScalarType, class MV, class OP>
644 {
645  int i;
648  const LinearProblem<ScalarType,MV,OP>& lp = iSolver->getProblem();
649  // Compute scaling term (done once for each block that's being solved)
650  if (firstcallCheckStatus_) {
651  //
652  // Get some current solver information.
653  //
654  firstcallCheckStatus_ = false;
655 
656  if (scaletype_== NormOfRHS) {
657  Teuchos::RCP<const MV> rhs = lp.getRHS();
658  numrhs_ = MVT::GetNumberVecs( *rhs );
659  scalevector_.resize( numrhs_ );
660  MVT::MvNorm( *rhs, scalevector_, scalenormtype_ );
661  }
662  else if (scaletype_==NormOfInitRes || scaletype_==NormOfFullInitRes || scaletype_==NormOfFullScaledInitRes) {
663  Teuchos::RCP<const MV> init_res = lp.getInitResVec();
664  numrhs_ = MVT::GetNumberVecs( *init_res );
665  scalevector_.resize( numrhs_ );
666  MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
667  }
668  else if (scaletype_==NormOfPrecInitRes || scaletype_==NormOfFullPrecInitRes || scaletype_==NormOfFullScaledPrecInitRes) {
670  numrhs_ = MVT::GetNumberVecs( *init_res );
671  scalevector_.resize( numrhs_ );
672  MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
673  }
674  else {
675  numrhs_ = MVT::GetNumberVecs( *(lp.getRHS()) );
676  }
677 
678  resvector_.resize( numrhs_ );
679  testvector_.resize( numrhs_ );
680 
681  curLSNum_ = lp.getLSNumber();
682  curLSIdx_ = lp.getLSIndex();
683  curBlksz_ = (int)curLSIdx_.size();
684  int validLS = 0;
685  for (i=0; i<curBlksz_; ++i) {
686  if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
687  validLS++;
688  }
689  curNumRHS_ = validLS;
690  //
691  // Initialize the testvector.
692  for (i=0; i<numrhs_; i++) { testvector_[i] = one; }
693 
694  // Return an error if the scaling is zero.
695  if (scalevalue_ == zero) {
696  return Failed;
697  }
698  }
699  return Undefined;
700 }
701 
702 } // end namespace Belos
703 
704 #endif /* BELOS_STATUS_TEST_RESNORM_H */
virtual Teuchos::RCP< const MV > getNativeResiduals(std::vector< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > *norms) const =0
Get the residuals native to the solver.
bool firstcallDefineResForm_
Is this the first time DefineResForm is called?
bool showMaxResNormOnly_
Determines if the entries for all of the residuals are shown or just the max.
std::string resFormStr() const
Description of current residual form.
int defineScaleForm(ScaleType TypeOfScaling, NormType TypeOfNorm, MagnitudeType ScaleValue=Teuchos::ScalarTraits< MagnitudeType >::one())
Define form of the scaling, its norm, its optional weighting std::vector, or, alternatively, define an explicit value.
MagnitudeType getTolerance() const
Returns the value of the tolerance, , set in the constructor.
const std::vector< MagnitudeType > * getTestValue() const
Returns the test value, , computed in most recent call to CheckStatus.
int curLSNum_
The current number of linear systems that have been loaded into the linear problem.
ScaleType
The type of scaling to use on the residual norm value.
Definition: BelosTypes.hpp:120
Exception thrown to signal error in a status test during Belos::StatusTest::checkStatus().
int setShowMaxResNormOnly(bool showMaxResNormOnly)
Set whether the only maximum residual norm is displayed when the print() method is called...
std::vector< int > curLSIdx_
The indices of the current number of right-hand sides being solved for.
StatusType firstCallCheckStatusSetup(Iteration< ScalarType, MV, OP > *iSolver)
Call to setup initial scaling std::vector.
Teuchos::RCP< MV > updateSolution(const Teuchos::RCP< MV > &update=Teuchos::null, bool updateLP=false, ScalarType scale=Teuchos::ScalarTraits< ScalarType >::one())
Compute the new solution to the linear system using the given update vector.
ScaleType scaletype_
Type of scaling to use (Norm of RHS, Norm of Initial Residual, None or User provided) ...
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual Teuchos::RCP< MV > getCurrentUpdate() const =0
Get the current update to the linear system.
Declaration of basic traits for the multivector type.
std::string description() const
Method to return description of the maximum iteration status test.
std::vector< MagnitudeType > testvector_
Test std::vector = resvector_ / scalevector_.
std::vector< int > ind_
Vector containing the indices for the vectors that passed the test.
Teuchos::RCP< const MV > getInitPrecResVec() const
A pointer to the preconditioned initial residual vector.
An abstract class of StatusTest for stopping criteria using residual norms.
void printStatus(std::ostream &os, StatusType type) const
Print message for each status specific to this stopping test.
An implementation of StatusTestResNorm using a family of residual norms.
int curBlksz_
The current blocksize of the linear system being solved.
Teuchos::RCP< const MV > getRHS() const
A pointer to the right-hand side B.
std::vector< MagnitudeType > resvector_
Residual norm std::vector.
StatusType
Whether the StatusTest wants iteration to stop.
Definition: BelosTypes.hpp:189
const std::vector< MagnitudeType > * getResNormValue() const
Returns the residual norm value, , computed in most recent call to CheckStatus.
int numrhs_
The total number of right-hand sides being solved for.
Traits class which defines basic operations on multivectors.
MultiVecTraits< ScalarType, MV > MVT
NormType resnormtype_
Type of norm to use on residual (OneNorm, TwoNorm, or InfNorm).
bool getShowMaxResNormOnly()
Returns whether the only maximum residual norm is displayed when the print() method is called...
int defineResForm(ResType TypeOfResidual, NormType TypeOfNorm)
Define form of the residual, its norm and optional weighting std::vector.
bool firstcallCheckStatus_
Is this the first time CheckStatus is called?
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.
int getLSNumber() const
The number of linear systems that have been set.
NormType scalenormtype_
Type of norm to use on the scaling (OneNorm, TwoNorm, or InfNorm)
int setQuorum(int quorum)
Sets the number of residuals that must pass the convergence test before Passed is returned...
std::vector< int > convIndices()
Returns the std::vector containing the indices of the residuals that passed the test.
virtual const LinearProblem< ScalarType, MV, OP > & getProblem() const =0
Get a constant reference to the linear problem.
Teuchos::RCP< const MV > getInitResVec() const
A pointer to the initial unpreconditioned residual vector.
Belos::StatusTest abstract class for specifying a residual norm stopping criteria.
StatusType getStatus() const
Return the result of the most recent CheckStatus call.
int curNumRHS_
The current number of right-hand sides being solved for.
int setTolerance(MagnitudeType tolerance)
Set the value of the tolerance.
NormType
The type of vector norm to compute.
Definition: BelosTypes.hpp:97
const std::vector< int > getLSIndex() const
(Zero-based) indices of the linear system(s) currently being solved.
StatusType checkStatus(Iteration< ScalarType, MV, OP > *iSolver)
Check convergence status: Passed, Failed, or Undefined.
std::vector< MagnitudeType > scalevector_
Scaling std::vector.
void computeCurrResVec(MV *R, const MV *X=0, const MV *B=0) const
Compute a residual R for this operator given a solution X, and right-hand side B. ...
ResType restype_
Type of residual to use (explicit or implicit)
void print(std::ostream &os, int indent=0) const
Output formatted description of stopping test to output stream.
ResType
Select how the residual std::vector is produced.
int quorum_
Number of residuals that must pass the convergence test before Passed is returned.
void reset()
Resets the internal configuration to the initial state.
bool firstcallDefineScaleForm_
Is this the first time DefineScaleForm is called?
int getQuorum() const
Returns the number of residuals that must pass the convergence test before Passed is returned...
MagnitudeType tolerance_
Tolerance used to determine convergence.
StatusTestGenResNorm(MagnitudeType Tolerance, int quorum=-1, bool showMaxResNormOnly=false)
Constructor.
MagnitudeType scalevalue_
Scaling value.
Teuchos::RCP< MV > curSoln_
Most recent solution vector used by this status test.
const std::vector< MagnitudeType > * getScaledNormValue() const
Returns the scaled norm value, .
Teuchos::RCP< MV > getSolution()
Returns the current solution estimate that was computed for the most recent residual test...
bool getLOADetected() const
Returns a boolean indicating a loss of accuracy has been detected in computing the residual...