Hey guys More weirdness related to facet integrals from IP viscosity that only appears with COFFEE optimisations. This time both on the bendy branches and on master (of last week). In the code below, the matrices M, Ms and M1+M2 should all be the same. For DG1 they are not equal on a recent (last week, before any sprint merges) install of firedrake/ffc/ufl/coffee on master. It is also incorrect on an older install I had, so it's not a recent thing. DG1 does however produce the right answer on the bendy branches: firedrake (bendy_changes) + ffc/ufl (fd_bendy) + coffee (master). For RT1 on the other hand it's the opposite: the matrices are correct and equal on master, but with bendy some of them are different from the correct answer. All of this only if use COFFEE optimisation, i.e. with parameters['coffee']={} the matrices are correct and all equal. I realize at this point the best thing to do is probably to wait if this all gets magically fixed when the sprint is over - will report back if not - but just thought to let you guys know already Cheers Stephan from firedrake import * mesh2d = UnitSquareMesh(1,1) U = VectorFunctionSpace(mesh2d, "DG", 1) v = TestFunction(U) u = Function(U) #parameters['coffee']={} n = FacetNormal(mesh2d) def outer_jump(v, n): return outer(v('+'), n('+'))+outer(v('-'), n('-')) F1 = inner(outer_jump(v, n), outer_jump(u,n))*dS F2 = inner(outer_jump(v, n), outer_jump(n,u))*dS F = inner(outer_jump(v, n), outer_jump(u,n)+outer_jump(n,u))*dS M1 = assemble(derivative(F1, u)).M.values M2 = assemble(derivative(F2, u)).M.values Ms = assemble(derivative(F1+F2, u)).M.values M = assemble(derivative(F, u)).M.values print abs(M-(M1+M2)).max() print abs(M-Ms).max()