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/types.h>
10 #include <dolfinx/mesh/cell_types.h>
11 #include <functional>
12 #include <memory>
13 #include <unsupported/Eigen/CXX11/Tensor>
14 #include <vector>
15 
16 struct ufc_coordinate_mapping;
17 struct ufc_finite_element;
18 
19 namespace dolfinx::fem
20 {
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;
49 
52  mesh::CellType cell_shape() const;
53 
56  int space_dimension() const;
57 
62  int block_size() const;
63 
66  int value_size() const;
67 
71  int reference_value_size() const;
72 
75  int value_rank() const;
76 
78  int value_dimension(int i) const;
79 
82  std::string family() const;
83 
85  // reference_values[num_points][num_dofs][reference_value_size]
87  Eigen::Tensor<double, 3, Eigen::RowMajor>& reference_values,
88  const Eigen::Ref<const Eigen::Array<
89  double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>>& X) const;
90 
93  // reference_value_derivatives[num_points][num_dofs][reference_value_size][num_derivatives]
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;
98 
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 Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic, 1>>& detJ,
107  const Eigen::Tensor<double, 3, Eigen::RowMajor>& K,
108  const std::uint32_t permutation_info) const;
109 
112  Eigen::Tensor<double, 4, Eigen::RowMajor>& values, std::size_t order,
113  const Eigen::Tensor<double, 4, Eigen::RowMajor>& reference_values,
114  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
115  Eigen::Dynamic, Eigen::RowMajor>>& X,
116  const Eigen::Tensor<double, 3, Eigen::RowMajor>& J,
117  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic, 1>>& detJ,
118  const Eigen::Tensor<double, 3, Eigen::RowMajor>& K,
119  const std::uint32_t permutation_info) const;
120 
125  const Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>&
127 
130  void transform_values(
131  ufc_scalar_t* reference_values,
132  const Eigen::Ref<const Eigen::Array<ufc_scalar_t, Eigen::Dynamic,
133  Eigen::Dynamic, Eigen::RowMajor>>&
134  physical_values,
135  const Eigen::Ref<const Eigen::Array<double, Eigen::Dynamic,
136  Eigen::Dynamic, Eigen::RowMajor>>&
137  coordinate_dofs) const;
138 
140  int num_sub_elements() const;
141 
143  std::size_t hash() const;
144 
146  std::shared_ptr<const FiniteElement>
147  extract_sub_element(const std::vector<int>& component) const;
148 
149 private:
150  std::string _signature, _family;
151 
152  mesh::CellType _cell_shape;
153 
154  int _tdim, _space_dim, _value_size, _reference_value_size;
155 
156  // Dof coordinates on the reference element
157  Eigen::Array<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> _refX;
158 
159  // List of sub-elements (if any)
160  std::vector<std::shared_ptr<const FiniteElement>> _sub_elements;
161 
162  // Recursively extract sub finite element
163  static std::shared_ptr<const FiniteElement>
164  extract_sub_element(const FiniteElement& finite_element,
165  const std::vector<int>& component);
166 
167  // Simple hash of the signature string
168  std::size_t _hash;
169 
170  // Dimension of each value space
171  std::vector<int> _value_dimension;
172 
173  // Functions for basis and derivatives evaluation
174  std::function<int(double*, int, const double*)> _evaluate_reference_basis;
175 
176  std::function<int(double*, int, int, const double*)>
177  _evaluate_reference_basis_derivatives;
178 
179  std::function<int(double*, int, int, const double*, const double*,
180  const double*, const double*, const double*,
181  const std::uint32_t)>
182  _transform_reference_basis_derivatives;
183 
184  std::function<int(ufc_scalar_t*, const ufc_scalar_t*, const double*,
185  const ufc_coordinate_mapping*)>
186  _transform_values;
187 
188  // Block size for VectorElements and TensorElements
189  // This gives the number of DOFs colocated at each point
190  int _block_size;
191 };
192 } // namespace dolfinx::fem
dolfinx::fem::FiniteElement::FiniteElement
FiniteElement(const ufc_finite_element &ufc_element)
Create finite element from UFC finite element.
Definition: FiniteElement.cpp:19
dolfinx::fem::FiniteElement::extract_sub_element
std::shared_ptr< const FiniteElement > extract_sub_element(const std::vector< int > &component) const
Extract sub finite element for component.
Definition: FiniteElement.cpp:216
dolfinx::fem::FiniteElement::value_dimension
int value_dimension(int i) const
Return the dimension of the value space for axis i.
Definition: FiniteElement.cpp:99
dolfinx::fem::FiniteElement::evaluate_reference_basis
void evaluate_reference_basis(Eigen::Tensor< double, 3, Eigen::RowMajor > &reference_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:108
dolfinx::fem::FiniteElement::num_sub_elements
int num_sub_elements() const
Return the number of sub elements (for a mixed element)
Definition: FiniteElement.cpp:211
dolfinx::fem::FiniteElement::evaluate_reference_basis_derivatives
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:124
dolfinx::fem::FiniteElement::FiniteElement
FiniteElement(FiniteElement &&element)=default
Move constructor.
dolfinx::fem::FiniteElement::dof_reference_coordinates
const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > & dof_reference_coordinates() const
Tabulate the reference coordinates of all dofs on an element.
Definition: FiniteElement.cpp:185
dolfinx::mesh::CellType
CellType
Cell type identifier.
Definition: cell_types.h:21
dolfinx::fem::FiniteElement::cell_shape
mesh::CellType cell_shape() const
Cell shape.
Definition: FiniteElement.cpp:84
dolfinx::fem::FiniteElement::operator=
FiniteElement & operator=(FiniteElement &&element)=default
Move assignment.
dolfinx::fem::FiniteElement
Finite Element, containing the dof layout on a reference element, and various methods for evaluating ...
Definition: FiniteElement.h:25
dolfinx::fem::FiniteElement::transform_reference_basis
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 Eigen::Ref< const Eigen::Array< double, Eigen::Dynamic, 1 >> &detJ, const Eigen::Tensor< double, 3, Eigen::RowMajor > &K, const std::uint32_t permutation_info) const
Push basis functions forward to physical element.
Definition: FiniteElement.cpp:140
dolfinx::fem::FiniteElement::value_rank
int value_rank() const
Rank of the value space.
Definition: FiniteElement.cpp:95
dolfinx::fem::FiniteElement::hash
std::size_t hash() const
Return simple hash of the signature string.
Definition: FiniteElement.cpp:213
dolfinx::fem::FiniteElement::~FiniteElement
virtual ~FiniteElement()=default
Destructor.
dolfinx::fem::FiniteElement::transform_reference_basis_derivatives
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 Eigen::Ref< const Eigen::Array< double, Eigen::Dynamic, 1 >> &detJ, const Eigen::Tensor< double, 3, Eigen::RowMajor > &K, const std::uint32_t permutation_info) const
Push basis function (derivatives) forward to physical element.
Definition: FiniteElement.cpp:162
dolfinx::fem::FiniteElement::space_dimension
int space_dimension() const
Dimension of the finite element function space.
Definition: FiniteElement.cpp:86
dolfinx::fem::FiniteElement::operator=
FiniteElement & operator=(const FiniteElement &element)=default
Copy assignment.
dolfinx::fem::FiniteElement::transform_values
void transform_values(ufc_scalar_t *reference_values, const Eigen::Ref< const Eigen::Array< ufc_scalar_t, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &physical_values, const Eigen::Ref< const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &coordinate_dofs) const
Map values of field from physical to reference space which has been evaluated at points given by dof_...
Definition: FiniteElement.cpp:196
dolfinx::fem::FiniteElement::signature
std::string signature() const
String identifying the finite element.
Definition: FiniteElement.cpp:82
dolfinx::fem::FiniteElement::family
std::string family() const
The finite element family.
Definition: FiniteElement.cpp:106
dolfinx::fem::FiniteElement::block_size
int block_size() const
Block size of the finite element function space. For VectorElements and TensorElements,...
Definition: FiniteElement.cpp:97
dolfinx::fem
Finite element method functionality.
Definition: assemble_matrix_impl.h:23
dolfinx::fem::FiniteElement::value_size
int value_size() const
The value size, e.g. 1 for a scalar function, 2 for a 2D vector.
Definition: FiniteElement.cpp:88
dolfinx::fem::FiniteElement::reference_value_size
int reference_value_size() const
The value size, e.g. 1 for a scalar function, 2 for a 2D vector for the reference element.
Definition: FiniteElement.cpp:90
dolfinx::fem::FiniteElement::FiniteElement
FiniteElement(const FiniteElement &element)=default
Copy constructor.