45#ifndef THYRA_GENERAL_SOLVE_CRITERIA_BELOS_STATUS_TEST_DEF_HPP
46#define THYRA_GENERAL_SOLVE_CRITERIA_BELOS_STATUS_TEST_DEF_HPP
48#include "Thyra_GeneralSolveCriteriaBelosStatusTest.hpp"
49#include "Thyra_VectorSpaceBase.hpp"
50#include "Thyra_MultiVectorBase.hpp"
51#include "Thyra_VectorBase.hpp"
52#include "Thyra_MultiVectorStdOps.hpp"
53#include "Thyra_VectorStdOps.hpp"
55#include "Teuchos_DebugDefaultAsserts.hpp"
56#include "Teuchos_VerbosityLevel.hpp"
57#include "Teuchos_as.hpp"
68 :convergenceTestFrequency_(-1),
80 const SolveCriteria<Scalar> &solveCriteria,
81 const int convergenceTestFrequency
86 typedef ScalarTraits<ScalarMag> SMT;
90 TEUCHOS_ASSERT_INEQUALITY(solveCriteria.requestedTol, >=, SMT::zero());
91 TEUCHOS_ASSERT_INEQUALITY(solveCriteria.solveMeasureType.numerator, !=, SOLVE_MEASURE_ONE);
92 TEUCHOS_ASSERT(nonnull(solveCriteria.numeratorReductionFunc) ||
93 nonnull(solveCriteria.denominatorReductionFunc) );
97 solveCriteria_ = solveCriteria;
98 convergenceTestFrequency_ = convergenceTestFrequency;
102 compute_r_ = solveCriteria.solveMeasureType.contains(SOLVE_MEASURE_NORM_RESIDUAL);
104 compute_x_ = (compute_r_ ||
105 solveCriteria.solveMeasureType.contains(SOLVE_MEASURE_NORM_SOLUTION));
110template<
class Scalar>
111ArrayView<const typename ScalarTraits<Scalar>::magnitudeType>
114 return lastAchievedTol_;
121template <
class Scalar>
131 TEUCHOS_ASSERT(iSolver);
136 if (currIter == 0 || currIter % convergenceTestFrequency_ != 0) {
141 const RCP<FancyOStream> out = this->getOStream();
142 const Teuchos::EVerbosityLevel verbLevel = this->getVerbLevel();
145 const int numRhs = lp.
getRHS()->domain()->dim();
148 if (solveCriteria_.solveMeasureType.contains(SOLVE_MEASURE_NORM_RHS)) {
149 TEUCHOS_TEST_FOR_EXCEPT_MSG(
true,
"ToDo: Handle ||b||");
153 if (solveCriteria_.solveMeasureType.contains(SOLVE_MEASURE_NORM_INIT_RESIDUAL)) {
154 TEUCHOS_TEST_FOR_EXCEPT_MSG(
true,
"ToDo: Handle ||r0||");
158 RCP<MultiVectorBase<Scalar> > X;
165 RCP<MultiVectorBase<Scalar> > R;
167 R = createMembers(lp.
getOperator()->range(), X->domain());
173 lastNumerator_.resize(numRhs);
174 lastDenominator_.resize(numRhs);
176 for (
int j = 0; j < numRhs; ++j) {
177 const RCP<const VectorBase<Scalar> > x_j = (nonnull(X) ? X->col(j) : null);
178 const RCP<const VectorBase<Scalar> > r_j = (nonnull(R) ? R->col(j) : null);
179 lastNumerator_[j] = computeReductionFunctional(
180 solveCriteria_.solveMeasureType.numerator,
181 solveCriteria_.numeratorReductionFunc.ptr(),
182 x_j.ptr(), r_j.ptr() );
183 lastDenominator_[j] = computeReductionFunctional(
184 solveCriteria_.solveMeasureType.denominator,
185 solveCriteria_.denominatorReductionFunc.ptr(),
186 x_j.ptr(), r_j.ptr() );
191 bool systemsAreConverged =
true;
192 lastAchievedTol_.resize(numRhs);
194 for (
int j = 0; j < numRhs; ++j) {
195 const ScalarMag convRatio = lastNumerator_[j] / lastDenominator_[j];
196 lastAchievedTol_[j] = convRatio;
197 const bool sys_converged_j = (convRatio <= solveCriteria_.requestedTol);
198 if (includesVerbLevel(verbLevel, Teuchos::VERB_MEDIUM)) {
199 printRhsStatus(currIter, j, *out);
201 if (!sys_converged_j) {
202 systemsAreConverged =
false;
207 lastCurrIter_ = currIter;
209 return lastRtnStatus_;
213template <
class Scalar>
217 return lastRtnStatus_;
221template <
class Scalar>
226 lastNumerator_.clear();
227 lastDenominator_.clear();
228 lastAchievedTol_.clear();
234template <
class Scalar>
236 std::ostream& os,
int indent
239 const int numRhs = lastNumerator_.size();
240 for (
int j = 0; j < numRhs; ++j) {
241 printRhsStatus(lastCurrIter_, j, os, indent);
249template <
class Scalar>
252 ESolveMeasureNormType measureType,
253 const Ptr<
const ReductionFunctional<Scalar> > &reductFunc,
254 const Ptr<
const VectorBase<Scalar> > &x,
255 const Ptr<
const VectorBase<Scalar> > &r
258 typedef ScalarTraits<ScalarMag> SMT;
259 ScalarMag rtn = -SMT::one();
260 Ptr<const VectorBase<Scalar> > v;
261 switch(measureType) {
262 case SOLVE_MEASURE_ONE:
265 case SOLVE_MEASURE_NORM_RESIDUAL:
268 case SOLVE_MEASURE_NORM_SOLUTION:
271 case SOLVE_MEASURE_NORM_INIT_RESIDUAL:
272 TEUCHOS_TEST_FOR_EXCEPT_MSG(
true,
"ToDo: Handle ||r0||!)")
273 case SOLVE_MEASURE_NORM_RHS:
274 TEUCHOS_TEST_FOR_EXCEPT_MSG(true, "ToDo: Handle ||b||!)");
275 TEUCHOS_SWITCH_DEFAULT_DEBUG_ASSERT();
277 if (rtn >= SMT::zero()) {
280 else if (nonnull(v) && rtn < SMT::zero()) {
281 if (nonnull(reductFunc)) {
282 rtn = reductFunc->reduce(*v);
288 TEUCHOS_IF_ELSE_DEBUG_ASSERT();
293template <
class Scalar>
295GeneralSolveCriteriaBelosStatusTest<Scalar>::printRhsStatus(
296 const int currIter,
const int j, std::ostream &out,
300 const ScalarMag convRatio = lastNumerator_[j] / lastDenominator_[j];
301 const bool sys_converged_j = (convRatio <= solveCriteria_.requestedTol);
302 for (
int i = 0; i < indent; ++i) { out <<
" "; }
304 <<
"["<<currIter<<
"] "
305 <<
"gN(vN("<<j<<
"))/gD(vD("<<j<<
")) = "
306 << lastNumerator_[j] <<
"/" << lastDenominator_[j] <<
" = "
307 << convRatio <<
" <= " << solveCriteria_.requestedTol <<
" : "
308 << (sys_converged_j ?
" true" :
"false")
Thyra specializations of MultiVecTraits and OperatorTraits.
virtual int getNumIters() const=0
virtual Teuchos::RCP< MV > getCurrentUpdate() const=0
virtual const LinearProblem< ScalarType, MV, OP > & getProblem() const=0
Teuchos::RCP< const OP > getOperator() const
Teuchos::RCP< const MV > getRHS() const
virtual void computeCurrResVec(MV *R, const MV *X=0, const MV *B=0) const
virtual Teuchos::RCP< MV > updateSolution(const Teuchos::RCP< MV > &update=Teuchos::null, bool updateLP=false, ScalarType scale=Teuchos::ScalarTraits< ScalarType >::one())
Subclass of Belos::StatusTest that implements every possible form of SolveCriteria that exists by for...
GeneralSolveCriteriaBelosStatusTest()
void setSolveCriteria(const SolveCriteria< Scalar > &solveCriteria, const int convergenceTestFrequency)
virtual Belos::StatusType getStatus() const
virtual Belos::StatusType checkStatus(Belos::Iteration< Scalar, MV, OP > *iSolver)
ScalarTraits< Scalar >::magnitudeType ScalarMag
ArrayView< const ScalarMag > achievedTol() const
virtual void print(std::ostream &os, int indent) const