Loading [MathJax]/extensions/tex2jax.js
ug4
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
ug::ReferenceMapping< TRefElem, TWorldDim > Class Template Reference

#include <reference_mapping.h>

Public Member Functions

void global_to_local (MathVector< dim > &locPos, const MathVector< worldDim > &globPos, const size_t maxIter=1000, const number tol=1e-10) const
 map global coordinate to local coordinate
 
void global_to_local (MathVector< dim > *vLocPos, const MathVector< worldDim > *vGlobPos, size_t n, const size_t maxIter=1000, const number tol=1e-10) const
 map global coordinate to local coordinate for n local positions
 
void global_to_local (std::vector< MathVector< dim > > &vLocPos, const std::vector< MathVector< worldDim > > &vGlobPos, const size_t maxIter=1000, const number tol=1e-10) const
 map global coordinate to local coordinate for a vector of local positions
 
bool is_linear () const
 returns if mapping is affine
 
void jacobian (MathMatrix< worldDim, dim > &J, const MathVector< dim > &locPos) const
 returns jacobian
 
void jacobian (MathMatrix< worldDim, dim > *vJ, const MathVector< dim > *vLocPos, size_t n) const
 returns jacobian for n local positions
 
void jacobian (std::vector< MathMatrix< worldDim, dim > > &J, const std::vector< MathVector< dim > > &vLocPos) const
 returns jacobian for a vector of local positions
 
void jacobian_transposed (MathMatrix< dim, worldDim > &JT, const MathVector< dim > &locPos) const
 returns transposed of jacobian
 
void jacobian_transposed (MathMatrix< dim, worldDim > *vJT, const MathVector< dim > *vLocPos, size_t n) const
 returns transposed of jacobian for n local positions
 
void jacobian_transposed (std::vector< MathMatrix< dim, worldDim > > &vJT, const std::vector< MathVector< dim > > &vLocPos) const
 returns transposed of jacobian for a vector of positions
 
number jacobian_transposed_inverse (MathMatrix< worldDim, dim > &JTInv, const MathVector< dim > &locPos) const
 returns transposed of the inverse of the jacobian and sqrt of gram determinante
 
void jacobian_transposed_inverse (MathMatrix< worldDim, dim > *vJTInv, const MathVector< dim > *vLocPos, size_t n) const
 returns transposed of the inverse of the jacobian for n local positions
 
void jacobian_transposed_inverse (MathMatrix< worldDim, dim > *vJTInv, number *vDet, const MathVector< dim > *vLocPos, size_t n) const
 returns transposed of the inverse of the jacobian for n local positions
 
void jacobian_transposed_inverse (std::vector< MathMatrix< worldDim, dim > > &vJTInv, const std::vector< MathVector< dim > > &vLocPos) const
 returns transposed of the inverse of the jacobian for a vector of positions
 
void jacobian_transposed_inverse (std::vector< MathMatrix< worldDim, dim > > &vJTInv, std::vector< number > &vDet, const std::vector< MathVector< dim > > &vLocPos) const
 returns transposed of the inverse of the jacobian for a vector of positions
 
