Finite element (dolfinx::fem)

Under development

Finite elements

class dolfinx::fem::FiniteElement

Finite Element, containing the dof layout on a reference element, and various methods for evaluating and transforming the basis.

Public Functions

explicit FiniteElement(const ufcx_finite_element &e)

Create finite element from UFC finite element.

Parameters

e[in] UFC finite element

FiniteElement(const basix::FiniteElement &element, int bs)

Create finite element from a Basix finite element.

Parameters
  • element[in] Basix finite element

  • bs[in] The block size

FiniteElement(const FiniteElement &element) = delete

Copy constructor.

FiniteElement(FiniteElement &&element) = default

Move constructor.

virtual ~FiniteElement() = default

Destructor.

FiniteElement &operator=(const FiniteElement &element) = delete

Copy assignment.

FiniteElement &operator=(FiniteElement &&element) = default

Move assignment.

bool operator==(const FiniteElement &e) const

Check if two elements are equivalent.

Note

Equality can be checked only for non-mixed elements. For a mixed element, this function will raise an exception.

Returns

True is the two elements are the same

bool operator!=(const FiniteElement &e) const

Check if two elements are not equivalent.

Note

Equality can be checked only for non-mixed elements. For a mixed element, this function will raise an exception.

Returns

True is the two elements are not the same

std::string signature() const noexcept

String identifying the finite element.

Warning

The function is provided for convenience, but it should not be relied upon for determining the element type. Use other functions, commonly returning enums, to determine element properties.

Returns

Element signature

mesh::CellType cell_shape() const noexcept

Cell shape.

Returns

Element cell shape

int space_dimension() const noexcept

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

Returns

Dimension of the finite element space

int block_size() const noexcept

Block size of the finite element function space. For VectorElements and TensorElements, this is the number of DOFs colocated at each DOF point. For other elements, this is always 1.

Returns

Block size of the finite element space

int value_size() const

The value size, e.g. 1 for a scalar function, 2 for a 2D vector, 9 for a second-order tensor in 3D.

Note

The return value of this function is equal to std::accumulate(value_shape().begin(), value_shape().end(), 1, std::multiplies<int>()).

Returns

The value size

int reference_value_size() const

The value size, e.g. 1 for a scalar function, 2 for a 2D vector, 9 for a second-order tensor in 3D, for the reference element.

Returns

The value size for the reference element

xtl::span<const int> value_shape() const noexcept

Shape of the value space. The rank is the size of the value_shape.

std::string family() const noexcept

The finite element family.

Returns

The string of the finite element family

void tabulate(xt::xtensor<double, 4> &values, const xt::xtensor<double, 2> &X, int order) const

Evaluate all derivatives of the basis functions up to given order at given points in reference cell.

Parameters
  • values[inout] Four dimensional xtensor that will be filled with the tabulated values. Should be of shape {num_derivatives, num_points, num_dofs, reference_value_size}

  • X[in] Two dimensional xtensor of shape [num_points, geometric dimension] containing the points at the reference element

  • order[in] The number of derivatives (up to and including this order) to tabulate for

template<typename O, typename P, typename Q, typename R>
inline std::function<void(O&, const P&, const Q&, double, const R&)> map_fn() const

Return a function that performs the appropriate push-forward (pull-back) for the element type.

For a pull-back the passed arguments should be:

  • U [out] The data on the reference cell after the pull-back, flattened with row-major layout, shape=(num_points, ref value_size)

  • u [in] The data on the physical cell that should be pulled back , flattened with row-major layout, shape=(num_points, value_size)

  • K [in] The inverse oif the Jacobian matrix of the map, shape=(tdim, gdim)

  • detJ_inv [in] 1/det(J)

  • J [in] The Jacobian matrix, shape=(gdim, tdim)

Template Parameters
  • O – The type that hold the computed pushed-forward (pulled-back) data (ndim==1)

  • P – The type that hold the data to be pulled back (pushed forwarded) (ndim==1)

  • Q – The type that holds the Jacobian (inverse Jacobian) matrix (ndim==2)

  • R – The type that holds the inverse Jacobian (Jacobian) matrix (ndim==2)

Returns

A function that for a push-forward takes arguments

  • u [out] The data on the physical cell after the push-forward flattened with row-major layout, shape=(num_points, value_size)

  • U [in] The data on the reference cell physical field to push forward, flattened with row-major layout, shape=(num_points, ref_value_size)

  • J [in] The Jacobian matrix of the map ,shape=(gdim, tdim)

  • detJ [in] det(J)

  • K [in] The inverse of the Jacobian matrix, shape=(tdim, gdim)

int num_sub_elements() const noexcept

Get the number of sub elements (for a mixed or blocked element)

Returns

The number of sub elements

bool is_mixed() const noexcept

Check if element is a mixed element, i.e. composed of two or more elements of different types. A block element, e.g. a Lagrange element with block size > 1 is not considered mixed.

Returns

True is element is mixed.

const std::vector<std::shared_ptr<const FiniteElement>> &sub_elements() const noexcept

Subelements (if any)

std::shared_ptr<const FiniteElement> extract_sub_element(const std::vector<int> &component) const

Extract sub finite element for component.

basix::maps::type map_type() const

Get the map type used by the element.

bool interpolation_ident() const noexcept

Check if interpolation into the finite element space is an identity operation given the evaluation on an expression at specific points, i.e. the degree-of-freedom are equal to point evaluations. The function will return true for Lagrange elements.

Returns

True if interpolation is an identity operation

bool map_ident() const noexcept

Check if the push forward/pull back map from the values on reference to the values on a physical cell for this element is the identity map.

Returns

True if the map is the identity

const xt::xtensor<double, 2> &interpolation_points() const

Points on the reference cell at which an expression need to be evaluated in order to interpolate the expression in the finite element space. For Lagrange elements the points will just be the nodal positions. For other elements the points will typically be the quadrature points used to evaluate moment degrees of freedom.

Returns

Points on the reference cell. Shape is (num_points, tdim).

const xt::xtensor<double, 2> &interpolation_operator() const

Interpolation operator (matrix) Pi that maps a function evaluated at the points provided by FiniteElement::interpolation_points to the element degrees of freedom, i.e. dofs = Pi f_x. See the Basix documentation for basix::FiniteElement::interpolation_matrix for how the data in f_x should be ordered.

Returns

The interpolation operator Pi. Shape is (num_dofs, num_points*value_size)

xt::xtensor<double, 2> create_interpolation_operator(const FiniteElement &from) const

Create a matrix that maps degrees of freedom from one element to this element (interpolation).

Note

The two elements must use the same mapping between the reference and physical cells

Note

Does not support mixed elements

Parameters

from[in] The element to interpolate from

Returns

Matrix operator that maps the from degrees-of-freedom to the degrees-of-freedom of this element. Shape is (num_dofs of this element, num_dofs of from).

bool needs_dof_transformations() const noexcept

Check if DOF transformations are needed for this element.

DOF transformations will be needed for elements which might not be continuous when two neighbouring cells disagree on the orientation of a shared subentity, and when this cannot be corrected for by permuting the DOF numbering in the dofmap.

For example, Raviart-Thomas elements will need DOF transformations, as the neighbouring cells may disagree on the orientation of a basis function, and this orientation cannot be corrected for by permuting the DOF numbers on each cell.

Returns

True if DOF transformations are required

bool needs_dof_permutations() const noexcept

Check if DOF permutations are needed for this element.

DOF permutations will be needed for elements which might not be continuous when two neighbouring cells disagree on the orientation of a shared subentity, and when this can be corrected for by permuting the DOF numbering in the dofmap.

