TensorElement(...,symmetry=True) is totally broken
Generated code for element = TensorElement("Discontinuous Lagrange", triangle, 0, symmetry=True)
shows that the implementation in FFC is totally wrong.
-
Element has dimension 3 and value size 4.
evaluate_basis
sets only first three values. Obviously, one of basis functions should contribute to two entries invalues
. -
evaluate_dof
is wrong as well.
It is a mystery why the element seems to work at leat a bit, see DOLFIN report when both evaluate_basis
and evaluate_dof
is wrong.
Possibly, this is another candidate for implementing in FIAT.
Comments (20)
-
-
reporter I don't see why, the tensor element is basically the same as vector and mixed elements, and the symmetry must be known to the form compiler to detect common subexpressions in forms.
The reason is that basis (and dual basis) tabulated by FIAT would already respect symmetry and that would avoid the need for working it around again and again on every place (like now, it's taken into account in form compilation but not in element compilation). But my opinion is not strong here.
-
For efficient evaluation, ffc should generate code that executes (sum_i u_i phi_i(x)) once and copies the result into output values that by symmetry are equal. That's not fiats job. (Maybe it will be finats job in the future.)
-
reporter Sure, but it seems to me that once basis (tabulated by FIAT) is symmetric-tensor-valued (and so dimension is adequately smaller than for non-symmetric tensor element), then there would be no extra evaluation. But maybe I'm wrong here.
-
reporter - marked as critical
Either this needs a fix or error from
evaluate_basis/dof
. Otherwise we give to users code which computes rubbish. -
reporter -
assigned issue to
-
assigned issue to
-
reporter Could be fixable by modifying
ffc.representation._generate_foo_offsets
to take into account for symmetry. -
I tried to fix some symmetry handling in uflacs and hit this bug while testing.
Here's a test snippet that may be useful:
from dolfin import * from ufl.classes import FloatValue parameters["form_compiler"]["representation"] = "uflacs" def test_symmetries(): mesh = UnitSquareMesh(1, 1) element = TensorElement("CG", triangle, 1, symmetry=True) V = FunctionSpace(mesh, element) f = Function(V) c = Constant(((1,2), (2,3))) f.interpolate(c) fv = f((0.5,0.5)) value1 = assemble(c**2*dx(mesh)) value2 = assemble(f**2*dx(mesh)) print value1 # 18 = 1 + 4 + 4 + 9 print value2 # 13 = ??? assert abs(value1 - float(sum(fvc**2 for fvc in [1,2,2,3]))) < 1e-12 assert abs(value2 - float(sum(fvc**2 for fvc in [1,2,2,3]))) < 1e-12 # this is wrong, 0 instead of 3 for last component... for i in range(len(fv)): assert abs(fv[i] - [1.0, 2.0, 2.0, 3.0][i]) < 1e-12
-
Also look for FIXME comments around physical vs reference value shape/size in representation.py, there are a couple of fixmes pointing to bugs there.
-
And note that the UFL TensorElement supports generic symmetries, i.e.
symmetry[(i,j)] -> (k,l)
and not just symmetric-tensor-valued. -
reporter My first attempt was to store
ufl.FiniteElement.symmetry()
into representation and take this into account in code generation. But that's probably not right. Dealing right with value dimensions and offsets at generation stage seems to need a cleanup, as you say, and that might fix the issue.UFL symmetries, yes!
Also note how quadrature and tensor deal with symmetry. They use
from ufl.permutation import build_component_numbering
:- https://bitbucket.org/fenics-project/ffc/src/32692902865a7623c57c13228caa028f207d8308/ffc/quadrature/quadraturetransformerbase.py?at=master&fileviewer=file-view-default#quadraturetransformerbase.py-880
- https://bitbucket.org/fenics-project/ffc/src/32692902865a7623c57c13228caa028f207d8308/ffc/tensor/monomialtransformation.py?at=master&fileviewer=file-view-default#monomialtransformation.py-469
-
Good, uflacs also uses
build_component_numbering
, inuflacs/analysis/modified_terminals.py
. This is the code I've tried to clean up a bit: when the transformation to e.g.ReferenceValue(Coefficient(tensor_element))
has been applied, the symmetry dict in the element is no longer relevant, as theelement.reference_value_shape()
is flattened with the symmetries subtracted from the dimension (seeufl/finiteelement/mixedelement.py
).The element generation should probably also use
build_component_numbering
so we have one definition of the numbering. -
reporter -
assigned issue to
I believe you're about to fix this, aren't you? Feel free to reassign back to me.
-
assigned issue to
-
- changed milestone to 2017.1
-
reporter - changed component to finite-elements
-
reporter - changed milestone to 2017.2
-
- removed responsible
-
reporter This is still an issue with the new finite element code.
-
reporter - changed milestone to 2018.1
-
reporter It would be good to fix this. It is a bug of type "return junk values without warning".
- Log in to comment
It's not such a mystery. It's handled correctly in the form compilation but not in the element compilation. Those are separate code paths.
"Possibly, this is another candidate for implementing in FIAT."
I don't see why, the tensor element is basically the same as vector and mixed elements, and the symmetry must be known to the form compiler to detect common subexpressions in forms.