Basix 0.8.0

Home     Installation     Demos     C++ docs     Python docs

Loading...
Searching...
No Matches
basix Namespace Reference

Basix: FEniCS runtime basis evaluation library. More...

Namespaces

namespace  cell
 Information about reference cells.
 
namespace  doftransforms
 Functions to transform DOFs in high degree Lagrange spaces. The functions in this namespace calculate the permutations that can be used to rotate and reflect DOF points in Lagrange spaces.
 
namespace  element
 Interfaces for creating finite elements.
 
namespace  indexing
 Indexing.
 
namespace  lattice
 Lattices of points.
 
namespace  maps
 Information about finite element maps.
 
namespace  math
 Mathematical functions.
 
namespace  moments
 Functions to create integral moment DOFs.
 
namespace  polynomials
 Polynomials.
 
namespace  polyset
 Polynomial expansion sets.
 
namespace  precompute
 Matrix and permutation pre-computation.
 
namespace  quadrature
 Quadrature rules.
 
namespace  sobolev
 Information about Sobolev spaces.
 

Classes

class  FiniteElement
 A finite element. More...
 

Functions

template<std::floating_point T>
FiniteElement< T > create_custom_element (cell::type cell_type, const std::vector< std::size_t > &value_shape, impl::mdspan_t< const T, 2 > wcoeffs, const std::array< std::vector< impl::mdspan_t< const T, 2 > >, 4 > &x, const std::array< std::vector< impl::mdspan_t< const T, 4 > >, 4 > &M, int interpolation_nderivs, maps::type map_type, sobolev::space sobolev_space, bool discontinuous, int embedded_subdegree, int embedded_superdegree, polyset::type poly_type)
 
template<std::floating_point T>
FiniteElement< T > create_element (element::family family, cell::type cell, int degree, element::lagrange_variant lvariant, element::dpc_variant dvariant, bool discontinuous, std::vector< int > dof_ordering={})
 
std::vector< int > tp_dof_ordering (element::family family, cell::type cell, int degree, element::lagrange_variant lvariant, element::dpc_variant dvariant, bool discontinuous)
 
template<std::floating_point T>
std::vector< std::vector< FiniteElement< T > > > tp_factors (element::family family, cell::type cell, int degree, element::lagrange_variant lvariant, element::dpc_variant dvariant, bool discontinuous, std::vector< int > dof_ordering)
 
template<std::floating_point T>
FiniteElement< T > create_tp_element (element::family family, cell::type cell, int degree, element::lagrange_variant lvariant, element::dpc_variant dvariant, bool discontinuous)
 
std::string version ()
 
template<std::floating_point T>
std::pair< std::vector< T >, std::array< std::size_t, 2 > > compute_interpolation_operator (const FiniteElement< T > &element_from, const FiniteElement< T > &element_to)
 Compute a matrix that represents the interpolation between two elements.
 

Detailed Description

Basix: FEniCS runtime basis evaluation library.

Function Documentation

◆ create_custom_element()

template<std::floating_point T>
FiniteElement< T > basix::create_custom_element ( cell::type  cell_type,
const std::vector< std::size_t > &  value_shape,
impl::mdspan_t< const T, 2 >  wcoeffs,
const std::array< std::vector< impl::mdspan_t< const T, 2 > >, 4 > &  x,
const std::array< std::vector< impl::mdspan_t< const T, 4 > >, 4 > &  M,
int  interpolation_nderivs,
maps::type  map_type,
sobolev::space  sobolev_space,
bool  discontinuous,
int  embedded_subdegree,
int  embedded_superdegree,
polyset::type  poly_type 
)

Create a custom finite element

Parameters
[in]cell_typeThe cell type
[in]value_shapeThe value shape of the element
[in]wcoeffsMatrices for the kth value index containing the expansion coefficients defining a polynomial basis spanning the polynomial space for this element. Shape is (dim(finite element polyset), dim(Legendre polynomials))
[in]xInterpolation points. Indices are (tdim, entity index, point index, dim)
[in]MThe interpolation matrices. Indices are (tdim, entity index, dof, vs, point_index, derivative)
[in]interpolation_nderivsThe number of derivatives that need to be used during interpolation
[in]map_typeThe type of map to be used to map values from the reference to a cell
[in]sobolev_spaceThe underlying Sobolev space for the element
[in]discontinuousIndicates whether or not this is the discontinuous version of the element
[in]embedded_subdegreeThe highest degree n such that a Lagrange (or vector Lagrange) element of degree n is a subspace of this element
[in]embedded_superdegreeThe degree of a polynomial in this element's polyset
[in]poly_typeThe type of polyset to use for this element
Returns
A custom finite element

