3 #ifndef DUNE_GEOMETRY_AFFINEGEOMETRY_HH
4 #define DUNE_GEOMETRY_AFFINEGEOMETRY_HH
11 #include <dune/common/fmatrix.hh>
12 #include <dune/common/fvector.hh>
24 template<
class ctype,
int dim >
27 template<
class ctype,
int dim >
37 template<
class ct,
int mydim,
int cdim>
76 : refElement_(&refElement), origin_(origin), jacobianTransposed_(jt)
78 integrationElement_ = MatrixHelper::template rightInvA< mydimension, coorddimension >( jacobianTransposed_, jacobianInverseTransposed_ );
84 : refElement_( &
ReferenceElements::general( gt ) ), origin_(origin), jacobianTransposed_( jt )
86 integrationElement_ = MatrixHelper::template rightInvA< mydimension, coorddimension >( jacobianTransposed_, jacobianInverseTransposed_ );
90 template<
class CoordVector >
92 : refElement_(&refElement), origin_(coordVector[0])
95 jacobianTransposed_[ i ] = coordVector[ i+1 ] - origin_;
96 integrationElement_ = MatrixHelper::template rightInvA< mydimension, coorddimension >( jacobianTransposed_, jacobianInverseTransposed_ );
100 template<
class CoordVector >
105 jacobianTransposed_[ i ] = coordVector[ i+1 ] - origin_;
106 integrationElement_ = MatrixHelper::template rightInvA< mydimension, coorddimension >( jacobianTransposed_, jacobianInverseTransposed_ );
136 jacobianTransposed_.umtv( local,
global );
156 jacobianInverseTransposed_.mtv(
global - origin_, local );
172 DUNE_UNUSED_PARAMETER(local);
173 return integrationElement_;
179 return integrationElement_ * refElement_->
volume();
190 DUNE_UNUSED_PARAMETER(local);
191 return jacobianTransposed_;
202 DUNE_UNUSED_PARAMETER(local);
203 return jacobianInverseTransposed_;
213 ctype integrationElement_;
A unique label for each type of element that can occur in a grid.
Definition: affinegeometry.hh:19
This class provides access to geometric and topological properties of a reference element.
Definition: referenceelements.hh:29
int size(int c) const
number of subentities of codimension c
Definition: referenceelements.hh:80
ctype volume() const
obtain the volume of the reference element
Definition: referenceelements.hh:186
const FieldVector< ctype, dim > & position(int i, int c) const
position of the barycenter of entity (i,c)
Definition: referenceelements.hh:148
const GeometryType & type(int i, int c) const
obtain the type of subentity (i,c)
Definition: referenceelements.hh:130
Class providing access to the singletons of the reference elements.
Definition: referenceelements.hh:455
Implementation of the Geometry interface for affine geometries.
Definition: affinegeometry.hh:39
AffineGeometry(const ReferenceElement &refElement, const CoordVector &coordVector)
Create affine geometry from reference element and a vector of vertex coordinates.
Definition: affinegeometry.hh:91
AffineGeometry(Dune::GeometryType gt, const GlobalCoordinate &origin, const JacobianTransposed &jt)
Create affine geometry from GeometryType, one vertex, and the Jacobian matrix.
Definition: affinegeometry.hh:82
const JacobianInverseTransposed & jacobianInverseTransposed(const LocalCoordinate &local) const
Obtain the transposed of the Jacobian's inverse.
Definition: affinegeometry.hh:200
FieldVector< ctype, mydimension > LocalCoordinate
Type for local coordinate vector.
Definition: affinegeometry.hh:52
Dune::GeometryType type() const
Obtain the type of the reference element.
Definition: affinegeometry.hh:113
static const int mydimension
Dimension of the geometry.
Definition: affinegeometry.hh:46
AffineGeometry(const ReferenceElement &refElement, const GlobalCoordinate &origin, const JacobianTransposed &jt)
Create affine geometry from reference element, one vertex, and the Jacobian matrix.
Definition: affinegeometry.hh:74
ctype volume() const
Obtain the volume of the element.
Definition: affinegeometry.hh:177
AffineGeometry(Dune::GeometryType gt, const CoordVector &coordVector)
Create affine geometry from GeometryType and a vector of vertex coordinates.
Definition: affinegeometry.hh:101
ctype integrationElement(const LocalCoordinate &local) const
Obtain the integration element.
Definition: affinegeometry.hh:170
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type for the transposed Jacobian matrix.
Definition: affinegeometry.hh:58
GlobalCoordinate corner(int i) const
Obtain coordinates of the i-th corner.
Definition: affinegeometry.hh:119
int corners() const
Obtain number of corners of the corresponding reference element.
Definition: affinegeometry.hh:116
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Type for the transposed inverse Jacobian matrix.
Definition: affinegeometry.hh:61
static const int coorddimension
Dimension of the world space.
Definition: affinegeometry.hh:49
GlobalCoordinate global(const LocalCoordinate &local) const
Evaluate the mapping.
Definition: affinegeometry.hh:133
GlobalCoordinate center() const
Obtain the centroid of the mapping's image.
Definition: affinegeometry.hh:125
ct ctype
Type used for coordinates.
Definition: affinegeometry.hh:43
FieldVector< ctype, coorddimension > GlobalCoordinate
Type for coordinate vector in world space.
Definition: affinegeometry.hh:55
bool affine() const
Always true: this is an affine geometry.
Definition: affinegeometry.hh:110
const JacobianTransposed & jacobianTransposed(const LocalCoordinate &local) const
Obtain the transposed of the Jacobian.
Definition: affinegeometry.hh:188
friend const ReferenceElement & referenceElement(const AffineGeometry &geometry)
Definition: affinegeometry.hh:206
Definition: matrixhelper.hh:36
Unique label for each type of entities that can occur in DUNE grids.
Definition: type.hh:25