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());
177 int setQuorum(
int quorum) {quorum_ = quorum;
return(0);}
193 StatusType checkStatus(Iteration<ScalarType,MV,OP>* iSolver);
211 void print(std::ostream& os,
int indent = 0)
const;
214 void printStatus(std::ostream& os, StatusType type)
const;
223 Teuchos::RCP<MV>
getSolution() {
if (restype_==Implicit) {
return Teuchos::null; }
else {
return curSoln_; } }
239 const std::vector<MagnitudeType>*
getTestValue()
const {
return(&testvector_);}
265 StatusType firstCallCheckStatusSetup(Iteration<ScalarType,MV,OP>* iSolver);
274 std::ostringstream oss;
275 oss <<
"Belos::StatusTestGenResNorm<>: " << resFormStr();
276 oss <<
", tol = " << tolerance_;
290 std::ostringstream oss;
292 oss << ((resnormtype_==OneNorm) ?
"1-Norm" : (resnormtype_==TwoNorm) ?
"2-Norm" :
"Inf-Norm");
293 oss << ((restype_==Explicit) ?
" Exp" :
" Imp");
297 if (scaletype_!=
None)
303 if (scaletype_==UserProvided)
304 oss <<
" (User Scale)";
307 oss << ((scalenormtype_==OneNorm) ?
"1-Norm" : (resnormtype_==TwoNorm) ?
"2-Norm" :
"Inf-Norm");
308 if (scaletype_==NormOfInitRes || scaletype_==NormOfFullInitRes || scaletype_==NormOfFullScaledInitRes)
310 else if (scaletype_==NormOfPrecInitRes || scaletype_==NormOfFullPrecInitRes || scaletype_==NormOfFullScaledPrecInitRes)
402 template <
class StorageType,
class MV,
class OP>
405 : tolerance_(Tolerance),
407 showMaxResNormOnly_(showMaxResNormOnly),
409 resnormtype_(TwoNorm),
410 scaletype_(NormOfInitRes),
411 scalenormtype_(TwoNorm),
412 scalevalue_(
Teuchos::ScalarTraits<MagnitudeType>::one ()),
418 firstcallCheckStatus_(
true),
419 firstcallDefineResForm_(
true),
420 firstcallDefineScaleForm_(
true),
421 ensemble_converged(StorageType::static_size, 0),
422 ensemble_iterations(StorageType::static_size, 0)
428 template <
class StorageType,
class MV,
class OP>
429 StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::~StatusTestGenResNorm()
432 template <
class StorageType,
class MV,
class OP>
433 void StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::reset()
441 firstcallCheckStatus_ =
true;
442 curSoln_ = Teuchos::null;
443 ensemble_converged = std::vector<int>(StorageType::static_size, 0);
444 ensemble_iterations = std::vector<int>(StorageType::static_size, 0);
447 template <
class StorageType,
class MV,
class OP>
448 int StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::defineResForm( ResType TypeOfResidual, NormType TypeOfNorm )
450 TEUCHOS_TEST_FOR_EXCEPTION(firstcallDefineResForm_==
false,StatusTestError,
451 "StatusTestGenResNorm::defineResForm(): The residual form has already been defined.");
452 firstcallDefineResForm_ =
false;
454 restype_ = TypeOfResidual;
455 resnormtype_ = TypeOfNorm;
460 template <
class StorageType,
class MV,
class OP>
461 int StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::defineScaleForm(ScaleType TypeOfScaling, NormType TypeOfNorm,
462 MagnitudeType ScaleValue )
464 TEUCHOS_TEST_FOR_EXCEPTION(firstcallDefineScaleForm_==
false,StatusTestError,
465 "StatusTestGenResNorm::defineScaleForm(): The scaling type has already been defined.");
466 firstcallDefineScaleForm_ =
false;
468 scaletype_ = TypeOfScaling;
469 scalenormtype_ = TypeOfNorm;
470 scalevalue_ = ScaleValue;
475 template <
class StorageType,
class MV,
class OP>
476 StatusType StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::checkStatus( Iteration<ScalarType,MV,OP>* iSolver )
478 MagnitudeType zero = Teuchos::ScalarTraits<MagnitudeType>::zero();
479 const LinearProblem<ScalarType,MV,OP>& lp = iSolver->getProblem();
481 if (firstcallCheckStatus_) {
482 StatusType status = firstCallCheckStatusSetup(iSolver);
491 if ( curLSNum_ != lp.getLSNumber() ) {
495 curLSNum_ = lp.getLSNumber();
496 curLSIdx_ = lp.getLSIndex();
497 curBlksz_ = (int)curLSIdx_.size();
499 for (
int i=0; i<curBlksz_; ++i) {
500 if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
503 curNumRHS_ = validLS;
504 curSoln_ = Teuchos::null;
510 if (status_==Passed) {
return status_; }
512 if (restype_==Implicit) {
518 std::vector<MagnitudeType> tmp_resvector( curBlksz_ );
519 Teuchos::RCP<const MV> residMV = iSolver->getNativeResiduals( &tmp_resvector );
520 if ( residMV != Teuchos::null ) {
521 tmp_resvector.resize( MVT::GetNumberVecs( *residMV ) );
522 MVT::MvNorm( *residMV, tmp_resvector, resnormtype_ );
523 typename std::vector<int>::iterator p = curLSIdx_.begin();
524 for (
int i=0; p<curLSIdx_.end(); ++p, ++i) {
527 resvector_[*p] = tmp_resvector[i];
530 typename std::vector<int>::iterator p = curLSIdx_.begin();
531 for (
int i=0; p<curLSIdx_.end(); ++p, ++i) {
534 resvector_[*p] = tmp_resvector[i];
538 else if (restype_==Explicit) {
542 Teuchos::RCP<MV> cur_update = iSolver->getCurrentUpdate();
543 curSoln_ = lp.updateSolution( cur_update );
544 Teuchos::RCP<MV> cur_res = MVT::Clone( *curSoln_, MVT::GetNumberVecs( *curSoln_ ) );
545 lp.computeCurrResVec( &*cur_res, &*curSoln_ );
546 std::vector<MagnitudeType> tmp_resvector( MVT::GetNumberVecs( *cur_res ) );
547 MVT::MvNorm( *cur_res, tmp_resvector, resnormtype_ );
548 typename std::vector<int>::iterator p = curLSIdx_.begin();
549 for (
int i=0; p<curLSIdx_.end(); ++p, ++i) {
552 resvector_[*p] = tmp_resvector[i];
559 if ( scalevector_.size() > 0 ) {
560 typename std::vector<int>::iterator p = curLSIdx_.begin();
561 for (; p<curLSIdx_.end(); ++p) {
565 if ( scalevector_[ *p ] != zero ) {
567 testvector_[ *p ] = resvector_[ *p ] / scalevector_[ *p ] / scalevalue_;
569 testvector_[ *p ] = resvector_[ *p ] / scalevalue_;
575 typename std::vector<int>::iterator p = curLSIdx_.begin();
576 for (; p<curLSIdx_.end(); ++p) {
579 testvector_[ *p ] = resvector_[ *p ] / scalevalue_;
585 ind_.resize( curLSIdx_.size() );
586 typename std::vector<int>::iterator p = curLSIdx_.begin();
587 for (; p<curLSIdx_.end(); ++p) {
591 const int ensemble_size = StorageType::static_size;
592 bool all_converged =
true;
593 for (
int i=0; i<ensemble_size; ++i) {
594 if (!ensemble_converged[i]) {
595 if (testvector_[ *p ].coeff(i) > tolerance_) {
596 ++ensemble_iterations[i];
597 all_converged =
false;
599 else if (testvector_[ *p ].coeff(i) <= tolerance_) {
600 ensemble_converged[i] = 1;
605 TEUCHOS_TEST_FOR_EXCEPTION(
true,StatusTestError,
"StatusTestGenResNorm::checkStatus(): NaN has been detected.");
616 int need = (quorum_ == -1) ? curNumRHS_: quorum_;
617 status_ = (have >= need) ? Passed : Failed;
623 template <
class StorageType,
class MV,
class OP>
624 void StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::print(std::ostream& os,
int indent)
const 626 for (
int j = 0;
j < indent;
j ++)
628 printStatus(os, status_);
630 if (status_==Undefined)
631 os <<
", tol = " << tolerance_ << std::endl;
634 if(showMaxResNormOnly_ && curBlksz_ > 1) {
635 const MagnitudeType maxRelRes = *std::max_element(
636 testvector_.begin()+curLSIdx_[0],testvector_.begin()+curLSIdx_[curBlksz_-1]
638 for (
int j = 0;
j < indent + 13;
j ++)
640 os <<
"max{residual["<<curLSIdx_[0]<<
"..."<<curLSIdx_[curBlksz_-1]<<
"]} = " << maxRelRes
641 << ( maxRelRes <= tolerance_ ?
" <= " :
" > " ) << tolerance_ << std::endl;
644 for (
int i=0; i<numrhs_; i++ ) {
645 for (
int j = 0;
j < indent + 13;
j ++)
647 os <<
"residual [ " << i <<
" ] = " << testvector_[ i ];
648 os << ((testvector_[i]<tolerance_) ?
" < " : (testvector_[i]==tolerance_) ?
" == " : (testvector_[i]>tolerance_) ?
" > " :
" " ) << tolerance_ << std::endl;
655 template <
class StorageType,
class MV,
class OP>
656 void StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::printStatus(std::ostream& os, StatusType type)
const 658 os << std::left << std::setw(13) << std::setfill(
'.');
671 os << std::left << std::setfill(
' ');
675 template <
class StorageType,
class MV,
class OP>
676 StatusType StatusTestGenResNorm<Sacado::MP::Vector<StorageType>,MV,OP>::firstCallCheckStatusSetup( Iteration<ScalarType,MV,OP>* iSolver )
679 MagnitudeType zero = Teuchos::ScalarTraits<MagnitudeType>::zero();
680 MagnitudeType one = Teuchos::ScalarTraits<MagnitudeType>::one();
681 const LinearProblem<ScalarType,MV,OP>& lp = iSolver->getProblem();
683 if (firstcallCheckStatus_) {
687 firstcallCheckStatus_ =
false;
689 if (scaletype_== NormOfRHS) {
690 Teuchos::RCP<const MV> rhs = lp.getRHS();
691 numrhs_ = MVT::GetNumberVecs( *rhs );
692 scalevector_.resize( numrhs_ );
693 MVT::MvNorm( *rhs, scalevector_, scalenormtype_ );
695 else if (scaletype_==NormOfInitRes || scaletype_==NormOfFullInitRes || scaletype_==NormOfFullScaledInitRes) {
696 Teuchos::RCP<const MV> init_res = lp.getInitResVec();
697 numrhs_ = MVT::GetNumberVecs( *init_res );
698 scalevector_.resize( numrhs_ );
699 MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
701 else if (scaletype_==NormOfPrecInitRes || scaletype_==NormOfFullPrecInitRes || scaletype_==NormOfFullScaledPrecInitRes) {
702 Teuchos::RCP<const MV> init_res = lp.getInitPrecResVec();
703 numrhs_ = MVT::GetNumberVecs( *init_res );
704 scalevector_.resize( numrhs_ );
705 MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
708 numrhs_ = MVT::GetNumberVecs( *(lp.getRHS()) );
711 resvector_.resize( numrhs_ );
712 testvector_.resize( numrhs_ );
714 curLSNum_ = lp.getLSNumber();
715 curLSIdx_ = lp.getLSIndex();
716 curBlksz_ = (int)curLSIdx_.size();
718 for (i=0; i<curBlksz_; ++i) {
719 if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
722 curNumRHS_ = validLS;
725 for (i=0; i<numrhs_; i++) { testvector_[i] = one; }
728 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.
NormType getResNormType()
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.