dune-pdelab  2.7-git
Public Member Functions | Static Public Attributes | List of all members
Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver > Class Template Reference

Local operator that can be used to create a fully matrix-free Jacobi preconditioner. More...

#include <dune/pdelab/backend/istl/matrixfree/iterativeblockjacobipreconditioner.hh>

Inheritance diagram for Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver >:
Inheritance graph

Public Types

Flags selective assembly
enum  { doSkipEntity }
 Whether to do selective assembly on the elements, i.e. whether or not skip_entity() should be called. More...
 
enum  { doSkipIntersection }
 Whether to do selective assembly on the intersections, i.e. whether or not skip_intersection() should be called. More...
 
Flags for the sparsity pattern
enum  { doPatternVolume }
 Whether to assemble the pattern on the elements, i.e. whether or not pattern_volume() should be called. More...
 
enum  { doPatternVolumePostSkeleton }
 Whether to assemble the pattern on the elements after the skeleton has been handled, i.e. whether or not pattern_volume_post_skeleton() should be called. More...
 
enum  { doPatternSkeleton }
 Whether to assemble the pattern on the interior intersections, i.e. whether or not pattern_skeleton() should be called. More...
 
enum  { doPatternBoundary }
 Whether to assemble the pattern on the boundary intersections, i.e. whether or not pattern_boundary() should be called. More...
 
Flags for the non-constant part of the residual and the jacobian
enum  { doAlphaVolume }
 Whether to call the local operator's alpha_volume(), jacobian_apply_volume() and jacobian_volume(). More...
 
enum  { doAlphaVolumePostSkeleton }
 Whether to call the local operator's alpha_volume_post_skeleton(), jacobian_apply_volume_post_skeleton() and jacobian_volume_post_skeleton(). More...
 
enum  { doAlphaSkeleton }
 Whether to call the local operator's alpha_skeleton(), jacobian_apply_skeleton() and jacobian_skeleton(). More...
 
enum  { doAlphaBoundary }
 Whether to call the local operator's alpha_boundary(), jacobian_apply_boundary() and jacobian_boundary(). More...
 
Flags for the constant part of the residual
enum  { doLambdaVolume }
 Whether to call the local operator's lambda_volume(). More...
 
enum  { doLambdaVolumePostSkeleton }
 Whether to call the local operator's lambda_volume_post_skeleton(). More...
 
enum  { doLambdaSkeleton }
 Whether to call the local operator's lambda_skeleton(). More...
 
enum  { doLambdaBoundary }
 Whether to call the local operator's lambda_boundary(). More...
 
Special flags
enum  { doSkeletonTwoSided }
 Whether to visit the skeleton methods from both sides. More...
 
enum  { isLinear }
 Wheter the local operator describes a linear problem. More...
 

Public Member Functions

 IterativeBlockJacobiPreconditionerLocalOperator (const BlockDiagonalLocalOperator &blockDiagonalLocalOperator, const PointDiagonalLocalOperator &pointDiagonalLocalOperator, const GridFunctionSpace &gridFunctionSpace, SolverStatistics< int > &solverStatiscits, BlockSolverOptions solveroptions, const bool verbose=0)
 Constructor. More...
 
bool requireSetup ()
 
void setRequireSetup (bool v)
 
template<typename EG , typename LFSU , typename X , typename LFSV , typename Y >
void alpha_volume (const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, Y &y) const
 Prepare fully matrix-free preconditioner. More...
 
template<typename EG , typename LFSU , typename Z , typename LFSV , typename Y >
void jacobian_apply_volume (const EG &eg, const LFSU &lfsu, const Z &z, const LFSV &lfsv, Y &y) const
 Apply fully matrix-free preconditioner - linear case. More...
 
template<typename EG , typename LFSU , typename X , typename Z , typename LFSV , typename Y >
void jacobian_apply_volume (const EG &eg, const LFSU &lfsu, const X &x, const Z &z, const LFSV &lfsv, Y &y) const
 Apply fully matrix-free preconditioner - nonlinear case. More...
 
template<typename LFSU , typename LFSV , typename LocalPattern >
void pattern_volume (const LFSU &lfsu, const LFSV &lfsv, LocalPattern &pattern) const
 

Static Public Attributes

static constexpr bool doPatternVolume = true
 
static constexpr bool doAlphaVolume = true
 
static constexpr bool isLinear = BlockDiagonalLocalOperator::isLinear
 

Detailed Description

template<typename BlockDiagonalLocalOperator, typename PointDiagonalLocalOperator, typename GridFunctionSpace, typename DomainField, template< typename > class IterativeSolver>
class Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver >

