From The FEniCS project
The Unified Form Language (UFL) is a domain specific language for declaration of finite element discretizations of variational forms. More precisely, it defines a flexible interface for choosing finite element spaces and defining expressions for weak forms in a notation close to mathematical notation. The form compilers FFC and SyFi use UFL as their end-user interface, producing UFC implementations as their output.
Contents |
Finite element definitions
Cell definitions
Some ways to define a cell:
- cell = Cell(type, degree)
- cell = Cell("triangle", 1)
- cell = triangle
- cell = tetrahedron
- etc.
The object triangle is predefined and equivalent with Cell("triangle", 1).
Supported cell types (same as UFC):
- "interval"
- "triangle"
- "tetrahedron"
- "quadrilateral"
- "hexahedron"
Predefined Cell objects exist for each of these.
Element families
Long name | short name:
- family = "Lagrange" | "CG"
- family = "Discontinuous Lagrange" | "DG"
- family = "Bubble" | "B"
- family = "Crouzeix-Raviart" | "CR"
- family = "Brezzi-Douglas-Marini" | "BDM"
- family = "Brezzi-Douglas-Fortin-Marini" | "BDFM"
- family = "Raviart-Thomas" | "RT"
- family = "Nedelec 1st kind H(div)" | "N1div"
- family = "Nedelec 2nd kind H(div)" | "N2div"
- family = "Nedelec 1st kind H(curl)" | "N1curl"
- family = "Nedelec 2nd kind H(curl)" | "N2curl"
- family = "Quadrature" | "Q"
- family = "Boundary Quadrature" | "BQ"
New elements can be added dynamically by the form compiler.
Element definitions
Basic elements:
- element = FiniteElement(family, cell, degree)
Compound elements:
- element = VectorElement(family, cell, degree, dim=None)
- element = TensorElement(family, cell, degree, shape=None, symmetry=None)
Mixed elements:
- mixed_element = MixedElement(element1, element2[, element3, ...])
- mixed_element = element1 + element2 # equivalent with MixedElement(element1, element2)
Terminal expression types
Atomic expressions that do not depend on other UFL expressions.
Basis functions
- u = BasisFunction(element)
- v = TrialFunction(element)
- u = TestFunction(element)
- u, p = BasisFunctions(mixed_element): basisfunctions over the subelements of a mixed element
- u, p = TrialFunctions(mixed_element): trialfunctions over the subelements of a mixed element
- v, q = TestFunctions(mixed_element): testfunctions over the subelements of a mixed element
Coefficient functions
- f = Function(element)
- c = Constant(cell)
- v = VectorConstant(cell)
- t = TensorConstant(cell)
Geometric quantities
- n = cell.n # Facet normal
- x = cell.x # Global spatial coordinates
- d = cell.d # Geometric dimension
MeshSize has been removed, use a Constant instead.
Literals
- I = Identity(cell.d)
- Python scalars
Expression operators
Compound tensor algebra operators
These are operators that can be defined in terms of the basic operators. None of the operands here may have free indices.
- dot(a, b): contraction over the last axis of a and the first axis of b, rank(dot(a,b)) == rank(a)+rank(b)-2
- outer(a, b): outer product of a and b, rank(outer(a,b)) == rank(a) + rank(b)
- inner(a, b): contraction over all axes of a and b, they must have the exact same dimensions, rank(inner(a,b)) == 0
- cross(a, b): cross product of 3D vectors a and b
- transpose(A): transpose of matrix A, transpose(A)[i,j] == A[j,i]
- tr(A): trace of matrix A
- det(A): determinant of matrix A
- dev(A): deviatoric part of matrix A
- inv(A): inverse of matrix A (limited to 1x1, 2x2 and 3x3 matrices)
- cofac(A): cofactor of matrix A, cofactor(A) == det(A)*inverse(A) (limited to 1x1, 2x2 and 3x3 matrices)
- skew(A): skew symmetric part of A, skew(A)[i,j] = (A[i,j] - A[j,i]) / 2
Index notation
Defining indices:
- i = Index()
- i, j, k = indices(3)
- Predefined indices i, j, k, l and p, q, r, s
In the following, i-l and p-s are assumed distinct indices, the name v refers to some vector valued expression, and the name A refers to some matrix valued expression. The name C refers to a tensor expression of arbitrary rank.
Basic fixed indexing of a vector valued expression v or matrix valued expression A:
- v[0]: component access, representing the scalar value of the first component of v
- A[0,1]: component access, representing the scalar value of the first row, second column of A
Basic indexing:
- v[i]: component access, representing the scalar value of some component of v
- A[i,j]: component access, representing the scalar value of some component i,j of A
More advanced indexing:
- A[i,0]: component access, representing the scalar value of some component i of the first column of A
- A[i,:]: row access, representing some row i of A, i.e. rank(A[i,:]) == 1
- A[:,j]: column access, representing some column j of A, i.e. rank(A[:,j]) == 1
- C[...,0]: subtensor access, representing the subtensor of A with the last axis fixed, e.g., A[...,0] == A[:,0]
- C[j,...]: subtensor access, representing the subtensor of A with the last axis fixed, e.g., A[j,...] == A[j,:]
Implicit summation can occur in only a few situations. A product of two terms that shares the same free index is implicitly treated as a sum over that free index:
- v[i]*v[i]: sum_i v_i v_i (in LaTeX notation)
- A[i,j]*v[i]*v[j]: sum_j (sum_i A_{ij} v_i) v_j (in LaTeX notation)
A tensor valued expression indexed twice with the same free index is treated as a sum over that free index:
- A[i,i]: sum_i A_{ii} (in LaTeX notation)
- C[i,j,j,i]: sum_i sum_j C_{ijji} (in LaTeX notation)
The spatial derivative, in the direction of a free index, of an expression with the same free index, is treated as a sum over that free index:
- v[i].dx(i): sum_i v_i (in LaTeX notation)
- A[i,j].dx(i): sum_i d(A_{ij})/(dx_i) (in approximate LaTeX notation)
Note that these examples are some times written v_{i,i} and A_{ij,i} in pen-and-paper index notation.
Differential operators
Basic spatial derivatives:
- v.dx(i): partial derivative in coordinate direction i
- Dx(v, i): alternative notation for partial derivative
Compound differential operators:
- grad(f): gradient of expression f, defined as grad(f)[i,...] = f.dx(i)
- div(f): divergence of expression f, defined as div(f) = f[i,...].dx(i)
- curl(f): curl of expression f
- rot(f): rot of expression f
Derivatives w.r.t. user-defined variable:
- w = variable(w): annotate expression w as a variable that can be used in diff
- diff(f, w): the derivative of expression f w.r.t. the variable w
Restriction operators
- v('+'): restriction of v on positive side of facet
- v('-'): restriction of v on negative side of facet
- jump(v): jump of v across facet
- jump(v, n): jump of v across facet weighted by n
- avg(v): average of v across facet
Conditional operators
This is experimental and may need more work. In particular, we probably need and/or, and the syntax will be a bit clumsy.
- eq(a, b): condition that a == b
- ne(a, b): condition that a != b
- le(a, b): condition that a <= b
- ge(a, b): condition that a >= b
- lt(a, b): condition that a < b
- gt(a, b): condition that a > b
- conditional(condition, true_value, false_value): evaluate to true_value if condition is true, to false_value otherwise, similar to the C++ syntax (condition ? true_value: false_value)
- sign(f): the sign of f (+1 or -1)
Forms
A form is defined as the sum of integrals, like:
a = u*v*dx + f*v*ds
Integrals
Default integral operators like in FFC:
- (...) * dx: integrate (...) over all cells
- (...) * ds: integrate (...) over all exterior facets (the boundary)
- (...) * dS: integrate (...) over all interior facets (all facets not on the boundary)
Integrals over subdomains: dx(0), dx(1), ...
- (...) * dx(k) where k is an integer: integrate (...) over all cells in cell subdomain k (subset of domain)
- (...) * ds(k) where k is an integer: integrate (...) over all exterior facets in exterior facet subdomain k (subset of the boundary)
- (...) * dS(k) where k is an integer: integrate (...) over all interior facets in interior facet subdomain k (subset of facets not on the boundary)
Subdomains of cells, exterior facets, and interior facets have no correlation, ie subdomain k of the exterior facets does not have to be the boundary of subdomain k of the cells.
Integral over subdomain k with integration metadata:
- (...)*dx(k, metadata)
Transformations of entire forms
Automated computation of the derivative of a form with respect to a coefficient: If f is a rank 0 form that is nonlinear in coefficient w, F is the rank 1 derivative of f with respect to w and J is the rank 2 Jacobi of F with respect to w.
In the following examples, assume:
v = TestFunction(element) u = TrialFunction(element) w = Function(element)
Example:
f = (w**2)/2 * dx F = derivative(f, w, v) J = derivative(F, w, u)
is equivalent to:
f = (w**2)/2 * dx F = w*v*dx J = u*v*dx
Automated Action computation:
a = u*v*dx L = action(a, w)
is equivalent to:
a = u*v*dx L = w*v*dx
If a is a rank 2 form used to assemble the matrix A, L is a rank 1 form that can be used to assemble the vector b = Ax. Useful to define both the form of a matrix and the form of its action, and for the action of a Jacobi matrix computed using derivative.
Extraction of left- and right-hand sides:
F = v*(u - w)*dx + k*dot(grad(v), grad(0.5*(w + u)))*dx a, L = lhs(F), rhs(F)
should be equivalent to
a = v*u*dx + 0.5*k*dot(grad(v), grad(u))*dx L = v*w*dx - 0.5*k*dot(grad(v), grad(w))*dx
Main authors
Contributions
In no particular order:
Kristian B. Ølgaard, Garth Wells, Marie Rognes, Kent-Andre Mardal
License
UFL is licenced under GPLv3.
Dependencies and requirements
The UFL python code does not specify any additional dependencies.

