Error converting VTK to XDMF and importing into dolfinx

General
1 

        import gmsh
        import meshio
        import pyvista as pv
        
        # Initialize Gmsh
        gmsh.initialize()
        gmsh.logger.start()
        
        # Define the dimensions of the beam
        w = 4
        h = 8
        d = 120  # depth (assuming a 3D beam)
        
        # Create the 3D box (volume)
        box_tag = gmsh.model.occ.add_box(0, 0, 0, w, h, d)
        
        # Synchronize the CAD model with the Gmsh model
        gmsh.model.occ.synchronize()
        
        # Create a volume from the box
        volumes = gmsh.model.getEntities(dim=3)
        
        # Define transfinite lines
        gmsh.model.mesh.set_transfinite_automatic()
        gmsh.model.mesh.set_transfinite_surface(volumes[0][1], "Left")
        gmsh.model.mesh.set_transfinite_volume(volumes[0][1])
        
        # Set the recombination algorithm to ensure hexahedral elements
        gmsh.option.setNumber("Mesh.RecombineAll", 1)
        gmsh.model.mesh.generate(3)
        
        # Output to MSH file (supported format)
        gmsh.write("beam.msh")
        
        # Output to VTK file (supported by meshio)
        gmsh.write("beam.vtk")
        
        # Finalize Gmsh
        gmsh.finalize()
        
        # Check for errors
        #log = gmsh.logger.get()
        #print(log)
        
        # Read the VTK file using meshio
        mesh = meshio.read("beam.vtk")
        
        # Inspect the VTK file using PyVista
        vtk_mesh = pv.read("beam.vtk")
        
        # Get the cell types
        cell_types = vtk_mesh.celltypes
        
        # Check if the cell type is HEXAHEDRON
        if pv.CellType.HEXAHEDRON in cell_types:
            print("The mesh contains hexahedron cells.")
        else:
            print("The mesh does not contain hexahedron cells.")
        
        # Plot the mesh
        plotter = pv.Plotter()
        plotter.add_mesh(vtk_mesh, show_edges=True, edge_color="red", line_width=2)
        plotter.show()



        Info    : Meshing 1D...
        Info    : [  0%] Meshing curve 1 (Line)
        Info    : [ 10%] Meshing curve 2 (Line)
        Info    : [ 20%] Meshing curve 3 (Line)
        Info    : [ 30%] Meshing curve 4 (Line)
        Info    : [ 40%] Meshing curve 5 (Line)
        Info    : [ 50%] Meshing curve 6 (Line)
        Info    : [ 60%] Meshing curve 7 (Line)
        Info    : [ 60%] Meshing curve 8 (Line)
        Info    : [ 70%] Meshing curve 9 (Line)
        Info    : [ 80%] Meshing curve 10 (Line)
        Info    : [ 90%] Meshing curve 11 (Line)
        Info    : [100%] Meshing curve 12 (Line)
        Info    : Done meshing 1D (Wall 0.0169652s, CPU 0s)
        Info    : Meshing 2D...
        Info    : [  0%] Meshing surface 1 (Transfinite)
        Info    : [ 20%] Meshing surface 2 (Transfinite)
        Info    : [ 40%] Meshing surface 3 (Transfinite)
        Info    : [ 60%] Meshing surface 4 (Transfinite)
        Info    : [ 70%] Meshing surface 5 (Transfinite)
        Info    : [ 90%] Meshing surface 6 (Transfinite)
        Info    : Done meshing 2D (Wall 0.0114789s, CPU 0s)
        Info    : Meshing 3D...
        Info    : Meshing volume 1 (Transfinite)
        Info    : Done meshing 3D (Wall 0.00202203s, CPU 0s)
        Info    : Optimizing mesh...
        Info    : Done optimizing mesh (Wall 7.15256e-06s, CPU 0s)
        Info    : 48 nodes 117 elements
        Info    : Writing 'beam.msh'...
        Info    : Done writing 'beam.msh'
        Info    : Writing 'beam.vtk'...
        Info    : Done writing 'beam.vtk'
        The mesh contains hexahedron cells.

``

