How to create identifying the facets for hexahedral elements

dolfinx
,
1

Hi,
I used hexahedral elements to mesh a spherical shell (I uploaded the mesh file to github). There was a problem adding MeshTags to boundaries as suggested by Dokken. Some boundary facets are not marked.

The code and result are as follows:

boundaries = [(1, lambda x: np.isclose(x[0]**2+x[1]**2+x[2]**2, 1)),
              (2, lambda x: np.isclose(x[0]**2+x[1]**2+x[2]**2, 0.35**2))]

facet_indices, facet_markers = [], []
fdim = msh.topology.dim - 1
for (marker, locator) in boundaries:
    facets = mesh.locate_entities(msh, fdim, locator)
    facet_indices.append(facets)
    facet_markers.append(np.full_like(facets, marker))
facet_indices = np.hstack(facet_indices).astype(np.int32)
facet_markers = np.hstack(facet_markers).astype(np.int32)
sorted_facets = np.argsort(facet_indices)
facet_tag = mesh.meshtags(msh, fdim, facet_indices[sorted_facets], facet_markers[sorted_facets])
ds = ufl.Measure("ds", domain=msh, subdomain_data=facet_tag)

msh.topology.create_connectivity(msh.topology.dim-1, msh.topology.dim)
with io.XDMFFile(msh.comm, "facet_tags.xdmf", "w") as xdmf:
    xdmf.write_mesh(msh)
    xdmf.write_meshtags(facet_tag, msh.geometry)

image
image1053×964 51.4 KB

There is no problem with tetrahedral elements. So I want to know how to create identifying the facets for hexahedral elements correctly?
Thanks!

2

Try passing a tolerance (larger than the default one) as trailing arguments to np.isclose, see the numpy documentation.

If you are worried that doing so may end up locating entities which are not on the boundary, you may use mesh.locate_entities_boundary instead of mesh.locate_entities.

3

Thanks francesco-ballarin.
Visualizations using both methods have not changed. But I find that the size of facets is correct.
Could this be a display issue with Paraview?

4

Hi Zhaoyx, I tried your mesh file, all the boundary facets are marked correct with your code. 2 for inner ball, 1 for out ball surface. I use paraview version 5.11.2.

Screenshot from 2024-04-05 17-22-48

2 Likes
5

Please ensure that you usr the extract block filter and the xdmf3t reader when opening files with multiple tags in paraview.

6

Thank you all!
It worked after using Threshold.