void local_to_global (MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
 map local coordinate to global coordinate
 
void local_to_global (MathVector< worldDim > *vGlobPos, const MathVector< dim > *vLocPos, size_t n) const
 map local coordinate to global coordinate for n local positions
 
void local_to_global (std::vector< MathVector< worldDim > > &vGlobPos, const std::vector< MathVector< dim > > &vLocPos) const
 map local coordinate to global coordinate for a vector of local positions
 
 ReferenceMapping ()
 Default Constructor.
 
 ReferenceMapping (const MathVector< worldDim > *vCornerCoord)
 Constructor setting the corners of the element.
 
 ReferenceMapping (const std::vector< MathVector< worldDim > > &vCornerCoord)
 Constructor setting the corners of the element.
 
number sqrt_gram_det (const MathVector< dim > &locPos) const
 returns the determinate of the jacobian
 
void sqrt_gram_det (number *vDet, const MathVector< dim > *vLocPos, size_t n) const
 returns the determinate of the jacobian for n local positions
 
void sqrt_gram_det (std::vector< number > vDet, const std::vector< MathVector< dim > > &vLocPos) const
 returns the determinate of the jacobian for a vector of local positions
 
void update (const MathVector< worldDim > *vCornerCoord)
 refresh mapping for new set of corners
 
void update (const std::vector< MathVector< worldDim > > &vCornerCoord)
 refresh mapping for new set of corners
 

Static Public Attributes

static const int dim = TRefElem::dim
 reference dimension (domain space dimension)
 
static const bool isLinear = false
 flag if mapping is linear (i.e. Jacobian does not depend on x)
 
static const int worldDim = TWorldDim
 world dimension (range space dimension)
 

Detailed Description

template<typename TRefElem, int TWorldDim>
class ug::ReferenceMapping< TRefElem, TWorldDim >

This class describes the mapping from a reference element into the real (physical) world. The mapping is initialized by the physical positions of the vertices of the real world element. The order of those points must be given as indicated by the corresponding reference element.

Let \(R\) be the reference element and \(T\) be the element. Then, the reference mapping is a mapping:

\[ \phi: R \mapsto T \]

Template Parameters
TRefElemreference element
TWorldDimworld dimension

Constructor & Destructor Documentation

◆ ReferenceMapping() [1/3]

template<typename TRefElem , int TWorldDim>
ug::ReferenceMapping< TRefElem, TWorldDim >::ReferenceMapping ( )

Default Constructor.

◆ ReferenceMapping() [2/3]

template<typename TRefElem , int TWorldDim>
ug::ReferenceMapping< TRefElem, TWorldDim >::ReferenceMapping ( const MathVector< worldDim > *  vCornerCoord)

Constructor setting the corners of the element.

◆ ReferenceMapping() [3/3]

template<typename TRefElem , int TWorldDim>
ug::ReferenceMapping< TRefElem, TWorldDim >::ReferenceMapping ( const std::vector< MathVector< worldDim > > &  vCornerCoord)

Constructor setting the corners of the element.

Member Function Documentation

◆ global_to_local() [1/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::global_to_local ( MathVector< dim > &  locPos,
const MathVector< worldDim > &  globPos,
const size_t  maxIter = 1000,
const number  tol = 1e-10 
) const

map global coordinate to local coordinate

◆ global_to_local() [2/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::global_to_local ( MathVector< dim > *  vLocPos,
const MathVector< worldDim > *  vGlobPos,
size_t  n,
const size_t  maxIter = 1000,
const number  tol = 1e-10 
) const

map global coordinate to local coordinate for n local positions

◆ global_to_local() [3/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::global_to_local ( std::vector< MathVector< dim > > &  vLocPos,
const std::vector< MathVector< worldDim > > &  vGlobPos,
const size_t  maxIter = 1000,
const number  tol = 1e-10 
) const

map global coordinate to local coordinate for a vector of local positions

◆ is_linear()

template<typename TRefElem , int TWorldDim>
bool ug::ReferenceMapping< TRefElem, TWorldDim >::is_linear ( ) const

returns if mapping is affine

◆ jacobian() [1/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::jacobian ( MathMatrix< worldDim, dim > &  J,
const MathVector< dim > &  locPos 
) const

returns jacobian

◆ jacobian() [2/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::jacobian ( MathMatrix< worldDim, dim > *  vJ,
const MathVector< dim > *  vLocPos,
size_t  n 
) const

returns jacobian for n local positions

◆ jacobian() [3/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::jacobian ( std::vector< MathMatrix< worldDim, dim > > &  J,
const std::vector< MathVector< dim > > &  vLocPos 
) const

returns jacobian for a vector of local positions

◆ jacobian_transposed() [1/3]

◆ jacobian_transposed() [2/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::jacobian_transposed ( MathMatrix< dim, worldDim > *  vJT,
const MathVector< dim > *  vLocPos,
size_t  n 
) const

returns transposed of jacobian for n local positions

◆ jacobian_transposed() [3/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::jacobian_transposed ( std::vector< MathMatrix< dim, worldDim > > &  vJT,
const std::vector< MathVector< dim > > &  vLocPos 
) const

returns transposed of jacobian for a vector of positions

◆ jacobian_transposed_inverse() [1/5]

template<typename TRefElem , int TWorldDim>
number ug::ReferenceMapping< TRefElem, TWorldDim >::jacobian_transposed_inverse ( MathMatrix< worldDim, dim > &  JTInv,
const MathVector< dim > &  locPos 
) const

returns transposed of the inverse of the jacobian and sqrt of gram determinante

Referenced by ug::StdLinConsistentGravityX< refDim >::prepare_simplex().

◆ jacobian_transposed_inverse() [2/5]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::jacobian_transposed_inverse ( MathMatrix< worldDim, dim > *  vJTInv,
const MathVector< dim > *  vLocPos,
size_t  n 
) const

returns transposed of the inverse of the jacobian for n local positions

◆ jacobian_transposed_inverse() [3/5]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::jacobian_transposed_inverse ( MathMatrix< worldDim, dim > *  vJTInv,
number vDet,
const MathVector< dim > *  vLocPos,
size_t  n 
) const

returns transposed of the inverse of the jacobian for n local positions

◆ jacobian_transposed_inverse() [4/5]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::jacobian_transposed_inverse ( std::vector< MathMatrix< worldDim, dim > > &  vJTInv,
const std::vector< MathVector< dim > > &  vLocPos 
) const

returns transposed of the inverse of the jacobian for a vector of positions

◆ jacobian_transposed_inverse() [5/5]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::jacobian_transposed_inverse ( std::vector< MathMatrix< worldDim, dim > > &  vJTInv,
std::vector< number > &  vDet,
const std::vector< MathVector< dim > > &  vLocPos 
) const

returns transposed of the inverse of the jacobian for a vector of positions

◆ local_to_global() [1/3]

◆ local_to_global() [2/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::local_to_global ( MathVector< worldDim > *  vGlobPos,
const MathVector< dim > *  vLocPos,
size_t  n 
) const

map local coordinate to global coordinate for n local positions

◆ local_to_global() [3/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::local_to_global ( std::vector< MathVector< worldDim > > &  vGlobPos,
const std::vector< MathVector< dim > > &  vLocPos 
) const

map local coordinate to global coordinate for a vector of local positions

◆ sqrt_gram_det() [1/3]

template<typename TRefElem , int TWorldDim>
number ug::ReferenceMapping< TRefElem, TWorldDim >::sqrt_gram_det ( const MathVector< dim > &  locPos) const

◆ sqrt_gram_det() [2/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::sqrt_gram_det ( number vDet,
const MathVector< dim > *  vLocPos,
size_t  n 
) const

returns the determinate of the jacobian for n local positions

◆ sqrt_gram_det() [3/3]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::sqrt_gram_det ( std::vector< number vDet,
const std::vector< MathVector< dim > > &  vLocPos 
) const

returns the determinate of the jacobian for a vector of local positions

◆ update() [1/2]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::update ( const MathVector< worldDim > *  vCornerCoord)

refresh mapping for new set of corners

Referenced by ug::InterpolateOnDiffElements(), ug::MaxErrorOnElements(), ug::StdLinConsistentGravity< refDim >::prepare_edge(), ug::StdLinConsistentGravity< refDim >::prepare_hexahedron(), ug::StdLinConsistentGravity< refDim >::prepare_prism(), ug::StdLinConsistentGravity< refDim >::prepare_pyramid(), ug::StdLinConsistentGravity< refDim >::prepare_quadrilateral(), ug::StdLinConsistentGravityX< refDim >::prepare_simplex(), ug::StdLinConsistentGravity< refDim >::prepare_tetrahedron(), ug::StdLinConsistentGravity< refDim >::prepare_triangle(), ug::ReferenceMapping< ReferenceVertex, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceEdge, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceTriangle, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceQuadrilateral, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceTetrahedron, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferencePyramid, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferencePrism, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceHexahedron, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceOctahedron, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceVertex, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceEdge, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceTriangle, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceQuadrilateral, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceTetrahedron, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferencePyramid, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferencePrism, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceHexahedron, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceOctahedron, TWorldDim >::ReferenceMapping(), ug::ReferenceMapping< ReferenceVertex, TWorldDim >::update(), ug::ReferenceMapping< ReferenceEdge, TWorldDim >::update(), ug::ReferenceMapping< ReferenceTriangle, TWorldDim >::update(), ug::ReferenceMapping< ReferenceQuadrilateral, TWorldDim >::update(), ug::ReferenceMapping< ReferenceTetrahedron, TWorldDim >::update(), ug::ReferenceMapping< ReferencePyramid, TWorldDim >::update(), ug::ReferenceMapping< ReferencePrism, TWorldDim >::update(), ug::ReferenceMapping< ReferenceHexahedron, TWorldDim >::update(), and ug::ReferenceMapping< ReferenceOctahedron, TWorldDim >::update().

◆ update() [2/2]

template<typename TRefElem , int TWorldDim>
void ug::ReferenceMapping< TRefElem, TWorldDim >::update ( const std::vector< MathVector< worldDim > > &  vCornerCoord)

refresh mapping for new set of corners

Member Data Documentation

◆ dim

template<typename TRefElem , int TWorldDim>
const int ug::ReferenceMapping< TRefElem, TWorldDim >::dim = TRefElem::dim
static

reference dimension (domain space dimension)

◆ isLinear

template<typename TRefElem , int TWorldDim>
const bool ug::ReferenceMapping< TRefElem, TWorldDim >::isLinear = false
static

flag if mapping is linear (i.e. Jacobian does not depend on x)

◆ worldDim


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