Possible bug on ufl.conditional

Hello!

I’am trying to run the following code, but it seems impossible to create the matrix related to the second derivative when the function to derive is built performing the multiplication between this function by ufl.conditional.

from mpi4py import MPI
import dolfinx
import ufl
mesh = dolfinx.generation.RectangleMesh(
    MPI.COMM_WORLD, [[0.0, 0.0, 0.0], [1, 1, 0.0]], [2, 2]
)
element = ufl.FiniteElement("Lagrange", mesh.ufl_cell(), degree=1)
V = dolfinx.FunctionSpace(mesh, element)
u = dolfinx.Function(V)
Id = ufl.Identity(2)

print("----- c_ok works")
c_ok = ufl.conditional(ufl.ge(1.0, 0.0), u**2 * Id, u**2 * Id)
F = ufl.inner(c_ok, c_ok) * ufl.dx
F_u = ufl.derivative(F, u, ufl.TestFunction(V))
F_u_u = ufl.derivative(F_u, u, ufl.TrialFunction(V))
dolfinx.fem.assemble_matrix(F_u_u)

print("----- c_ko does not work")
c_ko = u * ufl.conditional(ufl.ge(1.0, 0.0), u * Id, u * Id)
F = ufl.inner(c_ko, c_ko) * ufl.dx
F_u = ufl.derivative(F, u, ufl.TestFunction(V))
F_u_u = ufl.derivative(F_u, u, ufl.TrialFunction(V))
dolfinx.fem.assemble_matrix(F_u_u)

Therefore, for the second case I get the following error:

Traceback (most recent call last):
  File "/home/try.py", line 25, in <module>
    dolfinx.fem.assemble_matrix(F_u_u)
  File "/usr/lib/python3.9/functools.py", line 877, in wrapper
    return dispatch(args[0].__class__)(*args, **kw)
  File "/usr/local/lib/python3.9/dist-packages/dolfinx/fem/assemble.py", line 212, in assemble_matrix
    A = cpp.fem.create_matrix(_create_cpp_form(a))
  File "/usr/local/lib/python3.9/dist-packages/dolfinx/fem/assemble.py", line 27, in _create_cpp_form
    return Form(form)._cpp_object
  File "/usr/local/lib/python3.9/dist-packages/dolfinx/fem/form.py", line 42, in __init__
    self._ufc_form, module, self._code = jit.ffcx_jit(
  File "/usr/local/lib/python3.9/dist-packages/dolfinx/jit.py", line 61, in mpi_jit
    return local_jit(*args, **kwargs)
  File "/usr/local/lib/python3.9/dist-packages/dolfinx/jit.py", line 215, in ffcx_jit
    r = ffcx.codegeneration.jit.compile_forms([ufl_object], parameters=p_ffcx, **p_jit)
  File "/usr/local/lib/python3.9/dist-packages/ffcx/codegeneration/jit.py", line 167, in compile_forms
    impl = _compile_objects(decl, forms, form_names, module_name, p, cache_dir,
  File "/usr/local/lib/python3.9/dist-packages/ffcx/codegeneration/jit.py", line 232, in _compile_objects
    _, code_body = ffcx.compiler.compile_ufl_objects(ufl_objects, prefix=module_name, parameters=parameters)
  File "/usr/local/lib/python3.9/dist-packages/ffcx/compiler.py", line 103, in compile_ufl_objects
    ir = compute_ir(analysis, object_names, prefix, parameters, visualise)
  File "/usr/local/lib/python3.9/dist-packages/ffcx/ir/representation.py", line 94, in compute_ir
    irs = [
  File "/usr/local/lib/python3.9/dist-packages/ffcx/ir/representation.py", line 95, in <listcomp>
    _compute_integral_ir(fd, i, analysis.element_numbers, integral_names, parameters, visualise)
  File "/usr/local/lib/python3.9/dist-packages/ffcx/ir/representation.py", line 357, in _compute_integral_ir
    integral_ir = compute_integral_ir(itg_data.domain.ufl_cell(), itg_data.integral_type,
  File "/usr/local/lib/python3.9/dist-packages/ffcx/ir/integral.py", line 130, in compute_integral_ir
    F = compute_argument_factorization(S, rank)
  File "/usr/local/lib/python3.9/dist-packages/ffcx/ir/analysis/factorization.py", line 296, in compute_argument_factorization
    factors = handler(v, fac, sf, F)
  File "/usr/lib/python3.9/functools.py", line 877, in wrapper
    return dispatch(args[0].__class__)(*args, **kw)
  File "/usr/local/lib/python3.9/dist-packages/ffcx/ir/analysis/factorization.py", line 82, in handle_sum
    raise RuntimeError("Expecting equal argument rank terms among summands.")
RuntimeError: Expecting equal argument rank terms among summands.

This appears to be related to elimination of zero terms by the form compiler. You can see that it works fine with

Id = dolfinx.Constant(mesh,((1,0),
                            (0,1)))

and

z = dolfinx.Constant(mesh,0)
Id = ufl.as_matrix(((1,z),
                    (z,1)))

but fails with

Id = ufl.as_matrix(((1,0),
                    (0,1)))

In older versions of FEniCS there was a global parameter to switch off this optimization, but I don’t know how to do that in DOLFIN-X.

1 Like

I am getting the same error with FEniCS 2019. Do you remember the name of the global parameter?