Basix
Public Member Functions | List of all members
basix::FiniteElement Class Reference

#include <finite-element.h>

Public Member Functions

 FiniteElement (element::family family, 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_transformations, const Eigen::ArrayXXd &points, const Eigen::MatrixXd M={}, mapping::type map_type=mapping::type::identity)
 
 FiniteElement (const FiniteElement &element)=default
 Copy constructor.
 
 FiniteElement (FiniteElement &&element)=default
 Move constructor.
 
 ~FiniteElement ()=default
 Destructor.
 
FiniteElementoperator= (const FiniteElement &element)=default
 Assignment operator.
 
FiniteElementoperator= (FiniteElement &&element)=default
 Move assignment operator.
 
std::vector< Eigen::ArrayXXd > tabulate (int nd, const Eigen::ArrayXXd &x) const
 
void tabulate (int nd, const Eigen::ArrayXXd &x, double *basis_data) const
 
cell::type cell_type () const
 
int degree () const
 
int value_size () const
 
const std::vector< int > & value_shape () const
 
int dim () const
 
element::family family () const
 
const mapping::type mapping_type () const
 
Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > map_push_forward (const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &reference_data, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &J, const tcb::span< const double > &detJ, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &K) const
 
template<typename T >
void map_push_forward_m (const Eigen::Array< T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &reference_data, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &J, const tcb::span< const double > &detJ, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &K, T *physical_data) const
 
Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > map_pull_back (const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &physical_data, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &J, const tcb::span< const double > &detJ, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &K) const
 
template<typename T >
void map_pull_back_m (const Eigen::Array< T, Eigen::Dynamic, Eigen::Dynamic > &physical_data, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &J, const tcb::span< const double > &detJ, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &K, T *reference_data) const
 
const std::vector< std::vector< int > > & entity_dofs () const
 
std::vector< Eigen::MatrixXd > base_permutations () const
 
const Eigen::ArrayXXd & points () const
 Return a set of interpolation points.
 
int num_points () const
 Return the number of interpolation points.
 
const Eigen::MatrixXd & interpolation_matrix () const
 

Detailed Description

Finite Element The basis is stored as a set of coefficients, which are applied to the underlying expansion set for that cell type, when tabulating.

Constructor & Destructor Documentation

◆ FiniteElement()

FiniteElement::FiniteElement ( element::family  family,
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_transformations,
const Eigen::ArrayXXd &  points,
const Eigen::MatrixXd  M = {},
mapping::type  map_type = mapping::type::identity 
)
Todo:
Document A finite element
Parameters
[in]family
[in]cell_type
[in]degree
[in]value_shape
[in]coeffs
[in]entity_dofs
[in]base_transformationsBase transformations
[in]points
[in]MThe interpolation matrix
[in]map_type

Member Function Documentation

◆ base_permutations()

std::vector< Eigen::MatrixXd > FiniteElement::base_permutations ( ) const

Get the base transformations The base transformations represent the effect of rotating or reflecting a subentity of the cell on the numbering and orientation of the DOFs. This returns a list of matrices with one matrix for each subentity permutation in the following order: Reversing edge 0, reversing edge 1, ... Rotate face 0, reflect face 0, rotate face 1, reflect face 1, ...

Example: Order 3 Lagrange on a triangle

This space has 10 dofs arranged like:

2
|\
6 4
| \
5 9 3
| \
0-7-8-1

For this element, the base transformations are: [Matrix swapping 3 and 4, Matrix swapping 5 and 6, Matrix swapping 7 and 8] The first row shows the effect of reversing the diagonal edge. The second row shows the effect of reversing the vertical edge. The third row shows the effect of reversing the horizontal edge.

Example: Order 1 Raviart-Thomas on a triangle

This space has 3 dofs arranged like:

|\
| \
| \
<-1 0
| / \
| L ^ \
| | \
---2---

These DOFs are integrals of normal components over the edges: DOFs 0 and 2 are oriented inward, DOF 1 is oriented outwards. For this element, the base transformation matrices are:

0: [[-1, 0, 0],
[ 0, 1, 0],
[ 0, 0, 1]]
1: [[1, 0, 0],
[0, -1, 0],
[0, 0, 1]]
2: [[1, 0, 0],
[0, 1, 0],
[0, 0, -1]]

The first matrix reverses DOF 0 (as this is on the first edge). The second matrix reverses DOF 1 (as this is on the second edge). The third matrix reverses DOF 2 (as this is on the third edge).

Example: DOFs on the face of Order 2 Nedelec first kind on a tetrahedron

On a face of this tetrahedron, this space has two face tangent DOFs:

|\ |\
| \ | \
| \ | ^\
| \ | | \
| 0->\ | 1 \
| \ | \
------ ------

For these DOFs, the subblocks of the base transformation matrices are:

rotation: [[-1, 1],
[ 1, 0]]
reflection: [[0, 1],
[1, 0]]

