42 #ifndef BELOS_STATUS_TEST_GEN_RESSUBNORM_H
43 #define BELOS_STATUS_TEST_GEN_RESSUBNORM_H
55 #ifdef HAVE_BELOS_THYRA
56 #include <Thyra_MultiVectorBase.hpp>
57 #include <Thyra_MultiVectorStdOps.hpp>
58 #include <Thyra_ProductMultiVectorBase.hpp>
71 template <
class ScalarType,
class MV,
class OP>
98 "StatusTestGenResSubNorm::StatusTestGenResSubNorm(): StatusTestGenResSubNorm only available for blocked operators (e.g., Thyra).");
117 "StatusTestGenResSubNorm::defineResForm(): StatusTestGenResSubNorm only available for blocked operators (e.g., Thyra).");
144 "StatusTestGenResSubNorm::defineScaleForm(): StatusTestGenResSubNorm only available for blocked operators (e.g., Thyra).");
157 int setSubIdx (
size_t subIdx ) {
return 0;}
195 void print(std::ostream& ,
int = 0)
const { }
214 size_t getSubIdx()
const {
return 0; }
220 std::vector<int>
convIndices() {
return std::vector<int>(0); }
226 const std::vector<MagnitudeType>*
getTestValue()
const {
return NULL;};
229 const std::vector<MagnitudeType>* getResNormValue()
const {
return NULL;};
232 const std::vector<MagnitudeType>* getScaledNormValue()
const {
return NULL;};
249 StatusType firstCallCheckStatusSetup(Iteration<ScalarType,MV,OP>* iSolver) {
259 {
return std::string(
""); }
263 #ifdef HAVE_BELOS_THYRA
266 template <
class ScalarType>
267 class StatusTestGenResSubNorm<ScalarType,Thyra::MultiVectorBase<ScalarType>,Thyra::LinearOpBase<ScalarType> >
268 :
public StatusTestResNorm<ScalarType,Thyra::MultiVectorBase<ScalarType>,Thyra::LinearOpBase<ScalarType> > {
272 typedef Thyra::MultiVectorBase<ScalarType> MV;
273 typedef Thyra::LinearOpBase<ScalarType> OP;
277 typedef MultiVecTraits<ScalarType,MV>
MVT;
278 typedef OperatorTraits<ScalarType,MV,OP> OT;
295 StatusTestGenResSubNorm(
MagnitudeType Tolerance,
size_t subIdx,
int quorum = -1,
bool showMaxResNormOnly =
false )
296 : tolerance_(Tolerance),
299 showMaxResNormOnly_(showMaxResNormOnly),
309 firstcallCheckStatus_(true),
310 firstcallDefineResForm_(true),
311 firstcallDefineScaleForm_(true) { }
314 virtual ~StatusTestGenResSubNorm() { };
327 int defineResForm(
NormType TypeOfNorm) {
329 "StatusTestGenResSubNorm::defineResForm(): The residual form has already been defined.");
330 firstcallDefineResForm_ =
false;
332 resnormtype_ = TypeOfNorm;
360 "StatusTestGenResSubNorm::defineScaleForm(): The scaling type has already been defined.");
361 firstcallDefineScaleForm_ =
false;
363 scaletype_ = TypeOfScaling;
364 scalenormtype_ = TypeOfNorm;
365 scalevalue_ = ScaleValue;
379 int setSubIdx (
size_t subIdx ) { subIdx_ = subIdx;
return(0);}
383 int setQuorum(
int quorum) {quorum_ = quorum;
return(0);}
386 int setShowMaxResNormOnly(
bool showMaxResNormOnly) {showMaxResNormOnly_ = showMaxResNormOnly;
return(0);}
401 const LinearProblem<ScalarType,MV,OP>& lp = iSolver->getProblem();
403 if (firstcallCheckStatus_) {
404 StatusType status = firstCallCheckStatusSetup(iSolver);
414 if ( curLSNum_ != lp.getLSNumber() ) {
418 curLSNum_ = lp.getLSNumber();
419 curLSIdx_ = lp.getLSIndex();
420 curBlksz_ = (int)curLSIdx_.size();
422 for (
int i=0; i<curBlksz_; ++i) {
423 if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
426 curNumRHS_ = validLS;
427 curSoln_ = Teuchos::null;
433 if (status_==
Passed) {
return status_; }
440 curSoln_ = lp.updateSolution( cur_update );
442 lp.computeCurrResVec( &*cur_res, &*curSoln_ );
444 MvSubNorm( *cur_res, subIdx_, tmp_resvector, resnormtype_ );
446 typename std::vector<int>::iterator pp = curLSIdx_.begin();
447 for (
int i=0; pp<curLSIdx_.end(); ++pp, ++i) {
450 resvector_[*pp] = tmp_resvector[i];
457 if ( scalevector_.size() > 0 ) {
458 typename std::vector<int>::iterator p = curLSIdx_.begin();
459 for (; p<curLSIdx_.end(); ++p) {
463 if ( scalevector_[ *p ] != zero ) {
465 testvector_[ *p ] = resvector_[ *p ] / scalevector_[ *p ] / scalevalue_;
467 testvector_[ *p ] = resvector_[ *p ] / scalevalue_;
473 typename std::vector<int>::iterator ppp = curLSIdx_.begin();
474 for (; ppp<curLSIdx_.end(); ++ppp) {
477 testvector_[ *ppp ] = resvector_[ *ppp ] / scalevalue_;
482 ind_.resize( curLSIdx_.size() );
483 typename std::vector<int>::iterator p2 = curLSIdx_.begin();
484 for (; p2<curLSIdx_.end(); ++p2) {
488 if (testvector_[ *p2 ] > tolerance_) {
490 }
else if (testvector_[ *p2 ] <= tolerance_) {
501 int need = (quorum_ == -1) ? curNumRHS_: quorum_;
522 firstcallCheckStatus_ =
true;
523 curSoln_ = Teuchos::null;
532 void print(std::ostream& os,
int indent = 0)
const {
533 os.setf(std::ios_base::scientific);
534 for (
int j = 0; j < indent; j ++)
539 os <<
", tol = " << tolerance_ << std::endl;
542 if(showMaxResNormOnly_ && curBlksz_ > 1) {
544 testvector_.begin()+curLSIdx_[0],testvector_.begin()+curLSIdx_[curBlksz_-1]
546 for (
int j = 0; j < indent + 13; j ++)
548 os <<
"max{residual["<<curLSIdx_[0]<<
"..."<<curLSIdx_[curBlksz_-1]<<
"]} = " << maxRelRes
549 << ( maxRelRes <= tolerance_ ?
" <= " :
" > " ) << tolerance_ << std::endl;
552 for (
int i=0; i<numrhs_; i++ ) {
553 for (
int j = 0; j < indent + 13; j ++)
555 os <<
"residual [ " << i <<
" ] = " << testvector_[ i ];
556 os << ((testvector_[i]<tolerance_) ?
" < " : (testvector_[i]==tolerance_) ?
" == " : (testvector_[i]>tolerance_) ?
" > " :
" " ) << tolerance_ << std::endl;
565 os << std::left << std::setw(13) << std::setfill(
'.');
578 os << std::left << std::setfill(
' ');
591 int getQuorum()
const {
return quorum_; }
594 size_t getSubIdx()
const {
return subIdx_; }
606 const std::vector<MagnitudeType>*
getTestValue()
const {
return(&testvector_);};
609 const std::vector<MagnitudeType>* getResNormValue()
const {
return(&resvector_);};
612 const std::vector<MagnitudeType>* getScaledNormValue()
const {
return(&scalevector_);};
629 StatusType firstCallCheckStatusSetup(Iteration<ScalarType,MV,OP>* iSolver) {
633 const LinearProblem<ScalarType,MV,OP>& lp = iSolver->getProblem();
635 if (firstcallCheckStatus_) {
639 firstcallCheckStatus_ =
false;
644 scalevector_.resize( numrhs_ );
645 MvSubNorm( *rhs, subIdx_, scalevector_, scalenormtype_ );
650 scalevector_.resize( numrhs_ );
651 MvSubNorm( *init_res, subIdx_, scalevector_, scalenormtype_ );
656 scalevector_.resize( numrhs_ );
657 MvSubNorm( *init_res, subIdx_, scalevector_, scalenormtype_ );
662 scalevector_.resize( numrhs_ );
663 MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
669 scalevector_.resize( numrhs_ );
670 MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
676 scalevector_.resize( numrhs_ );
677 MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
678 MvScalingRatio( *init_res, subIdx_, scalevalue_ );
683 scalevector_.resize( numrhs_ );
684 MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
685 MvScalingRatio( *init_res, subIdx_, scalevalue_ );
691 resvector_.resize( numrhs_ );
692 testvector_.resize( numrhs_ );
694 curLSNum_ = lp.getLSNumber();
695 curLSIdx_ = lp.getLSIndex();
696 curBlksz_ = (int)curLSIdx_.size();
698 for (i=0; i<curBlksz_; ++i) {
699 if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
702 curNumRHS_ = validLS;
705 for (i=0; i<numrhs_; i++) { testvector_[i] = one; }
708 if (scalevalue_ == zero) {
722 std::ostringstream oss;
723 oss <<
"Belos::StatusTestGenResSubNorm<>: " << resFormStr();
724 oss <<
", tol = " << tolerance_;
736 std::string resFormStr()
const
738 std::ostringstream oss;
740 oss << ((resnormtype_==
OneNorm) ?
"1-Norm" : (resnormtype_==
TwoNorm) ?
"2-Norm" :
"Inf-Norm");
742 oss <<
" Res Vec [" << subIdx_ <<
"]) ";
745 if (scaletype_!=
None)
752 oss <<
" (User Scale)";
755 oss << ((scalenormtype_==
OneNorm) ?
"1-Norm" : (resnormtype_==
TwoNorm) ?
"2-Norm" :
"Inf-Norm");
757 oss <<
" Res0 [" << subIdx_ <<
"]";
759 oss <<
" Prec Res0 [" << subIdx_ <<
"]";
761 oss <<
" Full Res0 [" << subIdx_ <<
"]";
763 oss <<
" Full Prec Res0 [" << subIdx_ <<
"]";
765 oss <<
" scaled Full Res0 [" << subIdx_ <<
"]";
767 oss <<
" scaled Full Prec Res0 [" << subIdx_ <<
"]";
769 oss <<
" RHS [" << subIdx_ <<
"]";
788 typedef typename Thyra::ProductMultiVectorBase<ScalarType> TPMVB;
792 "Belos::StatusTestGenResSubNorm::MvSubNorm (Thyra specialization): "
793 "mv must be a Thyra::ProductMultiVector, but is of type " << thyProdVec);
797 typedef MultiVecTraits<ScalarType, MV>
MVT;
802 void MvScalingRatio(
const MV& mv,
size_t block,
MagnitudeType& lengthRatio) {
805 typedef typename Thyra::ProductMultiVectorBase<ScalarType> TPMVB;
809 "Belos::StatusTestGenResSubNorm::MvScalingRatio (Thyra specialization): "
810 "mv must be a Thyra::ProductMultiVector, but is of type " << thyProdVec);
814 lengthRatio = Teuchos::as<MagnitudeType>(thySubVec->range()->dim()) / Teuchos::as<MagnitudeType>(thyProdVec->range()->dim());
832 bool showMaxResNormOnly_;
847 std::vector<MagnitudeType> scalevector_;
850 std::vector<MagnitudeType> resvector_;
853 std::vector<MagnitudeType> testvector_;
856 std::vector<int> ind_;
871 std::vector<int> curLSIdx_;
880 bool firstcallCheckStatus_;
883 bool firstcallDefineResForm_;
886 bool firstcallDefineScaleForm_;
892 #endif // HAVE_BELOS_THYRA
virtual std::vector< int > convIndices()=0
Returns the std::vector containing the indices of the residuals that passed the test.
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().
virtual int setQuorum(int quorum)=0
Sets the number of residuals that must pass the convergence test before Passed is returned...
virtual int getQuorum() const =0
Returns the number of residuals that must pass the convergence test before Passed is returned...
An implementation of StatusTestResNorm using a family of norms of subvectors of the residual vectors...
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Declaration of basic traits for the multivector type.
An abstract class of StatusTest for stopping criteria using residual norms.
virtual MagnitudeType getTolerance() const =0
Returns the value of the tolerance, , set in the constructor.
virtual StatusType checkStatus(Iteration< ScalarType, MV, OP > *iSolver)=0
Check convergence status: Unconverged, Converged, Failed.
Class which defines basic traits for the operator type.
StatusType
Whether the StatusTest wants iteration to stop.
static int GetNumberVecs(const MV &mv)
Obtain the number of vectors in mv.
virtual void print(std::ostream &os, int indent=0) const =0
Output formatted description of stopping test to output stream.
virtual int defineScaleForm(ScaleType TypeOfScaling, NormType TypeOfNorm, MagnitudeType ScaleValue=Teuchos::ScalarTraits< MagnitudeType >::one())=0
Define the form of the scaling for the residual.
MultiVecTraits< ScalarType, MV > MVT
Traits class which defines basic operations on multivectors.
static Teuchos::RCP< MV > Clone(const MV &mv, const int numvecs)
Creates a new empty MV containing numvecs columns.
virtual Teuchos::RCP< MV > getSolution()=0
Returns the current solution estimate that was computed for the most recent residual test...
virtual int setShowMaxResNormOnly(bool showMaxResNormOnly)=0
Set whether the only maximum residual norm is displayed when the print() method is called...
virtual bool getShowMaxResNormOnly()=0
Returns whether the only maximum residual norm is displayed when the print() method is called...
Class which describes the linear problem to be solved by the iterative solver.
virtual std::string description() const
SCT::magnitudeType MagnitudeType
virtual const std::vector< MagnitudeType > * getTestValue() const =0
Returns the test value, , computed in most recent call to CheckStatus.
Belos::StatusTest abstract class for specifying a residual norm stopping criteria.
Teuchos::ScalarTraits< ScalarType > SCT
static magnitudeType magnitude(ScalarTypea)
NormType
The type of vector norm to compute.
#define TEUCHOS_UNREACHABLE_RETURN(dummyReturnVal)
static void MvNorm(const MV &mv, std::vector< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > &normvec, NormType type=TwoNorm)
Compute the norm of each individual vector of mv. Upon return, normvec[i] holds the value of ...
virtual void reset()=0
Informs the convergence test that it should reset its internal configuration to the initialized state...
virtual bool getLOADetected() const =0
Returns a boolean indicating a loss of accuracy has been detected in computing the residual...
virtual StatusType getStatus() const =0
Return the result of the most recent CheckStatus call.
virtual int setTolerance(MagnitudeType tolerance)=0
Set the value of the tolerance.
virtual void printStatus(std::ostream &os, StatusType type) const
Output the result of the most recent CheckStatus call.