42#ifndef BELOS_STATUS_TEST_IMPRESNORM_H
43#define BELOS_STATUS_TEST_IMPRESNORM_H
53#include "Teuchos_as.hpp"
103template <
class ScalarType,
class MV,
class OP>
107 typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType
MagnitudeType;
112 typedef Teuchos::ScalarTraits<ScalarType>
STS;
113 typedef Teuchos::ScalarTraits<MagnitudeType>
STM;
134 bool showMaxResNormOnly =
false);
226 void print(std::ostream& os,
int indent = 0)
const;
301 std::ostringstream oss;
302 oss <<
"Belos::StatusTestImpResNorm<>: " <<
resFormStr();
317 std::ostringstream oss;
330 oss <<
" (User Scale)";
419template <
class ScalarType,
class MV,
class OP>
422 : tolerance_(Tolerance),
423 currTolerance_(Tolerance),
425 showMaxResNormOnly_(showMaxResNormOnly),
435 firstcallCheckStatus_(true),
436 firstcallDefineResForm_(true),
437 firstcallDefineScaleForm_(true),
444template <
class ScalarType,
class MV,
class OP>
448template <
class ScalarType,
class MV,
class OP>
457 currTolerance_ = tolerance_;
458 firstcallCheckStatus_ =
true;
459 lossDetected_ =
false;
460 curSoln_ = Teuchos::null;
463template <
class ScalarType,
class MV,
class OP>
466 TEUCHOS_TEST_FOR_EXCEPTION(firstcallDefineResForm_==
false,
StatusTestError,
467 "StatusTestResNorm::defineResForm(): The residual form has already been defined.");
468 firstcallDefineResForm_ =
false;
470 resnormtype_ = TypeOfNorm;
475template <
class ScalarType,
class MV,
class OP>
479 TEUCHOS_TEST_FOR_EXCEPTION(firstcallDefineScaleForm_==
false,
StatusTestError,
480 "StatusTestResNorm::defineScaleForm(): The scaling type has already been defined.");
481 firstcallDefineScaleForm_ =
false;
483 scaletype_ = TypeOfScaling;
484 scalenormtype_ = TypeOfNorm;
485 scalevalue_ = ScaleValue;
490template <
class ScalarType,
class MV,
class OP>
501 if (firstcallCheckStatus_) {
502 StatusType status = firstCallCheckStatusSetup (iSolver);
518 curBlksz_ = (int)curLSIdx_.size();
520 for (
int i=0; i<curBlksz_; ++i) {
521 if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
524 curNumRHS_ = validLS;
525 curSoln_ = Teuchos::null;
554 std::vector<MagnitudeType> tmp_resvector( curBlksz_ );
556 if (! residMV.is_null ()) {
558 tmp_resvector.resize (MVT::GetNumberVecs (*residMV));
559 MVT::MvNorm (*residMV, tmp_resvector, resnormtype_);
560 typename std::vector<int>::iterator p = curLSIdx_.begin();
561 for (
int i=0; p<curLSIdx_.end(); ++p, ++i) {
564 resvector_[*p] = tmp_resvector[i];
568 typename std::vector<int>::iterator p = curLSIdx_.begin();
569 for (
int i=0; p<curLSIdx_.end(); ++p, ++i) {
572 resvector_[*p] = tmp_resvector[i];
579 if (scalevector_.size () > 0) {
581 typename std::vector<int>::iterator p = curLSIdx_.begin();
582 for (; p<curLSIdx_.end(); ++p) {
586 if ( scalevector_[ *p ] != zero ) {
588 testvector_[ *p ] = resvector_[ *p ] / scalevector_[ *p ] / scalevalue_;
590 testvector_[ *p ] = resvector_[ *p ] / scalevalue_;
596 typename std::vector<int>::iterator p = curLSIdx_.begin();
597 for (; p<curLSIdx_.end(); ++p) {
600 testvector_[ *p ] = resvector_[ *p ] / scalevalue_;
613 ind_.resize( curLSIdx_.size() );
614 std::vector<int> lclInd( curLSIdx_.size() );
615 typename std::vector<int>::iterator p = curLSIdx_.begin();
616 for (
int i=0; p<curLSIdx_.end(); ++p, ++i) {
619 if (testvector_[ *p ] > currTolerance_) {
621 }
else if (testvector_[ *p ] <= currTolerance_) {
634 "StatusTestImpResNorm::checkStatus(): One or more of the current "
635 "implicit residual norms is NaN.");
650 RCP<MV> cur_res = MVT::Clone (*curSoln_, MVT::GetNumberVecs (*curSoln_));
652 tmp_resvector.resize (MVT::GetNumberVecs (*cur_res));
653 std::vector<MagnitudeType> tmp_testvector (have);
654 MVT::MvNorm (*cur_res, tmp_resvector, resnormtype_);
657 if ( scalevector_.size() > 0 ) {
658 for (
int i=0; i<have; ++i) {
660 if ( scalevector_[ ind_[i] ] != zero ) {
662 tmp_testvector[ i ] = tmp_resvector[ lclInd[i] ] / scalevector_[ ind_[i] ] / scalevalue_;
664 tmp_testvector[ i ] = tmp_resvector[ lclInd[i] ] / scalevalue_;
669 for (
int i=0; i<have; ++i) {
670 tmp_testvector[ i ] = tmp_resvector[ lclInd[i] ] / scalevalue_;
681 for (
int i = 0; i < have; ++i) {
691 if (tmp_testvector[i] <= tolerance_) {
692 ind_[have2] = ind_[i];
698 const MagnitudeType diff = STM::magnitude (testvector_[ind_[i]] - tmp_testvector[i]);
699 if (diff > currTolerance_) {
707 lossDetected_ =
true;
708 ind_[have2] = ind_[i];
728 const MagnitudeType onePointFive = as<MagnitudeType>(3) / as<MagnitudeType> (2);
729 const MagnitudeType oneTenth = STM::one () / as<MagnitudeType> (10);
731 currTolerance_ = currTolerance_ - onePointFive * diff;
732 while (currTolerance_ < STM::zero ()) {
733 currTolerance_ += oneTenth * diff;
744 int need = (quorum_ == -1) ? curNumRHS_: quorum_;
751template <
class ScalarType,
class MV,
class OP>
754 for (
int j = 0; j < indent; j ++)
756 printStatus(os, status_);
759 os <<
", tol = " << tolerance_ << std::endl;
762 if(showMaxResNormOnly_ && curBlksz_ > 1) {
764 testvector_.begin()+curLSIdx_[0],testvector_.begin()+curLSIdx_[curBlksz_-1]
766 for (
int j = 0; j < indent + 13; j ++)
768 os <<
"max{residual["<<curLSIdx_[0]<<
"..."<<curLSIdx_[curBlksz_-1]<<
"]} = " << maxRelRes
769 << ( maxRelRes <= tolerance_ ?
" <= " :
" > " ) << tolerance_ << std::endl;
772 for (
int i=0; i<numrhs_; i++ ) {
773 for (
int j = 0; j < indent + 13; j ++)
775 os <<
"residual [ " << i <<
" ] = " << testvector_[ i ];
776 os << ((testvector_[i]<tolerance_) ?
" < " : (testvector_[i]==tolerance_) ?
" == " : (testvector_[i]>tolerance_) ?
" > " :
" " ) << tolerance_ << std::endl;
783template <
class ScalarType,
class MV,
class OP>
786 os << std::left << std::setw(13) << std::setfill(
'.');
793 os <<
"Unconverged (LoA)";
802 os << std::left << std::setfill(
' ');
806template <
class ScalarType,
class MV,
class OP>
815 if (firstcallCheckStatus_) {
819 firstcallCheckStatus_ =
false;
822 Teuchos::RCP<const MV> rhs = lp.
getRHS();
823 numrhs_ = MVT::GetNumberVecs( *rhs );
824 scalevector_.resize( numrhs_ );
825 MVT::MvNorm( *rhs, scalevector_, scalenormtype_ );
829 numrhs_ = MVT::GetNumberVecs( *init_res );
830 scalevector_.resize( numrhs_ );
831 MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
835 numrhs_ = MVT::GetNumberVecs( *init_res );
836 scalevector_.resize( numrhs_ );
837 MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
840 numrhs_ = MVT::GetNumberVecs( *(lp.
getRHS()) );
843 resvector_.resize( numrhs_ );
844 testvector_.resize( numrhs_ );
848 curBlksz_ = (int)curLSIdx_.size();
850 for (i=0; i<curBlksz_; ++i) {
851 if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
854 curNumRHS_ = validLS;
857 for (i=0; i<numrhs_; i++) { testvector_[i] = one; }
860 if (scalevalue_ == zero) {
Class which describes the linear problem to be solved by the iterative solver.
Declaration of basic traits for the multivector type.
Belos::StatusTest abstract class for specifying a residual norm stopping criteria.
virtual Teuchos::RCP< const MV > getNativeResiduals(std::vector< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > *norms) const =0
Get the residuals native to the solver.
virtual Teuchos::RCP< MV > getCurrentUpdate() const =0
Get the current update to the linear system.
virtual const LinearProblem< ScalarType, MV, OP > & getProblem() const =0
Get a constant reference to the linear problem.
A linear system to solve, and its associated information.
int getLSNumber() const
The number of linear systems that have been set.
Teuchos::RCP< const MV > getInitPrecResVec() const
A pointer to the preconditioned initial residual vector.
const std::vector< int > & getLSIndex() const
(Zero-based) indices of the linear system(s) currently being solved.
Teuchos::RCP< const MV > getRHS() const
A pointer to the right-hand side B.
Teuchos::RCP< const MV > getInitResVec() const
A pointer to the initial unpreconditioned residual vector.
virtual 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.
virtual 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.
Traits class which defines basic operations on multivectors.
Exception thrown to signal error in a status test during Belos::StatusTest::checkStatus().
Convergence test using the implicit residual norm(s), with an explicit residual norm(s) check for los...
StatusTestImpResNorm(MagnitudeType Tolerance, int quorum=-1, bool showMaxResNormOnly=false)
Constructor.
Teuchos::ScalarTraits< ScalarType > STS
std::vector< MagnitudeType > scalevector_
Scaling vector.
std::vector< MagnitudeType > resvector_
Residual norm vector.
StatusType status_
Status.
const std::vector< MagnitudeType > * getScaledNormValue() const
Returns the scaled norm value, .
bool firstcallCheckStatus_
Is this the first time CheckStatus is called?
bool lossDetected_
Has a loss in accuracy been detected?
virtual ~StatusTestImpResNorm()
Destructor (virtual for memory safety).
Teuchos::ScalarTraits< ScalarType >::magnitudeType MagnitudeType
The type of the magnitude (absolute value) of a ScalarType.
ScaleType scaletype_
Type of scaling to use (Norm of RHS, Norm of Initial Residual, None or User provided)
MagnitudeType getTolerance() const
"Original" convergence tolerance as set by user.
int numrhs_
The total number of right-hand sides being solved for.
MagnitudeType currTolerance_
MagnitudeType getCurrTolerance() const
Current convergence tolerance; may be changed to prevent loss of accuracy.
bool showMaxResNormOnly_
Determines if the entries for all of the residuals are shown or just the max.
Teuchos::ScalarTraits< MagnitudeType > STM
std::string description() const
Method to return description of the maximum iteration status test
std::vector< int > convIndices()
Returns the vector containing the indices of the residuals that passed the test.
void printStatus(std::ostream &os, StatusType type) const
Print message for each status specific to this stopping test.
bool firstcallDefineScaleForm_
Is this the first time DefineScaleForm is called?
int defineScaleForm(ScaleType TypeOfScaling, NormType TypeOfNorm, MagnitudeType ScaleValue=Teuchos::ScalarTraits< MagnitudeType >::one())
Define form of the scaling, its norm, its optional weighting vector, or, alternatively,...
std::vector< MagnitudeType > testvector_
Test vector = resvector_ / scalevector_.
std::string resFormStr() const
Description of current residual form.
const std::vector< MagnitudeType > * getResNormValue() const
Returns the residual norm value, , computed in most recent call to CheckStatus.
Teuchos::RCP< MV > getSolution()
Returns the current solution estimate that was computed for the most recent residual test.
void print(std::ostream &os, int indent=0) const
Output formatted description of stopping test to output stream.
void reset()
Resets the internal configuration to the initial state.
MagnitudeType scalevalue_
Scaling value.
bool getShowMaxResNormOnly()
Returns whether the only maximum residual norm is displayed when the print() method is called.
NormType scalenormtype_
Type of norm to use on the scaling (OneNorm, TwoNorm, or InfNorm)
bool getLOADetected() const
Returns a boolean indicating a loss of accuracy has been detected in computing the residual.
MagnitudeType tolerance_
Current tolerance used to determine convergence and the default tolerance set by user.
int getQuorum() const
Returns the number of residuals that must pass the convergence test before Passed is returned.
StatusType getStatus() const
Return the result of the most recent CheckStatus call.
std::vector< int > ind_
Vector containing the indices for the vectors that passed the test.
MultiVecTraits< ScalarType, MV > MVT
StatusType checkStatus(Iteration< ScalarType, MV, OP > *iSolver)
Check convergence status: Passed, Failed, or Undefined.
int setTolerance(MagnitudeType tolerance)
Set the value of the tolerance.
int quorum_
Number of residuals that must pass the convergence test before Passed is returned.
std::vector< int > curLSIdx_
The indices of the current number of right-hand sides being solved for.
bool firstcallDefineResForm_
Is this the first time DefineResForm is called?
int curNumRHS_
The current number of right-hand sides being solved for.
int curBlksz_
The current blocksize of the linear system being solved.
const std::vector< MagnitudeType > * getTestValue() const
Returns the test value, , computed in most recent call to CheckStatus.
NormType resnormtype_
Type of norm to use on residual (OneNorm, TwoNorm, or InfNorm).
Teuchos::RCP< MV > curSoln_
Current solution vector.
int defineResForm(NormType TypeOfNorm)
Define form of the residual, its norm and optional weighting vector.
int setShowMaxResNormOnly(bool showMaxResNormOnly)
Set whether the only maximum residual norm is displayed when the print() method is called.
StatusType firstCallCheckStatusSetup(Iteration< ScalarType, MV, OP > *iSolver)
Call to setup initial scaling vector.
int curLSNum_
The current number of linear systems that have been loaded into the linear problem.
int setQuorum(int quorum)
Sets the number of residuals that must pass the convergence test before Passed is returned.
An abstract class of StatusTest for stopping criteria using residual norms.
NormType
The type of vector norm to compute.
StatusType
Whether the StatusTest wants iteration to stop.
ScaleType
The type of scaling to use on the residual norm value.