◆ cell_type()

cell::type FiniteElement::cell_type ( ) const

Get the element cell type

Returns
The cell type

◆ degree()

int FiniteElement::degree ( ) const

Get the element polynomial degree

Returns
Polynomial degree

◆ dim()

int FiniteElement::dim ( ) const

Dimension of the finite element space (number of degrees of freedom for the element)

Returns
Number of degrees of freedom

◆ entity_dofs()

const std::vector< std::vector< int > > & FiniteElement::entity_dofs ( ) const

Get the number of dofs on each topological entity: (vertices, edges, faces, cell) in that order. For example, Lagrange degree 2 on a triangle has vertices: [1, 1, 1], edges: [1, 1, 1], cell: [0] The sum of the entity dofs must match the total number of dofs reported by FiniteElement::dim,

Returns
List of entity dof counts on each dimension

◆ family()

element::family FiniteElement::family ( ) const

Get the finite element family

Returns
The family

◆ interpolation_matrix()

const Eigen::MatrixXd & FiniteElement::interpolation_matrix ( ) const

Return a matrix of weights interpolation To interpolate a function in this finite element, the functions should be evaluated at each point given by FiniteElement::points(). These function values should then be multiplied by the weight matrix to give the coefficients of the interpolated function.

◆ map_pull_back()

Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > FiniteElement::map_pull_back ( const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  physical_data,
const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  J,
const tcb::span< const double > &  detJ,
const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  K 
) const

Map function values from a physical cell to the reference

Parameters
physical_dataThe function values on the cell
JThe Jacobian of the mapping
detJThe determinant of the Jacobian of the mapping
KThe inverse of the Jacobian of the mapping
Returns
The function values on the reference

◆ map_pull_back_m()

template<typename T >
void basix::FiniteElement::map_pull_back_m ( const Eigen::Array< T, Eigen::Dynamic, Eigen::Dynamic > &  physical_data,
const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  J,
const tcb::span< const double > &  detJ,
const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  K,
T *  reference_data 
) const

Map function values from a physical cell to the reference

Parameters
physical_dataThe function values on the cell
JThe Jacobian of the mapping
detJThe determinant of the Jacobian of the mapping
KThe inverse of the Jacobian of the mapping
reference_dataMemory location to fill

◆ map_push_forward()

Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > FiniteElement::map_push_forward ( const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  reference_data,
const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  J,
const tcb::span< const double > &  detJ,
const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  K 
) const

Map function values from the reference to a physical cell

Parameters
reference_dataThe function values on the reference
JThe Jacobian of the mapping
detJThe determinant of the Jacobian of the mapping
KThe inverse of the Jacobian of the mapping
Returns
The function values on the cell

◆ map_push_forward_m()

template<typename T >
void basix::FiniteElement::map_push_forward_m ( const Eigen::Array< T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  reference_data,
const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  J,
const tcb::span< const double > &  detJ,
const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &  K,
T *  physical_data 
) const

Direct to memory push forward

Parameters
reference_dataThe function values on the reference
JThe Jacobian of the mapping
detJThe determinant of the Jacobian of the mapping
KThe inverse of the Jacobian of the mapping
physical_dataMemory location to fill

◆ mapping_type()

const mapping::type FiniteElement::mapping_type ( ) const

Get the mapping type used for this element

Returns
The mapping

◆ tabulate() [1/2]

std::vector< Eigen::ArrayXXd > FiniteElement::tabulate ( int  nd,
const Eigen::ArrayXXd &  x 
) const

Compute basis values and derivatives at set of points.

Parameters
[in]ndThe order of derivatives, up to and including, to compute. Use 0 for the basis functions only.
[in]xThe points at which to compute the basis functions. The shape of x is (number of points, geometric dimension).
Returns
The basis functions (and derivatives). The first entry in the list is the basis function. Higher derivatives are stored in triangular (2D) or tetrahedral (3D) ordering, i.e. for the (x,y) derivatives in 2D: (0,0),(1,0),(0,1),(2,0),(1,1),(0,2),(3,0)... The function basix::idx can be used to find the appropriate derivative. If a vector result is expected, it will be stacked with all x values, followed by all y-values (and then z, if any), likewise tensor-valued results will be stacked in index order.

◆ tabulate() [2/2]

void FiniteElement::tabulate ( int  nd,
const Eigen::ArrayXXd &  x,
double *  basis_data 
) const

Direct to memory block tabulation

Parameters
ndNumber of derivatives
xPoints
basis_dataMemory location to fill

◆ value_shape()

const std::vector< int > & FiniteElement::value_shape ( ) const

Get the element value tensor shape, e.g. returning [1] for scalars.

Returns
Value shape

◆ value_size()

int FiniteElement::value_size ( ) const

Get the element value size This is just a convenience function returning product(value_shape)

Returns
Value size

The documentation for this class was generated from the following files: