27 #ifndef EWOMS_ISTL_WEIGHTED_RESIDUAL_REDUCTION_CRITERION_HH
28 #define EWOMS_ISTL_WEIGHTED_RESIDUAL_REDUCTION_CRITERION_HH
54 template <
class Vector,
class CollectiveCommunication>
57 typedef typename Vector::field_type Scalar;
58 typedef typename Vector::block_type BlockType;
68 Scalar maxResidual = 0.0)
70 residualReductionTolerance_(residualReductionTolerance),
72 maxResidual_(maxResidual)
79 { residualReductionTolerance_ = tol; }
85 {
return residualReductionTolerance_; }
92 {
return residualError_/std::max<Scalar>(1e-20, initialResidualError_); }
98 { absResidualTolerance_ = tol; }
105 {
return absResidualTolerance_; }
111 {
return residualError_; }
116 void setInitial(
const Vector& curSol,
const Vector& curResid)
override
118 updateErrors_(curSol, curSol, curResid);
122 residualError_ = std::max<Scalar>(residualError_,
123 std::numeric_limits<Scalar>::min()*1e10);
124 initialResidualError_ = residualError_;
125 lastResidualError_ = residualError_;
131 void update(
const Vector& curSol,
const Vector& changeIndicator,
const Vector& curResid)
override
132 { updateErrors_(curSol, changeIndicator, curResid); }
150 {
return !
converged() && (stagnates_ || residualError_ > maxResidual_); }
158 {
return residualError_/initialResidualError_; }
165 os << std::setw(20) <<
"iteration ";
166 os << std::setw(20) <<
"residual ";
167 os << std::setw(20) <<
"reduction ";
168 os << std::setw(20) <<
"rate ";
175 void print(Scalar iter, std::ostream& os = std::cout)
const override
177 const Scalar eps = std::numeric_limits<Scalar>::min()*1e10;
179 os << std::setw(20) << iter <<
" ";
181 os << std::setw(20) <<
accuracy() <<
" ";
182 os << std::setw(20) << lastResidualError_/std::max<Scalar>(residualError_, eps) <<
" ";
183 os << std::endl << std::flush;
188 void updateErrors_(
const Vector& curSol OPM_UNUSED,
const Vector& changeIndicator,
const Vector& curResid)
190 lastResidualError_ = residualError_;
191 residualError_ = 0.0;
193 for (
size_t i = 0; i < curResid.size(); ++i) {
194 for (
unsigned j = 0; j < BlockType::dimension; ++j) {
196 std::max<Scalar>(residualError_,
197 std::abs(curResid[i][j]));
199 if (stagnates_ && changeIndicator[i][j] != 0.0)
205 residualError_ = comm_.max(residualError_);
208 stagnates_ = comm_.min(stagnates_);
211 const CollectiveCommunication& comm_;
214 Scalar lastResidualError_;
217 Scalar residualError_;
220 Scalar initialResidualError_;
224 Scalar residualReductionTolerance_;
228 Scalar absResidualTolerance_;
Base class for all convergence criteria which only defines an virtual API.
Definition: convergencecriterion.hh:57
bool failed() const override
Returns true if the convergence criterion cannot be met anymore because the solver has broken down...
Definition: combinedcriterion.hh:149
void update(const Vector &curSol, const Vector &changeIndicator, const Vector &curResid) override
Update the internal members of the convergence criterion with the current solution.
Definition: combinedcriterion.hh:131
void setInitial(const Vector &curSol, const Vector &curResid) override
Set the initial solution of the linear system of equations.
Definition: combinedcriterion.hh:116
Scalar absResidual() const
Returns the infinity norm of the absolute residual.
Definition: combinedcriterion.hh:110
Convergence criterion which looks at the absolute value of the residual and fails if the linear solve...
Definition: combinedcriterion.hh:55
Scalar accuracy() const override
Returns the accuracy of the solution at the last update.
Definition: combinedcriterion.hh:157
void setAbsResidualTolerance(Scalar tol)
Sets the maximum absolute tolerated residual.
Definition: combinedcriterion.hh:97
void setResidualReductionTolerance(Scalar tol)
Sets the residual reduction tolerance.
Definition: combinedcriterion.hh:78
Scalar residualReductionTolerance() const
Returns the tolerance of the residual reduction of the solution.
Definition: combinedcriterion.hh:84
Scalar residualReduction() const
Returns the reduction of the maximum of the residual compared to the initial solution.
Definition: combinedcriterion.hh:91
Scalar absResidualTolerance() const
Returns the tolerated maximum of the the infinity norm of the absolute residual.
Definition: combinedcriterion.hh:104
void print(Scalar iter, std::ostream &os=std::cout) const override
Prints the information about the convergence behaviour for the current iteration. ...
Definition: combinedcriterion.hh:175
bool converged() const override
Returns true if and only if the convergence criterion is met.
Definition: combinedcriterion.hh:137
void printInitial(std::ostream &os=std::cout) const override
Prints the initial information about the convergence behaviour.
Definition: combinedcriterion.hh:163