DOLFIN-X
DOLFIN-X C++ interface
FiniteElement.h
1 // Copyright (C) 2008-2020 Anders Logg and Garth N. Wells
2 //
3 // This file is part of DOLFINX (https://www.fenicsproject.org)
4 //
5 // SPDX-License-Identifier: LGPL-3.0-or-later
6 
7 #pragma once
8 
9 #include <dolfinx/common/span.hpp>
10 #include <dolfinx/common/types.h>
11 #include <dolfinx/mesh/cell_types.h>
12 #include <functional>
13 #include <memory>
14 #include <numeric>
15 #include <vector>
16 
17 struct ufc_coordinate_mapping;
18 struct ufc_finite_element;
19 
20 namespace dolfinx::fem
21 {
25 {
26 public:
29  explicit FiniteElement(const ufc_finite_element& ufc_element);
30 
32  FiniteElement(const FiniteElement& element) = default;
33 
35  FiniteElement(FiniteElement&& element) = default;
36 
38  virtual ~FiniteElement() = default;
39 
41  FiniteElement& operator=(const FiniteElement& element) = default;
42 
44  FiniteElement& operator=(FiniteElement&& element) = default;
45 
48  std::string signature() const noexcept;
49 
52  mesh::CellType cell_shape() const noexcept;
53 
56  int space_dimension() const noexcept;
57 
62  int block_size() const noexcept;
63 
66  int value_size() const noexcept;
67 
71  int reference_value_size() const noexcept;
72 
75  int value_rank() const noexcept;
76 
78  int value_dimension(int i) const;
79 
82  std::string family() const noexcept;
83 
85  // reference_values[num_points][num_dofs][reference_value_size]
86  void evaluate_reference_basis(std::vector<double>& values,
87  const array2d<double>& X) const;
88 
91  // reference_value_derivatives[num_points][num_dofs][reference_value_size][num_derivatives]
92  void
93  evaluate_reference_basis_derivatives(std::vector<double>& reference_values,
94  int order,
95  const array2d<double>& X) const;
96 
98  void transform_reference_basis(std::vector<double>& values,
99  const std::vector<double>& reference_values,
100  const array2d<double>& X,
101  const std::vector<double>& J,
102  const tcb::span<const double>& detJ,
103  const std::vector<double>& K) const;
104 
107  std::vector<double>& values, std::size_t order,
108  const std::vector<double>& reference_values, const array2d<double>& X,
109  const std::vector<double>& J, const tcb::span<const double>& detJ,
110  const std::vector<double>& K) const;
111 
114  int num_sub_elements() const noexcept;
115 
117  std::size_t hash() const noexcept;
118 
120  std::shared_ptr<const FiniteElement>
121  extract_sub_element(const std::vector<int>& component) const;
122 
129  bool interpolation_ident() const noexcept;
130 
137  array2d<double> interpolation_points() const;
138 
156  template <typename T>
157  constexpr void interpolate(const array2d<T>& values,
158  tcb::span<T> dofs) const
159  {
160  const std::size_t rows = _space_dim / _bs;
161  assert(_space_dim % _bs == 0);
162  assert(dofs.size() == rows);
163 
164  // Compute dofs = Pi * x (matrix-vector multiply)
165  assert(_interpolation_matrix.size() % rows == 0);
166  const std::size_t cols = _interpolation_matrix.size() / rows;
167  for (std::size_t i = 0; i < rows; ++i)
168  {
169  // Dot product between row i of the matrix and 'values'
170  dofs[i] = std::transform_reduce(
171  std::next(_interpolation_matrix.begin(), i * cols),
172  std::next(_interpolation_matrix.begin(), i * cols + cols),
173  values.data(), T(0.0));
174  }
175  }
176 
181  bool needs_permutation_data() const noexcept;
182 
188  template <typename T>
189  void apply_dof_transformation(T* data, std::uint32_t cell_permutation,
190  int block_size) const
191  {
192  if constexpr (std::is_same<T, double>::value)
193  _apply_dof_transformation(data, cell_permutation, block_size);
194  else
195  _apply_dof_transformation_to_scalar(data, cell_permutation, block_size);
196  }
197 
203  template <typename T>
205  T* data, std::uint32_t cell_permutation, int block_size) const
206  {
207  if constexpr (std::is_same<T, double>::value)
208  _apply_inverse_transpose_dof_transformation(data, cell_permutation,
209  block_size);
210  else
211  _apply_inverse_transpose_dof_transformation_to_scalar(
212  data, cell_permutation, block_size);
213  }
214 
217  void map_pull_back(double* physical_data, const double* reference_data,
218  const double* J, const double* detJ, const double* K,
219  const int physical_dim, const int physical_value_size,
220  const int nresults, const int npoints) const;
221 
223  void map_pull_back(std::complex<double>* physical_data,
224  const std::complex<double>* reference_data,
225  const double* J, const double* detJ, const double* K,
226  const int physical_dim, const int physical_value_size,
227  const int nresults, const int npoints) const;
228 
229 private:
230  std::string _signature, _family;
231 
232  mesh::CellType _cell_shape;
233 
234  int _tdim, _space_dim, _value_size, _reference_value_size;
235 
236  // List of sub-elements (if any)
237  std::vector<std::shared_ptr<const FiniteElement>> _sub_elements;
238 
239  // Recursively extract sub finite element
240  static std::shared_ptr<const FiniteElement>
241  extract_sub_element(const FiniteElement& finite_element,
242  const std::vector<int>& component);
243 
244  // Simple hash of the signature string
245  std::size_t _hash;
246 
247  // Dimension of each value space
248  std::vector<int> _value_dimension;
249 
250  std::function<int(double*, const std::uint32_t, const int)>
251  _apply_dof_transformation;
252 
253  std::function<int(ufc_scalar_t*, const std::uint32_t, const int)>
254  _apply_dof_transformation_to_scalar;
255 
256  std::function<int(double*, const std::uint32_t, const int)>
257  _apply_inverse_transpose_dof_transformation;
258 
259  std::function<int(ufc_scalar_t*, const std::uint32_t, const int)>
260  _apply_inverse_transpose_dof_transformation_to_scalar;
261 
262  // Block size for VectorElements and TensorElements. This gives the
263  // number of DOFs colocated at each point.
264  int _bs;
265 
266  // True if interpolation is indentity, i.e. call to
267  // _interpolate_into_cell is not required
268  bool _interpolation_is_ident;
269 
270  // True if element needs dof permutation
271  bool _needs_permutation_data;
272 
273  // The basix element identifier
274  int _basix_element_handle;
275 
276  // The interpolation matrix. It has shape (dim,
277  // num_interpolation_points*basix_value_size). Note that _space_dim =
278  // _bs * dim. Storage is row-major
279  std::vector<double> _interpolation_matrix;
280 };
281 } // namespace dolfinx::fem
This class provides a dynamic 2-dimensional row-wise array data structure.
Definition: array2d.h:21
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:116
FiniteElement & operator=(FiniteElement &&element)=default
Move assignment.
array2d< double > interpolation_points() const
Points on the reference cell at which an expression need to be evaluated in order to interpolate the ...
Definition: FiniteElement.cpp:289
constexpr void interpolate(const array2d< T > &values, tcb::span< T > dofs) const
Definition: FiniteElement.h:157
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:120
std::string family() const noexcept
The finite element family.
Definition: FiniteElement.cpp:140
std::size_t hash() const noexcept
Return simple hash of the signature string.
Definition: FiniteElement.cpp:230
int value_dimension(int i) const
Return the dimension of the value space for axis i.
Definition: FiniteElement.cpp:132
FiniteElement(const FiniteElement &element)=default
Copy constructor.
void apply_inverse_transpose_dof_transformation(T *data, std::uint32_t cell_permutation, int block_size) const
Apply inverse transpose permutation to some data.
Definition: FiniteElement.h:204
void evaluate_reference_basis(std::vector< double > &values, const array2d< double > &X) const
Evaluate all basis functions at given points in reference cell.
Definition: FiniteElement.cpp:142
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:284
std::string signature() const noexcept
String identifying the finite element.
Definition: FiniteElement.cpp:109
std::shared_ptr< const FiniteElement > extract_sub_element(const std::vector< int > &component) const
Extract sub finite element for component.
Definition: FiniteElement.cpp:233
void map_pull_back(double *physical_data, const double *reference_data, const double *J, const double *detJ, const double *K, const int physical_dim, const int physical_value_size, const int nresults, const int npoints) const
Pull physical data back to the reference element. This passes the inputs directly into Basix's map_pu...
Definition: FiniteElement.cpp:309
FiniteElement(const ufc_finite_element &ufc_element)
Create finite element from UFC finite element.
Definition: FiniteElement.cpp:17
void transform_reference_basis(std::vector< double > &values, const std::vector< double > &reference_values, const array2d< double > &X, const std::vector< double > &J, const tcb::span< const double > &detJ, const std::vector< double > &K) const
Push basis functions forward to physical element.
Definition: FiniteElement.cpp:208
int value_size() const noexcept
The value size, e.g. 1 for a scalar function, 2 for a 2D vector.
Definition: FiniteElement.cpp:118
int block_size() const noexcept
Block size of the finite element function space. For VectorElements and TensorElements,...
Definition: FiniteElement.cpp:130
int value_rank() const noexcept
Rank of the value space.
Definition: FiniteElement.cpp:125
void evaluate_reference_basis_derivatives(std::vector< double > &reference_values, int order, const array2d< double > &X) const
Evaluate all basis function derivatives of given order at given points in reference cell.
Definition: FiniteElement.cpp:171
bool needs_permutation_data() const noexcept
Definition: FiniteElement.cpp:304
mesh::CellType cell_shape() const noexcept
Cell shape.
Definition: FiniteElement.cpp:111
void transform_reference_basis_derivatives(std::vector< double > &values, std::size_t order, const std::vector< double > &reference_values, const array2d< double > &X, const std::vector< double > &J, const tcb::span< const double > &detJ, const std::vector< double > &K) const
Push basis function (derivatives) forward to physical element.
virtual ~FiniteElement()=default
Destructor.
int num_sub_elements() const noexcept
Get the number of sub elements (for a mixed element)
Definition: FiniteElement.cpp:225
FiniteElement(FiniteElement &&element)=default
Move constructor.
void apply_dof_transformation(T *data, std::uint32_t cell_permutation, int block_size) const
Apply permutation to some data.
Definition: FiniteElement.h:189
Finite element method functionality.
Definition: assemble_matrix_impl.h:23
CellType
Cell type identifier.
Definition: cell_types.h:22