Problem with dolfinx.io.gmshio.model_to_mesh

Hi, I am trying to create a simple disk mesh with an inner disk in it that is cut out, but I am getting an annoying bug.

This is the code:


import gmsh
from dolfinx import mesh
import numpy as np
from mpi4py import MPI
import dolfinx


class DiskwithHoleMesh():
	def __init__(self, 
              radius_inner = 0.3, 
              radius = 1.0, 
			  center = (0.0, 0.0), 
			  resolution: float = 0.005,
			  output_dir: str = './meshes', 
			  ):
		
		self.radius = radius
		self.center = center
		self.resolution = resolution
		self.name = 'disk' 
		self.radius_inner = radius_inner
		self.output_dir = output_dir

		if not os.path.exists(self.output_dir):
			os.makedirs(self.output_dir)


	def create_mesh(self, 
				 order: int = 1):

		domain = self.create_disk(self.name, order)

		return domain


	def create_disk(self, name: str = 'disk_with_hole', order: int = 1) -> gmsh.model:

		gmsh.initialize()
		model = gmsh.model()

		model.add(name)
		model.setCurrent(name)

		outer_sphere = model.occ.addDisk(self.center[0], self.center[1], 0, self.radius, self.radius)
	
		inner_sphere = model.occ.addDisk(self.center[0], self.center[1], 0, self.radius_inner, 
									  self.radius_inner)

		model_dim_tags = model.occ.cut([(2, outer_sphere)], [(2, inner_sphere)], 
								 removeObject=True, removeTool=True)

		model.occ.synchronize()

		gmsh.option.setNumber("Mesh.CharacteristicLengthMin", self.resolution)
		gmsh.option.setNumber("Mesh.CharacteristicLengthMax", self.resolution)

		model.mesh.generate(2)
		gmsh.write(f"{self.output_dir}/{self.name}.msh")

		domain, cell_marker, facet_marker = dolfinx.io.gmshio.model_to_mesh(
			model, MPI.COMM_WORLD, 0, gdim=2)
		gmsh.finalize()

		return domain
	


disk_mesh = DiskwithHoleMesh(radius_inner=0.3, radius=1.0, resolution=0.05)
domain = disk_mesh.create_mesh(order=2)
disk_mesh.plot_mesh()

I am getting the following error:
in model_to_mesh(model, comm, rank, gdim, partitioner, dtype) [258](https://file+.vscode-resource.vscode-cdn.net/home/sebastian/my_code/GPstuff/~/anaconda3/envs/fenicsx-env/lib/python3.13/site-packages/dolfinx/io/gmshio.py:258) perm_sort = np.argsort(cell_dimensions) [260](https://file+.vscode-resource.vscode-cdn.net/home/sebastian/my_code/GPstuff/~/anaconda3/envs/fenicsx-env/lib/python3.13/site-packages/dolfinx/io/gmshio.py:260) # Broadcast cell type data and geometric dimension --> [261](https://file+.vscode-resource.vscode-cdn.net/home/sebastian/my_code/GPstuff/~/anaconda3/envs/fenicsx-env/lib/python3.13/site-packages/dolfinx/io/gmshio.py:261) cell_id = cell_information[perm_sort[-1]]["id"] [262](https://file+.vscode-resource.vscode-cdn.net/home/sebastian/my_code/GPstuff/~/anaconda3/envs/fenicsx-env/lib/python3.13/site-packages/dolfinx/io/gmshio.py:262) tdim = cell_information[perm_sort[-1]]["dim"] [263](https://file+.vscode-resource.vscode-cdn.net/home/sebastian/my_code/GPstuff/~/anaconda3/envs/fenicsx-env/lib/python3.13/site-packages/dolfinx/io/gmshio.py:263) num_nodes = cell_information[perm_sort[-1]]["num_nodes"] IndexError: index -1 is out of bounds for axis 0 with size 0

Any help on this is more than appreciated!

You should add Physical Surfaces otherwise, elements without any physical id do not get loaded by dolfinx, hence the error.

like this?


		# Add physical groups for the surfaces
		outer_surface = model.add_physical_group(2, [outer_sphere], tag=1)
		inner_surface = model.add_physical_group(2, [inner_sphere], tag=2)
		model.set_physical_name(2, outer_surface, "Outer Surface")
		model.set_physical_name(2, inner_surface, "Inner Surface")