Local operator that can be used to create a fully matrix-free Jacobi preconditioner.

Similar to the partial matrix-free class AssembledBlockJacobiPreconditionerLocalOperator this implements a local operator that can be used to implement a matrix-free Jacobi preconditioner. In contrast to the other class this implementation will be fully matrix-free.

A matrix-free Jacobi preconditioner needs to be able to apply the inverse of the block diagonal D to a vector. In the partial matrix-free class mentioned above this was done by assembling the diagonal blocks and inverting this matrix. In order to get a fully matrix-free version we instead use a Krylow solver on the diagonal block. This can once again be done in a matrix-free way but we will need a preconditioner for iterative solver. For this purpose we use another Jacobi preconditioner that operates on a single block. This means we need to assemble the point diagonal of this block and apply the inverse.

For examples see dune-pdelab/dune/pdelab/test/matrixfree/.

Template Parameters
BlockDiagonalLocalOperatorA lop for local application of diagonal blocks
PointDiagonalLocalOperatorA lop for local assembly of point diagonal
GridFunctionSpaceA grid function space
DomainFieldDomain field
IterativeSolverSolver that will be used to 'invert' the diagonal blocks

Member Enumeration Documentation

◆ anonymous enum

anonymous enum
inherited

Whether to do selective assembly on the elements, i.e. whether or not skip_entity() should be called.

Enumerator
doSkipEntity 

◆ anonymous enum

anonymous enum
inherited

Whether to do selective assembly on the intersections, i.e. whether or not skip_intersection() should be called.

Enumerator
doSkipIntersection 

◆ anonymous enum

anonymous enum
inherited

Whether to assemble the pattern on the elements, i.e. whether or not pattern_volume() should be called.

Enumerator
doPatternVolume 

◆ anonymous enum

anonymous enum
inherited

Whether to assemble the pattern on the elements after the skeleton has been handled, i.e. whether or not pattern_volume_post_skeleton() should be called.

Enumerator
doPatternVolumePostSkeleton 

◆ anonymous enum

anonymous enum
inherited

Whether to assemble the pattern on the interior intersections, i.e. whether or not pattern_skeleton() should be called.

Enumerator
doPatternSkeleton 

◆ anonymous enum

anonymous enum
inherited

Whether to assemble the pattern on the boundary intersections, i.e. whether or not pattern_boundary() should be called.

Enumerator
doPatternBoundary 

◆ anonymous enum

anonymous enum
inherited

Whether to call the local operator's alpha_volume(), jacobian_apply_volume() and jacobian_volume().

Enumerator
doAlphaVolume 

◆ anonymous enum

anonymous enum
inherited

Whether to call the local operator's alpha_volume_post_skeleton(), jacobian_apply_volume_post_skeleton() and jacobian_volume_post_skeleton().

Enumerator
doAlphaVolumePostSkeleton 

◆ anonymous enum

anonymous enum
inherited

Whether to call the local operator's alpha_skeleton(), jacobian_apply_skeleton() and jacobian_skeleton().

Enumerator
doAlphaSkeleton 

◆ anonymous enum

anonymous enum
inherited

Whether to call the local operator's alpha_boundary(), jacobian_apply_boundary() and jacobian_boundary().

Enumerator
doAlphaBoundary 

◆ anonymous enum

anonymous enum
inherited

Whether to call the local operator's lambda_volume().

Enumerator
doLambdaVolume 

◆ anonymous enum

anonymous enum
inherited

Whether to call the local operator's lambda_volume_post_skeleton().

Enumerator
doLambdaVolumePostSkeleton 

◆ anonymous enum

anonymous enum
inherited

Whether to call the local operator's lambda_skeleton().

Enumerator
doLambdaSkeleton 

◆ anonymous enum

anonymous enum
inherited

Whether to call the local operator's lambda_boundary().

Enumerator
doLambdaBoundary 

◆ anonymous enum

anonymous enum
inherited

Whether to visit the skeleton methods from both sides.

Enumerator
doSkeletonTwoSided 

◆ anonymous enum

anonymous enum
inherited

Wheter the local operator describes a linear problem.

Enumerator
isLinear 

Constructor & Destructor Documentation

◆ IterativeBlockJacobiPreconditionerLocalOperator()

template<typename BlockDiagonalLocalOperator , typename PointDiagonalLocalOperator , typename GridFunctionSpace , typename DomainField , template< typename > class IterativeSolver>
Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver >::IterativeBlockJacobiPreconditionerLocalOperator ( const BlockDiagonalLocalOperator &  blockDiagonalLocalOperator,
const PointDiagonalLocalOperator &  pointDiagonalLocalOperator,
const GridFunctionSpace gridFunctionSpace,
SolverStatistics< int > &  solverStatiscits,
BlockSolverOptions  solveroptions,
const bool  verbose = 0 
)
inline