For example, higher order Lagrange elements will need DOF permutations, as the arrangement of DOFs on a shared subentity may be different from the point of view of neighbouring cells, and this can be corrected for by permuting the DOF numbers on each cell.

Returns

True if DOF transformations are required

template<typename T>
inline std::function<void(const xtl::span<T>&, const xtl::span<const std::uint32_t>&, std::int32_t, int)> get_dof_transformation_function(bool inverse = false, bool transpose = false, bool scalar_element = false) const

Return a function that applies DOF transformation to some data.

The returned function will take four inputs:

  • [in,out] data The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size)

  • [in] cell_info Permutation data for the cell. The size of this is num_cells. For elements where no transformations are required, an empty span can be passed in.

  • [in] cell The cell number

  • [in] block_size The block_size of the input data

Parameters
  • inverse[in] Indicates whether the inverse transformations should be returned

  • transpose[in] Indicates whether the transpose transformations should be returned

  • scalar_element[in] Indicates whether the scalar transformations should be returned for a vector element

template<typename T>
inline std::function<void(const xtl::span<T>&, const xtl::span<const std::uint32_t>&, std::int32_t, int)> get_dof_transformation_to_transpose_function(bool inverse = false, bool transpose = false, bool scalar_element = false) const

Return a function that applies DOF transformation to some transposed data.

The returned function will take three inputs:

  • [in,out] data The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size)

  • [in] cell_info Permutation data for the cell. The size of this is num_cells. For elements where no transformations are required, an empty span can be passed in.

  • [in] cell The cell number

  • [in] block_size The block_size of the input data

Parameters
  • inverse[in] Indicates whether the inverse transformations should be returned

  • transpose[in] Indicates whether the transpose transformations should be returned

  • scalar_element[in] Indicated whether the scalar transformations should be returned for a vector element

template<typename T>
inline void apply_dof_transformation(const xtl::span<T> &data, std::uint32_t cell_permutation, int block_size) const

Apply DOF transformation to some data.

Parameters
  • data[inout] The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size)

  • cell_permutation[in] Permutation data for the cell

  • block_size[in] The block_size of the input data

template<typename T>
inline void apply_inverse_transpose_dof_transformation(const xtl::span<T> &data, std::uint32_t cell_permutation, int block_size) const

Apply inverse transpose transformation to some data. For VectorElements, this applies the transformations for the scalar subelement.

Parameters
  • data[inout] The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size)

  • cell_permutation[in] Permutation data for the cell

  • block_size[in] The block_size of the input data

template<typename T>
inline void apply_transpose_dof_transformation(const xtl::span<T> &data, std::uint32_t cell_permutation, int block_size) const

Apply transpose transformation to some data. For VectorElements, this applies the transformations for the scalar subelement.

Parameters
  • data[inout] The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size)

  • cell_permutation[in] Permutation data for the cell

  • block_size[in] The block_size of the input data

template<typename T>
inline void apply_inverse_dof_transformation(const xtl::span<T> &data, std::uint32_t cell_permutation, int block_size) const

Apply inverse transformation to some data. For VectorElements, this applies the transformations for the scalar subelement.

Parameters
  • data[inout] The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size)

  • cell_permutation[in] Permutation data for the cell

  • block_size[in] The block_size of the input data

template<typename T>
inline void apply_dof_transformation_to_transpose(const xtl::span<T> &data, std::uint32_t cell_permutation, int block_size) const

Apply DOF transformation to some transposed data.

Parameters
  • data[inout] The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size)

  • cell_permutation[in] Permutation data for the cell

  • block_size[in] The block_size of the input data

template<typename T>
inline void apply_inverse_dof_transformation_to_transpose(const xtl::span<T> &data, std::uint32_t cell_permutation, int block_size) const

Apply inverse of DOF transformation to some transposed data.

Parameters
  • data[inout] The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size)

  • cell_permutation[in] Permutation data for the cell

  • block_size[in] The block_size of the input data

template<typename T>
inline void apply_transpose_dof_transformation_to_transpose(const xtl::span<T> &data, std::uint32_t cell_permutation, int block_size) const

Apply transpose of transformation to some transposed data.

Parameters
  • data[inout] The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size)

  • cell_permutation[in] Permutation data for the cell

  • block_size[in] The block_size of the input data

template<typename T>
inline void apply_inverse_transpose_dof_transformation_to_transpose(const xtl::span<T> &data, std::uint32_t cell_permutation, int block_size) const

Apply inverse transpose transformation to some transposed data.

Parameters
  • data[inout] The data to be transformed. This data is flattened with row-major layout, shape=(num_dofs, block_size)

  • cell_permutation[in] Permutation data for the cell

  • block_size[in] The block_size of the input data

void permute_dofs(const xtl::span<std::int32_t> &doflist, std::uint32_t cell_permutation) const

Permute the DOFs of the element.

Parameters
  • doflist[inout] The numbers of the DOFs, a span of length num_dofs

  • cell_permutation[in] Permutation data for the cell

void unpermute_dofs(const xtl::span<std::int32_t> &doflist, std::uint32_t cell_permutation) const

Unpermute the DOFs of the element.

Parameters
  • doflist[inout] The numbers of the DOFs, a span of length num_dofs

  • cell_permutation[in] Permutation data for the cell

std::function<void(const xtl::span<std::int32_t>&, std::uint32_t)> get_dof_permutation_function(bool inverse = false, bool scalar_element = false) const

Return a function that applies DOF permutation to some data.

The returned function will take three inputs:

  • [in,out] doflist The numbers of the DOFs, a span of length num_dofs

  • [in] cell_permutation Permutation data for the cell

  • [in] block_size The block_size of the input data

Parameters
  • inverse[in] Indicates whether the inverse transformations should be returned

  • scalar_element[in] Indicated whether the scalar transformations should be returned for a vector element

Function spaces and functions

Finite element functions, expressions and constants

Function spaces

class dolfinx::fem::FunctionSpace

This class represents a finite element function space defined by a mesh, a finite element, and a local-to-global map of the degrees of freedom (dofmap).

Public Functions

FunctionSpace(std::shared_ptr<const mesh::Mesh> mesh, std::shared_ptr<const FiniteElement> element, std::shared_ptr<const DofMap> dofmap)

Create function space for given mesh, element and dofmap.

Parameters
  • mesh[in] The mesh

  • element[in] The element

  • dofmap[in] The dofmap

FunctionSpace(FunctionSpace &&V) = default

Move constructor.

virtual ~FunctionSpace() = default

Destructor.

FunctionSpace &operator=(FunctionSpace &&V) = default

Move assignment operator.

std::shared_ptr<FunctionSpace> sub(const std::vector<int> &component) const

Extract subspace for component.

Parameters

component[in] The subspace component

Returns

The subspace

bool contains(const FunctionSpace &V) const

Check whether V is subspace of this, or this itself.

Parameters

V[in] The space to be tested for inclusion

Returns

True if V is contained in or is equal to this FunctionSpace

std::pair<FunctionSpace, std::vector<std::int32_t>> collapse() const

Collapse a subspace and return a new function space and a map from new to old dofs.

Returns

The new function space and a map from new to old dofs

std::vector<int> component() const

Get the component with respect to the root superspace.

Returns

The component with respect to the root superspace , i.e. W.sub(1).sub(0) == [1, 0]

xt::xtensor<double, 2> tabulate_dof_coordinates(bool transpose) const

Tabulate the physical coordinates of all dofs on this process.

Todo:

Remove - see function in interpolate.h

Parameters

transpose[in] If false the returned data has shape (num_points, gdim), otherwise it is transposed and has shape (gdim, num_points)

Returns

