5#ifndef DUNE_ALBERTA_GEOMETRY_HH
6#define DUNE_ALBERTA_GEOMETRY_HH
9#include <dune/grid/albertagrid/misc.hh>
20 template<
int dim,
int dimworld >
28 template<
int codim,
class Gr
idImp >
29 struct AlbertaGridCoordinateReader
31 typedef typename std::remove_const< GridImp >::type Grid;
34 static constexpr int codimension = codim;
35 static constexpr int mydimension = dimension - codimension;
38 typedef Alberta::Real ctype;
40 typedef Alberta::ElementInfo< dimension > ElementInfo;
43 AlbertaGridCoordinateReader (
const GridImp &grid,
44 const ElementInfo &elementInfo,
47 elementInfo_( elementInfo ),
48 subEntity_( subEntity )
51 const ElementInfo &elementInfo ()
const
56 void coordinate (
int i, Coordinate &x )
const
58 assert( !elementInfo_ ==
false );
59 assert( (i >= 0) && (i <= mydimension) );
61 const int k = mapVertices( subEntity_, i );
62 const Alberta::GlobalVector &coord = grid_.getCoord( elementInfo_, k );
63 for(
int j = 0; j < coorddimension; ++j )
67 bool hasDeterminant ()
const
72 ctype determinant ()
const
74 assert( hasDeterminant() );
79 static int mapVertices (
int subEntity,
int i )
81 return Alberta::MapVertices< dimension, codimension >::apply( subEntity, i );
85 const ElementInfo &elementInfo_;
106 template<
int mydim,
int cdim,
class Gr
idImp >
112 typedef GridImp Grid;
115 static constexpr int dimbary = mydim + 1;
122 static constexpr int mydimension = mydim;
123 static constexpr int codimension = dimension - mydimension;
124 static constexpr int coorddimension = cdim;
136 static constexpr int numCorners = mydimension + 1;
146 template<
class CoordReader >
149 build( coordReader );
170 assert( (i >= 0) && (i <
corners()) );
181 GlobalCoordinate
global (
const LocalCoordinate &local )
const;
193 assert( calcedDet_ );
217 const JacobianTransposed &
231 const JacobianInverseTransposed &
262 template<
class CoordReader >
263 void build (
const CoordReader &coordReader );
265 void print ( std::ostream &out )
const;
269 ctype elDeterminant ()
const
278 GlobalCoordinate centroid_;
281 mutable JacobianTransposed jT_;
284 mutable JacobianInverseTransposed jTInv_;
287 mutable bool builtJT_;
289 mutable bool builtJTInv_;
291 mutable bool calcedDet_;
292 mutable ctype elDet_;
300 template<
int mydim,
int cdim,
class Gr
idImp >
301 class AlbertaGridGlobalGeometry
302 :
public AlbertaGridGeometry< mydim, cdim, GridImp >
304 typedef AlbertaGridGlobalGeometry< mydim, cdim, GridImp > This;
305 typedef AlbertaGridGeometry< mydim, cdim, GridImp > Base;
308 AlbertaGridGlobalGeometry ()
312 template<
class CoordReader >
313 AlbertaGridGlobalGeometry (
const CoordReader &coordReader )
314 : Base( coordReader )
319#if !DUNE_ALBERTA_CACHE_COORDINATES
320 template<
int dim,
int cdim >
321 class AlbertaGridGlobalGeometry< dim, cdim, const AlbertaGrid< dim, cdim > >
323 typedef AlbertaGridGlobalGeometry< dim, cdim, const AlbertaGrid< dim, cdim > > This;
326 typedef AlbertaGrid< dim, cdim > Grid;
329 static constexpr int dimbary = dim + 1;
331 typedef Alberta::ElementInfo< dim > ElementInfo;
335 typedef Alberta::Real ctype;
338 static constexpr int mydimension = dimension;
339 static constexpr int codimension = dimension - mydimension;
340 static constexpr int coorddimension = cdim;
342 typedef FieldVector< ctype, mydimension > LocalCoordinate;
352 static constexpr int numCorners = mydimension + 1;
357 AlbertaGridGlobalGeometry ()
362 template<
class CoordReader >
363 AlbertaGridGlobalGeometry (
const CoordReader &coordReader )
365 build( coordReader );
381 GlobalCoordinate
corner (
const int i )
const
383 assert( (i >= 0) && (i <
corners()) );
384 const Alberta::GlobalCoordinate &x = elementInfo_.coordinate( i );
386 for(
int j = 0; j < coorddimension; ++j )
392 GlobalCoordinate
center ()
const
394 GlobalCoordinate centroid_ =
corner( 0 );
395 for(
int i = 1; i < numCorners; ++i )
397 centroid_ *= ctype( 1 ) / ctype( numCorners );
402 GlobalCoordinate
global (
const LocalCoordinate &local )
const;
405 LocalCoordinate local (
const GlobalCoordinate &global )
const;
414 return elementInfo_.geometryCache().integrationElement();
436 return elementInfo_.geometryCache().jacobianTransposed();
440 const JacobianTransposed &
453 return elementInfo_.geometryCache().jacobianInverseTransposed();
457 const JacobianInverseTransposed &
464 Jacobian
jacobian (
const LocalCoordinate &local )
const
482 elementInfo_ = ElementInfo();
486 template<
class CoordReader >
487 void build (
const CoordReader &coordReader )
489 elementInfo_ = coordReader.elementInfo();
493 ElementInfo elementInfo_;
502 template<
class Gr
id >
503 class AlbertaGridLocalGeometryProvider
505 typedef AlbertaGridLocalGeometryProvider< Grid > This;
512 template<
int codim >
515 typedef AlbertaGridGeometry< dimension-codim, dimension, Grid > LocalGeometry;
518 typedef typename Codim< 0 >::LocalGeometry LocalElementGeometry;
519 typedef typename Codim< 1 >::LocalGeometry LocalFaceGeometry;
521 static constexpr int numChildren = 2;
522 static constexpr int numFaces = dimension + 1;
524 static constexpr int minFaceTwist = Alberta::Twist< dimension, dimension-1 >::minTwist;
525 static constexpr int maxFaceTwist = Alberta::Twist< dimension, dimension-1 >::maxTwist;
526 static constexpr int numFaceTwists = maxFaceTwist - minFaceTwist + 1;
529 struct GeoInFatherCoordReader;
530 struct FaceCoordReader;
532 AlbertaGridLocalGeometryProvider ()
534 buildGeometryInFather();
538 ~AlbertaGridLocalGeometryProvider ()
542 delete geometryInFather_[
child ][ 0 ];
543 delete geometryInFather_[
child ][ 1 ];
546 for(
int i = 0; i < numFaces; ++i )
548 for(
int j = 0; j < numFaceTwists; ++j )
549 delete faceGeometry_[ i ][ j ];
553 void buildGeometryInFather();
554 void buildFaceGeometry();
557 const LocalElementGeometry &
558 geometryInFather (
int child,
const int orientation = 1 )
const
560 assert( (
child >= 0) && (
child < numChildren) );
561 assert( (orientation == 1) || (orientation == -1) );
562 return *geometryInFather_[
child ][ (orientation + 1) / 2 ];
565 const LocalFaceGeometry &
566 faceGeometry (
int face,
int twist = 0 )
const
568 assert( (face >= 0) && (face < numFaces) );
569 assert( (twist >= minFaceTwist) && (twist <= maxFaceTwist) );
570 return *faceGeometry_[ face ][ twist - minFaceTwist ];
573 static const This &instance ()
575 static This theInstance;
580 template<
int codim >
581 static int mapVertices (
int subEntity,
int i )
583 return Alberta::MapVertices< dimension, codim >::apply( subEntity, i );
586 const LocalElementGeometry *geometryInFather_[ numChildren ][ 2 ];
587 const LocalFaceGeometry *faceGeometry_[ numFaces ][ numFaceTwists ];
geometry implementation for AlbertaGrid
Definition: geometry.hh:108
GeometryType type() const
obtain the type of reference element
Definition: geometry.hh:153
GlobalCoordinate center() const
return center of geometry
Definition: geometry.hh:175
const JacobianTransposed & jacobianTransposed(const LocalCoordinate &local) const
transposed of the geometry mapping's Jacobian
Definition: geometry.hh:218
const JacobianTransposed & jacobianTransposed() const
transposed of the geometry mapping's Jacobian
Definition: geometry.cc:55
int corners() const
number of corner the geometry
Definition: geometry.hh:162
const JacobianInverseTransposed & jacobianInverseTransposed() const
transposed inverse of the geometry mapping's Jacobian
Definition: geometry.cc:73
GlobalCoordinate corner(const int i) const
obtain the i-th corner of this geometry
Definition: geometry.hh:168
ctype integrationElement() const
integration element of the geometry mapping
Definition: geometry.hh:191
ctype integrationElement(const LocalCoordinate &local) const
integration element of the geometry mapping
Definition: geometry.hh:198
GlobalCoordinate global(const LocalCoordinate &local) const
map a point from the reference element to the geometry
Definition: geometry.cc:34
Alberta::Real ctype
type of coordinates
Definition: geometry.hh:119
ctype volume() const
volume of geometry
Definition: geometry.hh:204
Jacobian jacobian(const LocalCoordinate &local) const
geometry mapping's Jacobian
Definition: geometry.hh:238
const JacobianInverseTransposed & jacobianInverseTransposed(const LocalCoordinate &local) const
transposed inverse of the geometry mapping's Jacobian
Definition: geometry.hh:232
void invalidate()
invalidate the geometry
Definition: geometry.hh:254
JacobianInverse jacobianInverse(const LocalCoordinate &local) const
inverse of the geometry mapping's Jacobian
Definition: geometry.hh:244
bool affine() const
returns always true since we only have simplices
Definition: geometry.hh:159
void build(const CoordReader &coordReader)
build the geometry from a coordinate reader
Definition: geometry.cc:90
FieldMatrix< K, COLS, ROWS > transposed() const
Return transposed of the matrix as FieldMatrix.
Definition: fmatrix.hh:171
Unique label for each type of entities that can occur in DUNE grids.
Definition: type.hh:114
static constexpr int dimension
The dimension of the grid.
Definition: grid.hh:387
static constexpr int dimensionworld
The dimension of the world the grid lives in.
Definition: grid.hh:390
ct ctype
Define type used for coordinates in grid module.
Definition: grid.hh:518
Wrapper and interface classes for element geometries.
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:132
provides a wrapper for ALBERTA's el_info structure
constexpr GeometryType simplex(unsigned int dim)
Returns a GeometryType representing a simplex of dimension dim.
Definition: type.hh:453
Dune namespace.
Definition: alignedallocator.hh:13
static constexpr T factorial(const T &n) noexcept
calculate the factorial of n as a constexpr
Definition: math.hh:111