9 #include <dolfinx/common/span.hpp>
10 #include <dolfinx/common/types.h>
11 #include <dolfinx/mesh/cell_types.h>
14 #include <unsupported/Eigen/CXX11/Tensor>
17 struct ufc_coordinate_mapping;
18 struct ufc_finite_element;
29 explicit FiniteElement(
const ufc_finite_element& ufc_element);
82 std::
string family() const noexcept;
87 Eigen::Tensor<
double, 3, Eigen::RowMajor>& values,
88 const Eigen::Ref<const Eigen::Array<
89 double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>>& X) const;
95 Eigen::Tensor<
double, 4, Eigen::RowMajor>& reference_values,
int order,
96 const Eigen::Ref<const Eigen::Array<
97 double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>>& X) const;
101 Eigen::Tensor<
double, 3, Eigen::RowMajor>& values,
102 const Eigen::Tensor<
double, 3, Eigen::RowMajor>& reference_values,
103 const Eigen::Ref<const Eigen::Array<
double, Eigen::Dynamic,
104 Eigen::Dynamic, Eigen::RowMajor>>& X,
105 const Eigen::Tensor<
double, 3, Eigen::RowMajor>& J,
106 const tcb::span<const
double>& detJ,
107 const Eigen::Tensor<
double, 3, Eigen::RowMajor>& K) const;
111 Eigen::Tensor<
double, 4, Eigen::RowMajor>& values, std::
size_t order,
112 const Eigen::Tensor<
double, 4, Eigen::RowMajor>& reference_values,
113 const Eigen::Ref<const Eigen::Array<
double, Eigen::Dynamic,
114 Eigen::Dynamic, Eigen::RowMajor>>& X,
115 const Eigen::Tensor<
double, 3, Eigen::RowMajor>& J,
116 const tcb::span<const
double>& detJ,
117 const Eigen::Tensor<
double, 3, Eigen::RowMajor>& K) const;
124 std::
size_t hash() const noexcept;
144 Eigen::Array<
double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
165 void interpolate(const Eigen::Array<ufc_scalar_t, Eigen::Dynamic,
166 Eigen::Dynamic, Eigen::RowMajor>& values,
167 std::uint32_t cell_permutation,
168 tcb::span<ufc_scalar_t> dofs) const;
190 std::uint32_t cell_permutation,
194 std::
string _signature, _family;
196 mesh::CellType _cell_shape;
198 int _tdim, _space_dim, _value_size, _reference_value_size;
201 std::vector<std::shared_ptr<const
FiniteElement>> _sub_elements;
206 const std::vector<
int>& component);
212 std::vector<
int> _value_dimension;
214 std::function<
int(
double*,
int,
int, const
double*, const
double*,
215 const
double*, const
double*, const
double*)>
216 _transform_reference_basis_derivatives;
218 std::function<
int(
double*, const std::uint32_t, const
int)>
219 _apply_dof_transformation;
221 std::function<
int(ufc_scalar_t*, const std::uint32_t, const
int)>
222 _apply_dof_transformation_to_scalar;
230 bool _interpolation_is_ident;
233 bool _needs_permutation_data;
236 int _basix_element_handle;
239 Eigen::MatrixXd _interpolation_matrix;
Finite Element, containing the dof layout on a reference element, and various methods for evaluating ...
Definition: FiniteElement.h:25
int space_dimension() const noexcept
Dimension of the finite element function space.
Definition: FiniteElement.cpp:104
void evaluate_reference_basis_derivatives(Eigen::Tensor< double, 4, Eigen::RowMajor > &reference_values, int order, const Eigen::Ref< const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &X) const
Evaluate all basis function derivatives of given order at given points in reference cell.
Definition: FiniteElement.cpp:155
FiniteElement & operator=(FiniteElement &&element)=default
Move assignment.
Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > interpolation_points() const noexcept
Points on the reference cell at which an expression need to be evaluated in order to interpolate the ...
Definition: FiniteElement.cpp:286
void transform_reference_basis_derivatives(Eigen::Tensor< double, 4, Eigen::RowMajor > &values, std::size_t order, const Eigen::Tensor< double, 4, Eigen::RowMajor > &reference_values, const Eigen::Ref< const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &X, const Eigen::Tensor< double, 3, Eigen::RowMajor > &J, const tcb::span< const double > &detJ, const Eigen::Tensor< double, 3, Eigen::RowMajor > &K) const
Push basis function (derivatives) forward to physical element.
Definition: FiniteElement.cpp:200
int reference_value_size() const noexcept
The value size, e.g. 1 for a scalar function, 2 for a 2D vector for the reference element.
Definition: FiniteElement.cpp:108
std::string family() const noexcept
The finite element family.
Definition: FiniteElement.cpp:128
std::size_t hash() const noexcept
Return simple hash of the signature string.
Definition: FiniteElement.cpp:226
int value_dimension(int i) const
Return the dimension of the value space for axis i.
Definition: FiniteElement.cpp:120
void transform_reference_basis(Eigen::Tensor< double, 3, Eigen::RowMajor > &values, const Eigen::Tensor< double, 3, Eigen::RowMajor > &reference_values, const Eigen::Ref< const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &X, const Eigen::Tensor< double, 3, Eigen::RowMajor > &J, const tcb::span< const double > &detJ, const Eigen::Tensor< double, 3, Eigen::RowMajor > &K) const
Push basis functions forward to physical element.
Definition: FiniteElement.cpp:178
FiniteElement(const FiniteElement &element)=default
Copy constructor.
FiniteElement & operator=(const FiniteElement &element)=default
Copy assignment.
bool interpolation_ident() const noexcept
Check if interpolation into the finite element space is an identity operation given the evaluation on...
Definition: FiniteElement.cpp:280
std::string signature() const noexcept
String identifying the finite element.
Definition: FiniteElement.cpp:97
std::shared_ptr< const FiniteElement > extract_sub_element(const std::vector< int > &component) const
Extract sub finite element for component.
Definition: FiniteElement.cpp:229
void interpolate(const Eigen::Array< ufc_scalar_t, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &values, std::uint32_t cell_permutation, tcb::span< ufc_scalar_t > dofs) const
Definition: FiniteElement.cpp:291
FiniteElement(const ufc_finite_element &ufc_element)
Create finite element from UFC finite element.
Definition: FiniteElement.cpp:18
void apply_dof_transformation_to_scalar(ufc_scalar_t *data, std::uint32_t cell_permutation, int block_size) const
Apply permutation to some data.
Definition: FiniteElement.cpp:320
int value_size() const noexcept
The value size, e.g. 1 for a scalar function, 2 for a 2D vector.
Definition: FiniteElement.cpp:106
int block_size() const noexcept
Block size of the finite element function space. For VectorElements and TensorElements,...
Definition: FiniteElement.cpp:118
int value_rank() const noexcept
Rank of the value space.
Definition: FiniteElement.cpp:113
void apply_dof_transformation(double *data, std::uint32_t cell_permutation, int block_size) const
Apply permutation to some data.
Definition: FiniteElement.cpp:313
bool needs_permutation_data() const noexcept
Definition: FiniteElement.cpp:308
void evaluate_reference_basis(Eigen::Tensor< double, 3, Eigen::RowMajor > &values, const Eigen::Ref< const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &X) const
Evaluate all basis functions at given points in reference cell.
Definition: FiniteElement.cpp:130
mesh::CellType cell_shape() const noexcept
Cell shape.
Definition: FiniteElement.cpp:99
virtual ~FiniteElement()=default
Destructor.
int num_sub_elements() const noexcept
Get the number of sub elements (for a mixed element)
Definition: FiniteElement.cpp:221
FiniteElement(FiniteElement &&element)=default
Move constructor.
Finite element method functionality.
Definition: assemble_matrix_impl.h:24