DUNE PDELab (git)

Dune::PDELab::SeqDGAMGPrec< DGMatrix, DGPrec, CGPrec, P > Class Template Reference

#include <dune/pdelab/backend/istl/seq_amg_dg_backend.hh>

Public Types

typedef DGPrec::domain_type domain_type
 The domain type of the preconditioner.
 
typedef DGPrec::range_type range_type
 The range type of the preconditioner.
 
typedef X::field_type field_type
 The field type of the preconditioner.
 

Public Member Functions

SolverCategory::Category category () const override
 Category of the preconditioner (see SolverCategory::Category)
 
 SeqDGAMGPrec (DGMatrix &dgmatrix_, DGPrec &dgprec_, CGPrec &cgprec_, P &p_, int n1_, int n2_)
 Constructor. More...
 
virtual void pre (X &x, Y &b) override
 Prepare the preconditioner. More...
 
virtual void apply (X &x, const Y &b) override
 Apply the precondioner. More...
 
virtual void post (X &x) override
 Clean up. More...
 

Detailed Description

template<class DGMatrix, class DGPrec, class CGPrec, class P>
class Dune::PDELab::SeqDGAMGPrec< DGMatrix, DGPrec, CGPrec, P >

An ISTL preconditioner for DG based on AMG applied to CG subspace

The template parameters are: DGMatrix BCRSMatrix assembled with DG DGPrec preconditioner to be used for DG CGPrec preconditioner to be used on CG subspace P BCRSMatrix for grid transfer

Constructor & Destructor Documentation

◆ SeqDGAMGPrec()

template<class DGMatrix , class DGPrec , class CGPrec , class P >
Dune::PDELab::SeqDGAMGPrec< DGMatrix, DGPrec, CGPrec, P >::SeqDGAMGPrec ( DGMatrix &  dgmatrix_,
DGPrec &  dgprec_,
CGPrec &  cgprec_,
P &  p_,
int  n1_,
int  n2_ 
)
inline

Constructor.

Constructor gets all parameters to operate the prec.

Parameters
AThe matrix to operate on.
nThe number of iterations to perform.
wThe relaxation factor.

Member Function Documentation

◆ apply()

template<class DGMatrix , class DGPrec , class CGPrec , class P >
virtual void Dune::PDELab::SeqDGAMGPrec< DGMatrix, DGPrec, CGPrec, P >::apply ( X &  x,
const Y &  b 
)
inlineoverridevirtual

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< DGPrec::domain_type, DGPrec::range_type >.

◆ post()

template<class DGMatrix , class DGPrec , class CGPrec , class P >
virtual void Dune::PDELab::SeqDGAMGPrec< DGMatrix, DGPrec, CGPrec, P >::post ( X &  x)
inlineoverridevirtual

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.

Implements Dune::Preconditioner< DGPrec::domain_type, DGPrec::range_type >.

◆ pre()

template<class DGMatrix , class DGPrec , class CGPrec , class P >
virtual void Dune::PDELab::SeqDGAMGPrec< DGMatrix, DGPrec, CGPrec, P >::pre ( X &  x,
Y &  b 
)
inlineoverridevirtual

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
prec.pre(x,b); // prepare the preconditioner
prec.apply(x,b); // can be called multiple times now...
prec.post(x); // cleanup internal state
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:33
Parameters
xThe left hand side of the equation.
bThe right hand side of the equation.

Implements Dune::Preconditioner< DGPrec::domain_type, DGPrec::range_type >.


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.111.3 (Nov 12, 23:30, 2024)