The dof coordinates [([x0, y0, z0], [x1, y1, z1], ...) if transpose is false, and otherwise the returned data is transposed.

std::shared_ptr<const mesh::Mesh> mesh() const

The mesh.

std::shared_ptr<const FiniteElement> element() const

The finite element.

std::shared_ptr<const DofMap> dofmap() const

The dofmap.

Functions

template<typename T>
class dolfinx::fem::Function

This class represents a function \( u_h \) in a finite element function space \( V_h \), given by.

\[ u_h = \sum_{i=1}^{n} U_i \phi_i \]
where \( \{\phi_i\}_{i=1}^{n} \) is a basis for \( V_h \), and \( U \) is a vector of expansion coefficients for \( u_h \).

Public Types

using value_type = T

Field type for the Function, e.g. double.

Public Functions

inline explicit Function(std::shared_ptr<const FunctionSpace> V)

Create function on given function space.

Parameters

V[in] The function space

inline Function(std::shared_ptr<const FunctionSpace> V, std::shared_ptr<la::Vector<T>> x)

Create function on given function space with a given vector.

Warning

This constructor is intended for internal library use only

Parameters
  • V[in] The function space

  • x[in] The vector

Function(Function &&v) = default

Move constructor.

~Function() = default

Destructor.

Function &operator=(Function &&v) = default

Move assignment.

inline Function sub(int i) const

Extract subfunction (view into the Function)

Parameters

i[in] Index of subfunction

Returns

The subfunction

inline Function collapse() const

Collapse a subfunction (view into a Function) to a stand-alone Function.

Returns

New collapsed Function

inline std::shared_ptr<const FunctionSpace> function_space() const

Access the function space.

Returns

The function space

inline std::shared_ptr<const la::Vector<T>> x() const

Underlying vector.

inline std::shared_ptr<la::Vector<T>> x()

Underlying vector.

inline void interpolate(const Function<T> &v, const xtl::span<const std::int32_t> &cells)

Interpolate a Function.

Parameters
  • v[in] The function to be interpolated

  • cells[in] The cells to interpolate on

inline void interpolate(const Function<T> &v)

Interpolate a Function.

Parameters

v[in] The function to be interpolated

inline void interpolate(const std::function<xt::xarray<T>(const xt::xtensor<double, 2>&)> &f, const xtl::span<const std::int32_t> &cells)

Interpolate an expression function on a list of cells.

Parameters
  • f[in] The expression function to be interpolated

  • cells[in] The cells to interpolate on

inline void interpolate(const std::function<xt::xarray<T>(const xt::xtensor<double, 2>&)> &f)

Interpolate an expression function on the whole domain.

Parameters

f[in] The expression to be interpolated

inline void interpolate(const Expression<T> &e, const xtl::span<const std::int32_t> &cells)

Interpolate an Expression (based on UFL)

Parameters
  • e[in] The Expression to be interpolated. The Expression must have been created using the reference coordinates FiniteElement::interpolation_points() for the element associated with u.

  • cells[in] The cells to interpolate on

inline void interpolate(const Expression<T> &e)

Interpolate an Expression (based on UFL) on all cells.

Parameters

e[in] The function to be interpolated

inline void eval(const xt::xtensor<double, 2> &x, const xtl::span<const std::int32_t> &cells, xt::xtensor<T, 2> &u) const

Evaluate the Function at points.

Parameters
  • x[in] The coordinates of the points. It has shape (num_points, 3).

  • cells[in] An array of cell indices. cells[i] is the index of the cell that contains the point x(i). Negative cell indices can be passed, and the corresponding point will be ignored.

  • u[inout] The values at the points. Values are not computed for points with a negative cell index. This argument must be passed with the correct size.

Public Members

std::string name = "u"

Name.

Constants

template<typename T>
class dolfinx::fem::Constant

Constant value which can be attached to a Form. Constants may be scalar (rank 0), vector (rank 1), or tensor valued.

Public Functions

template<typename = std::enable_if_t<std::is_arithmetic_v<T> || std::is_same_v<T, std::complex<float>> || std::is_same_v<T, std::complex<double>>>>
inline Constant(T c)

Create a rank-0 (scalar-valued) constant.

Parameters

c[in] Value of the constant

inline Constant(const xt::xarray<T> &c)

Create a rank-d constant.

Parameters

c[in] Value of the constant

Public Members

std::vector<int> shape

Shape.

std::vector<T> value

Values, stored as a row-major flattened array.

Forms

template<typename T>
class dolfinx::fem::Form

A representation of finite element variational forms.

A note on the order of trial and test spaces: FEniCS numbers argument spaces starting with the leading dimension of the corresponding tensor (matrix). In other words, the test space is numbered 0 and the trial space is numbered 1. However, in order to have a notation that agrees with most existing finite element literature, in particular

\[ a = a(u, v) \]

the spaces are numbered from right to left

\[ a: V_1 \times V_0 \rightarrow \mathbb{R} \]

This is reflected in the ordering of the spaces that should be supplied to generated subclasses. In particular, when a bilinear form is initialized, it should be initialized as a(V_1, V_0) = ..., where V_1 is the trial space and V_0 is the test space. However, when a form is initialized by a list of argument spaces (the variable function_spaces in the constructors below), the list of spaces should start with space number 0 (the test space) and then space number 1 (the trial space).

Public Types

using scalar_type = T

Scalar type (T)

Public Functions

inline Form(const std::vector<std::shared_ptr<const fem::FunctionSpace>> &function_spaces, const std::map<IntegralType, std::pair<std::vector<std::pair<int, std::function<void(T*, const T*, const T*, const double*, const int*, const std::uint8_t*)>>>, const mesh::MeshTags<int>*>> &integrals, const std::vector<std::shared_ptr<const fem::Function<T>>> &coefficients, const std::vector<std::shared_ptr<const fem::Constant<T>>> &constants, bool needs_facet_permutations, const std::shared_ptr<const mesh::Mesh> &mesh = nullptr)

Create a finite element form.

Note

User applications will normally call a fem::Form builder function rather using this interfcae directly.

Parameters
  • function_spaces[in] Function spaces for the form arguments

  • integrals[in] The integrals in the form. The first key is the domain type. For each key there is a pair (list[domain id, integration kernel], domain markers).

  • coefficients[in]

  • constants[in] Constants in the Form

  • needs_facet_permutations[in] Set to true is any of the integration kernels require cell permutation data

  • mesh[in] The mesh of the domain. This is required when there are not argument functions from which the mesh can be extracted, e.g. for functionals

Form(const Form &form) = delete

Copy constructor.

Form(Form &&form) = default

Move constructor.

virtual ~Form() = default

Destructor.

inline int rank() const

Rank of the form (bilinear form = 2, linear form = 1, functional = 0, etc)

Returns

The rank of the form

inline std::shared_ptr<const mesh::Mesh> mesh() const

Extract common mesh for the form.

Returns

The mesh

inline const std::vector<std::shared_ptr<const fem::FunctionSpace>> &function_spaces() const

Return function spaces for all arguments.

Returns

Function spaces

inline const std::function<void(T*, const T*, const T*, const double*, const int*, const std::uint8_t*)> &kernel(IntegralType type, int i) const

Get the function for ‘kernel’ for integral i of given type.

Parameters
  • type[in] Integral type

  • i[in] Domain index

Returns

Function to call for tabulate_tensor

inline std::set<IntegralType> integral_types() const

Get types of integrals in the form.

Returns

Integrals types

inline int num_integrals(IntegralType type) const

Number of integrals of given type.

Parameters

type[in] Integral type

Returns

Number of integrals

inline std::vector<int> integral_ids(IntegralType type) const

Get the IDs for integrals (kernels) for given integral type. The IDs correspond to the domain IDs which the integrals are defined for in the form. ID=-1 is the default integral over the whole domain.

Parameters

type[in] Integral type

Returns

List of IDs for given integral type

inline const std::vector<std::int32_t> &cell_domains(int i) const

Get the list of cell indices for the ith integral (kernel) for the cell domain type.

Parameters

i[in] Integral ID, i.e. (sub)domain index

Returns

List of active cell entities for the given integral (kernel)

inline const std::vector<std::pair<std::int32_t, int>> &exterior_facet_domains(int i) const

Get the list of (cell_index, local_facet_index) pairs for the ith integral (kernel) for the exterior facet domain type.

Parameters

i[in] Integral ID, i.e. (sub)domain index

Returns

List of (cell_index, local_facet_index) pairs

inline const std::vector<std::tuple<std::int32_t, int, std::int32_t, int>> &interior_facet_domains(int i) const

Get the list of (cell_index_0, local_facet_index_0, cell_index_1, local_facet_index_1) tuples for the ith integral (kernel) for the interior facet domain type,.

Parameters

i[in] Integral ID, i.e. (sub)domain index

Returns

List of tuples of the form (cell_index_0, local_facet_index_0, cell_index_1, local_facet_index_1)

inline const std::vector<std::shared_ptr<const fem::Function<T>>> &coefficients() const

Access coefficients.

inline bool needs_facet_permutations() const

Get bool indicating whether permutation data needs to be passed into these integrals.

Returns

True if cell permutation data is required

inline std::vector<int> coefficient_offsets() const

Offset for each coefficient expansion array on a cell. Used to pack data for multiple coefficients in a flat array. The last entry is the size required to store all coefficients.

inline const std::vector<std::shared_ptr<const fem::Constant<T>>> &constants() const

Access constants.

Dirichlet boundary conditions

template<typename T>
class dolfinx::fem::DirichletBC

Object for setting (strong) Dirichlet boundary conditions.

\(u = g \ \text{on} \ G\),

where \(u\) is the solution to be computed, \(g\) is a function and \(G\) is a sub domain of the mesh.

A DirichletBC is specified by the function \(g\), the function space (trial space) and degrees of freedom to which the boundary condition applies.

Public Functions

template<typename S, typename U, typename = std::enable_if_t<std::is_convertible_v<S, T> or std::is_convertible_v<S, xt::xarray<T>>>>
inline DirichletBC(const S &g, U &&dofs, const std::shared_ptr<const FunctionSpace> &V)

Create a representation of a Dirichlet boundary condition constrained by a scalar or tensor constant.

Note

Can be used only with point-evaluation elements.

Note

The indices in dofs are for blocks, e.g. a block index corresponds to 3 degrees-of-freedom if the dofmap associated with g has block size 3.

Note

The size of of g must be equal to the block size if V. Use the Function version if this is not the case, e.g. for some mixed spaces.

Parameters
  • g[in] The boundary condition value (T or xt::xarray<T>)

  • dofs[in] Degree-of-freedom block indices ( std::vector<std::int32_t>) to be constrained. The indices must be sorted.

  • V[in] The function space to be constrained

template<typename U>
inline DirichletBC(const std::shared_ptr<const Constant<T>> &g, U &&dofs, const std::shared_ptr<const FunctionSpace> &V)

Create a representation of a Dirichlet boundary condition constrained by a fem::Constant.

Note

Can be used only with point-evaluation elements.

Note

The indices in dofs are for blocks, e.g. a block index corresponds to 3 degrees-of-freedom if the dofmap associated with g has block size 3.

Note

The size of of g must be equal to the block size if V. Use the Function version if this is not the case, e.g. for some mixed spaces.

Parameters
  • g[in] The boundary condition value

  • dofs[in] Degree-of-freedom block indices (std::vector<std::int32_t>) to be constrained. The indices must be sorted.

  • V[in] The function space to be constrained

template<typename U>
inline DirichletBC(const std::shared_ptr<const Function<T>> &g, U &&dofs)

Create a representation of a Dirichlet boundary condition where the space being constrained is the same as the function that defines the constraint Function, i.e. share the same fem::FunctionSpace

Note

The indices in dofs are for blocks, e.g. a block index corresponds to 3 degrees-of-freedom if the dofmap associated with g has block size 3.

Parameters
  • g[in] The boundary condition value.

  • dofs[in] Degree-of-freedom block indices (std::vector<std::int32_t>) to be constrained. The indices must be sorted.

template<typename U>
inline DirichletBC(const std::shared_ptr<const Function<T>> &g, U &&V_g_dofs, const std::shared_ptr<const FunctionSpace> &V)

Create a representation of a Dirichlet boundary condition where the space being constrained and the function that defines the constraint values do not share the same FunctionSpace.

A typical examples is when applying a constraint on a subspace. The (sub)space and the constrain function must have the same finite element.

Note

The indices in dofs are unrolled and not for blocks.

Parameters
  • g[in] The boundary condition value

  • V_g_dofs[in] Two arrays of degree-of-freedom indices (std::array<std::vector<std::int32_t>, 2>). First array are indices in the space where boundary condition is applied (V), second array are indices in the space of the boundary condition value function g. The arrays must be sorted by the indices in the first array. The dof indices are unrolled, i.e. are not by dof block.

  • V[in] The function (sub)space on which the boundary condition is applied

DirichletBC(const DirichletBC &bc) = default

Copy constructor.

Parameters

bc[in] The object to be copied

DirichletBC(DirichletBC &&bc) = default

Move constructor.

Parameters

bc[in] The object to be moved

~DirichletBC() = default

Destructor.

DirichletBC &operator=(const DirichletBC &bc) = default

Assignment operator.

Parameters

bc[in] Another DirichletBC object

DirichletBC &operator=(DirichletBC &&bc) = default

Move assignment operator.

inline std::shared_ptr<const fem::FunctionSpace> function_space() const

The function space to which boundary conditions are applied.

Returns

The function space

inline std::variant<std::shared_ptr<const Function<T>>, std::shared_ptr<const Constant<T>>> value() const

Return boundary value function g.

Returns

The boundary values Function

inline std::pair<xtl::span<const std::int32_t>, std::int32_t> dof_indices() const

Access dof indices (local indices, unrolled), including ghosts, to which a Dirichlet condition is applied, and the index to the first non-owned (ghost) index. The array of indices is sorted.

Returns

Sorted array of dof indices (unrolled) and index to the first entry in the dof index array that is not owned. Entries dofs[:pos] are owned and entries dofs[pos:] are ghosts.

inline void set(xtl::span<T> x, double scale = 1.0) const

Set bc entries in x to scale * x_bc

Parameters
  • x[in] The array in which to set scale * x_bc[i], where x_bc[i] is the boundary value of x[i]. Entries in x that do not have a Dirichlet condition applied to them are unchanged. The length of x must be less than or equal to the index of the greatest boundary dof index. To set values only for degrees-of-freedom that are owned by the calling rank, the length of the array x should be equal to the number of dofs owned by this rank.

  • scale[in] The scaling value to apply

inline void set(xtl::span<T> x, const xtl::span<const T> &x0, double scale = 1.0) const

Set bc entries in x to scale * (x0 - x_bc)

Parameters
  • x[in] The array in which to set scale * (x0 - x_bc)

  • x0[in] The array used in compute the value to set

  • scale[in] The scaling value to apply

inline void dof_values(xtl::span<T> values) const

Todo:

Review this function - it is almost identical to the ‘DirichletBC::set’ function

Set boundary condition value for entries with an applied boundary condition. Other entries are not modified.

Parameters

values[out] The array in which to set the dof values. The array must be at least as long as the array associated with V1 (the space of the function that provides the dof values)

inline void mark_dofs(const xtl::span<std::int8_t> &markers) const

Set markers[i] = true if dof i has a boundary condition applied. Value of markers[i] is not changed otherwise.

Parameters

markers[inout] Entry makers[i] is set to true if dof i in V0 had a boundary condition applied, i.e. dofs which are fixed by a boundary condition. Other entries in markers are left unchanged.

Degree-of-freedom maps

graph::AdjacencyList<std::int32_t> dolfinx::fem::transpose_dofmap(const graph::AdjacencyList<std::int32_t> &dofmap, std::int32_t num_cells)

Create an adjacency list that maps a global index (process-wise) to the ‘unassembled’ cell-wise contributions.

It is built from the usual (cell, local index) -> global index dof map. An ‘unassembled’ vector is the stacked cell contributions, ordered by cell index. If the usual dof map is:

Cell: 0 1 2 3 Global index: [ [0, 3, 5], [3, 2, 4], [4, 3, 2], [2, 1, 0]]

the ‘transpose’ dof map will be:

Global index: 0 1 2 3 4 5 Unassembled index: [ [0, 11], [10], [4, 8, 9], [1, 3, 7], [5, 6], [2] ]

Parameters
  • dofmap[in] The standard dof map that for each cell (node) gives the global (process-wise) index of each local (cell-wise) index.

  • num_cells[in] The number of cells (nodes) in dofmap to consider. The first num_cells are used. This is argument is typically used to exclude ghost cell contributions.

Returns

Map from global (process-wise) index to positions in an unaassembled array. The links for each node are sorted.

class dolfinx::fem::DofMap

Degree-of-freedom map.

This class handles the mapping of degrees of freedom. It builds a dof map based on an ElementDofLayout on a specific mesh topology. It will reorder the dofs when running in parallel. Sub-dofmaps, both views and copies, are supported.

Public Functions

template<typename E, typename U, typename = std::enable_if_t<std::is_same<graph::AdjacencyList<std::int32_t>, std::decay_t<U>>::value>>
inline DofMap(E &&element, std::shared_ptr<const common::IndexMap> index_map, int index_map_bs, U &&dofmap, int bs)

Create a DofMap from the layout of dofs on a reference element, an IndexMap defining the distribution of dofs across processes and a vector of indices.

Parameters
  • element[in] The layout of the degrees of freedom on an element (fem::ElementDofLayout)

  • index_map[in] The map describing the parallel distribution of the degrees of freedom

  • index_map_bs[in] The block size associated with the index_map

  • dofmap[in] Adjacency list (graph::AdjacencyList<std::int32_t>) with the degrees-of-freedom for each cell

  • bs[in] The block size of the dofmap

DofMap(DofMap &&dofmap) = default

Move constructor.

DofMap &operator=(DofMap &&dofmap) = default

Move assignment.

bool operator==(const DofMap &map) const

Equality operator.

Returns

Returns true if the data for the two dofmaps is equal

inline xtl::span<const std::int32_t> cell_dofs(int cell) const

Local-to-global mapping of dofs on a cell.

Parameters

cell[in] The cell index

Returns

Local-global dof map for the cell (using process-local indices)

int bs() const noexcept

Return the block size for the dofmap.

DofMap extract_sub_dofmap(const std::vector<int> &component) const

Extract subdofmap component.

Parameters

component[in] The component indices

Returns

The dofmap for the component

std::pair<DofMap, std::vector<std::int32_t>> collapse(MPI_Comm comm, const mesh::Topology &topology, const std::function<std::vector<int>(const graph::AdjacencyList<std::int32_t>&)> &reorder_fn = [](const graph::AdjacencyList< std::int32_t > &g) { return graph::reorder_gps(g);}) const

Create a “collapsed” dofmap (collapses a sub-dofmap)

Parameters
  • comm[in] MPI Communicator

  • topology[in] The mesh topology that the dofmap is defined on

  • reorder_fn[in] The graph re-ordering function to apply to the dof data

Returns

The collapsed dofmap

const graph::AdjacencyList<std::int32_t> &list() const

Get dofmap data.

Returns

The adjacency list with dof indices for each cell

inline const ElementDofLayout &element_dof_layout() const

Layout of dofs on an element.

int index_map_bs() const

Block size associated with the index_map.

Public Members

std::shared_ptr<const common::IndexMap> index_map

Index map that describes the parallel distribution of the dofmap.

Assembly

Functions

template<typename T>
std::map<std::pair<dolfinx::fem::IntegralType, int>, std::pair<xtl::span<const T>, int>> make_coefficients_span(const std::map<std::pair<IntegralType, int>, std::pair<std::vector<T>, int>> &coeffs)

Create a map xtl::span from a map of std::vector.

template<typename T>
T assemble_scalar(const Form<T> &M, const xtl::span<const T> &constants, const std::map<std::pair<IntegralType, int>, std::pair<xtl::span<const T>, int>> &coefficients)

Assemble functional into scalar. The caller supplies the form constants and coefficients for this version, which has efficiency benefits if the data can be re-used for multiple calls.

Note

Caller is responsible for accumulation across processes.

Parameters
  • M[in] The form (functional) to assemble

  • constants[in] The constants that appear in M

  • coefficients[in] The coefficients that appear in M

Returns

The contribution to the form (functional) from the local process

template<typename T>
T assemble_scalar(const Form<T> &M)

Assemble functional into scalar.

Note

Caller is responsible for accumulation across processes.

Parameters

M[in] The form (functional) to assemble

Returns

The contribution to the form (functional) from the local process

template<typename T>
void assemble_vector(xtl::span<T> b, const Form<T> &L, const xtl::span<const T> &constants, const std::map<std::pair<IntegralType, int>, std::pair<xtl::span<const T>, int>> &coefficients)

Assemble linear form into a vector, The caller supplies the form constants and coefficients for this version, which has efficiency benefits if the data can be re-used for multiple calls.

Parameters
  • b[inout] The vector to be assembled. It will not be zeroed before assembly.

  • L[in] The linear forms to assemble into b

  • constants[in] The constants that appear in L

  • coefficients[in] The coefficients that appear in L

template<typename T>
void assemble_vector(xtl::span<T> b, const Form<T> &L)

Assemble linear form into a vector.

Parameters
  • b[inout] The vector to be assembled. It will not be zeroed before assembly.

  • L[in] The linear forms to assemble into b

template<typename T>
void apply_lifting(xtl::span<T> b, const std::vector<std::shared_ptr<const Form<T>>> &a, const std::vector<xtl::span<const T>> &constants, const std::vector<std::map<std::pair<IntegralType, int>, std::pair<xtl::span<const T>, int>>> &coeffs, const std::vector<std::vector<std::shared_ptr<const DirichletBC<T>>>> &bcs1, const std::vector<xtl::span<const T>> &x0, double scale)

Modify b such that:

b <- b - scale * A_j (g_j - x0_j)

where j is a block (nest) index. For a non-blocked problem j = 0. The boundary conditions bcs1 are on the trial spaces V_j. The forms in [a] must have the same test space as L (from which b was built), but the trial space may differ. If x0 is not supplied, then it is treated as zero.

Ghost contributions are not accumulated (not sent to owner). Caller is responsible for calling VecGhostUpdateBegin/End.

template<typename T>
void apply_lifting(xtl::span<T> b, const std::vector<std::shared_ptr<const Form<T>>> &a, const std::vector<std::vector<std::shared_ptr<const DirichletBC<T>>>> &bcs1, const std::vector<xtl::span<const T>> &x0, double scale)

Modify b such that:

b <- b - scale * A_j (g_j - x0_j)

where j is a block (nest) index. For a non-blocked problem j = 0. The boundary conditions bcs1 are on the trial spaces V_j. The forms in [a] must have the same test space as L (from which b was built), but the trial space may differ. If x0 is not supplied, then it is treated as zero.

Ghost contributions are not accumulated (not sent to owner). Caller is responsible for calling VecGhostUpdateBegin/End.

template<typename T, typename U>
void assemble_matrix(U mat_add, const Form<T> &a, const xtl::span<const T> &constants, const std::map<std::pair<IntegralType, int>, std::pair<xtl::span<const T>, int>> &coefficients, const std::vector<std::shared_ptr<const DirichletBC<T>>> &bcs)

Assemble bilinear form into a matrix.

Parameters
  • mat_add[in] The function for adding values into the matrix

  • a[in] The bilinear from to assemble

  • constants[in] Constants that appear in a

  • coefficients[in] Coefficients that appear in a

  • bcs[in] Boundary conditions to apply. For boundary condition dofs the row and column are zeroed. The diagonal entry is not set.

template<typename T, typename U>
void assemble_matrix(U mat_add, const Form<T> &a, const std::vector<std::shared_ptr<const DirichletBC<T>>> &bcs)

Assemble bilinear form into a matrix.

Parameters
  • mat_add[in] The function for adding values into the matrix

  • a[in] The bilinear from to assemble

  • bcs[in] Boundary conditions to apply. For boundary condition dofs the row and column are zeroed. The diagonal entry is not set.

template<typename T, typename U>
void assemble_matrix(U mat_add, const Form<T> &a, const xtl::span<const T> &constants, const std::map<std::pair<IntegralType, int>, std::pair<xtl::span<const T>, int>> &coefficients, const xtl::span<const std::int8_t> &dof_marker0, const xtl::span<const std::int8_t> &dof_marker1)

Assemble bilinear form into a matrix. Matrix must already be initialised. Does not zero or finalise the matrix.

Parameters
  • mat_add[in] The function for adding values into the matrix

  • a[in] The bilinear form to assemble

  • constants[in] Constants that appear in a

  • coefficients[in] Coefficients that appear in a

  • dof_marker0[in] Boundary condition markers for the rows. If bc[i] is true then rows i in A will be zeroed. The index i is a local index.

  • dof_marker1[in] Boundary condition markers for the columns. If bc[i] is true then rows i in A will be zeroed. The index i is a local index.

template<typename T, typename U>
void assemble_matrix(U mat_add, const Form<T> &a, const xtl::span<const std::int8_t> &dof_marker0, const xtl::span<const std::int8_t> &dof_marker1)

Assemble bilinear form into a matrix. Matrix must already be initialised. Does not zero or finalise the matrix.

Parameters
  • mat_add[in] The function for adding values into the matrix

  • a[in] The bilinear form to assemble

  • dof_marker0[in] Boundary condition markers for the rows. If bc[i] is true then rows i in A will be zeroed. The index i is a local index.

  • dof_marker1[in] Boundary condition markers for the columns. If bc[i] is true then rows i in A will be zeroed. The index i is a local index.

template<typename T, typename U>
void set_diagonal(U set_fn, const xtl::span<const std::int32_t> &rows, T diagonal = 1.0)

Sets a value to the diagonal of a matrix for specified rows. It is typically called after assembly. The assembly function zeroes Dirichlet rows and columns. For block matrices, this function should normally be called only on the diagonal blocks, i.e. blocks for which the test and trial spaces are the same.

Parameters
  • set_fn[in] The function for setting values to a matrix

  • rows[in] The row blocks, in local indices, for which to add a value to the diagonal

  • diagonal[in] The value to add to the diagonal for the specified rows

template<typename T, typename U>
void set_diagonal(U set_fn, const fem::FunctionSpace &V, const std::vector<std::shared_ptr<const DirichletBC<T>>> &bcs, T diagonal = 1.0)

Sets a value to the diagonal of the matrix for rows with a Dirichlet boundary conditions applied. This function is typically called after assembly. The assembly function zeroes Dirichlet rows and columns. This function adds the value only to rows that are locally owned, and therefore does not create a need for parallel communication. For block matrices, this function should normally be called only on the diagonal blocks, i.e. blocks for which the test and trial spaces are the same.

Parameters
  • set_fn[in] The function for setting values to a matrix

  • V[in] The function space for the rows and columns of the matrix. It is used to extract only the Dirichlet boundary conditions that are define on V or subspaces of V.

  • bcs[in] The Dirichlet boundary condtions

  • diagonal[in] The value to add to the diagonal for rows with a boundary condition applied

template<typename T>
void set_bc(xtl::span<T> b, const std::vector<std::shared_ptr<const DirichletBC<T>>> &bcs, const xtl::span<const T> &x0, double scale = 1.0)

Set bc values in owned (local) part of the vector, multiplied by ‘scale’. The vectors b and x0 must have the same local size. The bcs should be on (sub-)spaces of the form L that b represents.

template<typename T>
void set_bc(xtl::span<T> b, const std::vector<std::shared_ptr<const DirichletBC<T>>> &bcs, double scale = 1.0)

Set bc values in owned (local) part of the vector, multiplied by ‘scale’. The bcs should be on (sub-)spaces of the form L that b represents.

Interpolation

xt::xtensor<double, 2> dolfinx::fem::interpolation_coords(const fem::FiniteElement &element, const mesh::Mesh &mesh, const xtl::span<const std::int32_t> &cells)

Compute the evaluation points in the physical space at which an expression should be computed to interpolate it in a finite elemenet space.

Parameters
  • element[in] The element to be interpolated into

  • mesh[in] The domain

  • cells[in] Indices of the cells in the mesh to compute interpolation coordinates for

Returns

The coordinates in the physical space at which to evaluate an expression. The shape is (3, num_points).

template<typename T>
void dolfinx::fem::interpolate(Function<T> &u, const Function<T> &v, const xtl::span<const std::int32_t> &cells)

Interpolate from one finite element Function to another on the same mesh.

Parameters
  • u[out] The function to interpolate into

  • v[in] The function to be interpolated

  • cells[in] List of cell indices to interpolate on

template<typename T>
void dolfinx::fem::interpolate(Function<T> &u, const xt::xarray<T> &f, const xtl::span<const std::int32_t> &cells)

Interpolate an expression f(x) in a finite element space.

Parameters
  • u[out] The function to interpolate into

  • f[in] Evaluation of the function f(x) at the physical points x given by fem::interpolation_coords. The element used in fem::interpolation_coords should be the same element as associated with u. The shape of f should be (value_size, num_points), or if value_size=1 the shape can be (num_points,).

  • cells[in] Indices of the cells in the mesh on which to interpolate. Should be the same as the list used when calling fem::interpolation_coords.

Sparsity pattern construction

template<typename T>
la::SparsityPattern dolfinx::fem::create_sparsity_pattern(const Form<T> &a)

Create a sparsity pattern for a given form.

Note

The pattern is not finalised, i.e. the caller is responsible for calling SparsityPattern::assemble.

Parameters

a[in] A bilinear form

Returns

The corresponding sparsity pattern

la::SparsityPattern dolfinx::fem::create_sparsity_pattern(const mesh::Topology &topology, const std::array<std::reference_wrapper<const DofMap>, 2> &dofmaps, const std::set<IntegralType> &integrals)

Create a sparsity pattern for a given form.

Note

The pattern is not finalised, i.e. the caller is responsible for calling SparsityPattern::assemble.

PETSc helpers

Mat dolfinx::fem::petsc::create_matrix(const Form<PetscScalar> &a, const std::string &type = std::string())

Create a matrix.

Parameters
  • a[in] A bilinear form

  • type[in] The PETSc matrix type to create

Returns

A sparse matrix with a layout and sparsity that matches the bilinear form. The caller is responsible for destroying the Mat object.

Mat dolfinx::fem::petsc::create_matrix_block(const std::vector<std::vector<const fem::Form<PetscScalar>*>> &a, const std::string &type = std::string())

Initialise a monolithic matrix for an array of bilinear forms.

Parameters
  • a[in] Rectangular array of bilinear forms. The a(i, j) form will correspond to the (i, j) block in the returned matrix

  • type[in] The type of PETSc Mat. If empty the PETSc default is used.

Returns

A sparse matrix with a layout and sparsity that matches the bilinear forms. The caller is responsible for destroying the Mat object.

Mat dolfinx::fem::petsc::create_matrix_nest(const std::vector<std::vector<const fem::Form<PetscScalar>*>> &a, const std::vector<std::vector<std::string>> &types)

Create nested (MatNest) matrix.

The caller is responsible for destroying the Mat object

Vec dolfinx::fem::petsc::create_vector_block(const std::vector<std::pair<std::reference_wrapper<const common::IndexMap>, int>> &maps)

Initialise monolithic vector. Vector is not zeroed.

The caller is responsible for destroying the Mat object

Vec dolfinx::fem::petsc::create_vector_nest(const std::vector<std::pair<std::reference_wrapper<const common::IndexMap>, int>> &maps)

Create nested (VecNest) vector. Vector is not zeroed.

void dolfinx::fem::petsc::assemble_vector(Vec b, const Form<PetscScalar> &L, const xtl::span<const PetscScalar> &constants, const std::map<std::pair<IntegralType, int>, std::pair<xtl::span<const PetscScalar>, int>> &coeffs)

Assemble linear form into an already allocated PETSc vector. Ghost contributions are not accumulated (not sent to owner). Caller is responsible for calling VecGhostUpdateBegin/End.

Parameters
  • b[inout] The PETsc vector to assemble the form into. The vector must already be initialised with the correct size. The process-local contribution of the form is assembled into this vector. It is not zeroed before assembly.

  • L[in] The linear form to assemble

  • constants[in] The constants that appear in L

  • coeffs[in] The coefficients that appear in L

void dolfinx::fem::petsc::assemble_vector(Vec b, const Form<PetscScalar> &L)

Assemble linear form into an already allocated PETSc vector. Ghost contributions are not accumulated (not sent to owner). Caller is responsible for calling VecGhostUpdateBegin/End.

Parameters
  • b[inout] The PETsc vector to assemble the form into. The vector must already be initialised with the correct size. The process-local contribution of the form is assembled into this vector. It is not zeroed before assembly.

  • L[in] The linear form to assemble

void dolfinx::fem::petsc::apply_lifting(Vec b, const std::vector<std::shared_ptr<const Form<PetscScalar>>> &a, const std::vector<xtl::span<const PetscScalar>> &constants, const std::vector<std::map<std::pair<IntegralType, int>, std::pair<xtl::span<const PetscScalar>, int>>> &coeffs, const std::vector<std::vector<std::shared_ptr<const DirichletBC<PetscScalar>>>> &bcs1, const std::vector<Vec> &x0, double scale)

Modify b such that:

b <- b - scale * A_j (g_j - x0_j)

where j is a block (nest) index. For a non-blocked problem j = 0. The boundary conditions bcs1 are on the trial spaces V_j. The forms in [a] must have the same test space as L (from which b was built), but the trial space may differ. If x0 is not supplied, then it is treated as zero.

Ghost contributions are not accumulated (not sent to owner). Caller is responsible for calling VecGhostUpdateBegin/End.

void dolfinx::fem::petsc::apply_lifting(Vec b, const std::vector<std::shared_ptr<const Form<PetscScalar>>> &a, const std::vector<std::vector<std::shared_ptr<const DirichletBC<PetscScalar>>>> &bcs1, const std::vector<Vec> &x0, double scale)

Modify b such that:

b <- b - scale * A_j (g_j - x0_j)

where j is a block (nest) index. For a non-blocked problem j = 0. The boundary conditions bcs1 are on the trial spaces V_j. The forms in [a] must have the same test space as L (from which b was built), but the trial space may differ. If x0 is not supplied, then it is treated as zero.

Ghost contributions are not accumulated (not sent to owner). Caller is responsible for calling VecGhostUpdateBegin/End.

void dolfinx::fem::petsc::set_bc(Vec b, const std::vector<std::shared_ptr<const DirichletBC<PetscScalar>>> &bcs, const Vec x0, double scale = 1.0)

Set bc values in owned (local) part of the PETSc vector, multiplied by ‘scale’. The vectors b and x0 must have the same local size. The bcs should be on (sub-)spaces of the form L that b represents.

Misc

Functions

template<typename T>
std::vector<std::vector<std::array<std::shared_ptr<const FunctionSpace>, 2>>> extract_function_spaces(const std::vector<std::vector<const Form<T>*>> &a)

Extract test (0) and trial (1) function spaces pairs for each bilinear form for a rectangular array of forms.

Parameters

a[in] A rectangular block on bilinear forms

Returns

Rectangular array of the same shape as a with a pair of function spaces in each array entry. If a form is null, then the returned function space pair is (null, null).

la::SparsityPattern create_sparsity_pattern(const mesh::Topology &topology, const std::array<std::reference_wrapper<const DofMap>, 2> &dofmaps, const std::set<IntegralType> &integrals)

Create a sparsity pattern for a given form.

Note

The pattern is not finalised, i.e. the caller is responsible for calling SparsityPattern::assemble.

template<typename T>
la::SparsityPattern create_sparsity_pattern(const Form<T> &a)

Create a sparsity pattern for a given form.

Note

The pattern is not finalised, i.e. the caller is responsible for calling SparsityPattern::assemble.

Parameters

a[in] A bilinear form

Returns

The corresponding sparsity pattern

ElementDofLayout create_element_dof_layout(const ufcx_dofmap &dofmap, const mesh::CellType cell_type, const std::vector<int> &parent_map = {})

Create an ElementDofLayout from a ufcx_dofmap.

DofMap create_dofmap(MPI_Comm comm, const ElementDofLayout &layout, mesh::Topology &topology, const std::function<std::vector<int>(const graph::AdjacencyList<std::int32_t>&)> &reorder_fn, const FiniteElement &element)

Create a dof map on mesh.

Parameters
  • comm[in] MPI communicator

  • layout[in] The dof layout on an element

  • topology[in] The mesh topology

  • element[in] The finite element

  • reorder_fn[in] The graph reordering function called on the dofmap

Returns

A new dof map

std::vector<std::string> get_coefficient_names(const ufcx_form &ufcx_form)

Get the name of each coefficient in a UFC form.

Parameters

ufcx_form[in] The UFC form return The name of each coefficient

std::vector<std::string> get_constant_names(const ufcx_form &ufcx_form)

Get the name of each constant in a UFC form.

Parameters

ufcx_form[in] The UFC form

Returns

The name of each constant

template<typename T>
Form<T> create_form(const ufcx_form &ufcx_form, const std::vector<std::shared_ptr<const FunctionSpace>> &spaces, const std::vector<std::shared_ptr<const Function<T>>> &coefficients, const std::vector<std::shared_ptr<const Constant<T>>> &constants, const std::map<IntegralType, const mesh::MeshTags<int>*> &subdomains, const std::shared_ptr<const mesh::Mesh> &mesh = nullptr)

Create a Form from UFC input.

Parameters
  • ufcx_form[in] The UFC form

  • spaces[in] Vector of function spaces

  • coefficients[in] Coefficient fields in the form

  • constants[in] Spatial constants in the form

  • subdomains[in] Subdomain markers

  • mesh[in] The mesh of the domain

template<typename T>
Form<T> create_form(const ufcx_form &ufcx_form, const std::vector<std::shared_ptr<const FunctionSpace>> &spaces, const std::map<std::string, std::shared_ptr<const Function<T>>> &coefficients, const std::map<std::string, std::shared_ptr<const Constant<T>>> &constants, const std::map<IntegralType, const mesh::MeshTags<int>*> &subdomains, const std::shared_ptr<const mesh::Mesh> &mesh = nullptr)

Create a Form from UFC input.

Parameters
  • ufcx_form[in] The UFC form

  • spaces[in] The function spaces for the Form arguments

  • coefficients[in] Coefficient fields in the form (by name)

  • constants[in] Spatial constants in the form (by name)

  • subdomains[in] Subdomain makers

  • mesh[in] The mesh of the domain. This is required if the form has no arguments, e.g. a functional

Returns

A Form

template<typename T>
Form<T> create_form(ufcx_form *(*fptr)(), const std::vector<std::shared_ptr<const FunctionSpace>> &spaces, const std::map<std::string, std::shared_ptr<const Function<T>>> &coefficients, const std::map<std::string, std::shared_ptr<const Constant<T>>> &constants, const std::map<IntegralType, const mesh::MeshTags<int>*> &subdomains, const std::shared_ptr<const mesh::Mesh> &mesh = nullptr)

Create a Form using a factory function that returns a pointer to a ufcx_form.

Parameters
  • fptr[in] pointer to a function returning a pointer to ufcx_form

  • spaces[in] The function spaces for the Form arguments

  • coefficients[in] Coefficient fields in the form (by name)

  • constants[in] Spatial constants in the form (by name)

  • subdomains[in] Subdomain markers

  • mesh[in] The mesh of the domain. This is required if the form has no arguments, e.g. a functional.

Returns

A Form

FunctionSpace create_functionspace(const std::shared_ptr<mesh::Mesh> &mesh, const basix::FiniteElement &e, int bs, const std::function<std::vector<int>(const graph::AdjacencyList<std::int32_t>&)> &reorder_fn = nullptr)

Create a FunctionSpace from a Basix element.

Parameters
  • mesh[in] Mesh

  • e[in] Basix finite element

  • bs[in] The block size, e.g. 3 for a ‘vector’ Lagrange element in 3D

  • reorder_fn[in] The graph reordering function to call on the dofmap. If nullptr, the default re-ordering is used.

Returns

The created function space

FunctionSpace create_functionspace(ufcx_function_space *(*fptr)(const char*), const std::string &function_name, const std::shared_ptr<mesh::Mesh> &mesh, const std::function<std::vector<int>(const graph::AdjacencyList<std::int32_t>&)> &reorder_fn = nullptr)

Create a FunctionSpace from UFC data.

Parameters
  • fptr[in] Function Pointer to a ufcx_function_space_create function

  • function_name[in] Name of a function whose function space to create. Function name is the name of Python variable for ufl.Coefficient, ufl.TrialFunction or ufl.TestFunction as defined in the UFL file.

  • mesh[in] Mesh

  • reorder_fn[in] The graph reordering function to call on the dofmap. If nullptr, the default re-ordering is used.

Returns

The created function space

template<typename T>
xtl::span<const std::uint32_t> get_cell_orientation_info(const std::vector<std::shared_ptr<const Function<T>>> &coefficients)
template<typename T, int _bs, typename Functor>
static inline void pack(const xtl::span<T> &coeffs, std::int32_t cell, int bs, const xtl::span<const T> &v, const xtl::span<const std::uint32_t> &cell_info, const DofMap &dofmap, Functor transform)
template<typename T, typename E, typename Functor>
void pack_coefficient_entity(const xtl::span<T> &c, int cstride, const Function<T> &u, const xtl::span<const std::uint32_t> &cell_info, const E &entities, Functor fetch_cells, std::int32_t offset)

Pack a single coefficient for a set of active entities.

Parameters
  • c[out] The coefficient to be packed

  • cstride[in] The total number of coefficient values to pack for each entity

  • u[in] The function to extract data from

  • cell_info[in] Array of bytes describing which transformation has to be applied on the cell to map it to the reference element

  • entities[in] The set of active entities

  • fetch_cells[in] Function that fetches the cell index for an entity in active_entities (signature: std::function<std::int32_t(E::value_type)>)

  • offset[in] The offset for c

template<typename T>
std::pair<std::vector<T>, int> allocate_coefficient_storage(const Form<T> &form, IntegralType integral_type, int id)

Allocate storage for coefficients of a pair (integral_type, id) from a fem::Form form.

Parameters
  • form[in] The Form

  • integral_type[in] Type of integral

  • id[in] The id of the integration domain

Returns

A storage container and the column stride

template<typename T>
std::map<std::pair<IntegralType, int>, std::pair<std::vector<T>, int>> allocate_coefficient_storage(const Form<T> &form)

Allocate memory for packed coefficients of a Form.

Parameters

form[in] The Form

Returns

A map from a form (integral_type, domain_id) pair to a (coeffs, cstride) pair

template<typename T>
void pack_coefficients(const Form<T> &form, IntegralType integral_type, int id, const xtl::span<T> &c, int cstride)

Pack coefficients of a Form for a given integral type and domain id.

Parameters
  • form[in] The Form

  • integral_type[in] Type of integral

  • id[in] The id of the integration domain

  • c[in] The coefficient array

  • cstride[in] The coefficient stride

template<typename T>
fem::Expression<T> create_expression(const ufcx_expression &expression, const std::vector<std::shared_ptr<const fem::Function<T>>> &coefficients, const std::vector<std::shared_ptr<const fem::Constant<T>>> &constants, const std::shared_ptr<const mesh::Mesh> &mesh = nullptr, const std::shared_ptr<const fem::FunctionSpace> &argument_function_space = nullptr)

Create Expression from UFC.

template<typename T>
fem::Expression<T> create_expression(const ufcx_expression &expression, const std::map<std::string, std::shared_ptr<const fem::Function<T>>> &coefficients, const std::map<std::string, std::shared_ptr<const fem::Constant<T>>> &constants, const std::shared_ptr<const mesh::Mesh> &mesh = nullptr, const std::shared_ptr<const fem::FunctionSpace> &argument_function_space = nullptr)

Create Expression from UFC input (with named coefficients and constants)

template<typename T>
void pack_coefficients(const Form<T> &form, std::map<std::pair<IntegralType, int>, std::pair<std::vector<T>, int>> &coeffs)

Pack coefficients of a Form.

Warning

This is subject to change

Parameters
  • form[in] The Form

  • coeffs[in] A map from a (integral_type, domain_id) pair to a (coeffs, cstride) pair

template<typename T>
std::pair<std::vector<T>, int> pack_coefficients(const Expression<T> &u, const xtl::span<const std::int32_t> &cells)

Pack coefficients of a Expression u for a give list of active cells.

Parameters
  • u[in] The Expression

  • cells[in] A list of active cells

Returns

A pair of the form (coeffs, cstride)

template<typename U>
std::vector<typename U::scalar_type> pack_constants(const U &u)

Pack constants of u of generic type U ready for assembly.

Warning

This function is subject to change