basix.ufl

Functions to directly wrap Basix elements in UFL.

Functions

_compute_signature(element)

Compute a signature of a custom element.

_ufl_pullback_from_enum(m)

Convert an enum to a UFL pull back.

_ufl_sobolev_space_from_enum(s)

Convert a Basix Sobolev space enum to a UFL Sobolev space.

blocked_element(sub_element, shape[, symmetry])

Create a UFL compatible blocked element.

custom_element(cell_type, ...[, polyset_type])

Create a UFL compatible custom Basix element.

element(family, cell, degree[, ...])

Create a UFL compatible element using Basix.

enriched_element(elements[, map_type])

Create an UFL compatible enriched element from a list of elements.

mixed_element(elements)

Create a UFL compatible mixed element from a list of elements.

quadrature_element(cell[, value_shape, ...])

Create a quadrature element.

real_element(cell, value_shape)

Create a real element.

wrap_element(element)

Wrap a Basix element as a Basix UFL element.

Classes

_BasixElement(element)

A wrapper allowing Basix elements to be used directly with UFL.

_BlockedElement(sub_element, shape[, symmetry])

Element with a block size that contains multiple copies of a sub element.

_ComponentElement(element, component)

An element representing one component of a _BasixElement.

_ElementBase(repr, cellname, ...[, degree, ...])

A base wrapper to allow elements to be used with UFL.

_MixedElement(sub_elements)

A mixed element that combines two or more elements.

_QuadratureElement(cell, points, weights, ...)

A quadrature element.

_RealElement(cell, value_shape)

A real element.

basix.ufl.blocked_element(sub_element: _ElementBase, shape: tuple[int, ...], symmetry: bool | None = None) _ElementBase

Create a UFL compatible blocked element.

Parameters:
  • sub_element – Element used for each block.

  • shape – Value shape of the element. For scalar-valued families, this can be used to create vector and tensor elements.

  • symmetry – Set to True if the tensor is symmetric. Valid for rank 2 elements only.

Returns:

A blocked finite element.

basix.ufl.custom_element(cell_type: CellType, reference_value_shape: list[int] | tuple[int, ...], wcoeffs: ndarray[Any, dtype[float64]], x: list[list[numpy.ndarray[Any, numpy.dtype[numpy.float64]]]], M: list[list[numpy.ndarray[Any, numpy.dtype[numpy.float64]]]], interpolation_nderivs: int, map_type: MapType, sobolev_space: SobolevSpace, discontinuous: bool, embedded_subdegree: int, embedded_superdegree: int, polyset_type: PolysetType = PolysetType.standard) _ElementBase

Create a UFL compatible custom Basix element.

Parameters:
  • cell_type – The cell type

  • reference_value_shape – The reference value shape of the element

  • wcoeffs – Matrices 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(Legenre polynomials)).

  • x – Interpolation points. Indices are (tdim, entity index, point index, dim).

  • M – The interpolation matrices. Indices are (tdim, entity index, dof, vs, point_index, derivative).

  • interpolation_nderivs – The number of derivatives that need to be used during interpolation.

  • map_type – The type of map to be used to map values from the reference to a cell.

  • sobolev_space – Underlying Sobolev space for the element.

  • discontinuous – Indicates whether or not this is the discontinuous version of the element.

  • embedded_subdegree – The highest degree n such that a Lagrange (or vector Lagrange) element of degree n is a subspace of this element.

  • embedded_superdegree – The highest degree of a polynomial in this element’s polyset.

  • polyset_type – Polyset type for the element.

Returns:

A custom finite element.

basix.ufl.element(family: ~basix.finite_element.ElementFamily | str, cell: ~basix.cell.CellType | str, degree: int, lagrange_variant: ~basix.finite_element.LagrangeVariant = LagrangeVariant.unset, dpc_variant: ~basix.finite_element.DPCVariant = DPCVariant.unset, discontinuous: bool = False, shape: tuple[int, ...] | None = None, symmetry: bool | None = None, dtype: ~numpy.dtype[~typing.Any] | None | ~typing.Type[~typing.Any] | ~numpy._typing._dtype_like._SupportsDType[~numpy.dtype[~typing.Any]] | str | ~typing.Tuple[~typing.Any, int] | ~typing.Tuple[~typing.Any, ~typing.SupportsIndex | ~typing.Sequence[~typing.SupportsIndex]] | ~typing.List[~typing.Any] | ~numpy._typing._dtype_like._DTypeDict | ~typing.Tuple[~typing.Any, ~typing.Any] = <class 'numpy.float64'>) _ElementBase

Create a UFL compatible element using Basix.

Parameters:
  • family – Element family/type.

  • cell – Element cell type.

  • degree – Degree of the finite element.

  • lagrange_variant – Variant of Lagrange to be used.

  • dpc_variant – Variant of DPC to be used.

  • discontinuous – If True, the discontinuous version of the element is created.

  • shape – Value shape of the element. For scalar-valued families, this can be used to create vector and tensor elements.

  • symmetry – Set to True if the tensor is symmetric. Valid for rank 2 elements only.

  • dtype – Floating point data type.

Returns:

A finite element.

basix.ufl.enriched_element(elements: list[basix.ufl._ElementBase], map_type: MapType | None = None) _ElementBase

Create an UFL compatible enriched element from a list of elements.

Parameters:
  • elements – The list of elements

  • map_type – The map type for the enriched element.

Returns:

An enriched finite element.

basix.ufl.mixed_element(elements: list[basix.ufl._ElementBase]) _ElementBase

Create a UFL compatible mixed element from a list of elements.

Parameters:

elements – The list of elements

Returns:

A mixed finite element.

basix.ufl.quadrature_element(cell: str | CellType, value_shape: tuple[int, ...] = (), scheme: str | None = None, degree: int | None = None, points: ndarray[Any, dtype[float64]] | None = None, weights: ndarray[Any, dtype[float64]] | None = None, pullback: AbstractPullback = IdentityPullback()) _ElementBase

Create a quadrature element.

When creating this element, either the quadrature scheme and degree must be input or the quadrature points and weights must be.

Parameters:
  • cell – Cell to create the element on.

  • value_shape – Value shape of the element.

  • scheme – Quadrature scheme.

  • degree – Quadrature degree.

  • points – Quadrature points.

  • weights – Quadrature weights.

  • pullback – Map name.

Returns:

A ‘quadrature’ finite element.

basix.ufl.real_element(cell: CellType | str, value_shape: tuple[int, ...]) _ElementBase

Create a real element.

Parameters:
  • cell – Cell to create the element on.

  • value_shape – Value shape of the element.

Returns:

A ‘real’ finite element.

basix.ufl.wrap_element(element: FiniteElement) _ElementBase

Wrap a Basix element as a Basix UFL element.