Belos Package Browser (Single Doxygen Collection)
Development
|
Convergence test using the implicit residual norm(s), with an explicit residual norm(s) check for loss of accuracy if necessary. More...
#include <BelosStatusTestImpResNorm.hpp>
Public Types | |
typedef Teuchos::ScalarTraits < ScalarType >::magnitudeType | MagnitudeType |
The type of the magnitude (absolute value) of a ScalarType. More... | |
Public Types inherited from Belos::StatusTestResNorm< ScalarType, MV, OP > | |
typedef Teuchos::ScalarTraits < ScalarType > | SCT |
typedef SCT::magnitudeType | MagnitudeType |
typedef MultiVecTraits < ScalarType, MV > | MVT |
Abbreviations for method implementations | |
typedef Teuchos::ScalarTraits < ScalarType > | STS |
typedef Teuchos::ScalarTraits < MagnitudeType > | STM |
typedef MultiVecTraits < ScalarType, MV > | MVT |
Constructors and destructor. | |
StatusTestImpResNorm (MagnitudeType Tolerance, int quorum=-1, bool showMaxResNormOnly=false) | |
Constructor. More... | |
virtual | ~StatusTestImpResNorm () |
Destructor (virtual for memory safety). More... | |
Form and parameter definition methods. | |
int | defineResForm (NormType TypeOfNorm) |
Define form of the residual, its norm and optional weighting vector. More... | |
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. More... | |
int | setTolerance (MagnitudeType tolerance) |
Set the value of the tolerance. More... | |
int | setQuorum (int quorum) |
Sets the number of residuals that must pass the convergence test before Passed is returned. More... | |
int | setShowMaxResNormOnly (bool showMaxResNormOnly) |
Set whether the only maximum residual norm is displayed when the print() method is called. More... | |
Status methods | |
StatusType | checkStatus (Iteration< ScalarType, MV, OP > *iSolver) |
Check convergence status: Passed, Failed, or Undefined. More... | |
StatusType | getStatus () const |
Return the result of the most recent CheckStatus call. More... | |
Reset methods | |
void | reset () |
Resets the internal configuration to the initial state. More... | |
Print methods | |
void | print (std::ostream &os, int indent=0) const |
Output formatted description of stopping test to output stream. More... | |
void | printStatus (std::ostream &os, StatusType type) const |
Print message for each status specific to this stopping test. More... | |
Methods to access data members. | |
Teuchos::RCP< MV > | getSolution () |
Returns the current solution estimate that was computed for the most recent residual test. More... | |
int | getQuorum () const |
Returns the number of residuals that must pass the convergence test before Passed is returned. More... | |
bool | getShowMaxResNormOnly () |
Returns whether the only maximum residual norm is displayed when the print() method is called. More... | |
std::vector< int > | convIndices () |
Returns the vector containing the indices of the residuals that passed the test. More... | |
MagnitudeType | getTolerance () const |
"Original" convergence tolerance as set by user. More... | |
MagnitudeType | getCurrTolerance () const |
Current convergence tolerance; may be changed to prevent loss of accuracy. More... | |
const std::vector < MagnitudeType > * | getTestValue () const |
Returns the test value, , computed in most recent call to CheckStatus. More... | |
const std::vector < MagnitudeType > * | getResNormValue () const |
Returns the residual norm value, , computed in most recent call to CheckStatus. More... | |
const std::vector < MagnitudeType > * | getScaledNormValue () const |
Returns the scaled norm value, . More... | |
bool | getLOADetected () const |
Returns a boolean indicating a loss of accuracy has been detected in computing the residual. More... | |
Misc. | |
StatusType | firstCallCheckStatusSetup (Iteration< ScalarType, MV, OP > *iSolver) |
Call to setup initial scaling vector. More... | |
Overridden from Teuchos::Describable | |
std::string | description () const |
Method to return description of the maximum iteration status test. More... | |
Private methods. | |
std::string | resFormStr () const |
Description of current residual form. More... | |
Private data members. | |
MagnitudeType | tolerance_ |
Current tolerance used to determine convergence and the default tolerance set by user. More... | |
MagnitudeType | currTolerance_ |
int | quorum_ |
Number of residuals that must pass the convergence test before Passed is returned. More... | |
bool | showMaxResNormOnly_ |
Determines if the entries for all of the residuals are shown or just the max. More... | |
NormType | resnormtype_ |
Type of norm to use on residual (OneNorm, TwoNorm, or InfNorm). More... | |
ScaleType | scaletype_ |
Type of scaling to use (Norm of RHS, Norm of Initial Residual, None or User provided) More... | |
NormType | scalenormtype_ |
Type of norm to use on the scaling (OneNorm, TwoNorm, or InfNorm) More... | |
MagnitudeType | scalevalue_ |
Scaling value. More... | |
std::vector< MagnitudeType > | scalevector_ |
Scaling vector. More... | |
std::vector< MagnitudeType > | resvector_ |
Residual norm vector. More... | |
std::vector< MagnitudeType > | testvector_ |
Test vector = resvector_ / scalevector_. More... | |
Teuchos::RCP< MV > | curSoln_ |
Current solution vector. More... | |
std::vector< int > | ind_ |
Vector containing the indices for the vectors that passed the test. More... | |
StatusType | status_ |
Status. More... | |
int | curBlksz_ |
The current blocksize of the linear system being solved. More... | |
int | curNumRHS_ |
The current number of right-hand sides being solved for. More... | |
std::vector< int > | curLSIdx_ |
The indices of the current number of right-hand sides being solved for. More... | |
int | curLSNum_ |
The current number of linear systems that have been loaded into the linear problem. More... | |
int | numrhs_ |
The total number of right-hand sides being solved for. More... | |
bool | firstcallCheckStatus_ |
Is this the first time CheckStatus is called? More... | |
bool | firstcallDefineResForm_ |
Is this the first time DefineResForm is called? More... | |
bool | firstcallDefineScaleForm_ |
Is this the first time DefineScaleForm is called? More... | |
bool | lossDetected_ |
Has a loss in accuracy been detected? More... | |
Additional Inherited Members | |
Public Member Functions inherited from Belos::StatusTestResNorm< ScalarType, MV, OP > | |
virtual int | setTolerance (MagnitudeType tolerance)=0 |
Set the value of the tolerance. More... | |
virtual int | defineScaleForm (ScaleType TypeOfScaling, NormType TypeOfNorm, MagnitudeType ScaleValue=Teuchos::ScalarTraits< MagnitudeType >::one())=0 |
Define the form of the scaling for the residual. More... | |
Public Member Functions inherited from Belos::StatusTest< ScalarType, MV, OP > | |
StatusTest () | |
Constructor. More... | |
virtual | ~StatusTest () |
Destructor. More... | |
Public Member Functions inherited from Teuchos::Describable | |
DescribableStreamManipulatorState | describe (const Describable &describable, const EVerbosityLevel verbLevel=Describable::verbLevel_default) |
std::ostream & | operator<< (std::ostream &os, const DescribableStreamManipulatorState &d) |
virtual void | describe (FancyOStream &out, const EVerbosityLevel verbLevel=verbLevel_default) const |
void | describe (std::ostream &out, const EVerbosityLevel verbLevel=verbLevel_default) const |
virtual | ~Describable () |
LabeledObject () | |
virtual | ~LabeledObject () |
virtual void | setObjectLabel (const std::string &objectLabel) |
virtual std::string | getObjectLabel () const |
Static Public Attributes inherited from Teuchos::Describable | |
static const EVerbosityLevel | verbLevel_default |
Convergence test using the implicit residual norm(s), with an explicit residual norm(s) check for loss of accuracy if necessary.
This test passes when a quorum of residual vectors (for all right-hand sides) passes a residual norm test. For residual vector , the form of the test is , where is the convergence tolerance set in the constructor or by setTolerance(). The default quorum consists of all the vectors, but you can set the number of vectors that must pass before the whole test passes, either by the constructor argument or by calling setQuorum().
The default residual vector norm (the norm used in the numerator of the test, and applied to the current residual vectors) is the 2-norm, but you can change the norm type (1-norm, 2-norm, or infinity norm) using defineResForm(). The default scaling factor for residual vector i is the 2-norm of the initial residual vector, but you can change this using defineScaleForm(). The norm used for the scaling factors may be different than the norm used for the current residual vectors. Furthermore, each residual vector may use a different scaling factor .
This test starts by using the "implicit" residual norm, otherwise called the "recursive" or "native" residual norm. It comes from the iterative method's projected problem, unlike the "explicit" residual norm, which comes from explicitly computing . Once the implicit residual norm reaches the convergence tolerance, this test then checks the explicit residual norm to make sure that it has reached the convergence tolerance as well.
We say that there is a "potential loss of accuracy" when we first detect that the implicit residual norm(s) have met the desired original convergence tolerance, but the explicit residual norm(s) have not. We don't actually call this a "loss of accuracy" (in the sense of getLOADetected() returning true) unless we tried to remedy this situation, but couldn't fix it. Upon detecting a potential loss of accuracy, this test tells the solver to "iterate a few more steps" by making the "current" residual tolerance (the value of getCurrTolerance()) smaller, and forcing the implicit residual norm(s) to meet the current tolerance before moving on to test the explicit norm(s). If that doesn't work, this test declares a loss of accuracy.
Definition at line 72 of file BelosStatusTestImpResNorm.hpp.
typedef Teuchos::ScalarTraits<ScalarType>::magnitudeType Belos::StatusTestImpResNorm< ScalarType, MV, OP >::MagnitudeType |
The type of the magnitude (absolute value) of a ScalarType.
Definition at line 75 of file BelosStatusTestImpResNorm.hpp.
|
private |
Definition at line 80 of file BelosStatusTestImpResNorm.hpp.
|
private |
Definition at line 81 of file BelosStatusTestImpResNorm.hpp.
|
private |
Definition at line 82 of file BelosStatusTestImpResNorm.hpp.
Belos::StatusTestImpResNorm< ScalarType, MV, OP >::StatusTestImpResNorm | ( | MagnitudeType | Tolerance, |
int | quorum = -1 , |
||
bool | showMaxResNormOnly = false |
||
) |
Constructor.
Tolerance | [in] Convergence tolerance . |
quorum | [in] The number of vectors in the problem that must pass the convergence criteria before this StatusTest passes. -1 means that all residuals must pass. |
showMaxResNormOnly | [in] Whether only the maximum residual norm (of all vectors) is displayed when the print() method is called. |
Definition at line 389 of file BelosStatusTestImpResNorm.hpp.
|
virtual |
Destructor (virtual for memory safety).
Definition at line 413 of file BelosStatusTestImpResNorm.hpp.
int Belos::StatusTestImpResNorm< ScalarType, MV, OP >::defineResForm | ( | NormType | TypeOfNorm | ) |
Define form of the residual, its norm and optional weighting vector.
This method defines the form of . We specify:
Definition at line 432 of file BelosStatusTestImpResNorm.hpp.
int Belos::StatusTestImpResNorm< ScalarType, MV, OP >::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.
This method defines the form of how the residual is scaled (if at all). It operates in two modes:
User-provided scaling value:
Definition at line 444 of file BelosStatusTestImpResNorm.hpp.
|
inline |
Set the value of the tolerance.
We allow the tolerance to be reset for cases where, in the process of testing the residual, we find that the initial tolerance was too tight or too lax.
Definition at line 147 of file BelosStatusTestImpResNorm.hpp.
|
inlinevirtual |
Sets the number of residuals that must pass the convergence test before Passed is returned.
quorum=-1
then all residuals must pass the convergence test before Passed is returned. Implements Belos::StatusTestResNorm< ScalarType, MV, OP >.
Definition at line 154 of file BelosStatusTestImpResNorm.hpp.
|
inlinevirtual |
Set whether the only maximum residual norm is displayed when the print() method is called.
Implements Belos::StatusTestResNorm< ScalarType, MV, OP >.
Definition at line 160 of file BelosStatusTestImpResNorm.hpp.
|
virtual |
Check convergence status: Passed, Failed, or Undefined.
This method checks to see if the convergence criteria are met. Depending on how the residual test is constructed this method will return the appropriate status type.
Implements Belos::StatusTest< ScalarType, MV, OP >.
Definition at line 460 of file BelosStatusTestImpResNorm.hpp.
|
inlinevirtual |
Return the result of the most recent CheckStatus call.
Implements Belos::StatusTest< ScalarType, MV, OP >.
Definition at line 179 of file BelosStatusTestImpResNorm.hpp.
|
virtual |
Resets the internal configuration to the initial state.
Implements Belos::StatusTest< ScalarType, MV, OP >.
Definition at line 417 of file BelosStatusTestImpResNorm.hpp.
|
virtual |
Output formatted description of stopping test to output stream.
Implements Belos::StatusTest< ScalarType, MV, OP >.
Definition at line 720 of file BelosStatusTestImpResNorm.hpp.
|
virtual |
Print message for each status specific to this stopping test.
Reimplemented from Belos::StatusTest< ScalarType, MV, OP >.
Definition at line 752 of file BelosStatusTestImpResNorm.hpp.
|
inlinevirtual |
Returns the current solution estimate that was computed for the most recent residual test.
Implements Belos::StatusTestResNorm< ScalarType, MV, OP >.
Definition at line 204 of file BelosStatusTestImpResNorm.hpp.
|
inlinevirtual |
Returns the number of residuals that must pass the convergence test before Passed is returned.
quorum=-1
then all residuals must pass the convergence test before Passed is returned. Implements Belos::StatusTestResNorm< ScalarType, MV, OP >.
Definition at line 208 of file BelosStatusTestImpResNorm.hpp.
|
inlinevirtual |
Returns whether the only maximum residual norm is displayed when the print() method is called.
Implements Belos::StatusTestResNorm< ScalarType, MV, OP >.
Definition at line 211 of file BelosStatusTestImpResNorm.hpp.
|
inlinevirtual |
Returns the vector containing the indices of the residuals that passed the test.
Implements Belos::StatusTestResNorm< ScalarType, MV, OP >.
Definition at line 214 of file BelosStatusTestImpResNorm.hpp.
|
inlinevirtual |
"Original" convergence tolerance as set by user.
This value is the convergence tolerance as set by the user in the constructor, or by the setTolerance() method. See this class' main documentation and the documentation of getCurrTolerance() for an explanation of the difference between the "original" and "current" tolerances.
Implements Belos::StatusTestResNorm< ScalarType, MV, OP >.
Definition at line 223 of file BelosStatusTestImpResNorm.hpp.
|
inline |
Current convergence tolerance; may be changed to prevent loss of accuracy.
The difference between "original" tolerance (the value of getTolerance()) and "current" tolerance (the value of this method) relates to the idea of "loss of accuracy." See this class' main documentation for details. We do not allow users to set the "current" tolerance.
Definition at line 234 of file BelosStatusTestImpResNorm.hpp.
|
inlinevirtual |
Returns the test value, , computed in most recent call to CheckStatus.
Implements Belos::StatusTestResNorm< ScalarType, MV, OP >.
Definition at line 239 of file BelosStatusTestImpResNorm.hpp.
|
inline |
Returns the residual norm value, , computed in most recent call to CheckStatus.
Definition at line 242 of file BelosStatusTestImpResNorm.hpp.
|
inline |
Returns the scaled norm value, .
Definition at line 245 of file BelosStatusTestImpResNorm.hpp.
|
inlinevirtual |
Returns a boolean indicating a loss of accuracy has been detected in computing the residual.
Implements Belos::StatusTestResNorm< ScalarType, MV, OP >.
Definition at line 248 of file BelosStatusTestImpResNorm.hpp.
StatusType Belos::StatusTestImpResNorm< ScalarType, MV, OP >::firstCallCheckStatusSetup | ( | Iteration< ScalarType, MV, OP > * | iSolver | ) |
Call to setup initial scaling vector.
After this function is called getScaledNormValue()
can be called to get the scaling vector.
Definition at line 776 of file BelosStatusTestImpResNorm.hpp.
|
inlinevirtual |
Method to return description of the maximum iteration status test.
Reimplemented from Teuchos::Describable.
Definition at line 267 of file BelosStatusTestImpResNorm.hpp.
|
inlineprivate |
Description of current residual form.
Definition at line 283 of file BelosStatusTestImpResNorm.hpp.
|
private |
Current tolerance used to determine convergence and the default tolerance set by user.
Definition at line 319 of file BelosStatusTestImpResNorm.hpp.
|
private |
Definition at line 319 of file BelosStatusTestImpResNorm.hpp.
|
private |
Number of residuals that must pass the convergence test before Passed is returned.
Definition at line 322 of file BelosStatusTestImpResNorm.hpp.
|
private |
Determines if the entries for all of the residuals are shown or just the max.
Definition at line 325 of file BelosStatusTestImpResNorm.hpp.
|
private |
Type of norm to use on residual (OneNorm, TwoNorm, or InfNorm).
Definition at line 328 of file BelosStatusTestImpResNorm.hpp.
|
private |
Type of scaling to use (Norm of RHS, Norm of Initial Residual, None or User provided)
Definition at line 331 of file BelosStatusTestImpResNorm.hpp.
|
private |
Type of norm to use on the scaling (OneNorm, TwoNorm, or InfNorm)
Definition at line 334 of file BelosStatusTestImpResNorm.hpp.
|
private |
Scaling value.
Definition at line 337 of file BelosStatusTestImpResNorm.hpp.
|
private |
Scaling vector.
Definition at line 340 of file BelosStatusTestImpResNorm.hpp.
|
private |
Residual norm vector.
Definition at line 343 of file BelosStatusTestImpResNorm.hpp.
|
private |
Test vector = resvector_ / scalevector_.
Definition at line 346 of file BelosStatusTestImpResNorm.hpp.
|
private |
Current solution vector.
Definition at line 349 of file BelosStatusTestImpResNorm.hpp.
|
private |
Vector containing the indices for the vectors that passed the test.
Definition at line 352 of file BelosStatusTestImpResNorm.hpp.
|
private |
Status.
Definition at line 355 of file BelosStatusTestImpResNorm.hpp.
|
private |
The current blocksize of the linear system being solved.
Definition at line 358 of file BelosStatusTestImpResNorm.hpp.
|
private |
The current number of right-hand sides being solved for.
Definition at line 361 of file BelosStatusTestImpResNorm.hpp.
|
private |
The indices of the current number of right-hand sides being solved for.
Definition at line 364 of file BelosStatusTestImpResNorm.hpp.
|
private |
The current number of linear systems that have been loaded into the linear problem.
Definition at line 367 of file BelosStatusTestImpResNorm.hpp.
|
private |
The total number of right-hand sides being solved for.
Definition at line 370 of file BelosStatusTestImpResNorm.hpp.
|
private |
Is this the first time CheckStatus is called?
Definition at line 373 of file BelosStatusTestImpResNorm.hpp.
|
private |
Is this the first time DefineResForm is called?
Definition at line 376 of file BelosStatusTestImpResNorm.hpp.
|
private |
Is this the first time DefineScaleForm is called?
Definition at line 379 of file BelosStatusTestImpResNorm.hpp.
|
private |
Has a loss in accuracy been detected?
Definition at line 382 of file BelosStatusTestImpResNorm.hpp.