Dune Core Modules (unstable)

Richardson preconditioner. More...

#include <dune/istl/preconditioners.hh>

Public Types

typedef X domain_type
 The domain type of the preconditioner.
 
typedef Y range_type
 The range type of the preconditioner.
 
typedef X::field_type field_type
 The field type of the preconditioner.
 
typedef Simd::Scalar< field_typescalar_field_type
 scalar type underlying the field_type
 
typedef FieldTraits< scalar_field_type >::real_type real_field_type
 real scalar type underlying the field_type
 

Public Member Functions

 Richardson (real_field_type w=1.0)
 Constructor. More...
 
 Richardson (const ParameterTree &configuration)
 Constructor. More...
 
virtual void pre ([[maybe_unused]] X &x, [[maybe_unused]] Y &b)
 Prepare the preconditioner. More...
 
virtual void apply (X &v, const Y &d)
 Apply the precondioner. More...
 
virtual void post ([[maybe_unused]] X &x)
 Clean up. More...
 
virtual SolverCategory::Category category () const
 Category of the preconditioner (see SolverCategory::Category)
 
virtual void pre (X &x, Y &b)=0
 Prepare the preconditioner. More...
 
virtual void post (X &x)=0
 Clean up. More...
 

Detailed Description

template<class X, class Y>
class Dune::Richardson< X, Y >

Richardson preconditioner.

Multiply simply by a constant.

Template Parameters
XType of the update
YType of the defect

Constructor & Destructor Documentation

◆ Richardson() [1/2]

template<class X , class Y >
Dune::Richardson< X, Y >::Richardson ( real_field_type  w = 1.0)
inline

Constructor.

Constructor gets all parameters to operate the prec.

Parameters
wThe relaxation factor.

◆ Richardson() [2/2]

template<class X , class Y >
Dune::Richardson< X, Y >::Richardson ( const ParameterTree configuration)
inline

Constructor.

Parameters
configurationParameterTree containing preconditioner parameters.
ParameterTree Key Meaning
relaxation The relaxation factor. default=1.0

See ISTL_Factory for the ParameterTree layout and examples.

Member Function Documentation

◆ apply()

template<class X , class Y >
virtual void Dune::Richardson< X, Y >::apply ( X &  v,
const Y &  d 
)
inlinevirtual

Apply the precondioner.

Apply one step of the preconditioner to the system A(v)=d. On entry v=0 and d=b-A(x) (although this might not be computed in that way. On exit v contains the update, i.e one step computes \( v = M^{-1} d \) where \( M \) is the approximate inverse of the operator \( A \) characterizing the preconditioner.

Parameters
[out]vThe update to be computed
dThe current defect.

Implements Dune::Preconditioner< X, Y >.

◆ post() [1/2]

template<class X , class Y >
virtual void Dune::Richardson< X, Y >::post ( [[maybe_unused] ] X &  x)
inlinevirtual

Clean up.

Clean up. This method is called after the last apply call for the linear system to be solved. Memory may be deallocated safely here. x is the solution of the linear equation.

Parameters
xThe right hand side of the equation.

◆ post() [2/2]

template<class X , class Y >
virtual void Dune::Preconditioner< X, Y >::post ( X &  x)
pure virtualinherited

Clean up.

This method is called after the last apply call for the linear system to be solved. Memory may be deallocated safely here. x is the solution of the linear equation.

Parameters
xThe right hand side of the equation.

Implemented in Dune::BlockPreconditioner< X, Y, C, P >, Dune::NonoverlappingBlockPreconditioner< C, P >, Dune::InverseOperator2Preconditioner< O, c >, Dune::Amg::KAMG< M, X, S, PI, K, A >, Dune::Amg::FastAMG< M, X, PI, A >, Dune::Amg::AMG< M, X, S, PI, A >, Dune::Amg::AMG< M, X, S, SequentialInformation, std::allocator< X > >, and Dune::Amg::AMG< Operator, X, Smoother >.

◆ pre() [1/2]

template<class X , class Y >
virtual void Dune::Richardson< X, Y >::pre ( [[maybe_unused] ] X &  x,
[[maybe_unused] ] Y &  b 
)
inlinevirtual

Prepare the preconditioner.

Prepare the preconditioner. A solver solves a linear operator equation A(x)=b by applying one or several steps of the preconditioner. The method pre() is called before the first apply operation. b and x are right hand side and solution vector of the linear system respectively. It may. e.g., scale the system, allocate memory or compute a (I)LU decomposition. Note: The ILU decomposition could also be computed in the constructor or with a separate method of the derived method if several linear systems with the same matrix are to be solved.

Note
if a preconditioner is copied (e.g. for a second thread) again the pre() method has to be called to ensure proper memory management.
X x(0.0);
Y b = ...; // rhs
Preconditioner<X,Y> prec(...);
prec.pre(x,b); // prepare the preconditioner
prec.apply(x,b); // can be called multiple times now...
prec.post(x); // cleanup internal state
Parameters
xThe left hand side of the equation.
bThe right hand side of the equation.

◆ pre() [2/2]

template<class X , class Y >
virtual void Dune::Preconditioner< X, Y >::pre ( X &  x,
Y &  b 
)
pure virtualinherited

Prepare the preconditioner.

A solver solves a linear operator equation A(x)=b by applying one or several steps of the preconditioner. The method pre() is called before the first apply operation. b and x are right hand side and solution vector of the linear system respectively. It may. e.g., scale the system, allocate memory or compute a (I)LU decomposition. Note: The ILU decomposition could also be computed in the constructor or with a separate method of the derived method if several linear systems with the same matrix are to be solved.

Note
if a preconditioner is copied (e.g. for a second thread) again the pre() method has to be called to ensure proper memory management.
X x(0.0);
Y b = ...; // rhs
Preconditioner<X,Y> prec(...);
prec.pre(x,b); // prepare the preconditioner
prec.apply(x,b); // can be called multiple times now...
prec.post(x); // cleanup internal state
Parameters
xThe left hand side of the equation.
bThe right hand side of the equation.

Implemented in Dune::BlockPreconditioner< X, Y, C, P >, Dune::NonoverlappingBlockPreconditioner< C, P >, and Dune::InverseOperator2Preconditioner< O, c >.


The documentation for this class was generated from the following file:
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.80.0 (Mar 27, 23:31, 2024)