3 #ifndef DUNE_GEOMETRY_REFERENCEELEMENTS_HH 4 #define DUNE_GEOMETRY_REFERENCEELEMENTS_HH 9 #include <dune/common/array.hh> 10 #include <dune/common/forloop.hh> 11 #include <dune/common/nullptr.hh> 12 #include <dune/common/typetraits.hh> 13 #include <dune/common/visibility.hh> 14 #include <dune/common/unused.hh> 27 template<
class ctype,
int dim >
30 template<
class ctype,
int dim >
53 template<
class ctype,
int dim >
67 template<
int codim >
struct CreateGeometries;
84 assert( (c >= 0) && (c <= dim) );
85 return info_[ c ].size();
99 int size (
int i,
int c,
int cc )
const 101 assert( (i >= 0) && (i <
size( c )) );
102 return info_[ c ][ i ].size( cc );
120 assert( (i >= 0) && (i <
size( c )) );
121 return info_[ c ][ i ].number( ii, cc );
134 assert( (i >= 0) && (i <
size( c )) );
135 return info_[ c ][ i ].type();
150 const FieldVector< ctype, dim > &
position(
int i,
int c )
const 152 assert( (c >= 0) && (c <= dim) );
153 return baryCenters_[ c ][ i ];
165 const ctype tolerance = ctype( 64 ) * std::numeric_limits< ctype >::epsilon();
166 return GenericGeometry::template checkInside< ctype, dim >(
type().
id(), dim, local, tolerance );
182 template<
int codim >
184 DUNE_DEPRECATED_MSG(
"Method will be removed after dune-geometry-2.4")
185 checkInside ( const FieldVector< ctype, dim-codim > &local,
int i )
const 187 const ctype tolerance = ctype( 64 ) * std::numeric_limits< ctype >::epsilon();
188 return GenericGeometry::template
checkInside< ctype, dim-codim >(
type( i, codim ).
id(), dim-codim, local, tolerance );
202 template<
int codim >
205 integral_constant< int, codim > codimVariable;
206 return geometries_[ codimVariable ][ i ];
224 assert( (face >= 0) && (face <
int( integrationNormals_.size() )) );
225 return integrationNormals_[ face ];
229 void initialize (
unsigned int topologyId )
234 for(
int codim = 0; codim <= dim; ++codim )
237 info_[ codim ].resize( size );
238 for(
unsigned int i = 0; i <
size; ++i )
239 info_[ codim ][ i ].initialize( topologyId, codim, i );
243 const unsigned int numVertices =
size( dim );
244 baryCenters_[ dim ].resize( numVertices );
248 for(
int codim = 0; codim < dim; ++codim )
250 baryCenters_[ codim ].resize(
size(codim) );
251 for(
int i = 0; i <
size( codim ); ++i )
253 baryCenters_[ codim ][ i ] = FieldVector< ctype, dim >( ctype( 0 ) );
254 const unsigned int numCorners =
size( i, codim, dim );
255 for(
unsigned int j = 0; j < numCorners; ++j )
256 baryCenters_[ codim ][ i ] += baryCenters_[ dim ][
subEntity( i, codim, j, dim ) ];
257 baryCenters_[ codim ][ i ] *= ctype( 1 ) / ctype( numCorners );
262 volume_ = GenericGeometry::template referenceVolume< ctype >( topologyId, dim );
267 integrationNormals_.resize(
size( 1 ) );
272 Dune::ForLoop< CreateGeometries, 0, dim >::apply( *
this, geometries_ );
276 template<
int codim >
278 :
public std::vector< typename Codim< codim >::Geometry >
287 std::vector< FieldVector< ctype, dim > > baryCenters_[ dim+1 ];
288 std::vector< FieldVector< ctype, dim > > integrationNormals_;
291 GeometryTable geometries_;
293 std::vector< SubEntityInfo > info_[ dim+1 ];
297 template<
class ctype,
int dim >
301 : numbering_( nullptr )
303 std::fill( offset_.begin(), offset_.end(), 0 );
307 : offset_( other.offset_ ),
310 numbering_ = allocate();
311 std::copy( other.numbering_, other.numbering_ + capacity(), numbering_ );
319 offset_ = other.offset_;
321 deallocate( numbering_ );
322 numbering_ = allocate();
323 std::copy( other.numbering_, other.numbering_ + capacity(), numbering_ );
330 assert( (cc >= codim()) && (cc <= dim) );
331 return (offset_[ cc+1 ] - offset_[ cc ]);
336 assert( (ii >= 0) && (ii <
size( cc )) );
337 return numbering_[ offset_[ cc ] + ii ];
342 void initialize (
unsigned int topologyId,
int codim,
unsigned int i )
348 for(
int cc = 0; cc <= codim; ++cc )
350 for(
int cc = codim; cc <= dim; ++cc )
354 deallocate( numbering_ );
355 numbering_ = allocate();
356 for(
int cc = codim; cc <= dim; ++cc )
363 unsigned int *
allocate () {
return (capacity() != 0 ?
new unsigned int[ capacity() ] :
nullptr); }
365 unsigned int capacity ()
const {
return offset_[ dim+1 ]; }
368 unsigned int *numbering_;
369 array< unsigned int, dim+2 > offset_;
374 template<
class ctype,
int dim >
375 template<
int codim >
388 DUNE_UNUSED_PARAMETER(i);
394 const int size = refElement.
size( codim );
395 std::vector< FieldVector< ctype, dim > > origins( size );
396 std::vector< FieldMatrix< ctype, dim - codim, dim > > jacobianTransposeds( size );
399 integral_constant< int, codim > codimVariable;
400 geometries[ codimVariable ].reserve( size );
401 for(
int i = 0; i <
size; ++i )
404 geometries[ codimVariable ].push_back( geometry );
414 template<
class ctype,
int dim >
425 for(
unsigned int topologyId = 0; topologyId <
numTopologies; ++topologyId )
426 values_[ topologyId ].initialize( topologyId );
431 assert( type.
dim() == dim );
432 return values_[ type.
id() ];
440 const value_type &
cube ()
const 455 const_iterator
begin ()
const {
return values_; }
477 template<
class ctype,
int dim >
486 return container() (
type );
492 return container().simplex();
498 return container().cube();
501 static Iterator
begin () {
return container().begin(); }
502 static Iterator
end () {
return container().end(); }
514 #endif // #ifndef DUNE_GEOMETRY_REFERENCEELEMENTS_HH const value_type & pyramid() const
Definition: referenceelements.hh:445
unsigned int subTopologyId(unsigned int topologyId, int dim, int codim, unsigned int i)
Compute the topology id of a given subentity.
Definition: subtopologies.cc:47
Definition: topologytypes.hh:231
ReferenceElementContainer< ctype, dim >::const_iterator Iterator
Definition: referenceelements.hh:480
AffineGeometry< ctype, dim-codim, dim > Geometry
type of geometry embedding a subentity into the reference element
Definition: referenceelements.hh:75
~SubEntityInfo()
Definition: referenceelements.hh:314
Definition: referenceelements.hh:28
const value_type & prism() const
Definition: referenceelements.hh:450
const value_type & simplex() const
Definition: referenceelements.hh:435
Definition: topologytypes.hh:248
unsigned int capacity() const
Definition: referenceelements.hh:365
static const ReferenceElement< ctype, dim > & simplex()
get simplex reference elements
Definition: referenceelements.hh:490
int codim() const
Definition: referenceelements.hh:361
SubEntityInfo()
Definition: referenceelements.hh:300
unsigned int size(unsigned int topologyId, int dim, int codim)
Compute the number of subentities of a given codimension.
Definition: subtopologies.cc:16
static const ReferenceElement< ctype, dim > & general(const GeometryType &type)
get general reference elements
Definition: referenceelements.hh:484
static const ReferenceElement< ctype, dim > & cube()
get hypercube reference elements
Definition: referenceelements.hh:496
unsigned int id() const
Return the topology id the type.
Definition: type.hh:326
An implementation of the Geometry interface for affine geometries.
const_iterator begin() const
Definition: referenceelements.hh:455
Definition: topologytypes.hh:259
unsigned int referenceCorners(unsigned int topologyId, int dim, FieldVector< ct, cdim > *corners)
Definition: referencedomain.hh:308
int size(int i, int c, int cc) const
number of subentities of codimension cc of subentity (i,c)
Definition: referenceelements.hh:99
bool checkInside(const FieldVector< ctype, dim > &local) const
check if a coordinate is in the reference element
Definition: referenceelements.hh:163
const_iterator end() const
Definition: referenceelements.hh:456
int size(int c) const
number of subentities of codimension c
Definition: referenceelements.hh:82
SubEntityInfo(const SubEntityInfo &other)
Definition: referenceelements.hh:306
const FieldVector< ctype, dim > & integrationOuterNormal(int face) const
obtain the integration outer normal of the reference element
Definition: referenceelements.hh:222
unsigned int * allocate()
Definition: referenceelements.hh:363
void initialize(unsigned int topologyId, int codim, unsigned int i)
Definition: referenceelements.hh:342
const GeometryType & type() const
Definition: referenceelements.hh:340
const value_type & cube() const
Definition: referenceelements.hh:440
unsigned int numTopologies(int dim)
obtain the number of topologies of a given dimension
Definition: topologytypes.hh:134
ReferenceElementContainer()
Definition: referenceelements.hh:423
Codim< codim >::Geometry geometry(int i) const
obtain the embedding of subentity (i,codim) into the reference element
Definition: referenceelements.hh:203
int subEntity(int i, int c, int ii, int cc) const
obtain number of ii-th subentity with codim cc of (i,c)
Definition: referenceelements.hh:118
void subTopologyNumbering(unsigned int topologyId, int dim, int codim, unsigned int i, int subcodim, unsigned int *beginOut, unsigned int *endOut)
Definition: subtopologies.cc:85
void deallocate(unsigned int *ptr)
Definition: referenceelements.hh:364
ctype volume() const
obtain the volume of the reference element
Definition: referenceelements.hh:210
Unique label for each type of entities that can occur in DUNE grids.
Definition: type.hh:24
unsigned int dim() const
Return dimension of the type.
Definition: type.hh:321
static Iterator begin()
Definition: referenceelements.hh:501
Collection of types depending on the codimension.
Definition: referenceelements.hh:72
topological information about the subentities of a reference element
Definition: referenceelements.hh:298
This class provides access to geometric and topological properties of a reference element...
Definition: affinegeometry.hh:25
const value_type * const_iterator
Definition: referenceelements.hh:421
Definition: topologytypes.hh:214
unsigned int referenceIntegrationOuterNormals(unsigned int topologyId, int dim, const FieldVector< ct, cdim > *origins, FieldVector< ct, cdim > *normals)
Definition: referencedomain.hh:471
const FieldVector< ctype, dim > & position(int i, int c) const
position of the barycenter of entity (i,c)
Definition: referenceelements.hh:150
unsigned int referenceEmbeddings(unsigned int topologyId, int dim, int codim, FieldVector< ct, cdim > *origins, FieldMatrix< ct, mydim, cdim > *jacobianTransposeds)
Definition: referencedomain.hh:406
ReferenceElement< ctype, dim > value_type
Definition: referenceelements.hh:420
static Iterator end()
Definition: referenceelements.hh:502
const GeometryType & type() const
obtain the type of this reference element
Definition: referenceelements.hh:139
const GeometryType & type(int i, int c) const
obtain the type of subentity (i,c)
Definition: referenceelements.hh:132
Definition: affinegeometry.hh:18
Class providing access to the singletons of the reference elements.
Definition: affinegeometry.hh:28
int number(int ii, int cc) const
Definition: referenceelements.hh:334
int size(int cc) const
Definition: referenceelements.hh:328
Implementation of the Geometry interface for affine geometries.
Definition: affinegeometry.hh:38