![Screenshot 2024-07-19 090209|642x500](upload://3cz1rByo6BDi5FKnawcDgEEmnTu.png)




    import meshio
    
    # Read the VTK file
    
    vtk_mesh = meshio.read("beam.vtk")
    
    # Write to XDMF file
    
    meshio.write("beam.xdmf", vtk_mesh)
    
    print('done...')



    import dolfinx 
    
    from mpi4py import MPI
    
    print(dolfinx.__version__)
    
    with dolfinx.io.XDMFFile(MPI.COMM_WORLD, 'beam.xdmf', "r") as xdmf:
           mesh = xdmf.read_mesh(name="Grid")


Traceback:


    0.8.0
    Traceback (most recent call last):
      File "/home/prusso/dolfinx-demo/main.py", line 11, in <module>
        mesh = xdmf.read_mesh(name="Grid")
               ^^^^^^^^^^^^^^^^^^^^^^^^^^^
      File "/home/prusso/spack/var/spack/environments/fenicsx-env/.spack-env/view/lib/python3.11/site-packages/dolfinx/io/utils.py", line 258, in read_mesh
        cell_shape, cell_degree = super().read_cell_type(name, xpath)
                                  ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
    RuntimeError: Cannot recognise cell type. Unknown value: mixed




    I checked that the cell type of the mesh created with gmsh is hexahedron. Beyond that I really don't understand why the cell type can't be recognized or why it is being declared "mixed" since as far as I know there should be only hexahedrons.

Any way to get this resolved and have the mesh import into dolfinx?
2

Please consider using dolfinx.io.gmshio — DOLFINx 0.8.0 documentation rather than relying on meshio.

1 Like
3 
  import gmsh
  
  gmsh.initialize()
  gmsh.open("beam.msh")
  
  # Check for entities, nodes, and elements
  entities = gmsh.model.getEntities()
  nodes = gmsh.model.mesh.getNodes()
  elements = gmsh.model.mesh.getElements()
  
  print(f"Entities: {entities}")
  print(f"Nodes: {len(nodes[0])}")
  print(f"Elements: {len(elements[1][0])}")
  
  gmsh.finalize()



  Info    : Reading 'beam.msh'...
  Info    : 27 entities
  Info    : 48 nodes
  Info    : 117 elements
  Info    : Done reading 'beam.msh'
  Entities: [(0, 1), (0, 2), (0, 3), (0, 4), (0, 5), (0, 6), (0, 7), (0, 8), (1, 1), (1, 2), (1, 3), (1, 4), (1, 5), (1, 6), (1, 7), (1, 8), (1, 9), (1, 10), (1, 11), (1, 12), (2, 1), (2, 2), (2, 3), (2, 4), (2, 5), (2, 6), (3, 1)]
  Nodes: 48
  Elements: 52



  import dolfinx
  from dolfinx.io import gmshio
  from mpi4py import MPI
  from dolfinx.mesh import create_mesh
  
  # Initialize MPI
  comm = MPI.COMM_WORLD
  
  # Path to the .msh file
  msh_filename = "beam.msh"
  
  # Read the mesh from the .msh file
  gmsh_model, cell_type, geometric_dimension = gmshio.read_from_msh(msh_filename, comm)
  
  # Create the mesh
  mesh = create_mesh(comm, gmsh_model, cell_type, geometric_dimension)
  
  # Print some information about the mesh
  print(f"Mesh: {mesh}")
  print(f"Number of cells: {mesh.topology.index_map(dolfinx.mesh.CellType.tetrahedron).size_local}")
  print(f"Number of vertices: {mesh.topology.index_map(dolfinx.mesh.Vertex).size_local}")



Exception has occurred: IndexError
index -1 is out of bounds for axis 0 with size 0
  File "/home/prusso/dolfinx-demo/main.py", line 13, in <module>
    gmsh_model, cell_type, geometric_dimension = gmshio.read_from_msh(msh_filename, comm)
                                                 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
IndexError: index -1 is out of bounds for axis 0 with size 0

I’m not quire sure what it means yet…How can I go about solving this?

4

You must add physical groups to dimensions 3 and 2 to get gmshio.read_from_msh working. You’ve got the outputs of read_from_msh wrong, and the call to create_mesh is not needed.

Please in future do not unnecessarily indent the code.

5

You are correct about those… So far after correcting both thos problems, if I run Python in debug mode there is no difficulty loading the mesh from the beam.msh file:


from dolfinx.io import gmshio
from mpi4py import MPI
from dolfinx.mesh import create_mesh

# Initialize MPI
comm = MPI.COMM_WORLD

# Path to the .msh file
msh_filename = "beam.msh"


mesh = gmshio.read_from_msh(msh_filename, comm)

There does seem to be a difficulty if I use the regular run mode of Python I get this traceback:


Info    : Reading 'beam.msh'...
Info    : 27 entities
Info    : 48 nodes
Info    : 57 elements
Info    : Done reading 'beam.msh'
Entities: [(0, 1), (0, 2), (0, 3), (0, 4), (0, 5), (0, 6), (0, 7), (0, 8), (1, 1), (1, 2), (1, 3), (1, 4), (1, 5), (1, 6), (1, 7), (1, 8), (1, 9), (1, 10), (1, 11), (1, 12), (2, 1), (2, 2), (2, 3), (2, 4), (2, 5), (2, 6), (3, 1)]
Nodes: 48
Elements: 46
Traceback (most recent call last):
  File "/home/prusso/dolfinx-demo/main2.py", line 1, in <module>
    from dolfinx.io import gmshio
  File "/home/prusso/spack/var/spack/environments/fenicsx-env/.spack-env/view/lib/python3.11/site-packages/dolfinx/__init__.py", line 13, in <module>
    from petsc4py import PETSc as _PETSc
  File "/home/prusso/spack/var/spack/environments/fenicsx-env/.spack-env/view/lib/python3.11/site-packages/petsc4py/PETSc.py", line 4, in <module>
    PETSc = ImportPETSc(ARCH)
            ^^^^^^^^^^^^^^^^^
  File "/home/prusso/spack/var/spack/environments/fenicsx-env/.spack-env/view/lib/python3.11/site-packages/petsc4py/lib/__init__.py", line 30, in ImportPETSc
    return Import('petsc4py', 'PETSc', path, arch)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/prusso/spack/var/spack/environments/fenicsx-env/.spack-env/view/lib/python3.11/site-packages/petsc4py/lib/__init__.py", line 97, in Import
    module = import_module(pkg, name, path, arch)
             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/prusso/spack/var/spack/environments/fenicsx-env/.spack-env/view/lib/python3.11/site-packages/petsc4py/lib/__init__.py", line 76, in import_module
    spec.loader.exec_module(module)
  File "petsc4py/PETSc.pyx", line 1, in init petsc4py.PETSc
  File "/home/prusso/spack/var/spack/environments/fenicsx-env/.spack-env/view/lib/python3.11/site-packages/numpy/__init__.py", line 159, in <module>
    from . import ma
  File "/home/prusso/spack/var/spack/environments/fenicsx-env/.spack-env/view/lib/python3.11/site-packages/numpy/ma/__init__.py", line 42, in <module>
    from . import core
  File "/home/prusso/spack/var/spack/environments/fenicsx-env/.spack-env/view/lib/python3.11/site-packages/numpy/ma/core.py", line 7940, in <module>
    inner.__doc__ = doc_note(np.inner.__doc__,
                    ^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/prusso/spack/var/spack/environments/fenicsx-env/.spack-env/view/lib/python3.11/site-packages/numpy/ma/core.py", line 125, in doc_note
    notesplit = re.split(r'\n\s*?Notes\n\s*?-----', inspect.cleandoc(initialdoc))
                                                    ^^^^^^^^^^^^^^^^
AttributeError: module 'inspect' has no attribute 'cleandoc'

Is there a way to get this traceback to resolve?

6

The final question is unrelated to the original one, and very hard for anyone to reproduce. Double check the spack environments, and make sure that the inspect module isn’t overwritten by something else (e.g., a file called inspect.py in the PYTHONPATH).

1 Like