1#ifndef DUNE_FEM_INTERPOLATIONQUADRATURE_HH
2#define DUNE_FEM_INTERPOLATIONQUADRATURE_HH
5#include <dune/geometry/type.hh>
6#include <dune/geometry/quadraturerules.hh>
24 template<
typename FieldImp,
int dim,
template <
class,
unsigned int>
class PointSet >
25 class InterpolationQuadratureFactory
30 typedef PointSet< FieldType, dim > PointSetType;
35 typedef InterpolationQuadratureFactory< FieldType, dimension, PointSet > ThisType;
46 const GeometryType elementGeometry_;
47 mutable size_t numElementInterpolPoints_;
57 InterpolationQuadratureFactory(
const GeometryType &geometry,
61 elementGeometry_( geometry ),
62 numElementInterpolPoints_( 0 )
65 if( geometry.isCube() )
67 auto points = PointSetType::buildCubeQuadrature( order );
68 order_ = points.quadOrder();
70 assert(order_ >= order);
72 for(
unsigned int i=0; i<points.size(); ++i )
79 DUNE_THROW(InvalidStateException,
"InterpolationQuadratureFactory: Unsupported geometry type " << geometry);
95 if( reqDim == (dim+1) && dim < 3 )
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();
106 return std::vector< ElementCoordinateType >();
116 return numShapeFunctions == numElementInterpolPoints_;
123 return elementGeometry_;
128 static unsigned int maxOrder ()
130 return QuadratureRules<FieldType,dim>::
131 maxOrder( Dune::GeometryTypes::cube(dim), Dune::QuadratureType::Enum(PointSetType::pointSetId) );
135 template<
class FieldType,
int dim,
template <
class,
unsigned int>
class PointSet >
136 struct InterpolationQuadratureTraitsImpl
138 static const int pointSetId = PointSet<FieldType,dim>::pointSetId;
140 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > SimplexQuadratureType;
141 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > CubeQuadratureType;
145 typedef int QuadratureKeyType ;
148 template<
class FieldType,
template <
class,
unsigned int>
class PointSet >
149 struct InterpolationQuadratureTraitsImpl< FieldType, 0, PointSet >
151 static const int dim = 0;
153 static const int pointSetId = PointSet<FieldType,dim>::pointSetId;
155 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > PointQuadratureType;
159 typedef int QuadratureKeyType ;
162 template<
class FieldType,
template <
class,
unsigned int>
class PointSet >
163 struct InterpolationQuadratureTraitsImpl< FieldType, 1, PointSet >
165 static const int dim = 1;
167 static const int pointSetId = PointSet<FieldType,dim>::pointSetId;
169 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > LineQuadratureType;
173 typedef int QuadratureKeyType ;
176 template<
class FieldType,
template <
class,
unsigned int>
class PointSet >
177 struct InterpolationQuadratureTraitsImpl< FieldType, 3, PointSet >
179 static const int dim = 3;
181 static const int pointSetId = PointSet<FieldType,dim>::pointSetId;
183 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > SimplexQuadratureType;
184 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > CubeQuadratureType;
186 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > PrismQuadratureType;
187 typedef InterpolationQuadratureFactory< FieldType, dim, PointSet > PyramidQuadratureType;
191 typedef int QuadratureKeyType ;
195#if HAVE_DUNE_LOCALFUNCTIONS
202 template <
class FieldType,
int dim >
203 using GaussLobattoQuadratureTraits = detail::InterpolationQuadratureTraitsImpl< FieldType, dim, GaussLobattoPointSet >;
211 template <
class FieldType,
int dim >
212 using GaussLegendreQuadratureTraits = detail::InterpolationQuadratureTraitsImpl< FieldType, dim, GaussLegendrePointSet >;
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