dune-geometry
2.3.0
|
generic geometry implementation based on corner coordinates More...
#include <dune/geometry/multilineargeometry.hh>
Classes | |
class | JacobianInverseTransposed |
Public Types | |
typedef ct | ctype |
coordinate type More... | |
typedef FieldVector< ctype, mydimension > | LocalCoordinate |
type of user data More... | |
typedef FieldVector< ctype, coorddimension > | GlobalCoordinate |
type of global coordinates More... | |
typedef FieldMatrix< ctype, mydimension, coorddimension > | JacobianTransposed |
type of jacobian transposed More... | |
typedef JacobianInverseTransposed | Jacobian |
For backward-compatibility, export the type JacobianInverseTransposed as Jacobian. More... | |
typedef Dune::ReferenceElement < ctype, mydimension > | ReferenceElement |
type of reference element More... | |
Public Member Functions | |
template<class Corners > | |
MultiLinearGeometry (const ReferenceElement &refElement, const Corners &corners) | |
constructor More... | |
template<class Corners > | |
MultiLinearGeometry (Dune::GeometryType gt, const Corners &corners) | |
constructor More... | |
bool | affine () const |
is this mapping affine? More... | |
Dune::GeometryType | type () const |
obtain the name of the reference element More... | |
int | corners () const |
obtain number of corners of the corresponding reference element More... | |
GlobalCoordinate | corner (int i) const |
obtain coordinates of the i-th corner More... | |
GlobalCoordinate | center () const |
obtain the centroid of the mapping's image More... | |
GlobalCoordinate | global (const LocalCoordinate &local) const |
evaluate the mapping More... | |
LocalCoordinate | local (const GlobalCoordinate &global) const |
evaluate the inverse mapping More... | |
ctype | integrationElement (const LocalCoordinate &local) const |
obtain the integration element More... | |
ctype | volume () const |
obtain the volume of the mapping's image More... | |
const JacobianTransposed & | jacobianTransposed (const LocalCoordinate &local) const |
obtain the transposed of the Jacobian More... | |
const JacobianInverseTransposed & | jacobianInverseTransposed (const LocalCoordinate &local) const |
obtain the transposed of the Jacobian's inverse More... | |
Static Public Attributes | |
static const int | mydimension = mydim |
geometry dimension More... | |
static const int | coorddimension = cdim |
coordinate dimension More... | |
Protected Types | |
typedef Traits::MatrixHelper | MatrixHelper |
typedef conditional < hasSingleGeometryType, integral_constant< unsigned int, Traits::template hasSingleGeometryType < mydimension >::topologyId > , unsigned int >::type | TopologyId |
typedef Dune::ReferenceElements< ctype, mydimension > | ReferenceElements |
Protected Member Functions | |
const ReferenceElement & | refElement () const |
TopologyId | topologyId () const |
TopologyId | topologyId (integral_constant< bool, true >) const |
unsigned int | topologyId (integral_constant< bool, false >) const |
Static Protected Member Functions | |
template<bool add, int dim> | |
static void | global (TopologyId topologyId, integral_constant< int, dim >, CornerIterator &cit, const ctype &df, const LocalCoordinate &x, const ctype &rf, GlobalCoordinate &y) |
template<bool add> | |
static void | global (TopologyId topologyId, integral_constant< int, 0 >, CornerIterator &cit, const ctype &df, const LocalCoordinate &x, const ctype &rf, GlobalCoordinate &y) |
template<bool add, int rows, int dim> | |
static void | jacobianTransposed (TopologyId topologyId, integral_constant< int, dim >, CornerIterator &cit, const ctype &df, const LocalCoordinate &x, const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt) |
template<bool add, int rows> | |
static void | jacobianTransposed (TopologyId topologyId, integral_constant< int, 0 >, CornerIterator &cit, const ctype &df, const LocalCoordinate &x, const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt) |
template<int dim> | |
static bool | affine (TopologyId topologyId, integral_constant< int, dim >, CornerIterator &cit, JacobianTransposed &jt) |
static bool | affine (TopologyId topologyId, integral_constant< int, 0 >, CornerIterator &cit, JacobianTransposed &jt) |
Protected Attributes | |
JacobianTransposed | jacobianTransposed_ |
JacobianInverseTransposed | jacobianInverseTransposed_ |
generic geometry implementation based on corner coordinates
Based on the recursive definition of the reference elements, the MultiLinearGeometry provides a generic implementation of a geometry given the corner coordinates.
The geometric mapping is multilinear in the classical sense only in the case of cubes; for simplices it is linear. The name is still justified, because the mapping satisfies the important property of begin linear along edges.
ct | coordinate type |
mydim | geometry dimension |
cdim | coordinate dimension |
Traits | traits allowing to tweak some implementation details (optional) |
The requirements on the traits are documented along with their default, MultiLinearGeometryTraits.
As an additional feature, this class allows to attach arbitrary user data to each object. This is used in GeometryGrid to implement reference counting.
typedef ct Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::ctype |
coordinate type
typedef FieldVector< ctype, coorddimension > Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::GlobalCoordinate |
type of global coordinates
typedef JacobianInverseTransposed Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::Jacobian |
For backward-compatibility, export the type JacobianInverseTransposed as Jacobian.
typedef FieldMatrix< ctype, mydimension, coorddimension > Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::JacobianTransposed |
type of jacobian transposed
typedef FieldVector< ctype, mydimension > Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::LocalCoordinate |
type of user data
For example, GeometryGrid uses this to implement reference counting.type of local coordinates
|
protected |
typedef Dune::ReferenceElement< ctype, mydimension > Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::ReferenceElement |
type of reference element
|
protected |
|
protected |
|
inline |
constructor
[in] | refElement | reference element for the geometry |
[in] | corners | corners to store internally |
|
inline |
constructor
[in] | gt | geometry type |
[in] | corners | corners to store internally |
|
inline |
is this mapping affine?
References Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed_, and Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::topologyId().
|
inlinestaticprotected |
References Dune::GenericGeometry::isPrism().
|
inlinestaticprotected |
|
inline |
obtain the centroid of the mapping's image
References Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::global(), and Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::refElement().
|
inline |
obtain coordinates of the i-th corner
References Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::corners().
Referenced by Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::global(), and Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::local().
|
inline |
obtain number of corners of the corresponding reference element
References Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::mydimension, Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::refElement(), and Dune::ReferenceElement< ctype, dim >::size().
Referenced by Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::corner().
|
inline |
evaluate the mapping
[in] | local | local coordinate to map |
References Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::local(), and Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::topologyId().
Referenced by Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::center(), Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::global(), and Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::local().
|
inlinestaticprotected |
References Dune::GenericGeometry::isPrism(), and Dune::GenericGeometry::isPyramid().
|
inlinestaticprotected |
|
inline |
obtain the integration element
If the Jacobian of the mapping is denoted by $J(x)$, the integration integration element is given by
[in] | local | local coordinate to evaluate the integration element in |
References Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed().
Referenced by Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::integrationElement(), and Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::volume().
|
inline |
obtain the transposed of the Jacobian's inverse
The Jacobian's inverse is defined as a pseudo-inverse. If we denote the Jacobian by , the following condition holds:
References Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::JacobianInverseTransposed::setup().
Referenced by Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianInverseTransposed().
|
inline |
obtain the transposed of the Jacobian
[in] | local | local coordinate to evaluate Jacobian in |
References Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed_, Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::local(), and Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::topologyId().
Referenced by Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::integrationElement(), Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed(), and Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::local().
|
inlinestaticprotected |
References Dune::GenericGeometry::isPrism(), and Dune::GenericGeometry::isPyramid().
|
inlinestaticprotected |
|
inline |
evaluate the inverse mapping
[in] | global | global coordinate to map |
References Dune::GenericGeometry::checkInside(), Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::global(), Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed(), Dune::ReferenceElement< ctype, dim >::position(), and Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::refElement().
Referenced by Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::global(), Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed(), and Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::local().
|
inlineprotected |
Referenced by Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::center(), Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::center(), Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::corners(), Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::local(), Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::topologyId(), Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::volume(), and Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::volume().
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inline |
obtain the name of the reference element
References Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::mydimension, and Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::topologyId().
|
inline |
obtain the volume of the mapping's image
References Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::integrationElement(), Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::refElement(), and Dune::ReferenceElement< ctype, dim >::volume().
Referenced by Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::volume().
|
static |
coordinate dimension
|
mutableprotected |
|
mutableprotected |
Referenced by Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::affine(), Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::global(), Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::integrationElement(), Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianInverseTransposed(), Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed(), Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed(), and Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::local().
|
static |
geometry dimension
Referenced by Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::corners(), and Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::type().