dune-fem 2.8-git
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
5#include <dune/grid/common/capabilities.hh>
9
16namespace Dune
17{
18
19 namespace Fem
20 {
21
22 // DiscontinuousGalerkinLocalInterpolation
23 // ---------------------------------------
24
28 template< class DiscreteFunctionSpace >
30 {
32
33 public:
35 typedef typename DiscreteFunctionSpaceType::GridType GridType;
36 typedef typename DiscreteFunctionSpaceType::EntityType EntityType;
37
38 static const bool isAlwaysAffine = Dune::Capabilities::isCartesian< GridType >::v ||
39 ( Dune::Capabilities::hasSingleGeometryType< GridType >::v && ((Dune::Capabilities::hasSingleGeometryType< GridType >::topologyId >> 1) == 0)) ;
40 // always true for orthonormal spaces
41 //static const bool isAlwaysAffine = true;
42
43 private:
44 typedef typename DiscreteFunctionSpaceType::RangeType RangeType;
45 typedef typename RangeType::value_type RangeFieldType;
46
47 typedef typename DiscreteFunctionSpaceType::GridPartType GridPartType;
48
49 typedef typename DiscreteFunctionSpaceType :: LocalMassMatrixType
50 LocalMassMatrixType;
51
52 typedef typename LocalMassMatrixType::VolumeQuadratureType QuadratureType;
53
54 public:
56 : space_( space ),
57 massStorage_( space_.localMassMatrixStorage() )
58 {}
59
62
63 void bind( const EntityType& ) {}
64 void unbind() {}
65
66 ThisType &operator= ( const ThisType &other ) = delete;
67
68 template< class LocalFunction, class LocalDofVector >
69 void operator () ( const LocalFunction &localFunction, LocalDofVector &dofs ) const
70 {
71 // set all dofs to zero
72 std::fill( dofs.begin(), dofs.end(), typename LocalDofVector::value_type(0) );
73
74 // get entity and geometry
75 const EntityType &entity = localFunction.entity();
76
77 if( entity.type().isNone() )
78 {
80 ElementQuadratureType quadrature( entity, localFunction.order() + space_.order( entity ) );
81 bool isAffine = computeInterpolation( entity, quadrature, localFunction, dofs );
82 if( ! isAffine )
83 {
85 AggloMassMatrix massMat( space_, massMatrix().volumeQuadratureOrder( entity ) );
86 // apply inverse of mass matrix
87 auto basisFunctionSet = space_.basisFunctionSet(entity);
88 massMat.applyInverse( entity, basisFunctionSet, dofs );
89 }
90 }
91 else
92 {
93 QuadratureType quadrature( entity, localFunction.order() + space_.order( entity ) );
94 bool isAffine = computeInterpolation( entity, quadrature, localFunction, dofs );
95 if( ! isAffine )
96 {
97 // apply inverse of mass matrix
98 auto basisFunctionSet = space_.basisFunctionSet(entity);
99 massMatrix().applyInverse( entity, basisFunctionSet, dofs );
100 }
101 }
102
103 }
104
105 private:
106 template<class QuadImpl, class LocalFunction, class LocalDofVector >
107 bool computeInterpolation( const EntityType& entity,
108 const QuadImpl& quadrature,
110 LocalDofVector &dofs ) const
111 {
112 const int nop = quadrature.nop();
113
114 auto& values = massStorage_.second;
115 // adjust size of values
116 values.resize( nop );
117
118 // evaluate local function for all quadrature points
119 localFunction.evaluateQuadrature( quadrature, values );
120
121 bool isAffine = isAlwaysAffine ;
122 if( ! isAlwaysAffine )
123 {
124 const auto geometry = entity.geometry();
125 isAffine = geometry.affine();
126
127 if( ! isAffine )
128 {
129 // apply weight
130 for(auto qp : quadrature )
131 values[ qp.index() ] *= qp.weight() * geometry.integrationElement( qp.position() );
132 }
133 }
134
135 if( isAffine )
136 {
137 // apply weight only
138 for(auto qp : quadrature )
139 values[ qp.index() ] *= qp.weight();
140 }
141
142 // add values to local function
143 space_.basisFunctionSet(entity).axpy( quadrature, values, dofs );
144
145 return isAffine;
146 }
147
148 const LocalMassMatrixType &massMatrix () const
149 {
150 return massStorage_.first;
151 }
152
153 const DiscreteFunctionSpaceType& space_;
154 typename DiscreteFunctionSpaceType::LocalMassMatrixStorageType& massStorage_;
155 };
156
157 } // namespace Fem
158
159} // namespace Dune
160
161#endif // #ifndef DUNE_FEM_SPACE_DISCONTINUOUSGALERKIN_LOCALINTERPOLATION_HH
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
discrete function space