42 #ifndef BELOS_STATUS_TEST_IMPRESNORM_H
43 #define BELOS_STATUS_TEST_IMPRESNORM_H
103 template <
class ScalarType,
class MV,
class OP>
134 bool showMaxResNormOnly =
false);
180 tolerance_ = tolerance;
193 showMaxResNormOnly_ = showMaxResNormOnly;
226 void print(std::ostream& os,
int indent = 0)
const;
267 return currTolerance_;
271 const std::vector<MagnitudeType>*
getTestValue()
const {
return(&testvector_);};
301 std::ostringstream oss;
302 oss <<
"Belos::StatusTestImpResNorm<>: " << resFormStr();
303 oss <<
", tol = " << tolerance_;
315 std::string resFormStr()
const
317 std::ostringstream oss;
319 oss << ((resnormtype_==
OneNorm) ?
"1-Norm" : (resnormtype_==
TwoNorm) ?
"2-Norm" :
"Inf-Norm");
323 if (scaletype_!=
None)
330 oss <<
" (User Scale)";
333 oss << ((scalenormtype_==
OneNorm) ?
"1-Norm" : (resnormtype_==
TwoNorm) ?
"2-Norm" :
"Inf-Norm");
357 bool showMaxResNormOnly_;
372 std::vector<MagnitudeType> scalevector_;
375 std::vector<MagnitudeType> resvector_;
378 std::vector<MagnitudeType> testvector_;
384 std::vector<int> ind_;
396 std::vector<int> curLSIdx_;
405 bool firstcallCheckStatus_;
408 bool firstcallDefineResForm_;
411 bool firstcallDefineScaleForm_;
419 template <
class ScalarType,
class MV,
class OP>
422 : tolerance_(Tolerance),
423 currTolerance_(Tolerance),
425 showMaxResNormOnly_(showMaxResNormOnly),
435 firstcallCheckStatus_(true),
436 firstcallDefineResForm_(true),
437 firstcallDefineScaleForm_(true),
444 template <
class ScalarType,
class MV,
class OP>
448 template <
class ScalarType,
class MV,
class OP>
457 currTolerance_ = tolerance_;
458 firstcallCheckStatus_ =
true;
459 lossDetected_ =
false;
460 curSoln_ = Teuchos::null;
463 template <
class ScalarType,
class MV,
class OP>
467 "StatusTestResNorm::defineResForm(): The residual form has already been defined.");
468 firstcallDefineResForm_ =
false;
470 resnormtype_ = TypeOfNorm;
475 template <
class ScalarType,
class MV,
class OP>
480 "StatusTestResNorm::defineScaleForm(): The scaling type has already been defined.");
481 firstcallDefineScaleForm_ =
false;
483 scaletype_ = TypeOfScaling;
484 scalenormtype_ = TypeOfNorm;
485 scalevalue_ = ScaleValue;
490 template <
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_;
751 template <
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;
783 template <
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(
' ');
806 template <
class ScalarType,
class MV,
class OP>
815 if (firstcallCheckStatus_) {
819 firstcallCheckStatus_ =
false;
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) {
virtual Teuchos::RCP< const MV > getNativeResiduals(std::vector< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > *norms) const =0
Get the residuals native to the solver.
Teuchos::RCP< MV > getSolution()
Returns the current solution estimate that was computed for the most recent residual test...
virtual ~StatusTestImpResNorm()
Destructor (virtual for memory safety).
const std::vector< MagnitudeType > * getTestValue() const
Returns the test value, , computed in most recent call to CheckStatus.
ScaleType
The type of scaling to use on the residual norm value.
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...
StatusType checkStatus(Iteration< ScalarType, MV, OP > *iSolver)
Check convergence status: Passed, Failed, or Undefined.
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.
#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.
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.
Teuchos::RCP< const MV > getRHS() const
A pointer to the right-hand side B.
StatusType
Whether the StatusTest wants iteration to stop.
const std::vector< MagnitudeType > * getResNormValue() const
Returns the residual norm value, , computed in most recent call to CheckStatus.
Traits class which defines basic operations on multivectors.
int getQuorum() const
Returns the number of residuals that must pass the convergence test before Passed is returned...
Teuchos::ScalarTraits< ScalarType >::magnitudeType MagnitudeType
The type of the magnitude (absolute value) of a ScalarType.
int setShowMaxResNormOnly(bool showMaxResNormOnly)
Set whether the only maximum residual norm is displayed when the print() method is called...
const std::vector< MagnitudeType > * getScaledNormValue() const
Returns the scaled norm value, .
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.
int setQuorum(int quorum)
Sets the number of residuals that must pass the convergence test before Passed is returned...
StatusType firstCallCheckStatusSetup(Iteration< ScalarType, MV, OP > *iSolver)
Call to setup initial scaling vector.
MagnitudeType getCurrTolerance() const
Current convergence tolerance; may be changed to prevent loss of accuracy.
int defineResForm(NormType TypeOfNorm)
Define form of the residual, its norm and optional weighting vector.
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.
bool getLOADetected() const
Returns a boolean indicating a loss of accuracy has been detected in computing the residual...
NormType
The type of vector norm to compute.
const std::vector< int > getLSIndex() const
(Zero-based) indices of the linear system(s) currently being solved.
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. ...
StatusType getStatus() const
Return the result of the most recent CheckStatus call.
MagnitudeType getTolerance() const
"Original" convergence tolerance as set by user.
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, define an explicit value.
TypeTo as(const TypeFrom &t)
StatusTestImpResNorm(MagnitudeType Tolerance, int quorum=-1, bool showMaxResNormOnly=false)
Constructor.
int setTolerance(MagnitudeType tolerance)
Set the value of the tolerance.
std::string description() const
Method to return description of the maximum iteration status test.
void printStatus(std::ostream &os, StatusType type) const
Print message for each status specific to this stopping test.
void reset()
Resets the internal configuration to the initial state.
void print(std::ostream &os, int indent=0) const
Output formatted description of stopping test to output stream.
std::vector< int > convIndices()
Returns the vector containing the indices of the residuals that passed the test.
bool getShowMaxResNormOnly()
Returns whether the only maximum residual norm is displayed when the print() method is called...