dune-geometry  2.3.1
Classes | Public Types | Public Member Functions | Static Public Attributes | Protected Types | Protected Member Functions | Static Protected Member Functions | Protected Attributes | List of all members
Dune::MultiLinearGeometry< ct, mydim, cdim, Traits > Class Template Reference

generic geometry implementation based on corner coordinates More...

#include <dune/geometry/multilineargeometry.hh>

Inheritance diagram for Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >:
Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >

Classes

class  JacobianInverseTransposed

Public Types

typedef ct ctype
 coordinate type
typedef FieldVector< ctype,
mydimension
LocalCoordinate
 type of user data
typedef FieldVector< ctype,
coorddimension
GlobalCoordinate
 type of global coordinates
typedef FieldMatrix< ctype,
mydimension, coorddimension
JacobianTransposed
 type of jacobian transposed
typedef JacobianInverseTransposed Jacobian
 For backward-compatibility, export the type JacobianInverseTransposed as Jacobian.
typedef Dune::ReferenceElement
< ctype, mydimension
ReferenceElement
 type of reference element

Public Member Functions

template<class Corners >
 MultiLinearGeometry (const ReferenceElement &refElement, const Corners &corners)
 constructor
template<class Corners >
 MultiLinearGeometry (Dune::GeometryType gt, const Corners &corners)
 constructor
bool affine () const
 is this mapping affine?
Dune::GeometryType type () const
 obtain the name of the reference element
int corners () const
 obtain number of corners of the corresponding reference element
GlobalCoordinate corner (int i) const
 obtain coordinates of the i-th corner
GlobalCoordinate center () const
 obtain the centroid of the mapping's image
GlobalCoordinate global (const LocalCoordinate &local) const
 evaluate the mapping
LocalCoordinate local (const GlobalCoordinate &global) const
 evaluate the inverse mapping
ctype integrationElement (const LocalCoordinate &local) const
 obtain the integration element
ctype volume () const
 obtain the volume of the mapping's image
const JacobianTransposedjacobianTransposed (const LocalCoordinate &local) const
 obtain the transposed of the Jacobian
const JacobianInverseTransposedjacobianInverseTransposed (const LocalCoordinate &local) const
 obtain the transposed of the Jacobian's inverse

Static Public Attributes

static const int mydimension = mydim
 geometry dimension
static const int coorddimension = cdim
 coordinate dimension

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 ReferenceElementrefElement () 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_

Detailed Description

template<class ct, int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
class Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >

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.

Template Parameters
ctcoordinate type
mydimgeometry dimension
cdimcoordinate dimension
Traitstraits 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.

Member Typedef Documentation

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
typedef ct Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::ctype

coordinate type

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
typedef FieldVector< ctype, coorddimension > Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::GlobalCoordinate

type of global coordinates

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
typedef JacobianInverseTransposed Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::Jacobian

For backward-compatibility, export the type JacobianInverseTransposed as Jacobian.

Deprecated:
This typedef will be removed after the release of dune 2.3
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
typedef FieldMatrix< ctype, mydimension, coorddimension > Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::JacobianTransposed

type of jacobian transposed

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
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

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
typedef Traits::MatrixHelper Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::MatrixHelper
protected
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
typedef Dune::ReferenceElement< ctype, mydimension > Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::ReferenceElement

type of reference element

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
typedef Dune::ReferenceElements< ctype, mydimension > Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::ReferenceElements
protected
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
typedef conditional< hasSingleGeometryType, integral_constant< unsigned int, Traits::template hasSingleGeometryType< mydimension >::topologyId >, unsigned int >::type Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::TopologyId
protected

Constructor & Destructor Documentation

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
template<class Corners >
Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::MultiLinearGeometry ( const ReferenceElement refElement,
const Corners &  corners 
)
inline

constructor

Parameters
[in]refElementreference element for the geometry
[in]cornerscorners to store internally
Note
The type of corners is actually a template argument. It is only required that the internal corner storage can be constructed from this object.
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
template<class Corners >
Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::MultiLinearGeometry ( Dune::GeometryType  gt,
const Corners &  corners 
)
inline

constructor

Parameters
[in]gtgeometry type
[in]cornerscorners to store internally
Note
The type of corners is actually a template argument. It is only required that the internal corner storage can be constructed from this object.

