dune-geometry 2.9.1
Loading...
Searching...
No Matches
Classes | Public Types | Public Member Functions | Static Public Attributes | Protected Types | Protected Member Functions | Static Protected Member Functions | 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 >:
Inheritance graph

Classes

class  JacobianInverseTransposed
 

Public Types

typedef ct ctype
 coordinate type
 
typedef FieldVector< ctype, mydimensionLocalCoordinate
 type of local coordinates
 
typedef FieldVector< ctype, coorddimensionGlobalCoordinate
 type of global coordinates
 
typedef ctype Volume
 type of volume
 
typedef FieldMatrix< ctype, mydimension, coorddimensionJacobianTransposed
 type of jacobian transposed
 
typedef FieldMatrix< ctype, coorddimension, mydimensionJacobian
 Type for the Jacobian matrix.
 
typedef FieldMatrix< ctype, mydimension, coorddimensionJacobianInverse
 Type for the inverse Jacobian matrix.
 
typedef ReferenceElements::ReferenceElement 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 &globalCoord) const
 evaluate the inverse mapping
 
Volume integrationElement (const LocalCoordinate &local) const
 obtain the integration element
 
Volume volume () const
 obtain the volume of the mapping's image
 
JacobianTransposed jacobianTransposed (const LocalCoordinate &local) const
 obtain the transposed of the Jacobian
 
JacobianInverseTransposed jacobianInverseTransposed (const LocalCoordinate &local) const
 obtain the transposed of the Jacobian's inverse
 
Jacobian jacobian (const LocalCoordinate &local) const
 Obtain the Jacobian.
 
JacobianInverse jacobianInverse (const LocalCoordinate &local) const
 Obtain the Jacobian's inverse.
 

Static Public Attributes

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

Protected Types

typedef Dune::ReferenceElements< ctype, mydimensionReferenceElements
 
typedef Traits::MatrixHelper MatrixHelper
 
typedef std::conditional< hasSingleGeometryType, std::integral_constant< unsignedint, Traits::templatehasSingleGeometryType< mydimension >::topologyId >, unsignedint >::type TopologyId
 

Protected Member Functions

ReferenceElement refElement () const
 
TopologyId topologyId () const
 
bool affine (JacobianTransposed &jacobianT) const
 

Static Protected Member Functions

template<bool add, int dim, class CornerIterator >
static void global (TopologyId topologyId, std::integral_constant< int, dim >, CornerIterator &cit, const ctype &df, const LocalCoordinate &x, const ctype &rf, GlobalCoordinate &y)
 
template<bool add, class CornerIterator >
static void global (TopologyId topologyId, std::integral_constant< int, 0 >, CornerIterator &cit, const ctype &df, const LocalCoordinate &x, const ctype &rf, GlobalCoordinate &y)
 
template<bool add, int rows, int dim, class CornerIterator >
static void jacobianTransposed (TopologyId topologyId, std::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, class CornerIterator >
static void jacobianTransposed (TopologyId topologyId, std::integral_constant< int, 0 >, CornerIterator &cit, const ctype &df, const LocalCoordinate &x, const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt)
 
template<int dim, class CornerIterator >
static bool affine (TopologyId topologyId, std::integral_constant< int, dim >, CornerIterator &cit, JacobianTransposed &jt)
 
template<class CornerIterator >
static bool affine (TopologyId topologyId, std::integral_constant< int, 0 >, CornerIterator &cit, JacobianTransposed &jt)
 

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.

Member Typedef Documentation

◆ ctype

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

coordinate type

◆ GlobalCoordinate

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

◆ Jacobian

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

Type for the Jacobian matrix.

◆ JacobianInverse

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

Type for the inverse Jacobian matrix.

◆ JacobianTransposed

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

◆ LocalCoordinate

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

type of local coordinates

◆ MatrixHelper

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

◆ ReferenceElement

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

type of reference element

◆ ReferenceElements

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

◆ TopologyId

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

◆ Volume

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

type of volume

Constructor & Destructor Documentation

◆ MultiLinearGeometry() [1/2]

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.

◆ MultiLinearGeometry() [2/2]

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

◆ affine() [1/4]

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
bool Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::affine ( ) const
inline

is this mapping affine?

◆ affine() [2/4]

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
bool Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::affine ( JacobianTransposed & jacobianT) const
inlineprotected

◆ affine() [3/4]

template<class ct , int mydim, int cdim, class Traits >
template<class CornerIterator >
bool Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::affine ( TopologyId topologyId,
std::integral_constant< int, 0 > ,
CornerIterator & cit,
JacobianTransposed & jt )
inlinestaticprotected

◆ affine() [4/4]

template<class ct , int mydim, int cdim, class Traits >
template<int dim, class CornerIterator >
bool Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::affine ( TopologyId topologyId,
std::integral_constant< int, dim > ,
CornerIterator & cit,
JacobianTransposed & jt )
inlinestaticprotected

◆ center()

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
GlobalCoordinate Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::center ( ) const
inline

obtain the centroid of the mapping's image

◆ corner()

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
GlobalCoordinate Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::corner ( int i) const
inline

obtain coordinates of the i-th corner

◆ corners()

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

obtain number of corners of the corresponding reference element

◆ global() [1/3]

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

evaluate the mapping

Parameters
[in]locallocal coordinate to map
Returns
corresponding global coordinate

◆ global() [2/3]

template<class ct , int mydim, int cdim, class Traits >
template<bool add, class CornerIterator >
void Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::global ( TopologyId topologyId,
std::integral_constant< int, 0 > ,
CornerIterator & cit,
const ctype & df,
const LocalCoordinate & x,
const ctype & rf,
GlobalCoordinate & y )
inlinestaticprotected

◆ global() [3/3]

template<class ct , int mydim, int cdim, class Traits >
template<bool add, int dim, class CornerIterator >
void Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::global ( TopologyId topologyId,
std::integral_constant< int, dim > ,
CornerIterator & cit,
const ctype & df,
const LocalCoordinate & x,
const ctype & rf,
GlobalCoordinate & y )
inlinestaticprotected

◆ integrationElement()

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

◆ jacobian()

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

Obtain the Jacobian.

Parameters
[in]locallocal coordinate to evaluate Jacobian in
Returns
a copy of the transposed of the Jacobian

◆ jacobianInverse()

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

Obtain 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.\]

◆ jacobianInverseTransposed()

template<class ct , int mydim, int cdim, class Traits >
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.\]

◆ jacobianTransposed() [1/3]

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

◆ jacobianTransposed() [2/3]

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

◆ jacobianTransposed() [3/3]

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

◆ local()

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

evaluate the inverse mapping

Parameters
[in]globalCoordglobal 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()
GlobalCoordinate global(const LocalCoordinate &local) const
evaluate the mapping
Definition multilineargeometry.hh:290

◆ refElement()

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

◆ topologyId()

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
TopologyId Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::topologyId ( ) const
inlineprotected

◆ type()

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
Dune::GeometryType Dune::MultiLinearGeometry< ct, mydim, cdim, Traits >::type ( ) const
inline

obtain the name of the reference element

◆ volume()

template<class ct , int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct >>
Volume 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()
ReferenceElement refElement() const
Definition multilineargeometry.hh:425
Volume integrationElement(const LocalCoordinate &local) const
obtain the integration element
Definition multilineargeometry.hh:350
which is wrong for n-linear surface maps and other nonlinear maps.

Member Data Documentation

◆ coorddimension

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

coordinate dimension

◆ mydimension

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

geometry dimension


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