1#ifndef DUNE_FEM_SPACE_DISCONTINUOUSGALERKIN_LOCALINTERPOLATION_HH
2#define DUNE_FEM_SPACE_DISCONTINUOUSGALERKIN_LOCALINTERPOLATION_HH
6#include <dune/fem/operator/1order/localmassmatrix.hh>
7#include <dune/fem/quadrature/cachingquadrature.hh>
8#include <dune/fem/quadrature/agglomerationquadrature.hh>
28 template<
class DiscreteFunctionSpace >
36 typedef typename DiscreteFunctionSpaceType::EntityType EntityType;
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;
62 void bind(
const EntityType& ) {}
67 template<
class LocalFunction,
class LocalDofVector >
68 void operator () (
const LocalFunction &localFunction, LocalDofVector &dofs )
const
71 std::fill( dofs.begin(), dofs.end(),
typename LocalDofVector::value_type(0) );
74 const EntityType &entity = localFunction.
entity();
76 if( entity.type().isNone() )
78 typedef ElementQuadrature< GridPartType, EntityType::codimension > ElementQuadratureType;
79 ElementQuadratureType quadrature( entity, localFunction.
order() + space_.order( entity ) );
80 bool isAffine = computeInterpolation( entity, quadrature, localFunction, dofs );
84 AggloMassMatrix massMat( space_, massMatrix().volumeQuadratureOrder( entity ) );
86 auto basisFunctionSet = space_.basisFunctionSet(entity);
87 massMat.applyInverse( entity, basisFunctionSet, dofs );
92 QuadratureType quadrature( entity, localFunction.
order() + space_.order( entity ) );
93 bool isAffine = computeInterpolation( entity, quadrature, localFunction, dofs );
97 auto basisFunctionSet = space_.basisFunctionSet(entity);
98 massMatrix().
applyInverse( entity, basisFunctionSet, dofs );
105 template<
class QuadImpl,
class LocalFunction,
class LocalDofVector >
106 bool computeInterpolation(
const EntityType& entity,
107 const QuadImpl& quadrature,
109 LocalDofVector &dofs )
const
111 const int nop = quadrature.nop();
113 auto& values = space_.localMassMatrixStorage().second;
116 values.resize( nop );
121 bool isAffine = isAlwaysAffine ;
122 if( ! isAlwaysAffine )
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 );
150 return space_.localMassMatrixStorage().first;
Definition: localinterpolation.hh:30
GridPartType::GridType GridType
type of underlying dune grid
Definition: discretefunctionspace.hh:214
interface for local functions
Definition: localfunction.hh:77
void evaluateQuadrature(const QuadratureType &quad, Vectors &... vec) const
evaluate all basisfunctions for all quadrature points and store the results in the result vector
Definition: localfunction.hh:393
const EntityType & entity() const
obtain the entity, this local function lives on
Definition: localfunction.hh:305
int order() const
obtain the order of this local function
Definition: localfunction.hh:293
void applyInverse(MassCaller &caller, const EntityType &entity, const BasisFunctionSet &basisFunctionSet, LocalFunction &lf) const
Definition: localmassmatrix.hh:340
Local Mass Matrix for arbitrary spaces.
Definition: localmassmatrix.hh:967
A set of traits classes to store static information about grid implementation.
Dune namespace.
Definition: alignedallocator.hh:13
Specialize with 'true' for if the codimension 0 entity of the grid has only one possible geometry typ...
Definition: capabilities.hh:27
Specialize with 'true' if the grid is a Cartesian grid. Cartesian grids satisfy the following propert...
Definition: capabilities.hh:48