Expression

class dolfin.functions.expression.Expression(cppcode=None, element=None, cell=None, domain=None, degree=None, name=None, label=None, mpi_comm=None, **kwargs)

Bases: object

This class represents a user-defined expression.

Expressions can be used as coefficients in variational forms or interpolated into finite element spaces.

Arguments
cppcode
C++ argument, see below
element
Optional element argument
degree
Optional element degree when element is not given.
cell
Optional cell argument to used in code generation
domain
Optional argument to determine the geometric dimension
degree
Optional argument to determine the degree of the element
name
Optional argument to set name of the Variable
label
Optional argument to set label of Variable
mpi_comm
Optional argument to allow JIT compilation on mpi groups. The Expression is only available at ranks in the same group. The Expression is JIT compiled on the first rank of the group.
1. Simple user-defined JIT-compiled expressions

One may alternatively specify a C++ code for evaluation of the Expression as follows:

f0 = Expression('sin(x[0]) + cos(x[1])')
f1 = Expression(('cos(x[0])', 'sin(x[1])'), element = V.ufl_element())

Here, f0 is is scalar and f1 is vector-valued.

Tensor expressions of rank 2 (matrices) may also be created:

f2 = Expression((('exp(x[0])','sin(x[1])'),
                ('sin(x[0])','tan(x[1])')))

In general, a single string expression will be interpreted as a scalar, a tuple of strings as a tensor of rank 1 (a vector) and a tuple of tuples of strings as a tensor of rank 2 (a matrix).

The expressions may depend on x[0], x[1], and x[2] which carry information about the coordinates where the expression is evaluated. All math functions defined in <cmath> are available to the user.

User defined parameters can be included as follows:

f = Expression('A*sin(x[0]) + B*cos(x[1])', A=2.0, B=Constant(4.0))

The parameters can be scalars and any scalar valued GenericFunction and are all initialized to the passed default value. The user defined parameters can be accessed and set as attributes or via the dict-like user_parameters attribute:

f.A = 5.0
f.B = Expression("value", value=6.0)
f.user_parameters["A"] = 1.0
f.user_parameters["B"] = Constant(5.0)

A parameter can only be updated with its original value-type. So if a parameter is a float, it can only be updated with float.

2. Complex user-defined JIT-compiled Expressions

One may also define an Expression using more complicated logic with the ‘cppcode’ argument. This argument should be a string of C++ code that implements a class that inherits from dolfin::Expression.

The following code illustrates how to define an Expression that depends on material properties of the cells in a Mesh. A MeshFunction is used to mark cells with different properties.

Note the use of the ‘cell’ parameter.

code = '''
class MyFunc : public Expression
{
public:

  std::shared_ptr<MeshFunction<std::size_t> > cell_data;

  MyFunc() : Expression()
  {
  }

void eval(Array<double>& values, const Array<double>& x,
          const ufc::cell& c) const
  {
    assert(cell_data);
    const Cell cell(*cell_data->mesh(), c.index);
    switch ((*cell_data)[cell.index()])
    {
    case 0:
      values[0] = exp(-x[0]);
      break;
    case 1:
      values[0] = exp(-x[2]);
      break;
    default:
      values[0] = 0.0;
    }
  }
};'''

cell_data = CellFunction('uint', V.mesh())
f = Expression(code)
f.cell_data = cell_data

While JIT compiling an Expression the public interface is scanned for dependencies. A user therefore need to include header files (also dolfin header files) declaring types that is used inside a method. The following example illustrates this. Note the inclusion of the dolfin namespace.

... code-block:: python

code = ‘’’ #include “dolfin/fem/GenericDofMap.h” namespace dolfin {

class Delta : public Expression { public:

Delta() : Expression() {}

void eval(Array<double>& values, const Array<double>& data,
const ufc::cell& cell) const

{ }

void update(const std::shared_ptr<const Function> u,
double nu, double dt, double C1, double U_infty, double chord)
{

const std::shared_ptr<const Mesh> mesh = u->function_space()->mesh(); const std::shared_ptr<const GenericDofMap> dofmap = u->function_space()->dofmap(); const uint ncells = mesh->num_cells(); uint ndofs_per_cell; if (ncells > 0) {

CellIterator cell(*mesh); ndofs_per_cell = dofmap->cell_dimension(cell->index());

} else {

return;

}

}

};

}’‘’

e = Expression(code)

3. User-defined expressions by subclassing

The user can subclass Expression and overload the ‘eval’ function. The value_shape of such an Expression will default to 0. If a user wants a vector or tensor Expression, the value_shape method needs to be overloaded.

class MyExpression0(Expression):
    def eval(self, value, x):
        dx = x[0] - 0.5
        dy = x[1] - 0.5
        value[0] = 500.0*exp(-(dx*dx + dy*dy)/0.02)
        value[1] = 250.0*exp(-(dx*dx + dy*dy)/0.01)
    def value_shape(self):
        return (2,)
f0 = MyExpression0()

If a user wants to use the Expression in a UFL form and have more controll in which finite element should be used to interpolate the expression in, the user can pass this information using the element kwarg:

V = FunctionSpace(mesh, "BDM", 1)
f1 = MyExpression0(element=V.ufl_element())

The user can also subclass Expression by overloading the eval_cell function. By this the user gets access to more powerful data structures, such as cell, facet and normal information, during assembly.

class MyExpression1(Expression):
    def eval_cell(self, value, x, ufc_cell):
        if ufc_cell.index > 10:
            value[0] = 1.0
        else:
            value[0] = -1.0

f2 = MyExpression1()

The ufc_cell object can be queried for the following data:

ufc_cell.cell_shape
ufc_cell.index
ufc_cell.topological_dimension
ufc_cell.geometric_dimension
ufc_cell.local_facet # only available on boundaries, otherwise -1
ufc_cell.mesh_identifier

The user can customize initialization of derived Expressions. However, because of magic behind the scenes, a user needs to pass optional arguments to __init__ using **kwargs, and _not_ calling the base class __init__:

class MyExpression1(Expression):
    def __init__(self, mesh, domain):
        self._mesh = mesh
        self._domain = domain
    def eval(self, values, x):
        ...

f3 = MyExpression1(mesh=mesh, domain=domain)

Note that subclassing may be significantly slower than using JIT-compiled expressions. This is because a callback from C++ to Python will be involved each time a Expression needs to be evaluated during assembly.

str() <==> info(x)
ufl_element()

Return the ufl FiniteElement.

ufl_evaluate(x, component, derivatives)

Function used by ufl to evaluate the Expression

value_shape()

Returns the value shape of the expression