Ewoms::Linear::FixPointCriterion< Vector, CollectiveCommunication > Class Template Reference

Provides a convergence criterion for the linear solvers which looks at the weighted maximum of the difference between two iterations. More...

#include <fixpointcriterion.hh>

Inheritance diagram for Ewoms::Linear::FixPointCriterion< Vector, CollectiveCommunication >:
Ewoms::Linear::ConvergenceCriterion< Vector >

Public Member Functions

 FixPointCriterion (const CollectiveCommunication &comm)
 
 FixPointCriterion (const CollectiveCommunication &comm, const Vector &weightVec, Scalar reduction)
 
void setWeight (const Vector &weightVec)
 Sets the relative weight of a primary variable. More...
 
Scalar weight (int outerIdx, int innerIdx) const
 Return the relative weight of a primary variable. More...
 
void setTolerance (Scalar tol)
 Set the maximum allowed weighted maximum difference between two iterations. More...
 
Scalar tolerance () const
 Return the maximum allowed weighted difference between two iterations for the solution considered to be converged.
 
void setInitial (const Vector &curSol, const Vector &curResid OPM_UNUSED)
 Set the initial solution of the linear system of equations. More...
 
void update (const Vector &curSol, const Vector &changeIndicator OPM_UNUSED, const Vector &curResid OPM_UNUSED)
 Update the internal members of the convergence criterion with the current solution. More...
 
bool converged () const
 Returns true if and only if the convergence criterion is met. More...
 
Scalar accuracy () const
 Returns the accuracy of the solution at the last update. More...
 
- Public Member Functions inherited from Ewoms::Linear::ConvergenceCriterion< Vector >
virtual ~ConvergenceCriterion ()
 Destructor. More...
 
virtual void setInitial (const Vector &curSol, const Vector &curResid)=0
 Set the initial solution of the linear system of equations. More...
 
virtual void update (const Vector &curSol, const Vector &changeIndicator, const Vector &curResid)=0
 Update the internal members of the convergence criterion with the current solution. More...
 
virtual bool failed () const
 Returns true if the convergence criterion cannot be met anymore because the solver has broken down.
 
virtual void printInitial (std::ostream &os OPM_UNUSED=std::cout) const
 Prints the initial information about the convergence behaviour. More...
 
virtual void print (Scalar iter OPM_UNUSED, std::ostream &os OPM_UNUSED=std::cout) const
 Prints the information about the convergence behaviour for the current iteration. More...
 

Detailed Description

template<class Vector, class CollectiveCommunication>
class Ewoms::Linear::FixPointCriterion< Vector, CollectiveCommunication >

Provides a convergence criterion for the linear solvers which looks at the weighted maximum of the difference between two iterations.

For the FixPointCriterion, the error of the solution is defined as

\[ e^k = \max_i\{ \left| w_i \delta^k_i \right| \}\;, \]

where $\delta = x^k - x^{k + 1} $ is the difference between two consequtive iterative solution vectors $x^k$ and $x^{k + 1}$ and $w_i$ is the weight of the $i$-th degree of freedom.

This criterion requires that the block type of the vector is a Dune::FieldVector

Member Function Documentation

◆ accuracy()

template<class Vector , class CollectiveCommunication >
Scalar Ewoms::Linear::FixPointCriterion< Vector, CollectiveCommunication >::accuracy ( ) const
inlinevirtual

Returns the accuracy of the solution at the last update.

A value of zero means that the solution was exact.

Implements Ewoms::Linear::ConvergenceCriterion< Vector >.

◆ converged()

template<class Vector , class CollectiveCommunication >
bool Ewoms::Linear::FixPointCriterion< Vector, CollectiveCommunication >::converged ( ) const
inlinevirtual

Returns true if and only if the convergence criterion is met.

Implements Ewoms::Linear::ConvergenceCriterion< Vector >.

◆ setInitial()

template<class Vector , class CollectiveCommunication >
void Ewoms::Linear::FixPointCriterion< Vector, CollectiveCommunication >::setInitial ( const Vector &  curSol,
const Vector &curResid  OPM_UNUSED 
)
inline

Set the initial solution of the linear system of equations.

This version of the method does NOT take the two-norm of the residual as argument. If the two-norm of the defect is available for the linear solver, the version of the update() method with it should be called.

Parameters
curSolThe current iterative solution of the linear system of equations
curResidThe residual vector of the current iterative solution of the linear system of equations

◆ setTolerance()

template<class Vector , class CollectiveCommunication >
void Ewoms::Linear::FixPointCriterion< Vector, CollectiveCommunication >::setTolerance ( Scalar  tol)
inline

Set the maximum allowed weighted maximum difference between two iterations.

Set the maximum allowed maximum difference between two iterationsfor the solution considered to be converged.

◆ setWeight()

template<class Vector , class CollectiveCommunication >
void Ewoms::Linear::FixPointCriterion< Vector, CollectiveCommunication >::setWeight ( const Vector &  weightVec)
inline

Sets the relative weight of a primary variable.

For the FixPointCriterion, the error of the solution is defined as

\[ e^k = \max_i\{ \left| w_i \delta^k_i \right| \}\;, \]

where $\delta = x^k - x^{k + 1} $ is the difference between two consequtive iterative solution vectors $x^k$ and $x^{k + 1}$ and $w_i$ is the weight of the $i$-th degree of freedom.

This method is specific to the FixPointCriterion.

Parameters
weightVecA Dune::BlockVector<Dune::FieldVector<Scalar, n> > with the relative weights of the degrees of freedom

◆ update()

template<class Vector , class CollectiveCommunication >
void Ewoms::Linear::FixPointCriterion< Vector, CollectiveCommunication >::update ( const Vector &  curSol,
const Vector &changeIndicator  OPM_UNUSED,
const Vector &curResid  OPM_UNUSED 
)
inline

Update the internal members of the convergence criterion with the current solution.

This version of the method does NOT take the two-norm of the residual as argument. If the two-norm of the defect is available for the linear solver, the version of the update() method with it should be called.

Parameters
curSolThe current iterative solution of the linear system of equations
changeIndicatorA vector where all non-zero values indicate that the solution has changed since the last iteration.
curResidThe residual vector of the current iterative solution of the linear system of equations

◆ weight()

template<class Vector , class CollectiveCommunication >
Scalar Ewoms::Linear::FixPointCriterion< Vector, CollectiveCommunication >::weight ( int  outerIdx,
int  innerIdx 
) const
inline

Return the relative weight of a primary variable.

For the FixPointCriterion, the error of the solution is defined as

\[ e^k = \max_i\{ \left| w_i \delta^k_i \right| \}\;, \]

where $\delta = x^k - x^{k + 1} $ is the difference between two consequtive iterative solution vectors $x^k$ and $x^{k + 1}$ and $w_i$ is the weight of the $i$-th degree of freedom.

This method is specific to the FixPointCriterion.

Parameters
outerIdxThe index of the outer vector (i.e. Dune::BlockVector)
innerIdxThe index of the inner vector (i.e. Dune::FieldVector)

The documentation for this class was generated from the following file: