dune-fem 2.8.0
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
Dune::Fem::CGInverseOperator< DiscreteFunction, Op > Class Template Reference

Inverse operator base on CG method. Uses a runtime parameter fem.preconditioning which enables diagonal preconditioning if diagonal matrix entries are available, i.e., Op :: assembled is true. More...

#include <dune/fem/solver/cginverseoperator.hh>

Inheritance diagram for Dune::Fem::CGInverseOperator< DiscreteFunction, Op >:
Inheritance graph

Public Types

typedef SolverParameter SolverParameterType
 
typedef BaseType::DomainFunctionType DomainFunctionType
 
typedef BaseType::RangeFunctionType RangeFunctionType
 
typedef DomainFunctionType DestinationType
 
typedef Op OperatorType
 type of operator
 
typedef OperatorType::RangeFieldType RangeFieldType
 
typedef Dune::FieldTraits< RangeFieldType >::real_type RealType
 
typedef Fem::Operator< RangeFunctionType, DomainFunctionTypePreconditioningType
 
typedef Fem::Operator< RangeFunctionType, DomainFunctionTypePreconditionerType
 
typedef DomainFunction::RangeFieldType DomainFieldType
 field type of the operator's domain
 

Public Member Functions

 CGInverseOperator (const SolverParameter &param=SolverParameter(Parameter::container()))
 
 CGInverseOperator (RealType redEps, RealType absLimit, unsigned int maxIter, bool verbose, const ParameterReader &parameter=Parameter::container())
 constructor of CGInverseOperator
 
 CGInverseOperator (RealType redEps, RealType absLimit, unsigned int maxIter, const ParameterReader &parameter=Parameter::container())
 constructor of CGInverseOperator
 
 CGInverseOperator (RealType redEps, RealType absLimit, const ParameterReader &parameter=Parameter::container())
 
template<class LinearOperator , std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0>
 CGInverseOperator (const LinearOperator &op, RealType redEps, RealType absLimit, unsigned int maxIter, bool verbose, const ParameterReader &parameter=Parameter::container())
 constructor of CGInverseOperator
 
template<class LinearOperator , std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0>
 CGInverseOperator (const LinearOperator &op, RealType redEps, RealType absLimit, unsigned int maxIter, const ParameterReader &parameter=Parameter::container())
 constructor of CGInverseOperator
 
template<class LinearOperator , std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0>
 CGInverseOperator (const LinearOperator &op, RealType redEps, RealType absLimit, const ParameterReader &parameter=Parameter::container())
 
 CGInverseOperator (const OperatorType &op, const PreconditioningType &precond, RealType redEps, RealType absLimit, unsigned int maxIter, bool verbose, const ParameterReader &parameter=Parameter::container())
 constructor of CGInverseOperator
 
 CGInverseOperator (const OperatorType &op, const PreconditioningType &precond, RealType redEps, RealType absLimit, unsigned int maxIter, const ParameterReader &parameter=Parameter::container())
 constructor of CGInverseOperator
 
 CGInverseOperator (const OperatorType &op, const PreconditioningType &precond, RealType redEps, RealType absLimit, const ParameterReader &parameter=Parameter::container())
 
template<class LinearOperator , std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0>
void bind (const LinearOperator &op)
 
void unbind ()
 
void bind (const OperatorType &op)
 
void bind (const OperatorType &op, const PreconditionerType &precond)
 
virtual void operator() (const DomainFunctionType &arg, RangeFunctionType &dest) const
 application operator
 
template<typename... A>
void prepare (A...) const
 
virtual void apply (const DomainFunctionType &arg, RangeFunctionType &dest) const
 application operator
 
unsigned int iterations () const
 number of iterations needed for last solve
 
void setMaxIterations (unsigned int maxIter)
 
double averageCommTime () const
 return average communication time during last solve
 
virtual void finalize ()
 finalization of operator
 

Protected Member Functions

template<class LinearOperator >
void checkPreconditioning (const LinearOperator &linearOp)
 

Protected Attributes

std::unique_ptr< PreconditioningTypeprecondObj_
 
const PreconditionerTypepreconditioner_
 
SolverParameter parameter_
 
const OperatorTypeoperator_ = nullptr
 
SolverType solver_
 

Detailed Description

template<class DiscreteFunction, class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
class Dune::Fem::CGInverseOperator< DiscreteFunction, Op >

Inverse operator base on CG method. Uses a runtime parameter fem.preconditioning which enables diagonal preconditioning if diagonal matrix entries are available, i.e., Op :: assembled is true.

Member Typedef Documentation

◆ DestinationType

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
typedef DomainFunctionType Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::DestinationType

◆ DomainFieldType

typedef DomainFunction::RangeFieldType Dune::Fem::Operator< DiscreteFunction , DiscreteFunction >::DomainFieldType
inherited

field type of the operator's domain

◆ DomainFunctionType

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
typedef BaseType::DomainFunctionType Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::DomainFunctionType

◆ OperatorType

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
typedef Op Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::OperatorType

type of operator

◆ PreconditionerType

template<class DiscreteFunction >
typedef Fem::Operator< RangeFunctionType, DomainFunctionType > Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::PreconditionerType
inherited

◆ PreconditioningType

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
typedef Fem::Operator< RangeFunctionType,DomainFunctionType > Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::PreconditioningType

◆ RangeFieldType

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
typedef OperatorType::RangeFieldType Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::RangeFieldType

◆ RangeFunctionType

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
typedef BaseType::RangeFunctionType Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::RangeFunctionType

◆ RealType

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
typedef Dune::FieldTraits<RangeFieldType>::real_type Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::RealType

◆ SolverParameterType

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
typedef SolverParameter Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::SolverParameterType

Constructor & Destructor Documentation

◆ CGInverseOperator() [1/10]

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::CGInverseOperator ( const SolverParameter param = SolverParameter(Parameter::container()))
inline

◆ CGInverseOperator() [2/10]

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::CGInverseOperator ( RealType  redEps,
RealType  absLimit,
unsigned int  maxIter,
bool  verbose,
const ParameterReader parameter = Parameter::container() 
)
inline

constructor of CGInverseOperator

Parameters
[in]redEpsreduction epsilon
[in]absLimitabsolut limit of residual
[in]maxItermaximum number of iteration steps
[in]verboseverbosity

◆ CGInverseOperator() [3/10]

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::CGInverseOperator ( RealType  redEps,
RealType  absLimit,
unsigned int  maxIter,
const ParameterReader parameter = Parameter::container() 
)
inline

constructor of CGInverseOperator

Parameters
[in]redEpsreduction epsilon
[in]absLimitabsolut limit of residual
[in]maxItermaximum number of iteration steps

◆ CGInverseOperator() [4/10]

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::CGInverseOperator ( RealType  redEps,
RealType  absLimit,
const ParameterReader parameter = Parameter::container() 
)
inline

◆ CGInverseOperator() [5/10]

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
template<class LinearOperator , std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0>
Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::CGInverseOperator ( const LinearOperator op,
RealType  redEps,
RealType  absLimit,
unsigned int  maxIter,
bool  verbose,
const ParameterReader parameter = Parameter::container() 
)
inline

constructor of CGInverseOperator

Parameters
[in]opoperator to invert
[in]redEpsreduction epsilon
[in]absLimitabsolut limit of residual
[in]maxItermaximum number of iteration steps
[in]verboseverbosity

◆ CGInverseOperator() [6/10]

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
template<class LinearOperator , std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0>
Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::CGInverseOperator ( const LinearOperator op,
RealType  redEps,
RealType  absLimit,
unsigned int  maxIter,
const ParameterReader parameter = Parameter::container() 
)
inline

constructor of CGInverseOperator

Parameters
[in]opoperator to invert
[in]redEpsreduction epsilon
[in]absLimitabsolut limit of residual
[in]maxItermaximum number of iteration steps

◆ CGInverseOperator() [7/10]

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
template<class LinearOperator , std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0>
Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::CGInverseOperator ( const LinearOperator op,
RealType  redEps,
RealType  absLimit,
const ParameterReader parameter = Parameter::container() 
)
inline

