Representation choice logic is flawed
If representation is chosen by integral metadata and parameter system choice is non-matching, wrong code for tabulate_tensor
is emitted. The code combines results of two different representations resulting in a rubbish.
The problem is reproducible by the following MWE
from dolfin import *
import numpy as np
import dijitso
#dijitso.set_log_level("debug")
mesh = UnitSquareMesh(1, 1)
V = FunctionSpace(mesh, 'N1curl', 1)
u = TestFunction(V)
# Works if parameter matches metadata
parameters['form_compiler']['representation'] = 'uflacs'
L0 = u[0]*dx(metadata={"representation": "uflacs"})
x0 = assemble(L0)
parameters['form_compiler']['representation'] = 'quadrature'
L1 = u[0]*dx(metadata={"representation": "quadrature"})
x1 = assemble(L1)
print(x0.array(), x1.array())
x1 -= x0
assert np.isclose(x1.norm('l2'), 0.0)
# Fails otherwise
parameters['form_compiler']['representation'] = 'auto'
L0 = u[0]*dx(metadata={"representation": "uflacs"})
x0 = assemble(L0)
L1 = u[0]*dx(metadata={"representation": "quadrature"})
x1 = assemble(L1)
print(x0.array(), x1.array())
x1 -= x0
assert np.isclose(x1.norm('l2'), 0.0)
This might be related to #112 as the solution might require some clean-up in representation choice. Helpful tool for resolving this is pull request #40.
Comments (6)
-
-
reporter Thanks for pointing this out. Here could be a cause of the bug.
-
Basically, if any integral metadata has representation set to uflacs or tsfc, that must be used for all integrals or an error raised in case of conflicting specification.
-
reporter This will be resolved by pull request #45.
compute_form_data
is now called correctly to match an effective representation. Using conflicting parameters/metadata is now error. -
reporter - changed milestone to 2017.1
-
reporter - changed status to resolved
Fix in master since c0b1920.
- Log in to comment
The problem is that the choice is not really per integral. The call to compute_form_data must be performed differently for uflacs vs quadrature/tensor, and is only done once for each form. I think it's fine if an erroneous mix of representations within a form throws an error; obviously rubbish code is not fine.