ufl team mailing list archive

[Kristian Oelgaard <k.b.oelgaard@xxxxxxxxxx>]

Hi,

Garth, Anders and I have discussed the possibility of defining control
statements in FFC/UFL. This can be useful for many things but if we
take DG advection diffusion as an example we have the following on an
exterior facet:

scalar = FiniteElement("Discontinuous Lagrange", "triangle", 1)
vector = VectorElement("Lagrange", "triangle", 2)

v = TestFunction(scalar)
u = TrialFunction(scalar)

b   = Function(vector)
n   = FacetNormal("triangle")
# outflow facet switch, computed in DOLFIN (has values 1.0 or 0.0)
of  = Function(constant)

a = dot(mult(v, n), mult(b, of*u))*ds

This works, but means that one has to compute the value of 'of' in DOLFIN.
This is currently handled by OutflowFacet in SpecialFunctions.h but to me
it seems silly to add an additional argument to a form the value of which
depends on the value of other arguments to the form. Wouldn't it be much
nicer to let FFC/UFL produce code that can decide which values to use in
tabulate_tensor()?

So the question is how we can do this. I think the following could work
but I don't know for sure and I definitely do not know if it's the best
solution:


A new FFC/UFL type with possible interface like this:

Condition( left, logical_operator, right, return_value, else (optional) )

that can work as a placeholder for an expression or alternatively just
represent a real number (a switch/bool).

The syntax in case of a real number for advection-diffusion:

condition = Condition(dot(n, b), ">", 0.0, 1.0)
The logical operator ">" should be in Python syntax, then we can always
map it in ufcformat.py (or any other format)

Then

a = dot(mult(v, n), mult(b, condition*u))*ds


and in tabulate_tensor() we would have code that looks something like:

cond0 = 0.0;
if (w[0][0]*w[1][0] + w[0][1]*w[1][1] > 0.0)
   cond0 = 1.0

A[0] = .... cond0* .....


Alternatively, condition could work like a placeholder

condition = Condition(dot(n, b), ">", 0.0, mult(b,u))

a = dot(mult(v, n), condition))*ds

The code would now be:

cond0[0] = 0.0;
cond0[1] = 0.0;

if (w[0][0]*w[1][0] + w[0][1]*w[1][1] > 0.0)
{
   cond0[0] = w[1][0]*...
   cond0[1] = w[1][1]*...
}


For a given tensor we can keep track of the number of conditions through
its own index and check if any conditions are present like we do for
coefficients, transforms, derivatives etc.

How would this fit into the UFL framework? Could it even work, or does
anyone see a cleaner way of doing it?


Kristian