◆ create_element()

template<std::floating_point T>
template basix::FiniteElement< double > basix::create_element ( element::family  family,
cell::type  cell,
int  degree,
element::lagrange_variant  lvariant,
element::dpc_variant  dvariant,
bool  discontinuous,
std::vector< int >  dof_ordering = {} 
)

Create an element using a given Lagrange variant and a given DPC variant

Parameters
[in]familyThe element family
[in]cellThe reference cell type that the element is defined on
[in]degreeThe degree of the element
[in]lvariantThe variant of Lagrange to use
[in]dvariantThe variant of DPC to use
[in]discontinuousIndicates whether the element is discontinuous between cells points of the element. The discontinuous element will have the same DOFs, but they will all be associated with the interior of the cell.
[in]dof_orderingOrdering of dofs for ElementDofLayout
Returns
A finite element

◆ tp_dof_ordering()

std::vector< int > basix::tp_dof_ordering ( element::family  family,
cell::type  cell,
int  degree,
element::lagrange_variant  lvariant,
element::dpc_variant  dvariant,
bool  discontinuous 
)

Get the tensor product DOF ordering for an element

Parameters
[in]familyThe element family
[in]cellThe reference cell type that the element is defined on. Currently limited to quadrilateral or hexahedron.
[in]degreeThe degree of the element
[in]lvariantThe variant of Lagrange to use
[in]dvariantThe variant of DPC to use
[in]discontinuousIndicates whether the element is discontinuous between cells points of the element. The discontinuous element will have the same DOFs, but they will all be associated with the interior of the cell.
Returns
A vector containing the dof ordering

◆ tp_factors()

template<std::floating_point T>
template std::vector< std::vector< basix::FiniteElement< double > > > basix::tp_factors ( element::family  family,
cell::type  cell,
int  degree,
element::lagrange_variant  lvariant,
element::dpc_variant  dvariant,
bool  discontinuous,
std::vector< int >  dof_ordering 
)

Get the tensor factors of an element

Parameters
[in]familyThe element family
[in]cellThe reference cell type that the element is defined on. Currently limited to quadrilateral or hexahedron.
[in]degreeThe degree of the element
[in]lvariantThe variant of Lagrange to use
[in]dvariantThe variant of DPC to use
[in]discontinuousIndicates whether the element is discontinuous between cells points of the element. The discontinuous element will have the same DOFs, but they will all be associated with the interior of the cell.
[in]dof_orderingThe ordering of the DOFs
Returns
A list of lists of finite element factors

◆ create_tp_element()

template<std::floating_point T>
template basix::FiniteElement< double > basix::create_tp_element ( element::family  family,
cell::type  cell,
int  degree,
element::lagrange_variant  lvariant,
element::dpc_variant  dvariant,
bool  discontinuous 
)

Create an element with Tensor Product dof ordering

Parameters
[in]familyThe element family
[in]cellThe reference cell type that the element is defined on. Currently limited to quadrilateral or hexahedron.
[in]degreeThe degree of the element
[in]lvariantThe variant of Lagrange to use
[in]dvariantThe variant of DPC to use
[in]discontinuousIndicates whether the element is discontinuous between cells points of the element. The discontinuous element will have the same DOFs, but they will all be associated with the interior of the cell.
Returns
A finite element

◆ version()

std::string basix::version ( )

Return the Basix version number

Returns
version string

◆ compute_interpolation_operator()

template<std::floating_point T>
std::pair< std::vector< T >, std::array< std::size_t, 2 > > basix::compute_interpolation_operator ( const FiniteElement< T > &  element_from,
const FiniteElement< T > &  element_to 
)

Compute a matrix that represents the interpolation between two elements.

If the two elements have the same value size, this function returns the interpolation between them.

If element_from has value size 1 and element_to has value size > 1, then this function returns a matrix to interpolate from a blocked element_from (ie multiple copies of element_from) into element_to.

If element_to has value size 1 and element_from has value size > 1, then this function returns a matrix that interpolates the components of element_from into copies of element_to.

Note
If the elements have different value sizes and both are greater than 1, this function throws a runtime error

In order to interpolate functions between finite element spaces on arbitrary cells, the functions must be pulled back to the reference element (this pull back includes applying DOF transformations). The matrix that this function returns can then be applied, then the result pushed forward to the cell. If element_from and element_to have the same map type, then only the DOF transformations need to be applied, as the pull back and push forward cancel each other out.

Parameters
[in]element_fromThe element to interpolate from
[in]element_toThe element to interpolate to
Returns
Matrix operator that maps the 'from' degrees-of-freedom to the 'to' degrees-of-freedom. Shape is (ndofs(element_to), ndofs(element_from))