ufl team mailing list archive

["Martin Sandve Alnæs" <martinal@xxxxxxxxx>]

2008/8/20 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) )


I'm not a fan of the "else (optional)". I didn't understand
the example below intuitively before I read this line.

Maybe we can have named conditional operators (mirroring fortran):

  eq(left, right) # ==
  ne(left, right) # !=
  le(left, right) # <=
  lt(left, right) # <
  ge(left, right) # >=
  gt(left, right) # >

And the conditional can then be more "structured" (shorter argument
list, in a way):

  c = Conditional(a, lt(left, right), b)

mirroring "c = a if (left<right) else b" in python,
or:

  c = Conditional(lt(left, right), a, b)

mirroring "c = (left<right) ? a : b" in C++.

Which do you prefer?


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

a and b can easily be any UFL expression, while left and right must of
course be scalar.


> The syntax in case of a real number for advection-diffusion:
>
> condition = Condition(dot(n, b), ">", 0.0, 1.0)


I my syntax suggestions this will be:
  condition = Conditional(gt(dot(n, b), 0.0), 1.0, 0.0)
or
  condition = Conditional(1.0, gt(dot(n, b), 0.0), 0.0)

or written differently:

  cond = gt(dot(n, b), 0.0)
  condition = Conditional(cond, 1.0, 0.0)
or
  condition = Conditional(1.0, cond, 0.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


In UFL:

  a = dot(v*n, 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))

What's the difference?

> 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]*...
> }

Ok.

> 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.

What do you mean? Sounds like FFC details, maybe I don't need to know.

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


I think this should be ok in UFL, although of course it causes
some mathematical issues w.r.t. differentiation. I guess that
often the derivative of a Conditinal can simply be another
Conditinal with the same condition and the derivatives of
the other arguments. We just have to be careful.

I'll implement this in UFL, and see if I spot any difficulties.

--
Martin