1#ifndef DUNE_FEM_SPACE_DISCONTINUOUSGALERKIN_LOCALINTERPOLATION_HH
2#define DUNE_FEM_SPACE_DISCONTINUOUSGALERKIN_LOCALINTERPOLATION_HH
5#include <dune/grid/common/capabilities.hh>
28 template<
class DiscreteFunctionSpace >
35 typedef typename DiscreteFunctionSpaceType::GridType
GridType;
36 typedef typename DiscreteFunctionSpaceType::EntityType
EntityType;
38 static const bool isAlwaysAffine = Dune::Capabilities::isCartesian< GridType >::v ||
39 ( Dune::Capabilities::hasSingleGeometryType< GridType >::v && ((Dune::Capabilities::hasSingleGeometryType< GridType >::topologyId >> 1) == 0)) ;
44 typedef typename DiscreteFunctionSpaceType::RangeType RangeType;
45 typedef typename RangeType::value_type RangeFieldType;
47 typedef typename DiscreteFunctionSpaceType::GridPartType GridPartType;
49 typedef typename DiscreteFunctionSpaceType :: LocalMassMatrixType
52 typedef typename LocalMassMatrixType::VolumeQuadratureType QuadratureType;
57 massStorage_( space_.localMassMatrixStorage() )
68 template<
class LocalFunction,
class LocalDofVector >
72 std::fill(
dofs.begin(),
dofs.end(),
typename LocalDofVector::value_type(0) );
77 if( entity.type().isNone() )
80 ElementQuadratureType quadrature( entity,
localFunction.order() + space_.order( entity ) );
85 AggloMassMatrix massMat( space_, massMatrix().volumeQuadratureOrder( entity ) );
87 auto basisFunctionSet = space_.basisFunctionSet(entity);
88 massMat.applyInverse( entity, basisFunctionSet,
dofs );
93 QuadratureType quadrature( entity,
localFunction.order() + space_.order( entity ) );
98 auto basisFunctionSet = space_.basisFunctionSet(entity);
99 massMatrix().applyInverse( entity, basisFunctionSet,
dofs );
106 template<
class QuadImpl,
class LocalFunction,
class LocalDofVector >
107 bool computeInterpolation(
const EntityType& entity,
108 const QuadImpl& quadrature,
110 LocalDofVector &
dofs )
const
112 const int nop = quadrature.nop();
114 auto& values = massStorage_.second;
116 values.resize( nop );
124 const auto geometry = entity.geometry();
125 isAffine = geometry.affine();
130 for(
auto qp : quadrature )
131 values[ qp.index() ] *= qp.weight() * geometry.integrationElement( qp.position() );
138 for(
auto qp : quadrature )
139 values[ qp.index() ] *= qp.weight();
143 space_.basisFunctionSet(entity).axpy( quadrature, values,
dofs );
148 const LocalMassMatrixType &massMatrix ()
const
150 return massStorage_.first;
154 typename DiscreteFunctionSpaceType::LocalMassMatrixStorageType& massStorage_;
Definition: bindguard.hh:11
static GridFunctionView< GF > localFunction(const GF &gf)
Definition: gridfunctionview.hh:118
IteratorRange< typename DF::DofIteratorType > dofs(DF &df)
Iterates over all DOFs.
Definition: rangegenerators.hh:76
interface for local functions
Definition: localfunction.hh:77
Local Mass Matrix for arbitrary spaces.
Definition: localmassmatrix.hh:909
quadrature on the codim-0 reference element
Definition: elementquadrature.hh:18
Definition: discontinuousgalerkin/localinterpolation.hh:30
void operator()(const LocalFunction &localFunction, LocalDofVector &dofs) const
Definition: discontinuousgalerkin/localinterpolation.hh:69
void unbind()
Definition: discontinuousgalerkin/localinterpolation.hh:64
void bind(const EntityType &)
Definition: discontinuousgalerkin/localinterpolation.hh:63
DiscontinuousGalerkinLocalInterpolation(const ThisType &other)=default
DiscreteFunctionSpaceType::EntityType EntityType
Definition: discontinuousgalerkin/localinterpolation.hh:36
DiscreteFunctionSpace DiscreteFunctionSpaceType
Definition: discontinuousgalerkin/localinterpolation.hh:34
static const bool isAlwaysAffine
Definition: discontinuousgalerkin/localinterpolation.hh:38
DiscreteFunctionSpaceType::GridType GridType
Definition: discontinuousgalerkin/localinterpolation.hh:35
DiscontinuousGalerkinLocalInterpolation(const DiscreteFunctionSpaceType &space)
Definition: discontinuousgalerkin/localinterpolation.hh:55
ThisType & operator=(const ThisType &other)=delete
DiscontinuousGalerkinLocalInterpolation(ThisType &&other)=default