dune-fem 2.8.0
Loading...
Searching...
No Matches
interpolationquadrature.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_INTERPOLATIONQUADRATURE_HH
2#define DUNE_FEM_INTERPOLATIONQUADRATURE_HH
3
4//- Dune includes
5#include <dune/geometry/type.hh>
6#include <dune/geometry/quadraturerules.hh>
7
10
11namespace Dune
12{
13
14 namespace Fem
15 {
16
17 namespace detail
18 {
24 template< typename FieldImp, int dim, template <class,unsigned int> class PointSet >
25 class InterpolationQuadratureFactory
26 : public Dune::Fem::QuadratureImp< FieldImp, dim >
27 {
28 public:
29 typedef FieldImp FieldType;
30 typedef PointSet< FieldType, dim > PointSetType;
31
32 enum { dimension = dim };
33
34 private:
35 typedef InterpolationQuadratureFactory< FieldType, dimension, PointSet > ThisType;
37
38 protected:
39 typedef typename BaseType :: ElementCoordinateType ElementCoordinateType;
40 using BaseType :: addQuadraturePoint;
41
42 public:
43 typedef typename BaseType :: CoordinateType CoordinateType;
44
45 protected:
46 const GeometryType elementGeometry_;
47 mutable size_t numElementInterpolPoints_;
48 int order_;
49
50 public:
57 InterpolationQuadratureFactory( const GeometryType &geometry,
58 const int order,
59 const size_t id )
60 : BaseType( id ),
61 elementGeometry_( geometry ),
62 numElementInterpolPoints_( 0 ) // this is only set when interpolationPoints are requested
63 {
64 // revert quadrature order to polynomial order
65 if( geometry.isCube() )
66 {
67 auto points = PointSetType::buildCubeQuadrature( order );
68 order_ = points.quadOrder();
69
70 assert(order_ >= order);
71
72 for( unsigned int i=0; i<points.size(); ++i )
73 {
74 addQuadraturePoint( points[ i ].point(), points[ i ].weight() );
75 }
76 }
77 else
78 {
79 DUNE_THROW(InvalidStateException,"InterpolationQuadratureFactory: Unsupported geometry type " << geometry);
80 }
81 }
82
85 int order () const
86 {
87 return order_;
88 }
89
92 virtual std::vector< ElementCoordinateType > interpolationPoints(const int reqDim ) const
93 {
94 // if requested dimension matches potential element dimension
95 if( reqDim == (dim+1) && dim < 3 )
96 {
97 typedef PointSet< FieldType, dim+1 > ElementPointSet;
98 auto points = ElementPointSet::buildCubeQuadrature( order_ );
99 numElementInterpolPoints_ = points.size();
100 std::vector< ElementCoordinateType > pts( numElementInterpolPoints_ );
101 for( size_t i=0; i<numElementInterpolPoints_; ++i )
102 pts[ i ] = points[ i ].point();
103 return pts;
104 }
105 else
106 return std::vector< ElementCoordinateType >();
107 }
108
110 virtual bool isFaceInterpolationQuadrature( const size_t numShapeFunctions ) const
111 {
112 // when numElementInterpolPoints_ is set this means we have a face
113 // quadrature and then this is an interpolation quadrature if the
114 // number of shape functions matches the number of element
115 // interpolation points
116 return numShapeFunctions == numElementInterpolPoints_;
117 }
118
121 GeometryType geometryType () const
122 {
123 return elementGeometry_;
124 }
125
128 static unsigned int maxOrder ()
129 {
130 return QuadratureRules<FieldType,dim>::
131 maxOrder( Dune::GeometryTypes::cube(dim), Dune::QuadratureType::Enum(PointSetType::pointSetId) );
132 }
133 };
134
135 template< class FieldType, int dim, template <class,unsigned int> class PointSet >
136 struct InterpolationQuadratureTraitsImpl
137 {
138 static const int pointSetId = PointSet<FieldType,dim>::pointSetId;
139
140 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > SimplexQuadratureType;
141 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > CubeQuadratureType;
142
143 typedef Dune::Fem::QuadratureImp< FieldType, dim > IntegrationPointListType;
144
145 typedef int QuadratureKeyType ;
146 };
147
148 template< class FieldType, template <class,unsigned int> class PointSet >
149 struct InterpolationQuadratureTraitsImpl< FieldType, 0, PointSet >
150 {
151 static const int dim = 0;
152
153 static const int pointSetId = PointSet<FieldType,dim>::pointSetId;
154
155 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > PointQuadratureType;
156
157 typedef Dune::Fem::QuadratureImp< FieldType, dim > IntegrationPointListType;
158
159 typedef int QuadratureKeyType ;
160 };
161
162 template< class FieldType, template <class,unsigned int> class PointSet >
163 struct InterpolationQuadratureTraitsImpl< FieldType, 1, PointSet >
164 {
165 static const int dim = 1;
166
167 static const int pointSetId = PointSet<FieldType,dim>::pointSetId;
168
169 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > LineQuadratureType;
170
171 typedef Dune::Fem::QuadratureImp< FieldType, dim > IntegrationPointListType;
172
173 typedef int QuadratureKeyType ;
174 };
175
176 template< class FieldType, template <class,unsigned int> class PointSet >
177 struct InterpolationQuadratureTraitsImpl< FieldType, 3, PointSet >
178 {
179 static const int dim = 3;
180
181 static const int pointSetId = PointSet<FieldType,dim>::pointSetId;
182
183 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > SimplexQuadratureType;
184 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > CubeQuadratureType;
185
186 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > PrismQuadratureType;
187 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > PyramidQuadratureType;
188
189 typedef Dune::Fem::QuadratureImp< FieldType, dim > IntegrationPointListType;
190
191 typedef int QuadratureKeyType ;
192 };
193 } // end namespace detail
194
195#if HAVE_DUNE_LOCALFUNCTIONS
197 //
198 // GaussLobatto point set (same quadrature as in dune-geometry
199 // but different point ordering)
200 //
202 template <class FieldType, int dim >
203 using GaussLobattoQuadratureTraits = detail::InterpolationQuadratureTraitsImpl< FieldType, dim, GaussLobattoPointSet >;
204
206 //
207 // GaussLegendre point set (same quadrature as in dune-geometry
208 // but different point ordering)
209 //
211 template <class FieldType, int dim >
212 using GaussLegendreQuadratureTraits = detail::InterpolationQuadratureTraitsImpl< FieldType, dim, GaussLegendrePointSet >;
213#endif
214
215 } // namespace Fem
216
217} // namespace Dune
218
219#endif // #ifndef DUNE_FEM_LAGRANGEPOINTS_HH
Definition: bindguard.hh:11
virtual std::vector< ElementCoordinateType > interpolationPoints(const int reqDim) const
returns list of element interpolation points for a given face quadrature
Definition: quadratureimp.hh:152
virtual bool isFaceInterpolationQuadrature(const size_t numShapeFunctions) const
return true if quadrature is also a set of interpolation points for a given number of shape functions
Definition: quadratureimp.hh:160
size_t id() const
obtain the identifier of the integration point list
Definition: quadratureimp.hh:122
const CoordinateType & point(size_t i) const
obtain coordinates of i-th integration point
Definition: quadratureimp.hh:96
static const int dimension
dimension of quadrature
Definition: quadratureimp.hh:52
virtual GeometryType geometryType() const =0
obtain GeometryType for this integration point list
FieldVector< FieldType, dim+1 > ElementCoordinateType
Definition: quadratureimp.hh:43
Generic implementation of a Dune quadrature.
Definition: quadratureimp.hh:196
void addQuadraturePoint(const CoordinateType &point, const FieldType weight)
Adds a point-weight pair to the quadrature.
Definition: quadratureimp.hh:270
BaseType::CoordinateType CoordinateType
type of local coordinates
Definition: quadratureimp.hh:207
FieldImp FieldType
field type
Definition: quadratureimp.hh:199
const FieldType & weight(size_t i) const
obtain weight of i-th integration point
Definition: quadratureimp.hh:251