1#ifndef DUNE_FEM_AGGLOMERATIONQUADRATURE_HH
2#define DUNE_FEM_AGGLOMERATIONQUADRATURE_HH
5#include <dune/geometry/axisalignedcubegeometry.hh>
20 template<
typename Gr
idPartImp,
class IntegrationPo
intList >
37 typedef typename GridPartType :: ctype
RealType;
46 typedef typename GridPartType::template Codim< 0 >::EntityType
EntityType;
83 if( storage.size() < 20 )
85 storage.emplace( polyQuad );
106 template <
int dimw = dimensionworld >
107 static std::enable_if_t<dimension == 2 && dimw == 2, IntegrationPointListType>
111 Dune::GeometryType simplexType = Dune::GeometryTypes::simplex(
dimension );
113 typedef AxisAlignedCubeGeometry< RealType, dimension, dimension> CubeGeometryType;
115 const auto& elemGeo = entity.geometry();
120 const int corners = elemGeo.corners();
121 for(
int c = 1; c < corners; ++c )
123 const auto& corner = elemGeo.corner( c );
126 lower[ d ] =
std::min( lower[ d ], corner[ d ]);
127 upper[ d ] =
std::max( upper[ d ], corner[ d ]);
133 CubeGeometryType cubeGeom( lower, upper );
135 QuadratureType quad( simplexType, quadKey );
137 const int subEntities = entity.subEntities( 1 );
138 const int quadNop = quad.nop();
139 int order = quad.order();
143 quadImp.
reset( order, subEntities * quadNop );
145 Dune::FieldMatrix<double,dimension,dimension> A( 0 ) ;
146 const auto& center = entity.geometry().center();
148 for(
int i = 0; i<subEntities; ++i )
150 const auto subEntity = entity.template subEntity< 1 >( i );
151 const auto& geom = subEntity.geometry();
156 A[ c ] = geom.corner( c );
163 double vol =
std::abs(A.determinant()) / elemGeo.volume();
166 for(
int qp = 0; qp < quadNop; ++qp )
169 A.mtv( quad.point( qp ), point );
172 point = cubeGeom.local( point );
175 double weight = quad.weight( qp ) * vol;
182 std::shared_ptr< const IntegrationPointListImpl > quadPtr( &quadImp,
Deleter() );
186 template <
int dimw = dimensionworld >
187 static std::enable_if_t<dimension != 2 || dimw != 2, IntegrationPointListType>
Dune::Fem::Double abs(const Dune::Fem::Double &a)
Definition: double.hh:942
double max(const Dune::Fem::Double &v, const double p)
Definition: double.hh:965
double min(const Dune::Fem::Double &v, const double p)
Definition: double.hh:953
Definition: bindguard.hh:11
ThreadSafeValue realizes thread safety for a given variable by creating an instance of this variable ...
Definition: threadsafevalue.hh:18
Agglomeration is a simple quadrature for polyhedral cells based on sub triangulation
Definition: agglomerationquadrature.hh:22
std::stack< std::unique_ptr< PolyhedronQuadratureType > > PolyhedronQuadratureStorageType
Definition: agglomerationquadrature.hh:52
PolyhedronQuadrature< RealType, dimension > PolyhedronQuadratureType
Definition: agglomerationquadrature.hh:49
static std::enable_if_t< dimension !=2||dimw !=2, IntegrationPointListType > computeQuadrature(const EntityType &entity, const QuadratureKeyType &quadKey)
Definition: agglomerationquadrature.hh:188
IntegrationPointListType::CoordinateType CoordinateType
Definition: agglomerationquadrature.hh:41
IntegrationPointListType::IntegrationPointListType IntegrationPointListImpl
Definition: agglomerationquadrature.hh:50
static void pushObject(const PolyhedronQuadratureType *quad)
Definition: agglomerationquadrature.hh:79
@ dimensionworld
Definition: agglomerationquadrature.hh:34
GridPartType::ctype RealType
type for reals (usually double)
Definition: agglomerationquadrature.hh:37
GridPartImp GridPartType
type of the grid partition
Definition: agglomerationquadrature.hh:25
@ codimension
Definition: agglomerationquadrature.hh:28
static std::enable_if_t< dimension==2 &&dimw==2, IntegrationPointListType > computeQuadrature(const EntityType &entity, const QuadratureKeyType &quadKey)
returns quadrature points for polyhedral cells
Definition: agglomerationquadrature.hh:108
static PolyhedronQuadratureStorageType & quadStorage()
Definition: agglomerationquadrature.hh:55
static PolyhedronQuadratureType * getObject(const GeometryType &type)
Definition: agglomerationquadrature.hh:62
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: agglomerationquadrature.hh:46
@ dimension
Definition: agglomerationquadrature.hh:31
int QuadratureKeyType
Definition: agglomerationquadrature.hh:43
IntegrationPointList IntegrationPointListType
Definition: agglomerationquadrature.hh:39
Definition: agglomerationquadrature.hh:95
void operator()(const IntegrationPointListImpl *quad)
Definition: agglomerationquadrature.hh:96
quadrature on the codim-0 reference element
Definition: elementquadrature.hh:58
Definition: femquadratures.hh:287
void reset(const int order, const int nop)
Definition: femquadratures.hh:338
void addQuadraturePoint(const CoordinateType &point, const FieldType weight)
Adds a point-weight pair to the quadrature.
Definition: quadratureimp.hh:270
actual interface class for integration point lists
Definition: quadrature.hh:161