◆ CGInverseOperator() [8/10]

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::CGInverseOperator ( const OperatorType op,
const PreconditioningType precond,
RealType  redEps,
RealType  absLimit,
unsigned int  maxIter,
bool  verbose,
const ParameterReader parameter = Parameter::container() 
)
inline

constructor of CGInverseOperator

Parameters
[in]opoperator to invert
[in]precondprecondition operator
[in]redEpsreduction epsilon
[in]absLimitabsolut limit of residual
[in]maxItermaximum number of iteration steps

◆ CGInverseOperator() [9/10]

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::CGInverseOperator ( const OperatorType op,
const PreconditioningType precond,
RealType  redEps,
RealType  absLimit,
unsigned int  maxIter,
const ParameterReader parameter = Parameter::container() 
)
inline

constructor of CGInverseOperator

Parameters
[in]opoperator to invert
[in]precondprecondition operator
[in]redEpsreduction epsilon
[in]absLimitabsolut limit of residual
[in]maxItermaximum number of iteration steps

◆ CGInverseOperator() [10/10]

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::CGInverseOperator ( const OperatorType op,
const PreconditioningType precond,
RealType  redEps,
RealType  absLimit,
const ParameterReader parameter = Parameter::container() 
)
inline

Member Function Documentation

◆ apply()

template<class DiscreteFunction >
virtual void Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::apply ( const DomainFunctionType arg,
RangeFunctionType dest 
) const
inlinevirtualinherited

application operator

The application operator actually solves the linear system $op(dest) = arg$ using the CG method.

Parameters
[in]argargument discrete function
[out]destdestination discrete function

◆ averageCommTime()

template<class DiscreteFunction >
double Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::averageCommTime ( ) const
inlineinherited

return average communication time during last solve

◆ bind() [1/3]

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
template<class LinearOperator , std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0>
void Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::bind ( const LinearOperator op)
inline

◆ bind() [2/3]

template<class DiscreteFunction >
void Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::bind ( const OperatorType op)
inlineinherited

◆ bind() [3/3]

template<class DiscreteFunction >
void Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::bind ( const OperatorType op,
const PreconditionerType precond 
)
inlineinherited

◆ checkPreconditioning()

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
template<class LinearOperator >
void Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::checkPreconditioning ( const LinearOperator linearOp)
inlineprotected

◆ finalize()

virtual void Dune::Fem::Operator< DiscreteFunction , DiscreteFunction >::finalize ( )
inlinevirtualinherited

finalization of operator

Note
The default implementation is empty.

◆ iterations()

template<class DiscreteFunction >
unsigned int Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::iterations ( ) const
inlineinherited

number of iterations needed for last solve

◆ operator()()

template<class DiscreteFunction >
virtual void Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::operator() ( const DomainFunctionType arg,
RangeFunctionType dest 
) const
inlinevirtualinherited

application operator

The application operator actually solves the linear system $op(dest) = arg$ using the CG method.

Parameters
[in]argargument discrete function
[out]destdestination discrete function

Implements Dune::Fem::Operator< DiscreteFunction, DiscreteFunction >.

◆ prepare()

template<class DiscreteFunction >
template<typename... A>
void Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::prepare ( A...  ) const
inlineinherited

◆ setMaxIterations()

template<class DiscreteFunction >
void Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::setMaxIterations ( unsigned int  maxIter)
inlineinherited

◆ unbind()

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
void Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::unbind ( )
inline

Member Data Documentation

◆ operator_

template<class DiscreteFunction >
const OperatorType* Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::operator_ = nullptr
protectedinherited

◆ parameter_

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
SolverParameter Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::parameter_
protected

◆ preconditioner_

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
const PreconditionerType* Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::preconditioner_
protected

◆ precondObj_

template<class DiscreteFunction , class Op = Fem::Operator< DiscreteFunction, DiscreteFunction >>
std::unique_ptr< PreconditioningType > Dune::Fem::CGInverseOperator< DiscreteFunction, Op >::precondObj_
protected

◆ solver_

template<class DiscreteFunction >
SolverType Dune::Fem::Solver::CGInverseOperator< DiscreteFunction >::solver_
protectedinherited

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