CellSize of higher order meshes

[APaganini]

I'm thinking about how to generalize CellSize to higher order meshes. I suspect that, given that the polynomial boundary approximation theory relies on the assumption that the true boundary is h^(p+1) away from the flat cell edge (for p-FEM approximations), we could just use the CellSize of the underlying simplicial mesh (although I'm clueless about how to do that). Any advice?

Bonus question: Why does CellSize call ufl.Circumradius instead of ufl.CellDiameter ? For shape-regular meshes these are probrably equivalent, but most of the FEM-theory I read in the past always employed the equivalent of ufl.CellDiameter.

[wence-]

I'm thinking about how to generalize CellSize to higher order meshes. I suspect that, given that the polynomial boundary approximation theory relies on the assumption that the true boundary is h^(p+1) away from the flat cell edge (for p-FEM approximations), we could just use the CellSize of the underlying simplicial mesh (although I'm clueless about how to do that). Any advice?

Here's some discussion that may be relevant https://github.com/firedrakeproject/tsfc/issues/224

My approach is complicated, @medinaeder's is probably more straightforward: you just replace CellDiameter with ReferenceCellDiameter() * detJ. The issue there is that detJ is not constant over the quad points so you get a slightly different value at each quad point. I would also be worried about the dimensions not being correct.

Your suggestion of using the flat mesh is doable too, it needs some thought about how to implement things in the form compiler/UFL. Effectively you would say something like:

mesh = FlatMesh(...)
cmesh = CurvedMesh(mesh)

h = CellDiameter(mesh)

V = FunctionSpace(cmesh, ...)

v = TestFunction(V)

L = h * v*dx # integrate over cmesh, but use h from flat mesh.

Bonus question: Why does CellSize call ufl.Circumradius instead of ufl.CellDiameter ? For shape-regular meshes these are probrably equivalent, but most of the FEM-theory I read in the past always employed the equivalent of ufl.CellDiameter.

Historical accident I think, a patch to replace CellSize with CellDiameter will probably be happily accepted.

[APaganini]

1. The solution suggested by @medinaeder is useful if one needs to compute the volume of the element (but one should perform local integration over it). I agree with your suspicion that it does not scale correctly wrt the geometric dimension.

2. Your code snippet is exactly what I meant (although maybe one should have a solution that starts from the CurvedMesh). I have to confess that I am not 100% sure this is what is needed by DG, though. Any advice by DG experts is welcome.

3. I'll write the patch.

[APaganini]

The patch is here: #2401