27#ifndef EWOMS_WEIGHTED_RESIDUAL_REDUCTION_CRITERION_HH
28#define EWOMS_WEIGHTED_RESIDUAL_REDUCTION_CRITERION_HH
57template <
class Vector,
class CollectiveCommunication>
60 using Scalar =
typename Vector::field_type;
61 using BlockType =
typename Vector::block_type;
109 return (residWeightVec_.size() == 0)
118 { residualReductionTolerance_ =
tol; }
124 {
return residualReductionTolerance_; }
130 { absResidualTolerance_ =
tol; }
136 {
return absResidualTolerance_; }
143 {
return residualError_/std::max<Scalar>(1
e-20, initialResidualError_); }
149 { fixPointTolerance_ =
tol; }
157 {
return fixPointTolerance_; }
164 {
return fixPointError_; }
171 lastResidualError_ = 1
e100;
176 fixPointError_ = 1
e100;
180 residualError_ = std::max<Scalar>(residualError_, 1
e-20);
181 initialResidualError_ = residualError_;
191 lastResidualError_ = residualError_;
204 residualError_ <= absResidualTolerance_;
213 (!
converged() && fixPointError_ <= fixPointTolerance_)
214 || residualError_ > maxError_;
223 {
return residualError_; }
230 os << std::setw(20) <<
" Iter ";
231 os << std::setw(20) <<
" Delta ";
232 os << std::setw(20) <<
" Residual ";
233 os << std::setw(20) <<
" ResidRed ";
234 os << std::setw(20) <<
" Rate ";
241 void print(Scalar
iter, std::ostream&
os = std::cout)
const override
243 static constexpr Scalar eps = std::numeric_limits<Scalar>::min()*1
e10;
245 os << std::setw(20) <<
iter <<
" ";
247 os << std::setw(20) << residualError_ <<
" ";
249 os << std::setw(20) << lastResidualError_ / std::max<Scalar>(residualError_, eps) <<
" ";
250 os << std::endl << std::flush;
257 residualError_ = 0.0;
258 fixPointError_ = 0.0;
259 for (
size_t i = 0; i <
curResid.size(); ++i) {
260 for (
unsigned j = 0; j < BlockType::dimension; ++j) {
262 std::max<Scalar>(residualError_,
265 std::max<Scalar>(fixPointError_,
266 std::abs(
curSol[i][j] - lastSolVec_[i][j])
267 /std::max<Scalar>(1.0,
curSol[i][j]));
272 residualError_ = comm_.max(residualError_);
273 fixPointError_ = comm_.max(fixPointError_);
276 const CollectiveCommunication& comm_;
279 Vector residWeightVec_;
286 Scalar fixPointError_;
290 Scalar fixPointTolerance_;
294 Scalar residualError_;
298 Scalar lastResidualError_;
302 Scalar initialResidualError_;
306 Scalar residualReductionTolerance_;
310 Scalar absResidualTolerance_;
Base class for all convergence criteria which only defines an virtual API.
Definition convergencecriterion.hh:56
Convergence criterion which looks at the weighted absolute value of the residual.
Definition weightedresidreductioncriterion.hh:59
Scalar residualReductionTolerance() const
Returns the tolerance of the residual reduction of the solution.
Definition weightedresidreductioncriterion.hh:123
Scalar accuracy() const override
Returns the accuracy of the solution at the last update.
Definition weightedresidreductioncriterion.hh:222
bool converged() const override
Returns true if and only if the convergence criterion is met.
Definition weightedresidreductioncriterion.hh:198
Scalar residualAccuracy() const
Returns the reduction of the weighted maximum of the residual compared to the initial solution.
Definition weightedresidreductioncriterion.hh:142
Scalar residualWeight(size_t outerIdx, unsigned innerIdx) const
Return the relative weight of a row of the residual.
Definition weightedresidreductioncriterion.hh:107
void setInitial(const Vector &curSol, const Vector &curResid) override
Set the initial solution of the linear system of equations.
Definition weightedresidreductioncriterion.hh:169
Scalar fixPointTolerance() const
Returns the maximum tolerated difference between two iterations to be met before a solution is consid...
Definition weightedresidreductioncriterion.hh:156
bool failed() const override
Returns true if the convergence criterion cannot be met anymore because the solver has broken down.
Definition weightedresidreductioncriterion.hh:210
void setResidualTolerance(Scalar tol)
Sets the maximum absolute tolerated residual.
Definition weightedresidreductioncriterion.hh:129
void setResidualWeight(const Vector &residWeightVec)
Sets the relative weight of each row of the residual.
Definition weightedresidreductioncriterion.hh:98
void setFixPointTolerance(Scalar tol)
Sets the fix-point tolerance.
Definition weightedresidreductioncriterion.hh:148
void update(const Vector &curSol, const Vector &, const Vector &curResid) override
Definition weightedresidreductioncriterion.hh:187
void printInitial(std::ostream &os=std::cout) const override
Prints the initial information about the convergence behaviour.
Definition weightedresidreductioncriterion.hh:228
Scalar fixPointAccuracy() const
Returns the weighted maximum of the difference between the last two iterative solutions.
Definition weightedresidreductioncriterion.hh:163
void print(Scalar iter, std::ostream &os=std::cout) const override
Prints the information about the convergence behaviour for the current iteration.
Definition weightedresidreductioncriterion.hh:241
void setResidualReductionTolerance(Scalar tol)
Sets the residual reduction tolerance.
Definition weightedresidreductioncriterion.hh:117
Scalar absResidualTolerance() const
Returns the maximum absolute tolerated residual.
Definition weightedresidreductioncriterion.hh:135
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242