DOLFINx 0.9.0
DOLFINx C++ interface
|
Model of a finite element. More...
#include <FiniteElement.h>
Public Types | |
using | geometry_type = T |
Geometry type of the Mesh that the FunctionSpace is defined on. | |
Public Member Functions | |
FiniteElement (const basix::FiniteElement< geometry_type > &element, std::size_t block_size, bool symmetric=false) | |
Create a finite element from a Basix finite element. | |
FiniteElement (const std::vector< std::shared_ptr< const FiniteElement< geometry_type > > > &elements) | |
Create mixed finite element from a list of finite elements. | |
FiniteElement (mesh::CellType cell_type, std::span< const geometry_type > points, std::array< std::size_t, 2 > pshape, std::size_t block_size, bool symmetric=false) | |
Create a quadrature element. | |
FiniteElement (const FiniteElement &element)=delete | |
Copy constructor. | |
FiniteElement (FiniteElement &&element)=default | |
Move constructor. | |
~FiniteElement ()=default | |
Destructor. | |
FiniteElement & | operator= (const FiniteElement &element)=delete |
Copy assignment. | |
FiniteElement & | operator= (FiniteElement &&element)=default |
Move assignment. | |
bool | operator== (const FiniteElement &e) const |
bool | operator!= (const FiniteElement &e) const |
std::string | signature () const noexcept |
int | space_dimension () const noexcept |
int | block_size () const noexcept |
int | reference_value_size () const |
std::span< const std::size_t > | reference_value_shape () const noexcept |
The reference value shape. | |
const std::vector< std::vector< std::vector< int > > > & | entity_dofs () const noexcept |
The local DOFs associated with each subentity of the cell. | |
const std::vector< std::vector< std::vector< int > > > & | entity_closure_dofs () const noexcept |
The local DOFs associated with the closure of each subentity of the cell. | |
bool | symmetric () const |
Does the element represent a symmetric 2-tensor? | |
void | tabulate (std::span< geometry_type > values, std::span< const geometry_type > X, std::array< std::size_t, 2 > shape, int order) const |
Evaluate derivatives of the basis functions up to given order at points in the reference cell. | |
std::pair< std::vector< geometry_type >, std::array< std::size_t, 4 > > | tabulate (std::span< const geometry_type > X, std::array< std::size_t, 2 > shape, int order) const |
int | num_sub_elements () const noexcept |
Number of sub elements (for a mixed or blocked element). | |
bool | is_mixed () const noexcept |
Check if element is a mixed element. | |
const std::vector< std::shared_ptr< const FiniteElement< geometry_type > > > & | sub_elements () const noexcept |
Get subelements (if any) | |
std::shared_ptr< const FiniteElement< geometry_type > > | extract_sub_element (const std::vector< int > &component) const |
Extract sub finite element for component. | |
const basix::FiniteElement< geometry_type > & | basix_element () const |
Return underlying Basix element (if it exists). | |
basix::maps::type | map_type () const |
Get the map type used by the element. | |
bool | interpolation_ident () const noexcept |
bool | map_ident () const noexcept |
std::pair< std::vector< geometry_type >, std::array< std::size_t, 2 > > | interpolation_points () const |
Points on the reference cell at which an expression needs to be evaluated in order to interpolate the expression in the finite element space. | |
std::pair< std::vector< geometry_type >, std::array< std::size_t, 2 > > | interpolation_operator () const |
std::pair< std::vector< geometry_type >, std::array< std::size_t, 2 > > | create_interpolation_operator (const FiniteElement &from) const |
Create a matrix that maps degrees of freedom from one element to this element (interpolation). | |
bool | needs_dof_transformations () const noexcept |
Check if DOF transformations are needed for this element. | |
bool | needs_dof_permutations () const noexcept |
Check if DOF permutations are needed for this element. | |
template<typename U > | |
std::function< void(std::span< U >, std::span< const std::uint32_t >, std::int32_t, int)> | dof_transformation_fn (doftransform ttype, bool scalar_element=false) const |
Return a function that applies a DOF transformation operator to some data (see T_apply()). | |
template<typename U > | |
std::function< void(std::span< U >, std::span< const std::uint32_t >, std::int32_t, int)> | dof_transformation_right_fn (doftransform ttype, bool scalar_element=false) const |
Return a function that applies DOF transformation to some transposed data (see T_apply_right()). | |
template<typename U > | |
void | T_apply (std::span< U > data, std::uint32_t cell_permutation, int n) const |
Transform basis functions from the reference element ordering and orientation to the globally consistent physical element ordering and orientation. | |
template<typename U > | |
void | Tt_inv_apply (std::span< U > data, std::uint32_t cell_permutation, int n) const |
Apply the inverse transpose of the operator applied by T_apply(). | |
template<typename U > | |
void | Tt_apply (std::span< U > data, std::uint32_t cell_permutation, int n) const |
Apply the transpose of the operator applied by T_apply(). | |
template<typename U > | |
void | Tinv_apply (std::span< U > data, std::uint32_t cell_permutation, int n) const |
Apply the inverse of the operator applied by T_apply(). | |
template<typename U > | |
void | T_apply_right (std::span< U > data, std::uint32_t cell_permutation, int n) const |
Right(post)-apply the operator applied by T_apply(). | |
template<typename U > | |
void | Tinv_apply_right (std::span< U > data, std::uint32_t cell_permutation, int n) const |
Right(post)-apply the inverse of the operator applied by T_apply(). | |
template<typename U > | |
void | Tt_apply_right (std::span< U > data, std::uint32_t cell_permutation, int n) const |
Right(post)-apply the transpose of the operator applied by T_apply(). | |
template<typename U > | |
void | Tt_inv_apply_right (std::span< U > data, std::uint32_t cell_permutation, int n) const |
Right(post)-apply the transpose inverse of the operator applied by T_apply(). | |
void | permute (std::span< std::int32_t > doflist, std::uint32_t cell_permutation) const |
Permute indices associated with degree-of-freedoms on the reference element ordering to the globally consistent physical element degree-of-freedom ordering. | |
void | permute_inv (std::span< std::int32_t > doflist, std::uint32_t cell_permutation) const |
Perform the inverse of the operation applied by permute(). | |
std::function< void(std::span< std::int32_t >, std::uint32_t)> | dof_permutation_fn (bool inverse=false, bool scalar_element=false) const |
Return a function that applies a degree-of-freedom permutation to some data. | |
Model of a finite element.
Provides the dof layout on a reference element, and various methods for evaluating and transforming the basis.
FiniteElement | ( | const basix::FiniteElement< geometry_type > & | element, |
std::size_t | block_size, | ||
bool | symmetric = false ) |
Create a finite element from a Basix finite element.
[in] | element | Basix finite element |
[in] | block_size | The block size for the element |
[in] | symmetric | Is the element a symmetric tensor? |
FiniteElement | ( | const std::vector< std::shared_ptr< const FiniteElement< geometry_type > > > & | elements | ) |
Create mixed finite element from a list of finite elements.
[in] | elements | Basix finite elements |
FiniteElement | ( | mesh::CellType | cell_type, |
std::span< const geometry_type > | points, | ||
std::array< std::size_t, 2 > | pshape, | ||
std::size_t | block_size, | ||
bool | symmetric = false ) |
Create a quadrature element.
[in] | cell_type | The cell type |
[in] | points | Quadrature points |
[in] | pshape | Shape of points array |
[in] | block_size | The block size for the element |
[in] | symmetric | Is the element a symmetric tensor? |
const basix::FiniteElement< T > & basix_element | ( | ) | const |
Return underlying Basix element (if it exists).
Throws | a std::runtime_error is there no Basix element. |
|
noexcept |
Block size of the finite element function space. For BlockedElements, this is the number of DOFs colocated at each DOF point. For other elements, this is always 1.
std::pair< std::vector< T >, std::array< std::size_t, 2 > > create_interpolation_operator | ( | const FiniteElement< T > & | from | ) | const |
Create a matrix that maps degrees of freedom from one element to this element (interpolation).
[in] | from | The element to interpolate from |
from
degrees-of-freedom to the degrees-of-freedom of this element. The (0) matrix data (row-major storage) and (1) the shape (num_dofs of this
element, num_dofs of from
) are returned.std::function< void(std::span< std::int32_t >, std::uint32_t)> dof_permutation_fn | ( | bool | inverse = false, |
bool | scalar_element = false ) const |
Return a function that applies a degree-of-freedom permutation to some data.
The returned function can apply permute() to mixed-elements.
The signature of the returned function has three arguments:
[in] | inverse | Indicates if the inverse transformation should be returned. |
[in] | scalar_element | Indicates is the scalar transformations should be returned for a vector element. |
|
inline |
Return a function that applies a DOF transformation operator to some data (see T_apply()).
The transformation is applied from the left-hand side, i.e.
\[ u \leftarrow T u. \]
If the transformation for the (sub)element is a permutation only, the returned function will do change the ordering for the (sub)element as it is assumed that permutations are incorporated into the degree-of-freedom map.
See the documentation for T_apply() for a description of the transformation for a single element type. This function generates a function that can apply the transformation to a mixed element.
The signature of the returned function has four arguments:
[in] | ttype | The transformation type. Typical usage is:
|
[in] | scalar_element | Indicates whether the scalar transformations should be returned for a vector element |
|
inline |
Return a function that applies DOF transformation to some transposed data (see T_apply_right()).
The transformation is applied from the right-hand side, i.e.
\[ u^{t} \leftarrow u^{t} T. \]
If the transformation for the (sub)element is a permutation only, the returned function will do change the ordering for the (sub)element as it is assumed that permutations are incorporated into the degree-of-freedom map.
The signature of the returned function has four arguments:
[in] | ttype | Transformation type. See dof_transformation_fn(). |
[in] | scalar_element | Indicate if the scalar transformations should be returned for a vector element. |
|
noexcept |
Check if interpolation into the finite element space is an identity operation given the evaluation on an expression at specific points, i.e. the degree-of-freedom are equal to point evaluations. The function will return true
for Lagrange elements.
std::pair< std::vector< T >, std::array< std::size_t, 2 > > interpolation_operator | ( | ) | const |
Interpolation operator (matrix) Pi
that maps a function evaluated at the points provided by FiniteElement::interpolation_points to the element degrees of freedom, i.e. dofs = Pi f_x. See the Basix documentation for basix::FiniteElement::interpolation_matrix for how the data in f_x
should be ordered.
Pi
, returning the data for Pi
(row-major storage) and the shape (num_dofs, num_points * / value_size)
std::pair< std::vector< T >, std::array< std::size_t, 2 > > interpolation_points | ( | ) | const |
Points on the reference cell at which an expression needs to be evaluated in order to interpolate the expression in the finite element space.
For Lagrange elements the points will just be the nodal positions. For other elements the points will typically be the quadrature points used to evaluate moment degrees of freedom.
(num_points, tdim)
.
|
noexcept |
Check if element is a mixed element.
A mixed element i composed of two or more elements of different types (a block element, e.g. a Lagrange element with block size >= 1 is not considered mixed).
|
noexcept |
Check if the push forward/pull back map from the values on reference to the values on a physical cell for this element is the identity map.
|
noexcept |
Check if DOF permutations are needed for this element.
DOF permutations will be needed for elements which might not be continuous when two neighbouring cells disagree on the orientation of a shared subentity, and when this can be corrected for by permuting the DOF numbering in the dofmap.
For example, higher order Lagrange elements will need DOF permutations, as the arrangement of DOFs on a shared sub-entity may be different from the point of view of neighbouring cells, and this can be corrected for by permuting the DOF numbers on each cell.
|
noexcept |
Check if DOF transformations are needed for this element.
DOF transformations will be needed for elements which might not be continuous when two neighbouring cells disagree on the orientation of a shared sub-entity, and when this cannot be corrected for by permuting the DOF numbering in the dofmap.
For example, Raviart-Thomas elements will need DOF transformations, as the neighbouring cells may disagree on the orientation of a basis function, and this orientation cannot be corrected for by permuting the DOF numbers on each cell.
|
noexcept |
Number of sub elements (for a mixed or blocked element).
bool operator!= | ( | const FiniteElement< T > & | e | ) | const |
Check if two elements are not equivalent
bool operator== | ( | const FiniteElement< T > & | e | ) | const |
Check if two elements are equivalent
void permute | ( | std::span< std::int32_t > | doflist, |
std::uint32_t | cell_permutation ) const |
Permute indices associated with degree-of-freedoms on the reference element ordering to the globally consistent physical element degree-of-freedom ordering.
Given an array \(\tilde{d}\) that holds an integer associated with each degree-of-freedom and following the reference element degree-of-freedom ordering, this function computes
\[ d = P \tilde{d},\]
where \(P\) is a permutation matrix and \(d\) holds the integers in \(\tilde{d}\) but permuted to follow the globally consistent physical element degree-of-freedom ordering. The permutation is computed in-place.
[in,out] | doflist | Indicies associated with the degrees-of-freedom. Size=num_dofs . |
[in] | cell_permutation | Permutation data for the cell. |
void permute_inv | ( | std::span< std::int32_t > | doflist, |
std::uint32_t | cell_permutation ) const |
Perform the inverse of the operation applied by permute().
Given an array \(d\) that holds an integer associated with each degree-of-freedom and following the globally consistent physical element degree-of-freedom ordering, this function computes
\[ \tilde{d} = P^{T} d, \]
where \(P^{T}\) is a permutation matrix and \(\tilde{d}\) holds the integers in \(d\) but permuted to follow the reference element degree-of-freedom ordering. The permutation is computed in-place.
[in,out] | doflist | Indicies associated with the degrees-of-freedom. Size=num_dofs . |
[in] | cell_permutation | Permutation data for the cell. |
int reference_value_size | ( | ) | const |
The value size, e.g. 1 for a scalar function, 2 for a 2D vector, 9 for a second-order tensor in 3D, for the reference element
|
noexcept |
String identifying the finite element
|
noexcept |
Dimension of the finite element function space (the number of degrees-of-freedom for the element)
|
inline |
Transform basis functions from the reference element ordering and orientation to the globally consistent physical element ordering and orientation.
Consider that the value of a finite element function \(f_{h}\) at a point is given by
\[ f_{h} = \phi^{T} c, \]
where \(f_{h}\) has shape \(r \times 1\), \(\phi\) has shape \(d \times r\) and holds the finite element basis functions, and \(c\) has shape \(d \times 1\) and holds the degrees-of-freedom. The basis functions and degree-of-freedom are with respect to the physical element orientation. If the degrees-of-freedom on the physical element orientation are given by
\[ \phi = T \tilde{\phi}, \]
where \(T\) is a \(d \times d\) matrix, it follows from \(f_{h} = \phi^{T} c = \tilde{\phi}^{T} T^{T} c\) that
\[ \tilde{c} = T^{T} c. \]
This function applies \(T\) to data. The transformation is performed in-place. The operator \(T\) is orthogonal for many elements, but not all.
This function calls the corresponding Basix function.
[in,out] | data | Data to transform. The shape is (m, n) , where m is the number of dgerees-of-freedom and the storage is row-major. |
[in] | cell_permutation | Permutation data for the cell |
[in] | n | Number of columns in data . |
|
inline |
Right(post)-apply the operator applied by T_apply().
Computes
\[ v^{T} = u^{T} T \]
in-place.
[in,out] | data | The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size) . |
[in] | cell_permutation | Permutation data for the cell |
[in] | n | Block size of the input data |
std::pair< std::vector< T >, std::array< std::size_t, 4 > > tabulate | ( | std::span< const geometry_type > | X, |
std::array< std::size_t, 2 > | shape, | ||
int | order ) const |
Evaluate all derivatives of the basis functions up to given order at given points in reference cell
[in] | X | The reference coordinates at which to evaluate the basis functions. Shape is (num_points, topological dimension) (row-major storage) |
[in] | shape | The shape of X |
[in] | order | The number of derivatives (up to and including this order) to tabulate for |
void tabulate | ( | std::span< geometry_type > | values, |
std::span< const geometry_type > | X, | ||
std::array< std::size_t, 2 > | shape, | ||
int | order ) const |
Evaluate derivatives of the basis functions up to given order at points in the reference cell.
[in,out] | values | Array that will be filled with the tabulated basis values. Must have shape (num_derivatives, num_points, / num_dofs, reference_value_size) (row-major storage) |
[in] | X | The reference coordinates at which to evaluate the basis functions. Shape is (num_points, topological dimension) (row-major storage) |
[in] | shape | The shape of X |
[in] | order | The number of derivatives (up to and including this order) to tabulate for |
|
inline |
Apply the inverse of the operator applied by T_apply().
The transformation
\[ v = T^{-1} u \]
is performed in-place.
[in,out] | data | The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size) . |
[in] | cell_permutation | Permutation data for the cell. |
[in] | n | Block size of the input data. |
|
inline |
Right(post)-apply the inverse of the operator applied by T_apply().
Computes
\[ v^{T} = u^{T} T^{-1} \]
in-place.
[in,out] | data | Data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size) . |
[in] | cell_permutation | Permutation data for the cell |
[in] | n | Block size of the input data |
|
inline |
Apply the transpose of the operator applied by T_apply().
The transformation
\[ u \leftarrow T^{T} u \]
is performed in-place.
[in,out] | data | The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size) . |
[in] | cell_permutation | Permutation data for the cell. |
[in] | n | The block size of the input data. |
|
inline |
Right(post)-apply the transpose of the operator applied by T_apply().
Computes
\[ v^{T} = u^{T} T^{T} \]
in-place.
[in,out] | data | Data to be transformed. The data is flattened with row-major layout, shape=(num_dofs, block_size) . |
[in] | cell_permutation | Permutation data for the cell |
[in] | n | Block size of the input data. |
|
inline |
Apply the inverse transpose of the operator applied by T_apply().
The transformation
\[ v = T^{-T} u \]
is performed in-place.
[in,out] | data | The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size) . |
[in] | cell_permutation | Permutation data for the cell. |
[in] | n | Block_size of the input data. |
|
inline |
Right(post)-apply the transpose inverse of the operator applied by T_apply().
Computes
\[ v^{T} = u^{T} T^{-T} \]
in-place.
[in,out] | data | Data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size) . |
[in] | cell_permutation | Permutation data for the cell. |
[in] | n | Block size of the input data. |