11#include <basix/finite-element.h>
14#include <dolfinx/mesh/cell_types.h>
21struct ufcx_finite_element;
27template <std::
floating_po
int T>
37 inverse_transpose = 3,
124 std::array<std::
size_t, 2> shape,
int order) const;
135 std::pair<std::vector<
geometry_type>, std::array<std::
size_t, 4>>
190 std::pair<std::vector<
geometry_type>, std::array<std::
size_t, 2>>
202 std::pair<std::vector<
geometry_type>, std::array<std::
size_t, 2>>
217 std::pair<std::vector<
geometry_type>, std::array<std::
size_t, 2>>
265 template <typename U>
266 std::function<
void(std::span<U>, std::span<const std::uint32_t>, std::int32_t,
270 bool scalar_element = false)
const
275 return [](std::span<U>, std::span<const std::uint32_t>, std::int32_t, int)
281 if (_sub_elements.size() != 0)
286 std::vector<std::function<void(
287 std::span<U>, std::span<const std::uint32_t>, std::int32_t,
int)>>
288 sub_element_functions;
289 std::vector<int> dims;
290 for (std::size_t i = 0; i < _sub_elements.size(); ++i)
292 sub_element_functions.push_back(
293 _sub_elements[i]->
template get_pre_dof_transformation_function<U>(
298 return [dims, sub_element_functions](
299 std::span<U> data, std::span<const std::uint32_t> cell_info,
302 std::size_t offset = 0;
303 for (std::size_t e = 0; e < sub_element_functions.size(); ++e)
305 const std::size_t width = dims[e] *
block_size;
306 sub_element_functions[e](data.subspan(offset, width), cell_info,
312 else if (!scalar_element)
315 const std::function<void(std::span<U>, std::span<const std::uint32_t>,
318 = _sub_elements[0]->template get_pre_dof_transformation_function<U>(
321 return [ebs, sub_function](std::span<U> data,
322 std::span<const std::uint32_t> cell_info,
323 std::int32_t
cell,
int data_block_size)
324 { sub_function(data, cell_info,
cell, ebs * data_block_size); };
329 case doftransform::inverse_transpose:
330 return [
this](std::span<U> data, std::span<const std::uint32_t> cell_info,
336 case doftransform::transpose:
337 return [
this](std::span<U> data, std::span<const std::uint32_t> cell_info,
342 case doftransform::inverse:
343 return [
this](std::span<U> data, std::span<const std::uint32_t> cell_info,
347 case doftransform::standard:
348 return [
this](std::span<U> data, std::span<const std::uint32_t> cell_info,
352 throw std::runtime_error(
"Unknown transformation type");
370 template <
typename U>
371 std::function<void(std::span<U>, std::span<const std::uint32_t>, std::int32_t,
374 = doftransform::standard,
375 bool scalar_element =
false)
const
380 return [](std::span<U>, std::span<const std::uint32_t>, std::int32_t, int)
385 else if (_sub_elements.size() != 0)
390 std::vector<std::function<void(
391 std::span<U>, std::span<const std::uint32_t>, std::int32_t,
int)>>
392 sub_element_functions;
393 for (std::size_t i = 0; i < _sub_elements.size(); ++i)
395 sub_element_functions.push_back(
397 ->
template get_post_dof_transformation_function<U>(ttype));
400 return [
this, sub_element_functions](
401 std::span<U> data, std::span<const std::uint32_t> cell_info,
404 std::size_t offset = 0;
405 for (std::size_t e = 0; e < sub_element_functions.size(); ++e)
407 sub_element_functions[e](data.subspan(offset, data.size() - offset),
409 offset += _sub_elements[e]->space_dimension();
413 else if (!scalar_element)
416 const std::function<void(std::span<U>, std::span<const std::uint32_t>,
419 = _sub_elements[0]->template get_pre_dof_transformation_function<U>(
421 return [
this, sub_function](std::span<U> data,
422 std::span<const std::uint32_t> cell_info,
423 std::int32_t
cell,
int data_block_size)
426 const std::size_t dof_count = data.size() / data_block_size;
427 for (
int block = 0; block < data_block_size; ++block)
429 sub_function(data.subspan(block * dof_count, dof_count), cell_info,
438 case doftransform::inverse_transpose:
439 return [
this](std::span<U> data, std::span<const std::uint32_t> cell_info,
445 case doftransform::transpose:
446 return [
this](std::span<U> data, std::span<const std::uint32_t> cell_info,
451 case doftransform::inverse:
452 return [
this](std::span<U> data, std::span<const std::uint32_t> cell_info,
457 case doftransform::standard:
458 return [
this](std::span<U> data, std::span<const std::uint32_t> cell_info,
462 throw std::runtime_error(
"Unknown transformation type");
472 template <
typename U>
474 std::uint32_t cell_permutation,
478 _element->pre_apply_dof_transformation(data,
block_size, cell_permutation);
489 template <
typename U>
491 std::span<U> data, std::uint32_t cell_permutation,
int block_size)
const
494 _element->pre_apply_inverse_transpose_dof_transformation(data,
block_size,
505 template <
typename U>
507 std::uint32_t cell_permutation,
511 _element->pre_apply_transpose_dof_transformation(data,
block_size,
522 template <
typename U>
524 std::uint32_t cell_permutation,
528 _element->pre_apply_inverse_dof_transformation(data,
block_size,
538 template <
typename U>
540 std::uint32_t cell_permutation,
544 _element->post_apply_dof_transformation(data,
block_size, cell_permutation);
553 template <
typename U>
555 std::uint32_t cell_permutation,
559 _element->post_apply_inverse_dof_transformation(data,
block_size,
569 template <
typename U>
571 std::uint32_t cell_permutation,
575 _element->post_apply_transpose_dof_transformation(data,
block_size,
585 template <
typename U>
587 std::span<U> data, std::uint32_t cell_permutation,
int block_size)
const
590 _element->post_apply_inverse_transpose_dof_transformation(data,
block_size,
599 std::uint32_t cell_permutation)
const;
606 std::uint32_t cell_permutation)
const;
619 std::function<void(std::span<std::int32_t>, std::uint32_t)>
621 bool scalar_element =
false)
const;
624 std::string _signature;
631 std::vector<std::shared_ptr<const FiniteElement<geometry_type>>>
635 std::vector<std::size_t> _reference_value_shape;
645 bool _needs_dof_permutations;
646 bool _needs_dof_transformations;
649 std::unique_ptr<basix::FiniteElement<geometry_type>> _element;
653 std::pair<std::vector<geometry_type>, std::array<std::size_t, 2>> _points;
Definition CoordinateElement.h:26
Finite Element, containing the dof layout on a reference element, and various methods for evaluating ...
Definition FiniteElement.h:29
FiniteElement & operator=(const FiniteElement &element)=delete
Copy assignment.
const std::vector< std::shared_ptr< const FiniteElement< geometry_type > > > & sub_elements() const noexcept
Get subelements (if any)
Definition FiniteElement.cpp:424
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).
Definition FiniteElement.cpp:514
basix::maps::type map_type() const
Get the map type used by the element.
Definition FiniteElement.cpp:453
void post_apply_inverse_dof_transformation(std::span< U > data, std::uint32_t cell_permutation, int block_size) const
Apply inverse of DOF transformation to some transposed data.
Definition FiniteElement.h:554
int space_dimension() const noexcept
Dimension of the finite element function space (the number of degrees-of-freedom for the element)
Definition FiniteElement.cpp:368
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.
Definition FiniteElement.cpp:393
void post_apply_dof_transformation(std::span< U > data, std::uint32_t cell_permutation, int block_size) const
Apply DOF transformation to some transposed data.
Definition FiniteElement.h:539
FiniteElement & operator=(FiniteElement &&element)=default
Move assignment.
bool needs_dof_transformations() const noexcept
Check if DOF transformations are needed for this element.
Definition FiniteElement.cpp:557
FiniteElement(const FiniteElement &element)=delete
Copy constructor.
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...
Definition FiniteElement.cpp:485
std::shared_ptr< const FiniteElement< geometry_type > > extract_sub_element(const std::vector< int > &component) const
Extract sub finite element for component.
Definition FiniteElement.cpp:431
std::string signature() const noexcept
String identifying the finite element.
Definition FiniteElement.cpp:356
std::function< void(std::span< std::int32_t >, std::uint32_t)> get_dof_permutation_function(bool inverse=false, bool scalar_element=false) const
Return a function that applies DOF permutation to some data.
bool needs_dof_permutations() const noexcept
Check if DOF permutations are needed for this element.
Definition FiniteElement.cpp:563
void pre_apply_transpose_dof_transformation(std::span< U > data, std::uint32_t cell_permutation, int block_size) const
Apply transpose transformation to some data. For VectorElements, this applies the transformations for...
Definition FiniteElement.h:506
bool is_mixed() const noexcept
Check if element is a mixed element.
Definition FiniteElement.cpp:417
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,...
Definition FiniteElement.cpp:380
doftransform
DOF transformation type.
Definition FiniteElement.h:33
void pre_apply_dof_transformation(std::span< U > data, std::uint32_t cell_permutation, int block_size) const
Apply DOF transformation to some data.
Definition FiniteElement.h:473
void post_apply_transpose_dof_transformation(std::span< U > data, std::uint32_t cell_permutation, int block_size) const
Apply transpose of transformation to some transposed data.
Definition FiniteElement.h:570
bool operator!=(const FiniteElement &e) const
Check if two elements are not equivalent.
Definition FiniteElement.cpp:350
~FiniteElement()=default
Destructor.
T geometry_type
Geometry type of the Mesh that the FunctionSpace is defined on.
Definition FiniteElement.h:41
const basix::FiniteElement< geometry_type > & basix_element() const
Return underlying basix element (if it exists)
Definition FiniteElement.cpp:441
void pre_apply_inverse_dof_transformation(std::span< U > data, std::uint32_t cell_permutation, int block_size) const
Apply inverse transformation to some data. For VectorElements, this applies the transformations for t...
Definition FiniteElement.h:523
FiniteElement(FiniteElement &&element)=default
Move constructor.
int num_sub_elements() const noexcept
Number of sub elements (for a mixed or blocked element).
Definition FiniteElement.cpp:411
std::function< void(std::span< U >, std::span< const std::uint32_t >, std::int32_t, int)> get_post_dof_transformation_function(doftransform ttype=doftransform::standard, bool scalar_element=false) const
Return a function that applies DOF transformation to some transposed data.
Definition FiniteElement.h:373
std::pair< std::vector< geometry_type >, std::array< std::size_t, 2 > > interpolation_operator() const
Interpolation operator (matrix) Pi that maps a function evaluated at the points provided by FiniteEle...
Definition FiniteElement.cpp:501
int block_size() const noexcept
Block size of the finite element function space. For BlockedElements, this is the number of DOFs colo...
Definition FiniteElement.cpp:387
std::span< const std::size_t > reference_value_shape() const
The reference value shape.
Definition FiniteElement.cpp:374
void pre_apply_inverse_transpose_dof_transformation(std::span< U > data, std::uint32_t cell_permutation, int block_size) const
Apply inverse transpose transformation to some data. For VectorElements, this applies the transformat...
Definition FiniteElement.h:490
bool interpolation_ident() const noexcept
Check if interpolation into the finite element space is an identity operation given the evaluation on...
Definition FiniteElement.cpp:475
void unpermute_dofs(std::span< std::int32_t > doflist, std::uint32_t cell_permutation) const
Unpermute the DOFs of the element.
Definition FiniteElement.cpp:576
void post_apply_inverse_transpose_dof_transformation(std::span< U > data, std::uint32_t cell_permutation, int block_size) const
Apply inverse transpose transformation to some transposed data.
Definition FiniteElement.h:586
void permute_dofs(std::span< std::int32_t > doflist, std::uint32_t cell_permutation) const
Permute the DOFs of the element.
Definition FiniteElement.cpp:569
std::function< void(std::span< U >, std::span< const std::uint32_t >, std::int32_t, int)> get_pre_dof_transformation_function(doftransform ttype=doftransform::standard, bool scalar_element=false) const
Return a function that applies DOF transformation to some data.
Definition FiniteElement.h:268
bool operator==(const FiniteElement &e) const
Check if two elements are equivalent.
Definition FiniteElement.cpp:338
mesh::CellType cell_shape() const noexcept
Cell shape.
Definition FiniteElement.cpp:362
bool map_ident() const noexcept
Check if the push forward/pull back map from the values on reference to the values on a physical cell...
Definition FiniteElement.cpp:465
Finite element method functionality.
Definition assemble_matrix_impl.h:26
CellType
Cell type identifier.
Definition cell_types.h:22