dune-fem  2.4.1-rc
discontinuousgalerkin/localinterpolation.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_SPACE_DISCONTINUOUSGALERKIN_LOCALINTERPOLATION_HH
2 #define DUNE_FEM_SPACE_DISCONTINUOUSGALERKIN_LOCALINTERPOLATION_HH
3 
4 // dune-fem includes
7 
14 namespace Dune
15 {
16 
17  namespace Fem
18  {
19 
20  // DiscontinuousGalerkinLocalInterpolation
21  // ---------------------------------------
22 
26  template< class DiscreteFunctionSpace, template< class, int > class Quadrature = CachingQuadrature >
28  {
30 
31  public:
33  typedef typename DiscreteFunctionSpaceType::EntityType EntityType;
34 
35  private:
36  typedef typename DiscreteFunctionSpaceType::RangeType RangeType;
37  typedef typename RangeType::value_type RangeFieldType;
38 
39  typedef typename DiscreteFunctionSpaceType::GridPartType GridPartType;
41 
43 
44  public:
45  DiscontinuousGalerkinLocalInterpolation ( const DiscreteFunctionSpaceType &space, const int order = -1 )
46  : order_( space.order() ),
47  massMatrix_( space, (order < 0 ? 2*space.order() : order) )
48  {}
49 
50  private:
51  // forbid copying
52  DiscontinuousGalerkinLocalInterpolation ( const ThisType &other );
53  // forbid assignment
54  ThisType &operator= ( const ThisType &other );
55 
56  public:
57  template< class LocalFunction, class LocalDofVector >
58  void operator () ( const LocalFunction &localFunction, LocalDofVector &dofs ) const
59  {
60  // set all dofs to zero
61  dofs.clear();
62 
63  // get entity and geometry
64  const EntityType &entity = localFunction.entity();
65  const typename EntityType::Geometry geometry = entity.geometry();
66 
67  QuadratureType quadrature( entity, localFunction.order() + order_ );
68  const int nop = quadrature.nop();
69  for( int qp = 0; qp < nop; ++qp )
70  {
71  // evaluate local function
72  RangeType value;
73  localFunction.evaluate( quadrature[ qp ], value );
74 
75  // compute quadrature weight
76  const RangeFieldType intel = quadrature.weight( qp )*geometry.integrationElement( quadrature.point( qp ) );
77  value *= intel;
78 
79  dofs.axpy( quadrature[ qp ], value );
80  }
81 
82  // apply inverse of mass matrix
83  massMatrix().applyInverse( entity, dofs );
84  }
85 
86  private:
87  const LocalMassMatrixType &massMatrix () const { return massMatrix_; }
88 
89  const int order_;
90  LocalMassMatrixType massMatrix_;
91  };
92 
93  } // namespace Fem
94 
95 } // namespace Dune
96 
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
discrete function space
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