← Back to team overview

ufl team mailing list archive

Re: control statements in UFL

 

On Wed, Aug 20, 2008 at 10:02:37AM +0200, Martin Sandve Alnæs wrote:
> 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)

I think this second option looks better.

But wouldn't it be possible to just overload <, >, = etc for UFLObject
so that one may do

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

or even

  a = dot(v*n, (left < right)*u + (1 - left < right)*sin(u))*ds

We could have left < right return a Conditional object which would
return either 1 or 0 (in the generated code) depending on whether the
comparison is true or false.

-- 
Anders


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

Attachment: signature.asc
Description: Digital signature


Follow ups

References