1 #ifndef DUNE_FEM_SPACE_DISCONTINUOUSGALERKIN_LOCALINTERPOLATION_HH 2 #define DUNE_FEM_SPACE_DISCONTINUOUSGALERKIN_LOCALINTERPOLATION_HH 26 template<
class DiscreteFunctionSpace,
template<
class,
int >
class Quadrature = CachingQuadrature >
33 typedef typename DiscreteFunctionSpaceType::EntityType
EntityType;
36 typedef typename DiscreteFunctionSpaceType::RangeType RangeType;
37 typedef typename RangeType::value_type RangeFieldType;
39 typedef typename DiscreteFunctionSpaceType::GridPartType GridPartType;
46 : order_( space.order() ),
47 massMatrix_( space, (order < 0 ? 2*space.order() : order) )
54 ThisType &operator= (
const ThisType &other );
57 template<
class LocalFunction,
class LocalDofVector >
64 const EntityType &entity = localFunction.
entity();
65 const typename EntityType::Geometry geometry = entity.geometry();
67 QuadratureType quadrature( entity, localFunction.
order() + order_ );
68 const int nop = quadrature.nop();
69 for(
int qp = 0; qp < nop; ++qp )
73 localFunction.
evaluate( quadrature[ qp ], value );
76 const RangeFieldType intel = quadrature.weight( qp )*geometry.integrationElement( quadrature.point( qp ) );
79 dofs.axpy( quadrature[ qp ], value );
87 const LocalMassMatrixType &massMatrix ()
const {
return massMatrix_; }
90 LocalMassMatrixType massMatrix_;
97 #endif // #ifndef DUNE_FEM_SPACE_DISCONTINUOUSGALERKIN_LOCALINTERPOLATION_HH
void applyInverse(MassCaller &caller, const EntityType &entity, LocalFunction &lf) const
Definition: localmassmatrix.hh:231
interface for local functions
Definition: localfunction.hh:41
Definition: coordinate.hh:4
const EntityType & entity() const
obtain the entity, this local function lives on
Definition: localfunction.hh:285
void operator()(const LocalFunction &localFunction, LocalDofVector &dofs) const
Definition: discontinuousgalerkin/localinterpolation.hh:58
DiscreteFunctionSpaceType::EntityType EntityType
Definition: discontinuousgalerkin/localinterpolation.hh:33
void evaluate(const PointType &x, RangeType &ret) const
evaluate the local function
Definition: localfunction.hh:300
Definition: discontinuousgalerkin/localinterpolation.hh:27
DiscreteFunctionSpace DiscreteFunctionSpaceType
Definition: discontinuousgalerkin/localinterpolation.hh:32
DiscontinuousGalerkinLocalInterpolation(const DiscreteFunctionSpaceType &space, const int order=-1)
Definition: discontinuousgalerkin/localinterpolation.hh:45
int order() const
obtain the order of this local function
Definition: localfunction.hh:273
actual interface class for quadratures
Definition: quadrature.hh:320