Member Function Documentation

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
bool Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::affine ( ) const
inline
template<class ct , int mydim, int cdim, class Traits >
template<int dim>
bool Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::affine ( TopologyId  topologyId,
integral_constant< int, dim >  ,
CornerIterator &  cit,
JacobianTransposed jt 
)
inlinestaticprotected
template<class ct , int mydim, int cdim, class Traits >
bool Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::affine ( TopologyId  topologyId,
integral_constant< int, 0 >  ,
CornerIterator &  cit,
JacobianTransposed jt 
)
inlinestaticprotected
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
GlobalCoordinate Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::center ( ) const
inline
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
GlobalCoordinate Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::corner ( int  i) const
inline
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
int Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::corners ( ) const
inline
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
GlobalCoordinate Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::global ( const LocalCoordinate local) const
inline
template<class ct , int mydim, int cdim, class Traits >
template<bool add, int dim>
void Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::global ( TopologyId  topologyId,
integral_constant< int, dim >  ,
CornerIterator &  cit,
const ctype df,
const LocalCoordinate x,
const ctype rf,
GlobalCoordinate y 
)
inlinestaticprotected
template<class ct , int mydim, int cdim, class Traits >
template<bool add>
void Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::global ( TopologyId  topologyId,
integral_constant< int, 0 >  ,
CornerIterator &  cit,
const ctype df,
const LocalCoordinate x,
const ctype rf,
GlobalCoordinate y 
)
inlinestaticprotected
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
ctype Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::integrationElement ( const LocalCoordinate local) const
inline

obtain the integration element

If the Jacobian of the mapping is denoted by $J(x)$, the integration integration element $\mu(x)$ is given by

\[ \mu(x) = \sqrt{|\det (J^T(x) J(x))|}.\]

Parameters
[in]locallocal coordinate to evaluate the integration element in
Returns
the integration element $\mu(x)$.
Note
For affine mappings, it is more efficient to call jacobianInverseTransposed before integrationElement, if both are required.

Reimplemented in Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >.

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().

template<class ct , int mydim, int cdim, class Traits >
const MultiLinearGeometry< ct, mydim, cdim, Traits >::JacobianInverseTransposed & Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianInverseTransposed ( const LocalCoordinate local) const
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 $J(x)$, the following condition holds:

\[J^{-1}(x) J(x) = I.\]

Reimplemented in Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >.

References Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::JacobianInverseTransposed::setup().

Referenced by Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianInverseTransposed().

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
const JacobianTransposed& Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed ( const LocalCoordinate local) const
inline

obtain the transposed of the Jacobian

Parameters
[in]locallocal coordinate to evaluate Jacobian in
Returns
a reference to the transposed of the Jacobian
Note
The returned reference is reused on the next call to JacobianTransposed, destroying the previous value.

Reimplemented in Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >.

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().

template<class ct , int mydim, int cdim, class Traits >
template<bool add, int rows, int dim>
void Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed ( TopologyId  topologyId,
integral_constant< int, dim >  ,
CornerIterator &  cit,
const ctype df,
const LocalCoordinate x,
const ctype rf,
FieldMatrix< ctype, rows, cdim > &  jt 
)
inlinestaticprotected
template<class ct , int mydim, int cdim, class Traits >
template<bool add, int rows>
void Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed ( TopologyId  topologyId,
integral_constant< int, 0 >  ,
CornerIterator &  cit,
const ctype df,
const LocalCoordinate x,
const ctype rf,
FieldMatrix< ctype, rows, cdim > &  jt 
)
inlinestaticprotected
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
LocalCoordinate Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::local ( const GlobalCoordinate global) const
inline

evaluate the inverse mapping

Parameters
[in]globalglobal coordinate to map
Returns
corresponding local coordinate
Note
For given global coordinate y the returned local coordinate x that minimizes the following function over the local coordinate space spanned by the reference element.
(global( x ) - y).two_norm()

Reimplemented in Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >.

References 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().

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
const ReferenceElement& Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::refElement ( ) const
inlineprotected
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
TopologyId Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::topologyId ( ) const
inlineprotected
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
TopologyId Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::topologyId ( integral_constant< bool, true >  ) const
inlineprotected
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
unsigned int Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::topologyId ( integral_constant< bool, false >  ) const
inlineprotected
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
Dune::GeometryType Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::type ( ) const
inline
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
ctype Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::volume ( ) const
inline

obtain the volume of the mapping's image

Note
The current implementation just returns
integrationElement( refElement().position( 0, 0 ) ) * refElement().volume()
which is wrong for n-linear surface maps and other nonlinear maps.

Reimplemented in Dune::CachedMultiLinearGeometry< ct, mydim, cdim, Traits >.

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().

Member Data Documentation

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
const int Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::coorddimension = cdim
static

coordinate dimension

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
JacobianInverseTransposed Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianInverseTransposed_
mutableprotected
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
JacobianTransposed Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::jacobianTransposed_
mutableprotected
template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
const int Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::mydimension = mydim
static

The documentation for this class was generated from the following file: