Tabulation Project basix
finite-element.h
1 // FIXME: just include everything for now
2 // Need to define public API
3 
4 #pragma once
5 
6 #include "cell.h"
7 #include <Eigen/Dense>
8 #include <string>
9 #include <vector>
10 
11 namespace basix
12 {
13 
139 Eigen::MatrixXd
140 compute_expansion_coefficients(const Eigen::MatrixXd& span_coeffs,
141  const Eigen::MatrixXd& dual,
142  bool condition_check = false);
143 
148 {
149 
150 public:
153  const std::vector<int>& value_shape,
154  const Eigen::ArrayXXd& coeffs,
155  const std::vector<std::vector<int>>& entity_dofs,
156  const std::vector<Eigen::MatrixXd>& base_permutations,
157  const Eigen::ArrayXXd& points,
158  const Eigen::MatrixXd interpolation_matrix = {},
159  const std::string mapping_name = "affine");
160 
162  FiniteElement(const FiniteElement& element) = default;
163 
165  FiniteElement(FiniteElement&& element) = default;
166 
168  ~FiniteElement() = default;
169 
171  FiniteElement& operator=(const FiniteElement& element) = default;
172 
174  FiniteElement& operator=(FiniteElement&& element) = default;
175 
190  std::vector<Eigen::ArrayXXd> tabulate(int nd, const Eigen::ArrayXXd& x) const;
191 
194  cell::type cell_type() const;
195 
198  int degree() const;
199 
203  int value_size() const;
204 
207  const std::vector<int>& value_shape() const;
208 
212  int dim() const;
213 
216  const std::string& family_name() const;
217 
220  const std::string& mapping_name() const;
221 
228  const std::vector<std::vector<int>>& entity_dofs() const;
229 
308  std::vector<Eigen::MatrixXd> base_permutations() const;
309 
313  const Eigen::ArrayXXd& points() const;
314 
320  const Eigen::MatrixXd& interpolation_matrix() const;
321 
322 private:
323  // Cell type
324  cell::type _cell_type;
325 
326  // The name of the finite element family
327  std::string _family_name;
328 
329  // Degree
330  int _degree;
331 
332  // Value shape
333  std::vector<int> _value_shape;
334 
336  std::string _mapping_name;
337 
338  // Shape function coefficient of expansion sets on cell. If shape
339  // function is given by @f$\psi_i = \sum_{k} \phi_{k} \alpha^{i}_{k}@f$,
340  // then _coeffs(i, j) = @f$\alpha^i_k@f$. i.e., _coeffs.row(i) are the
341  // expansion coefficients for shape function i (@f$\psi_{i}@f$).
342  Eigen::MatrixXd _coeffs;
343 
344  // Number of dofs associated each subentity
345  // The dofs of an element are associated with entities of different
346  // topological dimension (vertices, edges, faces, cells). The dofs are listed
347  // in this order, with vertex dofs first. Each entry is the dof count on the
348  // associated entity, as listed by cell::topology.
349  std::vector<std::vector<int>> _entity_dofs;
350 
351  // Base permutations
352  std::vector<Eigen::MatrixXd> _base_permutations;
353 
354  // Set of points used for point evaluation
355  // Experimental - currently used for an implementation of
356  // "tabulate_dof_coordinates" Most useful for Lagrange. This may change or go
357  // away. For non-Lagrange elements, these points will be used in combination
358  // with _interpolation_matrix to perform interpolation
359  Eigen::ArrayXXd _points;
360 
362  Eigen::MatrixXd _interpolation_matrix;
363 };
364 
366 FiniteElement create_element(std::string family, std::string cell, int degree);
367 
370 const std::string& version();
371 
372 } // namespace basix
Definition: finite-element.h:148
FiniteElement & operator=(const FiniteElement &element)=default
Assignment operator.
cell::type cell_type() const
Definition: finite-element.cpp:102
FiniteElement & operator=(FiniteElement &&element)=default
Move assignment operator.
const Eigen::ArrayXXd & points() const
Definition: finite-element.cpp:168
int value_size() const
Definition: finite-element.cpp:106
FiniteElement(FiniteElement &&element)=default
Move constructor.
FiniteElement(std::string family_name, cell::type cell_type, int degree, const std::vector< int > &value_shape, const Eigen::ArrayXXd &coeffs, const std::vector< std::vector< int >> &entity_dofs, const std::vector< Eigen::MatrixXd > &base_permutations, const Eigen::ArrayXXd &points, const Eigen::MatrixXd interpolation_matrix={}, const std::string mapping_name="affine")
A finite element.
Definition: finite-element.cpp:78
int degree() const
Definition: finite-element.cpp:104
const std::vector< int > & value_shape() const
Definition: finite-element.cpp:114
~FiniteElement()=default
Destructor.
const std::vector< std::vector< int > > & entity_dofs() const
Definition: finite-element.cpp:130
int dim() const
Definition: finite-element.cpp:119
const std::string & family_name() const
Definition: finite-element.cpp:121
const Eigen::MatrixXd & interpolation_matrix() const
Definition: finite-element.cpp:125
FiniteElement(const FiniteElement &element)=default
Copy constructor.
std::vector< Eigen::ArrayXXd > tabulate(int nd, const Eigen::ArrayXXd &x) const
Definition: finite-element.cpp:136
const std::string & mapping_name() const
Definition: finite-element.cpp:123
std::vector< Eigen::MatrixXd > base_permutations() const
Definition: finite-element.cpp:163
type
Cell type.
Definition: cell.h:21
basix
Definition: basix.h:7
Eigen::MatrixXd compute_expansion_coefficients(const Eigen::MatrixXd &span_coeffs, const Eigen::MatrixXd &dual, bool condition_check=false)
Definition: finite-element.cpp:46
FiniteElement create_element(std::string family, std::string cell, int degree)
Create an element by name.
Definition: finite-element.cpp:22
const std::string & version()
Definition: finite-element.cpp:170
int degree(int handle)
Degree.
Definition: basix.cpp:50