27 #ifndef EWOMS_ISTL_WEIGHTED_RESIDUAL_REDUCTION_CRITERION_HH 28 #define EWOMS_ISTL_WEIGHTED_RESIDUAL_REDUCTION_CRITERION_HH 57 template <
class Vector,
class CollectiveCommunication>
60 typedef typename Vector::field_type Scalar;
61 typedef typename Vector::block_type BlockType;
69 const Vector& residWeights,
73 Scalar maxError = 0.0)
75 residWeightVec_(residWeights),
99 { residWeightVec_ = residWeightVec; }
109 return (residWeightVec_.size() == 0)
111 : residWeightVec_[outerIdx][innerIdx];
118 { residualReductionTolerance_ = tol; }
124 {
return residualReductionTolerance_; }
130 { absResidualTolerance_ = tol; }
136 {
return absResidualTolerance_; }
143 {
return residualError_/std::max<Scalar>(1e-20, initialResidualError_); }
149 { fixPointTolerance_ = tol; }
157 {
return fixPointTolerance_; }
164 {
return fixPointError_; }
169 void setInitial(
const Vector& curSol,
const Vector& curResid)
171 lastResidualError_ = 1e100;
173 lastSolVec_ = curSol;
174 updateErrors_(curSol, curResid);
176 fixPointError_ = 1e100;
180 residualError_ = std::max<Scalar>(residualError_, 1e-20);
181 initialResidualError_ = residualError_;
187 void update(
const Vector& curSol,
const Vector& curResid)
189 lastResidualError_ = residualError_;
190 updateErrors_(curSol, curResid);
202 residualError_ <= absResidualTolerance_;
211 (!
converged() && fixPointError_ <= fixPointTolerance_)
212 || residualError_ > maxError_;
221 {
return residualError_; }
228 os << std::setw(20) <<
" Iter ";
229 os << std::setw(20) <<
" Delta ";
230 os << std::setw(20) <<
" Residual ";
231 os << std::setw(20) <<
" ResidRed ";
232 os << std::setw(20) <<
" Rate ";
239 void print(Scalar iter, std::ostream& os = std::cout)
const 241 static constexpr Scalar eps = std::numeric_limits<Scalar>::min()*1e10;
243 os << std::setw(20) << iter <<
" ";
245 os << std::setw(20) << residualError_ <<
" ";
247 os << std::setw(20) << lastResidualError_ / std::max<Scalar>(residualError_, eps) <<
" ";
248 os << std::endl << std::flush;
253 void updateErrors_(
const Vector& curSol,
const Vector& curResid)
255 residualError_ = 0.0;
256 fixPointError_ = 0.0;
257 for (
size_t i = 0; i < curResid.size(); ++i) {
258 for (
unsigned j = 0; j < BlockType::dimension; ++j) {
260 std::max<Scalar>(residualError_,
263 std::max<Scalar>(fixPointError_,
264 std::abs(curSol[i][j] - lastSolVec_[i][j])
265 /std::max<Scalar>(1.0, curSol[i][j]));
268 lastSolVec_ = curSol;
270 residualError_ = comm_.max(residualError_);
271 fixPointError_ = comm_.max(fixPointError_);
274 const CollectiveCommunication& comm_;
277 Vector residWeightVec_;
284 Scalar fixPointError_;
288 Scalar fixPointTolerance_;
292 Scalar residualError_;
296 Scalar lastResidualError_;
300 Scalar initialResidualError_;
304 Scalar residualReductionTolerance_;
308 Scalar absResidualTolerance_;
Base class for all convergence criteria which only defines an virtual API.
Definition: convergencecriterion.hh:57
bool converged() const
Returns true if and only if the convergence criterion is met.
Definition: weightedresidreductioncriterion.hh:196
Scalar fixPointAccuracy() const
Returns the weighted maximum of the difference between the last two iterative solutions.
Definition: weightedresidreductioncriterion.hh:163
Definition: baseauxiliarymodule.hh:37
Scalar residualReductionTolerance() const
Returns the tolerance of the residual reduction of the solution.
Definition: weightedresidreductioncriterion.hh:123
void setFixPointTolerance(Scalar tol)
Sets the fix-point tolerance.
Definition: weightedresidreductioncriterion.hh:148
void setResidualReductionTolerance(Scalar tol)
Sets the residual reduction tolerance.
Definition: weightedresidreductioncriterion.hh:117
void print(Scalar iter, std::ostream &os=std::cout) const
Prints the information about the convergence behaviour for the current iteration. ...
Definition: weightedresidreductioncriterion.hh:239
Scalar residualAccuracy() const
Returns the reduction of the weighted maximum of the residual compared to the initial solution...
Definition: weightedresidreductioncriterion.hh:142
Scalar fixPointTolerance() const
Returns the maximum tolerated difference between two iterations to be met before a solution is consid...
Definition: weightedresidreductioncriterion.hh:156
void setInitial(const Vector &curSol, const Vector &curResid)
Set the initial solution of the linear system of equations.
Definition: weightedresidreductioncriterion.hh:169
Scalar accuracy() const
Returns the accuracy of the solution at the last update.
Definition: weightedresidreductioncriterion.hh:220
void setResidualWeight(const Vector &residWeightVec)
Sets the relative weight of each row of the residual.
Definition: weightedresidreductioncriterion.hh:98
void printInitial(std::ostream &os=std::cout) const
Prints the initial information about the convergence behaviour.
Definition: weightedresidreductioncriterion.hh:226
Scalar absResidualTolerance() const
Returns the maximum absolute tolerated residual.
Definition: weightedresidreductioncriterion.hh:135
bool failed() const
Returns true if the convergence criterion cannot be met anymore because the solver has broken down...
Definition: weightedresidreductioncriterion.hh:208
Convergence criterion which looks at the weighted absolute value of the residual. ...
Definition: weightedresidreductioncriterion.hh:58
void setResidualTolerance(Scalar tol)
Sets the maximum absolute tolerated residual.
Definition: weightedresidreductioncriterion.hh:129
void update(const Vector &curSol, const Vector &curResid)
Definition: weightedresidreductioncriterion.hh:187
Scalar residualWeight(size_t outerIdx, unsigned innerIdx) const
Return the relative weight of a row of the residual.
Definition: weightedresidreductioncriterion.hh:107