dolfinx.function

Collection of functions and function spaces

Functions

TensorFunctionSpace(mesh, element[, shape, …])

Create tensor finite element (composition of scalar elements) function space.

VectorFunctionSpace(mesh, element[, dim, …])

Create vector finite element (composition of scalar elements) function space.

Classes

Constant(domain, c, Sequence, float])

A constant with respect to a domain.

ElementMetaData(family, degree, form_degree)

Data for representing a finite element

Expression(ufl_expression, x, …)

Create dolfinx Expression.

Function(V, x, name)

A finite element function that is represented by a function space (domain, element and dofmap) and a vector holding the degrees-of-freedom

FunctionSpace(mesh, element, …)

A space on which Functions (fields) can be defined.

class dolfinx.function.Constant(domain, c: Union[numpy.ndarray, Sequence, float])[source]

Bases: ufl.constant.Constant

A constant with respect to a domain.

Parameters
  • domain (DOLFIN or UFL mesh) –

  • c – Value of the constant.

property value

Returns value of the constant.

class dolfinx.function.ElementMetaData(family: str, degree: int, form_degree: Optional[int] = None)[source]

Bases: tuple

Data for representing a finite element

Create new instance of ElementMetaData(family, degree, form_degree)

degree: int

Alias for field number 1

family: str

Alias for field number 0

form_degree: Optional[int]

Alias for field number 2

class dolfinx.function.Expression(ufl_expression: ufl.core.expr.Expr, x: numpy.ndarray, form_compiler_parameters: dict = {}, jit_parameters: dict = {})[source]

Bases: object

Create dolfinx Expression.

Represents a mathematical expression evaluated at a pre-defined set of points on the reference cell. This class closely follows the concept of a UFC Expression.

This functionality can be used to evaluate a gradient of a Function at the quadrature points in all cells. This evaluated gradient can then be used as input to a non-FEniCS function that calculates a material constitutive model.

Parameters
  • ufl_expression – Pure UFL expression

  • x – Array of points of shape (num_points, tdim) on the reference element.

  • form_compiler_parameters – Parameters used in FFCX compilation of this Expression. Run ffcx –help in the commandline to see all available options.

  • jit_parameters – Parameters controlling JIT compilation of C code.

Note

This wrapper is responsible for the FFCX compilation of the UFL Expr and attaching the correct data to the underlying C++ Expression.

eval(cells: numpy.ndarray, u: Optional[numpy.ndarray] = None) → numpy.ndarray[source]

Evaluate Expression in cells.

Parameters
  • cells – local indices of cells to evaluate expression.

  • u (optional) – array of shape (num_cells, num_points*value_size) to store result of expression evaluation.

Returns

u – The i-th row of u contains the expression evaluated on cells[i].

Return type

np.ndarray

Note

This function allocates u of the appropriate size if u is not passed.

property num_points

Return the number of evaluation points on the reference cell.

property ufl_expression

Return the original UFL Expression

property value_size

Return the value size of the expression

property x

Return the evaluation points on the reference cell

class dolfinx.function.Function(V: dolfinx.function.FunctionSpace, x: Optional[dolfinx.cpp.la.Vector] = None, name: Optional[str] = None)[source]

Bases: ufl.coefficient.Coefficient

A finite element function that is represented by a function space (domain, element and dofmap) and a vector holding the degrees-of-freedom

Initialize finite element Function.

collapse()[source]
compute_point_values()[source]
copy()[source]

Return a copy of the Function. The FunctionSpace is shared and the degree-of-freedom vector is copied.

eval(x: numpy.ndarray, cells: numpy.ndarray, u=None) → numpy.ndarray[source]

Evaluate Function at points x, where x has shape (num_points, 3), and cells has shape (num_points,) and cell[i] is the index of the cell containing point x[i]. If the cell index is negative the point is ignored.

property function_space

Return the FunctionSpace

property id

Return object id index.

interpolate(u) → None[source]

Interpolate an expression

property name

Name of the Function.

split()[source]

Extract any sub functions.

A sub function can be extracted from a discrete function that is in a mixed, vector, or tensor FunctionSpace. The sub function resides in the subspace of the mixed space.

sub(i: int)[source]

Return a sub function.

The sub functions are numbered from i = 0..N-1, where N is the total number of sub spaces.

ufl_evaluate(x, component, derivatives)[source]

Function used by ufl to evaluate the Expression

property vector

Return the vector holding Function degrees-of-freedom.

property x

Return the vector holding Function degrees-of-freedom.

class dolfinx.function.FunctionSpace(mesh: dolfinx.cpp.mesh.Mesh, element: Union[ufl.finiteelement.finiteelementbase.FiniteElementBase, dolfinx.function.ElementMetaData], cppV: Optional[dolfinx.cpp.function.FunctionSpace] = None, form_compiler_parameters: dict = {}, jit_parameters: dict = {})[source]

Bases: ufl.functionspace.FunctionSpace

A space on which Functions (fields) can be defined.

Create a finite element function space.

clone()dolfinx.function.FunctionSpace[source]

Return a new FunctionSpace \(W\) which shares data with this FunctionSpace \(V\), but with a different unique integer ID.

This function is helpful for defining mixed problems and using blocked linear algebra. For example, a matrix block defined on the spaces \(V \times W\) where, \(V\) and \(W\) are defined on the same finite element and mesh can be identified as an off-diagonal block whereas the \(V \times V\) and \(V \times V\) matrices can be identified as diagonal blocks. This is relevant for the handling of boundary conditions.

collapse(collapsed_dofs: bool = False)[source]
Collapse a subspace and return a new function space and a map from

new to old dofs.

Arguments
collapsed_dofs

Return the map from new to old dofs

Returns
FunctionSpace

The new function space.

dict

The map from new to old dofs (optional)

component()[source]

Return the component relative to the parent space.

contains(V) → bool[source]

Check whether a FunctionSpace is in this FunctionSpace, or is the same as this FunctionSpace.

property dim
property dofmap

Return the degree-of-freedom map associated with the function space.

dolfin_element()[source]

Return the DOLFIN element.

property element
property id

The unique identifier

property mesh

Return the mesh on which the function space is defined.

num_sub_spaces() → int[source]

Return the number of sub spaces.

sub(i: int)dolfinx.function.FunctionSpace[source]

Return the i-th sub space.

tabulate_dof_coordinates()[source]
ufl_cell()[source]
ufl_function_space() → ufl.functionspace.FunctionSpace[source]

Return the UFL function space

dolfinx.function.TensorFunctionSpace(mesh: dolfinx.cpp.mesh.Mesh, element: dolfinx.function.ElementMetaData, shape=None, symmetry: bool = None, restriction=None)dolfinx.function.FunctionSpace[source]

Create tensor finite element (composition of scalar elements) function space.

dolfinx.function.VectorFunctionSpace(mesh: dolfinx.cpp.mesh.Mesh, element: dolfinx.function.ElementMetaData, dim=None, restriction=None)dolfinx.function.FunctionSpace[source]

Create vector finite element (composition of scalar elements) function space.

dolfinx.function.singledispatch(func)[source]

Single-dispatch generic function decorator.

Transforms a function into a generic function, which can have different behaviours depending upon the type of its first argument. The decorated function acts as the default implementation, and additional implementations can be registered using the register() attribute of the generic function.