`ufl.extract_blocks` is broken for `MixedElement`

[jorgensd]

Minimal example with basix:

import basix.ufl
import ufl

cell = "triangle"
c_el = basix.ufl.element("Lagrange", cell,1, shape=(2, ))
mesh = ufl.Mesh(c_el)

element = basix.ufl.element("Lagrange", cell, 1)
me = basix.ufl.mixed_element([element, element])
W = ufl.FunctionSpace(mesh, me)
v,q = ufl.TestFunctions(W)

F = v*ufl.dx  + q*ufl.dx
F0 = ufl.extract_blocks(F,0)
F1 = ufl.extract_blocks(F,1)
print(F0)
print(F1)

First issue is that sub_elements is a property not a list. Fixing this, one get into the issue that an argument of a MixedElement has part=None, making extraction fail.

[jorgensd]

The manipulation of forms with mixed elements makes it hard to test. Here is an example from ufl 2019.1.0 (when the extract_block at least gives a result):

from ufl import *
from ufl.domain import default_domain
from ufl.algorithms.formsplitter import block_split
domain = default_domain(triangle)

el = FiniteElement("CG", triangle, 1)
el2 = FiniteElement("CG", triangle, 2)

me = MixedElement([el, el2])
W = FunctionSpace(domain, me)
v, q = TestFunctions(W)
F00 = v*dx 
F01 = q*dx
F = F01 + F00
F0 = block_split(F,0)
print(f"{str(F)=}")
print(f"{str(F0)=}")
print(f"{str(F00)=}")

yielding:

str(F)='{ v_0[1] } * dx(<Mesh #-1>[everywhere], {})\n  +  { v_0[0] } * dx(<Mesh #-1>[everywhere], {})'
str(F0)='{ ([v_0, 0])[1] } * dx(<Mesh #-1>[everywhere], {})\n  +  { ([v_0, 0])[0] } * dx(<Mesh #-1>[everywhere], {})'
str(F00)='{ v_0[0] } * dx(<Mesh #-1>[everywhere], {})'

I don't see any good way of checking equality between F0 and F00 (thus checking the capabilities of extract block). We could of course test this in FFCx, but I would like the test to be in ufl as the function can potentially be used by firedrake etc. My main issue is: '{ ([v_0, 0])[1] } * dx(<Mesh #-1>[everywhere], {})\n which I believe is 0 if I understand the brackets correctly

[wence-]

Run it through the form preprocessing and those lists should simplify out, I think.

FWIW, we don't use the UFL function in firedrake because we want to reconstruct our own Arguments (that have concrete information attached to them)

[jorgensd]

Run it through the form preprocessing and those lists should simplify out, I think.

FWIW, we don't use the UFL function in firedrake because we want to reconstruct our own Arguments (that have concrete information attached to them)

Unfortunately that doesn't work. Mwe using https://github.com/FEniCS/ufl/pull/285

import ufl
import ufl.algorithms
from ufl.finiteelement import FiniteElement, MixedElement



cell = ufl.triangle
domain = ufl.Mesh(FiniteElement("Lagrange", cell, 1, (2,), ufl.identity_pullback, ufl.H1))
fe_scalar = FiniteElement("Lagrange", cell, 1, (), ufl.identity_pullback, ufl.H1)
fe_vector = FiniteElement("Lagrange", cell, 2, (2, ), ufl.identity_pullback, ufl.H1)

me = MixedElement([fe_vector, fe_scalar])

W = ufl.FunctionSpace(domain, me)
wh =  ufl.Coefficient(W)
uh, ph = ufl.split(wh)
v, q = ufl.TestFunctions(W)
F_0 = ufl.inner(uh,v) * ufl.dx(domain=domain)
F_1 = ph * q *ufl. dx
F = F_0 + F_1
f_prec = ufl.algorithms.preprocess_form(ufl.extract_blocks(F, 0), False)
F_0_prec = ufl.algorithms.preprocess_form(F_0, False)
print(str(f_prec))
print()
print(str(F_0_prec))

yielding

{ sum_{i_8} ([([v_0[0], v_0[1], 0])[0], ([v_0[0], v_0[1], 0])[1]])[i_8] * ([w_0[0], w_0[1]])[i_8]  } * dx(<Mesh #0>[everywhere], {})
  +  { ([v_0[0], v_0[1], 0])[2] * w_0[2] } * dx(<Mesh #0>[everywhere], {})

{ sum_{i_9} ([v_0[0], v_0[1]])[i_9] * ([w_0[0], w_0[1]])[i_9]  } * dx(<Mesh #0>[everywhere], {})

Any suggestions as how to change the PR to restore the former capability and make a sensible check?

[wence-]

https://github.com/firedrakeproject/firedrake/blob/master/firedrake/formmanipulation.py#L18-L32

and https://github.com/firedrakeproject/firedrake/blob/master/firedrake/formmanipulation.py#L67-L72

[jorgensd]

How should I credit firedrake in the ufl file, as we would need these changes in form splitter? The appropriate thing to do would be to copy the whole firedrake license into ufl.

[dham]

Firedrake is LGPL and UFL is GPL so you can just lift the code under the licences. You are supposed to maintain copyright headers, but since git blame shows that Lawrence wrote that code, and Lawrence is already listed as an author of UFL, I think that's moot and you can just do it.

[jorgensd]

Firedrake is LGPL and UFL is GPL so you can just lift the code under the licences. You are supposed to maintain copyright headers, but since git blame shows that Lawrence wrote that code, and Lawrence is already listed as an author of UFL, I think that's moot and you can just do it.

Ok. That was what i wanted to clarify. @wence- let me know if this doesn’t suit you :)

[wence-]

That's fine.

[wence-]

Note that if you do that, the usage I implemented there has this performance antipattern which you may wish to try and address https://github.com/firedrakeproject/firedrake/issues/3560