Constructor.

Parameters
blockDiagonalLocalOperatorA local operator that can be used to (locally) apply a diagonal block through the <applyLocalDiagonalBlock> function. You can create such a local operator by wrapping your local operator into a <BlockDiagonalLocalOperatorWrapper>
pointDiagonalLocalOperatorA local operator that can be used to (locally) assemble the point diagonal of a diagonal block. You can create such a local operator by wrapping your local operator into a <PointDiagonalLocalOperatorWrapper>
gridFunctionSpaceA grid function space
solverStatA class for export solver statistics
solveroptionsA class for configuring your solver
verboseControls the amount of output

Member Function Documentation

◆ alpha_volume()

template<typename BlockDiagonalLocalOperator , typename PointDiagonalLocalOperator , typename GridFunctionSpace , typename DomainField , template< typename > class IterativeSolver>
template<typename EG , typename LFSU , typename X , typename LFSV , typename Y >
void Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver >::alpha_volume ( const EG &  eg,
const LFSU &  lfsu,
const X &  x,
const LFSV &  lfsv,
Y &  y 
) const
inline

Prepare fully matrix-free preconditioner.

This assembles the point diagonal of the diagonal block and stores its inverse. During the apply step this will be used as a preconditioner for the solver on the local block.

◆ jacobian_apply_volume() [1/2]

template<typename BlockDiagonalLocalOperator , typename PointDiagonalLocalOperator , typename GridFunctionSpace , typename DomainField , template< typename > class IterativeSolver>
template<typename EG , typename LFSU , typename X , typename Z , typename LFSV , typename Y >
void Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver >::jacobian_apply_volume ( const EG &  eg,
const LFSU &  lfsu,
const X &  x,
const Z &  z,
const LFSV &  lfsv,
Y &  y 
) const
inline

Apply fully matrix-free preconditioner - nonlinear case.

Compared to the linear case this needs to set the correct linearization point in the BlockDiagonalOperator

◆ jacobian_apply_volume() [2/2]

template<typename BlockDiagonalLocalOperator , typename PointDiagonalLocalOperator , typename GridFunctionSpace , typename DomainField , template< typename > class IterativeSolver>
template<typename EG , typename LFSU , typename Z , typename LFSV , typename Y >
void Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver >::jacobian_apply_volume ( const EG &  eg,
const LFSU &  lfsu,
const Z &  z,
const LFSV &  lfsv,
Y &  y 
) const
inline

Apply fully matrix-free preconditioner - linear case.

◆ pattern_volume()

template<typename LFSU , typename LFSV , typename LocalPattern >
void Dune::PDELab::FullVolumePattern::pattern_volume ( const LFSU &  lfsu,
const LFSV &  lfsv,
LocalPattern &  pattern 
) const
inlineinherited

◆ requireSetup()

template<typename BlockDiagonalLocalOperator , typename PointDiagonalLocalOperator , typename GridFunctionSpace , typename DomainField , template< typename > class IterativeSolver>
bool Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver >::requireSetup ( )
inline

◆ setRequireSetup()

template<typename BlockDiagonalLocalOperator , typename PointDiagonalLocalOperator , typename GridFunctionSpace , typename DomainField , template< typename > class IterativeSolver>
void Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver >::setRequireSetup ( bool  v)
inline

Member Data Documentation

◆ doAlphaVolume

template<typename BlockDiagonalLocalOperator , typename PointDiagonalLocalOperator , typename GridFunctionSpace , typename DomainField , template< typename > class IterativeSolver>
constexpr bool Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver >::doAlphaVolume = true
staticconstexpr

◆ doPatternVolume

template<typename BlockDiagonalLocalOperator , typename PointDiagonalLocalOperator , typename GridFunctionSpace , typename DomainField , template< typename > class IterativeSolver>
constexpr bool Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver >::doPatternVolume = true
staticconstexpr

◆ isLinear

template<typename BlockDiagonalLocalOperator , typename PointDiagonalLocalOperator , typename GridFunctionSpace , typename DomainField , template< typename > class IterativeSolver>
constexpr bool Dune::PDELab::IterativeBlockJacobiPreconditionerLocalOperator< BlockDiagonalLocalOperator, PointDiagonalLocalOperator, GridFunctionSpace, DomainField, IterativeSolver >::isLinear = BlockDiagonalLocalOperator::isLinear
staticconstexpr

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