dune-geometry 2.10
|
Geometry implementation based on local-basis function parametrization. More...
#include <dune/geometry/localfiniteelementgeometry.hh>
Public Types | |
using | ctype = typename LocalBasisTraits::DomainFieldType |
coordinate type | |
using | LocalCoordinate = FieldVector<ctype, mydimension> |
type of local coordinates | |
using | GlobalCoordinate = FieldVector<ctype, coorddimension> |
type of global coordinates | |
using | Volume = decltype(power(std::declval<ctype>(),mydimension)) |
type of volume | |
using | Jacobian = FieldMatrix<ctype, coorddimension, mydimension> |
type of jacobian | |
using | JacobianTransposed = FieldMatrix<ctype, mydimension, coorddimension> |
type of jacobian transposed | |
using | JacobianInverse = FieldMatrix<ctype, mydimension, coorddimension> |
type of jacobian inverse | |
using | JacobianInverseTransposed = FieldMatrix<ctype, coorddimension, mydimension> |
type of jacobian inverse transposed | |
using | ReferenceElements = Dune::ReferenceElements<ctype, mydimension> |
type of reference element | |
using | ReferenceElement = typename ReferenceElements::ReferenceElement |
Public Member Functions | |
LocalFiniteElementGeometry ()=default | |
Default constructed geometry results in an empty/invalid representation. | |
LocalFiniteElementGeometry (const ReferenceElement &refElement, const LocalFiniteElement &localFE, std::vector< GlobalCoordinate > vertices) | |
Constructor from a vector of coefficients of the LocalBasis parametrizing the geometry. | |
template<class Param , std::enable_if_t< std::is_invocable_r_v< GlobalCoordinate, Param, LocalCoordinate >, int > = 0> | |
LocalFiniteElementGeometry (const ReferenceElement &refElement, const LocalFiniteElement &localFE, Param &¶metrization) | |
Constructor from a local parametrization function, mapping local to (curved) global coordinates. | |
template<class... Args> | |
LocalFiniteElementGeometry (GeometryType gt, Args &&... args) | |
Constructor, forwarding to the other constructors that take a reference-element. | |
int | order () const |
Obtain the (highest) polynomial order of the parametrization. | |
bool | affine () const |
Is this mapping affine? Geometries of order 1 might be affine, but it needs to be checked on non-simplex geometries. | |
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 coordinate mapping. | |
LocalCoordinate | local (const GlobalCoordinate &y, Impl::GaussNewtonOptions< ctype > opts={}) const |
Evaluate the inverse coordinate mapping. | |
ctype | integrationElement (const LocalCoordinate &local) const |
Obtain the integration element. | |
Volume | volume (Impl::ConvergenceOptions< ctype > opts={}) const |
Obtain the volume of the mapping's image. | |
Volume | volume (const QuadratureRule< ctype, mydimension > &quadRule) const |
Obtain the volume of the mapping's image by given quadrature rules. | |
Jacobian | jacobian (const LocalCoordinate &local) const |
Obtain the Jacobian. | |
JacobianTransposed | jacobianTransposed (const LocalCoordinate &local) const |
Obtain the transposed of the Jacobian. | |
JacobianInverse | jacobianInverse (const LocalCoordinate &local) const |
Obtain the Jacobian's inverse. | |
JacobianInverseTransposed | jacobianInverseTransposed (const LocalCoordinate &local) const |
Obtain the transposed of the Jacobian's inverse. | |
const LocalFiniteElement & | finiteElement () const |
Obtain the local finite-element. | |
const std::vector< GlobalCoordinate > & | coefficients () const |
Obtain the coefficients of the parametrization. | |
const LocalBasis & | localBasis () const |
The local basis of the stored local finite-element. | |
Static Public Attributes | |
static const int | mydimension = LocalBasisTraits::dimDomain |
geometry dimension | |
static const int | coorddimension = cdim |
coordinate dimension | |
Protected Types | |
using | MatrixHelper = Impl::FieldMatrixHelper<ctype> |
Geometry implementation based on local-basis function parametrization.
Parametrization of the geometry by any localfunction interpolated into a local finite-element space.
LFE | Type of a local finite-element. |
cdim | Coordinate dimension. |
using Dune::LocalFiniteElementGeometry< LFE, cdim >::ctype = typename LocalBasisTraits::DomainFieldType |
coordinate type
using Dune::LocalFiniteElementGeometry< LFE, cdim >::GlobalCoordinate = FieldVector<ctype, coorddimension> |
type of global coordinates
using Dune::LocalFiniteElementGeometry< LFE, cdim >::Jacobian = FieldMatrix<ctype, coorddimension, mydimension> |
type of jacobian
using Dune::LocalFiniteElementGeometry< LFE, cdim >::JacobianInverse = FieldMatrix<ctype, mydimension, coorddimension> |
type of jacobian inverse
using Dune::LocalFiniteElementGeometry< LFE, cdim >::JacobianInverseTransposed = FieldMatrix<ctype, coorddimension, mydimension> |
type of jacobian inverse transposed
using Dune::LocalFiniteElementGeometry< LFE, cdim >::JacobianTransposed = FieldMatrix<ctype, mydimension, coorddimension> |
type of jacobian transposed
using Dune::LocalFiniteElementGeometry< LFE, cdim >::LocalCoordinate = FieldVector<ctype, mydimension> |
type of local coordinates
|
protected |
using Dune::LocalFiniteElementGeometry< LFE, cdim >::ReferenceElement = typename ReferenceElements::ReferenceElement |
using Dune::LocalFiniteElementGeometry< LFE, cdim >::ReferenceElements = Dune::ReferenceElements<ctype, mydimension> |
type of reference element
using Dune::LocalFiniteElementGeometry< LFE, cdim >::Volume = decltype(power(std::declval<ctype>(),mydimension)) |
type of volume
|
default |
Default constructed geometry results in an empty/invalid representation.
|
inline |
Constructor from a vector of coefficients of the LocalBasis parametrizing the geometry.
[in] | refElement | reference element for the geometry |
[in] | localFE | Local finite-element to use for the parametrization |
[in] | vertices | Coefficients of the local interpolation into the basis |
The vertices are the coefficients of a local interpolation in the local finite-element. For Lagrange local bases, these correspond to vertices on the curved geometry in the local Lagrange nodes.
|
inline |
Constructor from a local parametrization function, mapping local to (curved) global coordinates.
[in] | refElement | reference element for the geometry |
[in] | localFE | Local finite-element to use for the parametrization |
[in] | parametrization | parametrization function with signature GlobalCoordinate(LocalCoordinate) |
The parametrization function is not stored in the class, but interpolated into the local finite-element basis and the computed interpolation coefficients are stored.
|
inlineexplicit |
Constructor, forwarding to the other constructors that take a reference-element.
[in] | gt | geometry type |
[in] | args... | arguments passed to the other constructors |
|
inline |
Is this mapping affine? Geometries of order 1 might be affine, but it needs to be checked on non-simplex geometries.
|
inline |
Obtain the centroid of the mapping's image.
|
inline |
Obtain the coefficients of the parametrization.
|
inline |
Obtain coordinates of the i-th corner.
|
inline |
Obtain number of corners of the corresponding reference element.
|
inline |
Obtain the local finite-element.
|
inline |
Evaluate the coordinate mapping.
Implements a linear combination of local basis functions scaled by the vertices as coefficients.
[in] | local | local coordinate to map |
|
inline |
Obtain the integration element.
If the Jacobian of the geometry is denoted by
[in] | local | local coordinate to evaluate the integration element in. |
|
inline |
Obtain the Jacobian.
[in] | local | local coordinate to evaluate Jacobian in |
|
inline |
Obtain the Jacobian's inverse.
The Jacobian's inverse is defined as a pseudo-inverse. If we denote the Jacobian by
|
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
|
inline |
Obtain the transposed of the Jacobian.
[in] | local | local coordinate to evaluate Jacobian in |
|
inline |
Evaluate the inverse coordinate mapping.
[in] | y | Global coordinate to map |
[in] | opts | Parameters to control the behavior of the Gauss-Newton algorithm. |
y
the local coordinate x
that minimizes the function (global(x) - y).two_norm2()
over the local coordinate space spanned by the reference element.In | case of an error indicating that the local coordinate can not be obtained, a Dune::Exception is thrown, with an error code from `GaussNewtonErrorCode`. |
|
inline |
The local basis of the stored local finite-element.
|
inline |
Obtain the (highest) polynomial order of the parametrization.
|
inline |
Obtain the name of the reference element.
|
inline |
Obtain the volume of the mapping's image by given quadrature rules.
|
inline |
Obtain the volume of the mapping's image.
Calculates the volume of the entity by numerical integration. Since the polynomial order of the volume element is not known, iteratively compute numerical integrals with increasing order of the quadrature rules, until tolerance is reached.
opts | An optional control over the convergence, providing a break tolerance and a maximal iteration count. |
|
static |
coordinate dimension
|
static |
geometry dimension