![]() |
waves
Basic FE playground
|
An implementation of StatusTestResNorm using a family of weighted residual norms. More...
#include <BelosStatusTestWeightedGenResNorm_MP_Vector.hpp>
Public Types | |
typedef Teuchos::ScalarTraits< ScalarType > | SCT |
typedef SCT::magnitudeType | MagnitudeType |
typedef MultiVecTraits< ScalarType, MV > | MVT |
Enums. | |
enum | ResType { Implicit, Explicit } |
Select how the residual std::vector is produced. More... | |
Public Member Functions | |
Constructors/destructors. | |
StatusTestWeightedGenResNorm (MagnitudeType Tolerance, Teuchos::RCP< MV > weights, int quorum=-1, bool showMaxResNormOnly=false) | |
Constructor. More... | |
virtual | ~StatusTestWeightedGenResNorm () |
Destructor. More... | |
Form and parameter definition methods. | |
int | defineResForm (ResType TypeOfResidual, NormType TypeOfNorm) |
Define form of the residual, its norm and optional weighting std::vector. More... | |
int | defineScaleForm (ScaleType TypeOfScaling, NormType TypeOfNorm, MagnitudeType ScaleValue=Teuchos::ScalarTraits< MagnitudeType >::one()) |
Define form of the scaling, its norm, its optional weighting std::vector, or, alternatively, define an explicit value. More... | |
NormType | getResNormType () |
int | setTolerance (MagnitudeType tolerance) |
Set the value of the tolerance. More... | |
int | setQuorum (int quorum) |
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 () |
int | getQuorum () const |
bool | getShowMaxResNormOnly () |
Returns whether the only maximum residual norm is displayed when the print() method is called. More... | |
std::vector< int > | convIndices () |
Returns the std::vector containing the indices of the residuals that passed the test. More... | |
MagnitudeType | getTolerance () const |
Returns the value of the tolerance, \( \tau \), set in the constructor. More... | |
const std::vector< MagnitudeType > * | getTestValue () const |
Returns the test value, \( \frac{\|r\|}{\sigma} \), computed in most recent call to CheckStatus. More... | |
const std::vector< MagnitudeType > * | getResNormValue () const |
Returns the residual norm value, \( \|r\| \), computed in most recent call to CheckStatus. More... | |
const std::vector< MagnitudeType > * | getScaledNormValue () const |
Returns the scaled norm value, \( \sigma \). More... | |
bool | getLOADetected () const |
const std::vector< int > | getEnsembleIterations () const |
Returns number of ensemble iterations. More... | |
Misc. | |
StatusType | firstCallCheckStatusSetup (Iteration< ScalarType, MV, OP > *iSolver) |
Call to setup initial scaling std::vector. More... | |
Overridden from Teuchos::Describable | |
std::string | description () const |
Method to return description of the maximum iteration status test More... | |
Private Member Functions | |
Private methods. | |
std::string | resFormStr () const |
Description of current residual form. More... | |
Private helper functions | |
void | MvWNorm (const MV &mv, const MV &w, std::vector< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > &normVec, NormType type=TwoNorm) |
Private Attributes | |
Private data members. | |
MagnitudeType | tolerance_ |
Tolerance used to determine convergence. More... | |
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 std::vector. More... | |
std::vector< MagnitudeType > | resvector_ |
Residual norm std::vector. More... | |
std::vector< MagnitudeType > | testvector_ |
Test std::vector = resvector_ / scalevector_. More... | |
std::vector< int > | ind_ |
Vector containing the indices for the vectors that passed the test. More... | |
Teuchos::RCP< MV > | curSoln_ |
Most recent solution vector used by this status 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... | |
Teuchos::RCP< MV > | weights_ |
Weights of the used norm: More... | |
std::vector< int > | ensemble_iterations |
The number of iterations at which point each ensemble component converges. More... | |
An implementation of StatusTestResNorm using a family of weighted residual norms.
StatusTestWeightedGenResNorm is an implementation of StatusTestResNorm that allows a user to construct one of a family of residual tests for use as a status/convergence test for Belos. The form of the test is
\[ \frac{\|w\odot r_i\|}{\sigma_i} \le \tau \]
where
typedef SCT::magnitudeType Belos::StatusTestWeightedGenResNorm< ScalarType, MV, OP >::MagnitudeType |
typedef MultiVecTraits<ScalarType, MV> Belos::StatusTestWeightedGenResNorm< ScalarType, MV, OP >::MVT |
typedef Teuchos::ScalarTraits<ScalarType> Belos::StatusTestWeightedGenResNorm< ScalarType, MV, OP >::SCT |
enum Belos::StatusTestWeightedGenResNorm::ResType |
Belos::StatusTestWeightedGenResNorm< ScalarType, MV, OP >::StatusTestWeightedGenResNorm | ( | MagnitudeType | Tolerance, |
Teuchos::RCP< MV > | weights, | ||
int | quorum = -1 , |
||
bool | showMaxResNormOnly = false |
||
) |
Constructor.
The constructor takes a single argument specifying the tolerance ( \(\tau\)). If none of the form definition methods are called, we use \(\|r\|_2/\|r^{(0)}\|_2 \le \tau\) as the stopping criterion, where \(\|r\|_2\) uses the least costly form of the 2-norm of residual available from the iterative method and \(\|r^{(0)}\|_2\) is the corresponding norm of the initial residual. The least costly form of the 2-norm depends on the chosen iterative method. Most Krylov methods produce the preconditioned residual std::vector in a form that would be exact in infinite precision arithmetic. This std::vector may be different from the true residual either because left preconditioning was used, or because round-off error has introduced significant error, or both.
You can also state the number of vectors that must pass the convergence criteria before the status test passes by using the quorum
argument.
|
virtual |
Destructor.
StatusType Belos::StatusTestWeightedGenResNorm< ScalarType, MV, OP >::checkStatus | ( | Iteration< ScalarType, MV, OP > * | iSolver | ) |
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.
|
inline |
Returns the std::vector containing the indices of the residuals that passed the test.
int Belos::StatusTestWeightedGenResNorm< ScalarType, MV, OP >::defineResForm | ( | ResType | TypeOfResidual, |
NormType | TypeOfNorm | ||
) |
Define form of the residual, its norm and optional weighting std::vector.
This method defines the form of \(\|r\|\). We specify:
int Belos::StatusTestWeightedGenResNorm< ScalarType, MV, OP >::defineScaleForm | ( | ScaleType | TypeOfScaling, |
NormType | TypeOfNorm, | ||
MagnitudeType | ScaleValue = Teuchos::ScalarTraits<MagnitudeType>::one() |
||
) |
Define form of the scaling, its norm, its optional weighting std::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:
|
inline |
Method to return description of the maximum iteration status test
StatusType Belos::StatusTestWeightedGenResNorm< ScalarType, MV, OP >::firstCallCheckStatusSetup | ( | Iteration< ScalarType, MV, OP > * | iSolver | ) |
Call to setup initial scaling std::vector.
After this function is called getScaledNormValue()
can be called to get the scaling std::vector.
|
inline |
Returns number of ensemble iterations.
|
inline |
Returns a boolean indicating a loss of accuracy has been detected in computing the residual.
|
inline |
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.
|
inline |
|
inline |
Returns the residual norm value, \( \|r\| \), computed in most recent call to CheckStatus.
|
inline |
Returns the scaled norm value, \( \sigma \).
|
inline |
Returns whether the only maximum residual norm is displayed when the print() method is called.
|
inline |
Returns the current solution estimate that was computed for the most recent residual test.
|
inline |
Return the result of the most recent CheckStatus call.
|
inline |
Returns the test value, \( \frac{\|r\|}{\sigma} \), computed in most recent call to CheckStatus.
|
inline |
Returns the value of the tolerance, \( \tau \), set in the constructor.
|
inlineprivate |
void Belos::StatusTestWeightedGenResNorm< ScalarType, MV, OP >::print | ( | std::ostream & | os, |
int | indent = 0 |
||
) | const |
Output formatted description of stopping test to output stream.
void Belos::StatusTestWeightedGenResNorm< ScalarType, MV, OP >::printStatus | ( | std::ostream & | os, |
StatusType | type | ||
) | const |
Print message for each status specific to this stopping test.
void Belos::StatusTestWeightedGenResNorm< ScalarType, MV, OP >::reset |
Resets the internal configuration to the initial state.
|
inlineprivate |
Description of current residual form.
|
inline |
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.
|
inline |
Set whether the only maximum residual norm is displayed when the print() method is called.
|
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.
|
private |
The current blocksize of the linear system being solved.
|
private |
The indices of the current number of right-hand sides being solved for.
|
private |
The current number of linear systems that have been loaded into the linear problem.
|
private |
The current number of right-hand sides being solved for.
|
private |
Most recent solution vector used by this status test.
|
private |
The number of iterations at which point each ensemble component converges.
|
private |
Is this the first time CheckStatus is called?
|
private |
Is this the first time DefineResForm is called?
|
private |
Is this the first time DefineScaleForm is called?
|
private |
Vector containing the indices for the vectors that passed the test.
|
private |
The total number of right-hand sides being solved for.
|
private |
Number of residuals that must pass the convergence test before Passed is returned.
|
private |
Type of norm to use on residual (OneNorm, TwoNorm, or InfNorm).
|
private |
Residual norm std::vector.
|
private |
Type of norm to use on the scaling (OneNorm, TwoNorm, or InfNorm)
|
private |
Type of scaling to use (Norm of RHS, Norm of Initial Residual, None or User provided)
|
private |
Scaling value.
|
private |
Scaling std::vector.
|
private |
Determines if the entries for all of the residuals are shown or just the max.
|
private |
Status.
|
private |
Test std::vector = resvector_ / scalevector_.
|
private |
Tolerance used to determine convergence.
|
private |
Weights of the used norm: