DOF orderings for RT1 on quadrilaterals

Hi all,

I am using dolfinx and working with Raviart–Thomas elements (lowest order / order 1) on quadrilaterals on a rectangular mesh. For a given cell I need to access the dofs, which works

dofs = V.dofmap.cell_dofs(cell_idx)

and then as a next step I need to get the coordinates of the four DoFs. I understand that the tabulate_dof_coordinates does not work but the four DoFs are four flux values over the four edges (see DefElement: Q H(div)). So what I am looking for is the coordinate of the center of the facet for each of the DoFs as a two-dimensional coordinate. Or alternatively, speaking how are the four DoFs located with respect to the center of the cell (which one is at the top, at the bottom, to the left, to the right). Is there a straightforward way of extracting this information?

This seems to be related to DOF ordering H(div) elements but I could not make enough sense for me from that.

Thank you very much!

Hello,

DOFs on the RT element are no point evaluations, they represent integral moments. So there is no position on the facet at which they are located. If you take a look on the Basix definition you will find evaluation points which are used to calculated the DOF value. But these points are nothing else then the quadrature points, required to evaluate the integral moments on the facet.

I hope that helps.

Best wishes!

Hello mbr,

thanks for the quick reply. Unfortunately, I am not sure what to make of this. I am probably not familiar enough with the basics of dolfinx to understand your reply. Which “Basisx definition” do you mean? Could you provide an URL please?

For every facet of my (quadrilateral) cell I have one DoF. So the information I am seeking is: what is the location of the facet (or its coordinates) corresponding to the DoF? Is there some (straightforward) way to retrieve this information?

Thank you very much.

I believe that I have found a solution to my problem. If someone is interested, this seems to work for me.

facet_dim = mesh.topology.dim - 1
mesh.topology.create_connectivity(facet_dim, 0)
num_facets = mesh.topology.index_map(facet_dim).size_local
facets = dolfinx.cpp.mesh.entities_to_geometry(mesh, facet_dim, np.arange(num_facets, dtype=np.int32), False)
points = mesh.geometry.x
for idx, facet in enumerate(facets):
    dof = fem.locate_dofs_topological(V=RT0, entity_dim=1, entities=[idx])[0]
    dof_coordinates[dof,:] = .5 * (mesh.geometry.x[facet[0], :2] + mesh.geometry.x[facet[1], :2])