Mixed elements with higher value ranks cannot be compiled

[Jan Blechta]

from dolfin import *

mesh = UnitSquareMesh(3, 3)
RT = FiniteElement('RT', mesh, 3)
CG = FiniteElement('CG', mesh, 2)
RTCG = MixedElement((RT, CG))
vRTCG = MixedElement(3*(2*(RTCG,),), value_shape=(3, 2)+RTCG.value_shape())

W = FunctionSpaceBase(mesh, vRTCG)

Last line raises ufl.log.UFLException: Illegal component index '(0, 0)' (value rank 2)for element (value rank 1). Although it fails during compilation, it seems that problem is rather in UFL. Rank 1 or rank 2 version of vRTCG seems to be working.

Although complicated elements can always be flattened, the inconvenience is that one gets rid of tensor structure and may need to reconstruct component trial and test functions by unflattening.

EDIT: After discussion with Martin it showed up that value_shape argument of MixedElement was never intended for this usage. (It should be used only by Vector/TensorElement.) Martin also prefers not to enable construction of Vector/TensorElement from another arbitrary sub-element. Agreement is therefore to extend MixedElement constructor by shape keyword argument.

[Martin Sandve Alnæs]

Has this ever worked? I can't recall intending it to. I believe the value_shape has always been used only by the VectorElement and TensorElement.

[Jan Blechta]

No, it did not. It's a pity that one cannot utilize richness/generality of UFL for manipulation with arbitrarily high rank tensors.

[Jan Blechta]

@martinal On the other hand, TensorElements of higher ranks seem to be working:

from dolfin import *
mesh = UnitSquareMesh(3, 3) 
V = TensorFunctionSpace(mesh, 'CG', 1, shape=(3, 2, 2))

Is there a chance that TensorElement could be enhanced to allow construction from another element, not just from family? I mean

TensorElement.__init__(self, element, domain, degree, shape=None, symmetry=None, quad_scheme=None)

in addition to existing

TensorElement.__init__(self, family, domain, degree, shape=None, symmetry=None, quad_scheme=None)

?

[Martin Sandve Alnæs]

Maybe. But I think you have much bigger issues here than a simple ufl fix. Nested mixed elements don't work that well at the dolfin level already, consider e.g. splitting. It would be much easier to stick with flat structure at the dolfin level and implement some generic reshape functionality in ufl:

RT = FiniteElement('RT', mesh, 3)
CG = FiniteElement('CG', mesh, 2)
RTCG = MixedElement((RT, CG))
vRTCG = MixedElement( 3*2*(RTCG,)  )

vRTCG2 = reshape(vRTCG, (3, 2)+RTCG.value_shape())

[Jan Blechta]

Introducing TensorElement.__init__(self, family=None, domain=None, degree=None, sub_element=None, shape=None, symmetry=None, quad_scheme=None) where user supplies either family, domain, degree or sub_element is very easy.

But, as you say, other functionality is missing in UFL and DOLFIN. For instance, one may not create shallow/deep subfunctions. Maybe this is easy to fix, but large clean-up would be desired.

[Martin Sandve Alnæs]

Actually I wouldn't overload TensorElement with that possibility, it is explicitly defined as a special case of mixed element where the family of all subelements are the same. Extending MixedElement with shape would be a better enhancement.