27#ifndef EWOMS_RESID_REDUCTION_CRITERION_HH
28#define EWOMS_RESID_REDUCTION_CRITERION_HH
32#include <dune/istl/scalarproducts.hh>
50template <
class Vector>
53 using Scalar =
typename Vector::field_type;
72 {
return tolerance_; }
79 static constexpr Scalar eps = std::numeric_limits<Scalar>::min()*1
e10;
83 curDefect_ = scalarProduct_.norm(
curResid);
84 lastDefect_ = curDefect_;
85 initialDefect_ = std::max(curDefect_, eps);
95 lastDefect_ = curDefect_;
96 curDefect_ = scalarProduct_.norm(
curResid);
109 {
return curDefect_/initialDefect_; }
116 os << std::setw(20) <<
"iteration ";
117 os << std::setw(20) <<
"residual ";
118 os << std::setw(20) <<
"accuracy ";
119 os << std::setw(20) <<
"rate ";
126 void print(Scalar
iter, std::ostream&
os = std::cout)
const override
128 static constexpr Scalar eps = std::numeric_limits<Scalar>::min()*1
e10;
130 os << std::setw(20) <<
iter <<
" ";
131 os << std::setw(20) << curDefect_ <<
" ";
133 os << std::setw(20) << (lastDefect_/std::max(eps, curDefect_)) <<
" ";
138 Dune::ScalarProduct<Vector>& scalarProduct_;
141 Scalar initialDefect_;
Base class for all convergence criteria which only defines an virtual API.
Definition convergencecriterion.hh:56
Provides a convergence criterion which looks at the reduction of the two-norm of the residual for the...
Definition residreductioncriterion.hh:52
void printInitial(std::ostream &os=std::cout) const override
Prints the initial information about the convergence behaviour.
Definition residreductioncriterion.hh:114
Scalar accuracy() const override
Returns the accuracy of the solution at the last update.
Definition residreductioncriterion.hh:108
void print(Scalar iter, std::ostream &os=std::cout) const override
Prints the information about the convergence behaviour for the current iteration.
Definition residreductioncriterion.hh:126
bool converged() const override
Returns true if and only if the convergence criterion is met.
Definition residreductioncriterion.hh:102
void setInitial(const Vector &, const Vector &curResid) override
Set the initial solution of the linear system of equations.
Definition residreductioncriterion.hh:77
void setTolerance(Scalar tol)
Set the maximum allowed weighted maximum of the reduction of the linear residual.
Definition residreductioncriterion.hh:65
Scalar tolerance() const
Return the maximum allowed weighted maximum of the reduction of the linear residual.
Definition residreductioncriterion.hh:71
void update(const Vector &, const Vector &, const Vector &curResid) override
Definition residreductioncriterion.hh:91
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