Optimize assembly for jump terms

Hi all,

I have a form like:

I = \gamma * inner(jump(d2udn2) , jump(d2vdn2) )*dS

where \gamma is a function depending on u, d2udn2 and d2vdn2 are the second order normal derivatives for vector functions u and v, defined as:

d2udn2 = grad(grad(u)*n)*n

Each time I found the assembly for this term is extremely slow and memory-consuming, and in 3D problems, it will trigger the killed issue.

The parallelization is done by mpirun -n 24, with a 24-core Xeon CPU, and my cube mesh is constructed by 20x20x20 hexes. Thus I don’t think all memory will be eaten up.

parameters["ghost_mode"] = "shared_facet" is used for the integration on interior facets.

I’ve tried optimizers like:
parameters['form_compiler']['cpp_optimize'] = True
or
parameters["form_compiler"]['representation'] = 'uflacs'
But the killed issue still occurs.

I am really appreciated for your help!

Which function space and what order polynomial are you using?

It is a 2nd order Raviart-Thomas space using hex mesh:
RT_2 = P_{3,2,2} \times P_{2,3,2} \times P_{2,2,3}.

I guess my question is related to this issue posed on DOLFIN repo: https://bitbucket.org/fenics-project/dolfin/issues/1072/segmentation-fault-when-assembling-forms

Thanks!

Referring to: https://bitbucket.org/fenics-project/dolfin/issues/1072/segmentation-fault-when-assembling-forms

This is a MWE replicating this issue on the interior facets:

from dolfin import *
mesh = UnitCubeMesh.create(2,2,2,CellType.Type.hexahedron)
V = FunctionSpace(mesh,"Lagrange",3)
u = TrialFunction(V)
v = TestFunction(V)

# Working
Idx = inner(grad(grad(u)), grad(grad(v)))*dx
assemble(Idx)

# Working
IdS1 = inner(jump( grad(u) ), jump( grad(v) ) )*dS
assemble(IdS1)

# Not working
IdS2 = inner( jump( grad(grad(u)) ), jump( grad(grad(v)) ) )*dS
assemble(IdS2)

Thanks!

Hex support in old dolfin is really limited. I would suggest moving to dolfinx, Where @mscroggs had been working on extending the number of supported function spaces.

I’m unsure to what extent H(div) and H(curl) spaces are supported in the dolfin, but I’d use any functionality there with care as if these spaces are used on meshes whose cells are not correctly ordered, you could get wrong results (and there is no functionality to reorder non-simplex meshes).

H(div) and H(curl) spaces on quads and hexes are currently in development in dolfinx, but are not yet supported.