Illegal component index in UFL conditional()

Question

I have the following ULF code, where I try to make variable b conditional on 1< p, where p is a discrete Lagrange coefficient.

cell = triangle

V = VectorElement("Lagrange", cell, 2)
W = FiniteElement("DG", cell, 0)

p = Coefficient(W)

u = TestFunction(V)
v = TrialFunction(V)

b = conditional( lt(1,p), Identity(2), 2*Identity(2) )

a = inner(b,grad(u))*dx
L = inner(b,grad(v))*ds

but getting:

ufl.log.UFLException: Illegal component index '(0, 0)' (value rank 2)for element (value rank 0).

on compilation

Comments

Hi mdelmans,

I faced the same problem, and conditional() worked properly for me. so if at this point you haven't managed to solve it you should look for a bug previous to call conditional() because it works when it is used within a weak statement,

Regards,

Santiago

Answer 1

Hi, I think yours is an application of conditional which the authors of UFL/FFC did not have in mind. Consider asking them on the mailing list or registering an issue. In the mean time here's a fix - define your custom conditional which uses the scalar-valued one from UFL

from dolfin import *

def tensor_conditional(predicate, tvalue, fvalue):
    assert tvalue.ufl_shape == fvalue.ufl_shape

    # Scalar
    rank = len(tvalue.ufl_shape)
    if rank == 0: 
        return conditional(predicate, tvalue, fvalue)
    # Vector
    elif rank == 1:
        dim = tvalue.ufl_shape[0]
        conds = [conditional(predicate, tvalue[i], fvalue[i]) for i in range(dim)]
    # Matrices
    elif rank == 2:
        nrows, ncols = tvalue.ufl_shape
        conds = [[conditional(predicate, tvalue[i, j], fvalue[i, j]) for j in range(ncols)]
                 for i in range(nrows)]
    # Generalize later
    else:
        raise ValueError

    return as_tensor(conds)

# ----------------------------------------------------------------------------

mesh = UnitSquareMesh(10, 10)

V = VectorFunctionSpace(mesh, 'CG', 1)
W = FunctionSpace(mesh, 'DG', 0)

w = interpolate(Expression('x[0]*x[0]+x[1]*x[1]'), W)

x, y = SpatialCoordinate(mesh)

# This fails with your error
try:
    b = conditional(gt(w, 1), as_vector((x, y)), as_vector((-x, -y)))
    f = project(b, V)
except Exception as e:
    print '\033[1;37;31m%s\033[0m' % e
    print 'Try with tensor_conditional'

# Works!
b = tensor_conditional(gt(w, 1), as_vector((x, y)), as_vector((-x, -y)))
f = project(b, V)
plot(f)

V = FunctionSpace(mesh, 'CG', 1)
b = tensor_conditional(gt(w, 1),
                       as_tensor(((1, 0), (0, 2))),
                       Constant(1+1E-2)*as_tensor(((2, 0), (0, 1))))
f = project(det(b), V)
plot(f)

interactive()

You can put the tensor_conditional into your ufl file and use it in the other definitions.