Using XDMFFile.read_mesh to read XDMF files but the dimension is changed

General
1

I have stl files of the 3D-structures, and I am trying to make a stress-strain curve.
So first I use the following code to convert stl file to xdmf file:

stl_filename = '1.stl'
msh_filename = '1.msh'

stl_mesh = trimesh.load_mesh(stl_filename)
meshio_mesh = meshio.Mesh(
    points=stl_mesh.vertices,
    cells={"triangle": [list(face) for face in stl_mesh.faces]}
)
meshio.write("mesh.xdmf", meshio_mesh)

and mesh.xdmf can be opened by paraviwe, which shows a 3D structure.Here is the figure:

xdmf
xdmf805Ă—633 28.9 KB

Then I use the following code to read the xdmf file:

with io.XDMFFile(MPI.COMM_WORLD, mesh.xdmf, 'r') as infile:
    domain = infile.read_mesh(name="Grid")

but the domain.geometry.dim is 2, not 3. So when I try to apply force on the surface of the structure and the deformation shows zero.
I wonder what went wrong. When I use “mesh.create_box” to build a structure, the dimension is 3.
Looking forward to your help sincerely!

2

Please provide the XDMFFile (and preferrably the h5-file if you are using a binary format for the data).
I would also like to see the output of

nodes = infile.read_geometry_data(name="Grid")
print(nodes.shape)
printnodes)
3

Thanks for your reply.
The result of print('nodes_shape:',nodes.shape) is nodes_shape: (952, 3)
And I also find that

with io.XDMFFile(comm, xdmf_filename, 'r') as infile:
    domain = infile.read_mesh(name="Grid")
print('1:',domain.geometry.dim)
print('2:',domain.topology.dim)

the result is 3 and 2, respectively.
The stl, XDMFFile and the h5 -file has been uploaded to the URL below.
zrx-explorer/FEniCSx_files (github.com)

4

This is as expected then. The domain.geometry.dim should be 3, while the topology dimension should be 2 as you are reading in triangles.

Please note that you have now read in a shell structure. You would need to mesh the interior of the stl (for instance with GMSH or ftetwild) to get a tetrahedral or hexahedral volume mesh.

5

Thanks a lot!
I have tried using GMSH to open a stl file, the code is here:

def add_physical_for_gmsh(ifile, ofile = None):
    gmsh.initialize()
    gmsh.open(ifile)
    gmsh.model.geo.synchronize()    
    volumes = gmsh.model.getEntities(dim=3)
    if volumes:
        physical_group_tag = gmsh.model.addPhysicalGroup(3, [volume[1] for volume in volumes], name="AllVolumes")
    else:
        raise RuntimeError()
    if not ofile:
        ofile = ifile
    gmsh.write(ofile)
    gmsh.finalize()

def createGeometryAndMesh(stlfile, outfile, angle=40, mesh_size_factor=1.0, includeBoundary = True, forceParametrizablePatches=True):
    gmsh.clear()
    gmsh.initialize()
    gmsh.merge(stlfile)

    gmsh.model.occ.removeAllDuplicates()
    gmsh.model.mesh.classifySurfaces(angle * math.pi / 180., includeBoundary,
                                     forceParametrizablePatches,
                                     180 * math.pi / 180.)
    gmsh.model.mesh.createGeometry()
    gmsh.model.geo.synchronize()

    s = gmsh.model.getEntities(2)
    l = gmsh.model.geo.addSurfaceLoop([e[1] for e in s])
    gmsh.model.geo.addVolume([l])
    gmsh.model.geo.synchronize()
    volumes = gmsh.model.getEntities(dim=3)
    physical_group_tag = gmsh.model.addPhysicalGroup(3, [volume[1] for volume in volumes], name="AllVolumes")
    gmsh.model.geo.synchronize()
    gmsh.option.setNumber("Mesh.MeshSizeFactor", mesh_size_factor)
    gmsh.model.mesh.generate(3)
    gmsh.option.setNumber("Mesh.Optimize", 1)
    gmsh.model.mesh.optimize("Netgen")

    gmsh.model.geo.synchronize()
    gmsh.write(outfile)
    gmsh.finalize()

stl_filename = '1.stl'
msh_filename = '1.msh'

createGeometryAndMesh(stl_filename, msh_filename)

add_physical_for_gmsh(msh_filename, ofile = None)
print('add_finished')
domain, cell_markers, facet_markers = dfx.io.gmshio.read_from_msh(msh_filename, MPI.COMM_WORLD, gdim=3)
print('ok')

add_finished has been printed, while the ok is not been printed.The error info is shown as follows:

Info    : [ 90%] Processing parametrizations
Info    : Done reading '1.msh'
Segmentation fault

I wonder if it is the problem with my computer or with the structure 1.stl,which I have shown before.
And I have also tried other stl file using this code. And the error info went like this:

Traceback (most recent call last):
  File "/zjlab/create_mesh_createGeometryAndMesh.py", line 79, in <module>
    createGeometryAndMesh(stl_filename, msh_filename)
  File "/zjlab/create_mesh_createGeometryAndMesh.py", line 65, in createGeometryAndMesh
    gmsh.model.mesh.generate(3)
  File "/root/miniconda3/lib/python3.12/site-packages/gmsh.py", line 2060, in generate
    raise Exception(logger.getLastError())
Exception: Invalid boundary mesh (overlapping facets) on surface 74 surface 141

maybe this is related to the problems of the structure.
And I have also tried tetwild, which can implement format conversion successfully. However, I can’t find the python-api of the tetwild and I have thousands of stl files to deal with.
So could you plz help me with this question, or do you have any good ideas to batch convert stl files?

6

I would strongly suggest you start by

  1. Inspect your 1.msh file with meshio.
    There you can look at what cells have been stored in your mesh, and of what dimension.

One thing you could add after meshing is removing of duplicate nodes, as that is usually a place where things go wrong.

It would also be good if you could dig into dolfinx.io.gmshio.read_from_msh, as it is simply calling the GMSH Python API to extract data.
Figuring out which lines of: dolfinx/python/dolfinx/io/gmshio.py at ce09b5e8b56cbd69bcf5f17afb4c2b7217b21706 · FEniCS/dolfinx · GitHub throws the error is a start to shed some light on the issue.

1 Like
7

Thank you for your help. I’ll study the suggestions in detail.