Factorization failure for certain conditionals

[mrambausek]

While probably quite a corner case, I managed to stumble into that...

import ufl
from ffcx.codegeneration import jit

elmt = ufl.FiniteElement("CG", ufl.triangle, 1)

u1 = ufl.Coefficient(elmt)
u2 = ufl.Coefficient(elmt)
ctol = ufl.Constant(ufl.triangle)


def cond(tol):
    return ufl.conditional(u2 < tol, u1, u2) * ufl.dx
    

# works
jit.compile_forms([ufl.derivative(cond(0), (u1, u2))])

# works
jit.compile_forms([ufl.derivative(cond(ctol), u1)])

# works
jit.compile_forms([ufl.derivative(cond(ctol), u2)])

# fails
jit.compile_forms([ufl.derivative(cond(ctol), (u1, u2))])

[michalhabera]

@mrambausek Can you include also an error message? Thanks!

[mrambausek]

~/opt/fenics/2019.1.0-fedora31-workstation/fenics-py/lib64/python3.7/site-packages/ffcx/ir/uflacs/analysis/factorization.py in compute_argument_factorization(S, rank)
    340             for o in expr.ufl_operands:
--> 341                 F.add_edge(i, F.e2i[o])
    342 

KeyError: Zero((), (), ())

[mrambausek]

this fails as well:

jit.compile_forms([ufl.derivative(cond(1e-6), (u1, u2))]) 

which is quite surprising since the only change is from a float "0" to a float "1e-6"

[mrambausek]

~/opt/fenics/2019.1.0-fedora31-workstation/fenics-py/lib64/python3.7/site-packages/ffcx/ir/uflacs/analysis/factorization.py in compute_argument_factorization(S, rank)
    340             for o in expr.ufl_operands:
--> 341                 F.add_edge(i, F.e2i[o])
    342 

KeyError: Zero((), (), ())

What seems imporant, is that at least one of the derivatives wrt. (u1, u2) has to vanish. This conditional, for example, works:

def cond(tol):
    return ufl.conditional(u2 < tol, u1*u2, 2*u1*u2) * ufl.dx

[mrambausek]

At least for the MW(F)E, a possible workaround seems to be

v1, v2 = ufl.TestFunctions(u1.ufl_element() * u2.ufl_element())
jit.compile_forms([ufl.derivative(cond(ctol), u2, v2) + ufl.derivative(cond(ctol), u1, v1)])

[michalhabera]

Just note that taking a (!distributional) derivative of a conditional is in most use cases mathematical non-sense, unless you add an extra delta measure which comes from the differentiation of a discontinuous function. What UFL and FFC do is that they differentiate piecewise, which is not always correct.

But I'll have a better look at why this issue appeared. I rewrote some factorization parts a while ago when allowing tensor-valued expressions to be factorized. So the issue originated probably at that point.

[mrambausek]

Thanks a lot for looking into that -- if you find the time.

I would say, a valid use case for a distributional derivative of a conditional is when you have an energy density of which the derivative exists and is smooth, but has an unstable numerical evaluation (e.g., one needs to take limits). Antiderivatives of sigmoidal functions are candidates for such potentials. As a remedy, one can -- by means of a conditional -- replace the original energy density by a quadratic approximation around certain points. Now, in case of coupled problems, e.g., magneto-mechanics, problems as described here could appear in certain situations. Of course, this is still a corner case, given the special case that triggered the problem. I am happy with my current workaround, but I would highly appreciate a fix or at least a warning if it really is an unsupported operation.