42 #ifndef BELOS_STATUS_TEST_GEN_RESNORM_MP_VECTOR_H 43 #define BELOS_STATUS_TEST_GEN_RESNORM_MP_VECTOR_H 54 #include "BelosStatusTestGenResNorm.hpp" 84 template <
class Storage,
class MV,
class OP>
86 public StatusTestResNorm<Sacado::MP::Vector<Storage>,MV,OP> {
92 typedef Teuchos::ScalarTraits<ScalarType>
SCT;
94 typedef MultiVecTraits<ScalarType,MV>
MVT;
125 StatusTestGenResNorm( MagnitudeType Tolerance,
int quorum = -1,
bool showMaxResNormOnly =
false );
142 int defineResForm( ResType TypeOfResidual, NormType TypeOfNorm);
165 int defineScaleForm( ScaleType TypeOfScaling, NormType TypeOfNorm, MagnitudeType ScaleValue = Teuchos::ScalarTraits<MagnitudeType>::one());
175 int setQuorum(
int quorum) {quorum_ = quorum;
return(0);}
191 StatusType checkStatus(Iteration<ScalarType,MV,OP>* iSolver);
209 void print(std::ostream& os,
int indent = 0)
const;
212 void printStatus(std::ostream& os, StatusType type)
const;
221 Teuchos::RCP<MV>
getSolution() {
if (restype_==Implicit) {
return Teuchos::null; }
else {
return curSoln_; } }
237 const std::vector<MagnitudeType>*
getTestValue()
const {
return(&testvector_);};
263 StatusType firstCallCheckStatusSetup(Iteration<ScalarType,MV,OP>* iSolver);
272 std::ostringstream oss;
273 oss <<
"Belos::StatusTestGenResNorm<>: " << resFormStr();
274 oss <<
", tol = " << tolerance_;
288 std::ostringstream oss;
290 oss << ((resnormtype_==OneNorm) ?
"1-Norm" : (resnormtype_==TwoNorm) ?
"2-Norm" :
"Inf-Norm");
291 oss << ((restype_==Explicit) ?
" Exp" :
" Imp");
295 if (scaletype_!=
None)
301 if (scaletype_==UserProvided)
302 oss <<
" (User Scale)";
305 oss << ((scalenormtype_==OneNorm) ?
"1-Norm" : (resnormtype_==TwoNorm) ?
"2-Norm" :
"Inf-Norm");
306 if (scaletype_==NormOfInitRes || scaletype_==NormOfFullInitRes || scaletype_==NormOfFullScaledInitRes)
308 else if (scaletype_==NormOfPrecInitRes || scaletype_==NormOfFullPrecInitRes || scaletype_==NormOfFullScaledPrecInitRes)
400 template <
class StorageType,
class MV,
class OP>
403 : tolerance_(Tolerance),
405 showMaxResNormOnly_(showMaxResNormOnly),
407 resnormtype_(TwoNorm),
408 scaletype_(NormOfInitRes),
409 scalenormtype_(TwoNorm),
410 scalevalue_(
Teuchos::ScalarTraits<MagnitudeType>::one ()),
416 firstcallCheckStatus_(
true),
417 firstcallDefineResForm_(
true),
418 firstcallDefineScaleForm_(
true),
419 ensemble_converged(StorageType::static_size, 0),
420 ensemble_iterations(StorageType::static_size, 0)
426 template <
class StorageType,
class MV,
class OP>
427 StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::~StatusTestGenResNorm()
430 template <
class StorageType,
class MV,
class OP>
431 void StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::reset()
439 firstcallCheckStatus_ =
true;
440 curSoln_ = Teuchos::null;
441 ensemble_converged = std::vector<int>(StorageType::static_size, 0);
442 ensemble_iterations = std::vector<int>(StorageType::static_size, 0);
445 template <
class StorageType,
class MV,
class OP>
446 int StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::defineResForm( ResType TypeOfResidual, NormType TypeOfNorm )
448 TEUCHOS_TEST_FOR_EXCEPTION(firstcallDefineResForm_==
false,StatusTestError,
449 "StatusTestGenResNorm::defineResForm(): The residual form has already been defined.");
450 firstcallDefineResForm_ =
false;
452 restype_ = TypeOfResidual;
453 resnormtype_ = TypeOfNorm;
458 template <
class StorageType,
class MV,
class OP>
459 int StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::defineScaleForm(ScaleType TypeOfScaling, NormType TypeOfNorm,
460 MagnitudeType ScaleValue )
462 TEUCHOS_TEST_FOR_EXCEPTION(firstcallDefineScaleForm_==
false,StatusTestError,
463 "StatusTestGenResNorm::defineScaleForm(): The scaling type has already been defined.");
464 firstcallDefineScaleForm_ =
false;
466 scaletype_ = TypeOfScaling;
467 scalenormtype_ = TypeOfNorm;
468 scalevalue_ = ScaleValue;
473 template <
class StorageType,
class MV,
class OP>
474 StatusType StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::checkStatus( Iteration<ScalarType,MV,OP>* iSolver )
476 MagnitudeType zero = Teuchos::ScalarTraits<MagnitudeType>::zero();
477 const LinearProblem<ScalarType,MV,OP>& lp = iSolver->getProblem();
479 if (firstcallCheckStatus_) {
480 StatusType status = firstCallCheckStatusSetup(iSolver);
489 if ( curLSNum_ != lp.getLSNumber() ) {
493 curLSNum_ = lp.getLSNumber();
494 curLSIdx_ = lp.getLSIndex();
495 curBlksz_ = (
int)curLSIdx_.size();
497 for (
int i=0; i<curBlksz_; ++i) {
498 if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
501 curNumRHS_ = validLS;
502 curSoln_ = Teuchos::null;
508 if (status_==Passed) {
return status_; }
510 if (restype_==Implicit) {
516 std::vector<MagnitudeType> tmp_resvector( curBlksz_ );
517 Teuchos::RCP<const MV> residMV = iSolver->getNativeResiduals( &tmp_resvector );
518 if ( residMV != Teuchos::null ) {
519 tmp_resvector.resize( MVT::GetNumberVecs( *residMV ) );
520 MVT::MvNorm( *residMV, tmp_resvector, resnormtype_ );
521 typename std::vector<int>::iterator p = curLSIdx_.begin();
522 for (
int i=0; p<curLSIdx_.end(); ++p, ++i) {
525 resvector_[*p] = tmp_resvector[i];
528 typename std::vector<int>::iterator p = curLSIdx_.begin();
529 for (
int i=0; p<curLSIdx_.end(); ++p, ++i) {
532 resvector_[*p] = tmp_resvector[i];
536 else if (restype_==Explicit) {
540 Teuchos::RCP<MV> cur_update = iSolver->getCurrentUpdate();
541 curSoln_ = lp.updateSolution( cur_update );
542 Teuchos::RCP<MV> cur_res = MVT::Clone( *curSoln_, MVT::GetNumberVecs( *curSoln_ ) );
543 lp.computeCurrResVec( &*cur_res, &*curSoln_ );
544 std::vector<MagnitudeType> tmp_resvector( MVT::GetNumberVecs( *cur_res ) );
545 MVT::MvNorm( *cur_res, tmp_resvector, resnormtype_ );
546 typename std::vector<int>::iterator p = curLSIdx_.begin();
547 for (
int i=0; p<curLSIdx_.end(); ++p, ++i) {
550 resvector_[*p] = tmp_resvector[i];
557 if ( scalevector_.size() > 0 ) {
558 typename std::vector<int>::iterator p = curLSIdx_.begin();
559 for (; p<curLSIdx_.end(); ++p) {
563 if ( scalevector_[ *p ] != zero ) {
565 testvector_[ *p ] = resvector_[ *p ] / scalevector_[ *p ] / scalevalue_;
567 testvector_[ *p ] = resvector_[ *p ] / scalevalue_;
573 typename std::vector<int>::iterator p = curLSIdx_.begin();
574 for (; p<curLSIdx_.end(); ++p) {
577 testvector_[ *p ] = resvector_[ *p ] / scalevalue_;
583 ind_.resize( curLSIdx_.size() );
584 typename std::vector<int>::iterator p = curLSIdx_.begin();
585 for (; p<curLSIdx_.end(); ++p) {
589 const int ensemble_size = StorageType::static_size;
590 bool all_converged =
true;
591 for (
int i=0; i<ensemble_size; ++i) {
592 if (!ensemble_converged[i]) {
593 if (testvector_[ *p ].coeff(i) > tolerance_) {
594 ++ensemble_iterations[i];
595 all_converged =
false;
597 else if (testvector_[ *p ].coeff(i) <= tolerance_) {
598 ensemble_converged[i] = 1;
603 TEUCHOS_TEST_FOR_EXCEPTION(
true,StatusTestError,
"StatusTestGenResNorm::checkStatus(): NaN has been detected.");
614 int need = (quorum_ == -1) ? curNumRHS_: quorum_;
615 status_ = (have >= need) ? Passed : Failed;
621 template <
class StorageType,
class MV,
class OP>
622 void StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::print(std::ostream& os,
int indent)
const 624 for (
int j = 0;
j < indent;
j ++)
626 printStatus(os, status_);
628 if (status_==Undefined)
629 os <<
", tol = " << tolerance_ << std::endl;
632 if(showMaxResNormOnly_ && curBlksz_ > 1) {
633 const MagnitudeType maxRelRes = *std::max_element(
634 testvector_.begin()+curLSIdx_[0],testvector_.begin()+curLSIdx_[curBlksz_-1]
636 for (
int j = 0;
j < indent + 13;
j ++)
638 os <<
"max{residual["<<curLSIdx_[0]<<
"..."<<curLSIdx_[curBlksz_-1]<<
"]} = " << maxRelRes
639 << ( maxRelRes <= tolerance_ ?
" <= " :
" > " ) << tolerance_ << std::endl;
642 for (
int i=0; i<numrhs_; i++ ) {
643 for (
int j = 0;
j < indent + 13;
j ++)
645 os <<
"residual [ " << i <<
" ] = " << testvector_[ i ];
646 os << ((testvector_[i]<tolerance_) ?
" < " : (testvector_[i]==tolerance_) ?
" == " : (testvector_[i]>tolerance_) ?
" > " :
" " ) << tolerance_ << std::endl;
653 template <
class StorageType,
class MV,
class OP>
654 void StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::printStatus(std::ostream& os, StatusType type)
const 656 os << std::left << std::setw(13) << std::setfill(
'.');
669 os << std::left << std::setfill(
' ');
673 template <
class StorageType,
class MV,
class OP>
674 StatusType StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::firstCallCheckStatusSetup( Iteration<ScalarType,MV,OP>* iSolver )
677 MagnitudeType zero = Teuchos::ScalarTraits<MagnitudeType>::zero();
678 MagnitudeType one = Teuchos::ScalarTraits<MagnitudeType>::one();
679 const LinearProblem<ScalarType,MV,OP>& lp = iSolver->getProblem();
681 if (firstcallCheckStatus_) {
685 firstcallCheckStatus_ =
false;
687 if (scaletype_== NormOfRHS) {
688 Teuchos::RCP<const MV> rhs = lp.getRHS();
689 numrhs_ = MVT::GetNumberVecs( *rhs );
690 scalevector_.resize( numrhs_ );
691 MVT::MvNorm( *rhs, scalevector_, scalenormtype_ );
693 else if (scaletype_==NormOfInitRes || scaletype_==NormOfFullInitRes || scaletype_==NormOfFullScaledInitRes) {
694 Teuchos::RCP<const MV> init_res = lp.getInitResVec();
695 numrhs_ = MVT::GetNumberVecs( *init_res );
696 scalevector_.resize( numrhs_ );
697 MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
699 else if (scaletype_==NormOfPrecInitRes || scaletype_==NormOfFullPrecInitRes || scaletype_==NormOfFullScaledPrecInitRes) {
700 Teuchos::RCP<const MV> init_res = lp.getInitPrecResVec();
701 numrhs_ = MVT::GetNumberVecs( *init_res );
702 scalevector_.resize( numrhs_ );
703 MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
706 numrhs_ = MVT::GetNumberVecs( *(lp.getRHS()) );
709 resvector_.resize( numrhs_ );
710 testvector_.resize( numrhs_ );
712 curLSNum_ = lp.getLSNumber();
713 curLSIdx_ = lp.getLSIndex();
714 curBlksz_ = (
int)curLSIdx_.size();
716 for (i=0; i<curBlksz_; ++i) {
717 if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
720 curNumRHS_ = validLS;
723 for (i=0; i<numrhs_; i++) { testvector_[i] = one; }
726 if (scalevalue_ == zero) {
std::string description() const
Method to return description of the maximum iteration status test.
int curLSNum_
The current number of linear systems that have been loaded into the linear problem.
std::string resFormStr() const
Description of current residual form.
StatusType getStatus() const
Return the result of the most recent CheckStatus call.
NormType scalenormtype_
Type of norm to use on the scaling (OneNorm, TwoNorm, or InfNorm)
int curBlksz_
The current blocksize of the linear system being solved.
std::vector< int > ind_
Vector containing the indices for the vectors that passed the test.
int setTolerance(MagnitudeType tolerance)
Set the value of the tolerance.
Teuchos::RCP< MV > getSolution()
StatusType status_
Status.
const std::vector< MagnitudeType > * getResNormValue() const
Returns the residual norm value, , computed in most recent call to CheckStatus.
bool firstcallCheckStatus_
Is this the first time CheckStatus is called?
std::vector< MagnitudeType > resvector_
Residual norm std::vector.
bool getShowMaxResNormOnly()
Returns 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.
std::vector< MagnitudeType > testvector_
Test std::vector = resvector_ / scalevector_.
An implementation of StatusTestResNorm using a family of residual norms.
bool getLOADetected() const
int curNumRHS_
The current number of right-hand sides being solved for.
ResType restype_
Type of residual to use (explicit or implicit)
std::vector< int > ensemble_iterations
The number of iterations at which point each ensemble component converges.
bool showMaxResNormOnly_
Determines if the entries for all of the residuals are shown or just the max.
int setShowMaxResNormOnly(bool showMaxResNormOnly)
Set whether the only maximum residual norm is displayed when the print() method is called...
MultiVecTraits< ScalarType, MV > MVT
ResType
Select how the residual std::vector is produced.
const std::vector< MagnitudeType > * getScaledNormValue() const
Returns the scaled norm value, .
const std::vector< MagnitudeType > * getTestValue() const
Returns the test value, , computed in most recent call to CheckStatus.
int quorum_
Number of residuals that must pass the convergence test before Passed is returned.
MagnitudeType scalevalue_
Scaling value.
SCT::magnitudeType MagnitudeType
MagnitudeType getTolerance() const
Returns the value of the tolerance, , set in the constructor.
NormType resnormtype_
Type of norm to use on residual (OneNorm, TwoNorm, or InfNorm).
std::vector< int > convIndices()
Returns the std::vector containing the indices of the residuals that passed the test.
MagnitudeType tolerance_
Tolerance used to determine convergence.
ScaleType scaletype_
Type of scaling to use (Norm of RHS, Norm of Initial Residual, None or User provided) ...
const std::vector< int > getEnsembleIterations() const
Returns number of ensemble iterations.
bool firstcallDefineResForm_
Is this the first time DefineResForm is called?
Teuchos::ScalarTraits< ScalarType > SCT
int setQuorum(int quorum)
std::vector< MagnitudeType > scalevector_
Scaling std::vector.
Sacado::MP::Vector< Storage > ScalarType
int numrhs_
The total number of right-hand sides being solved for.
std::vector< int > ensemble_converged
Which ensemble components have converged.
bool firstcallDefineScaleForm_
Is this the first time DefineScaleForm is called?
Teuchos::RCP< MV > curSoln_
Most recent solution vector used by this status test.