Marking surfaces with pygmsh 3D

I would like to create a simple cube with a sphere cut out from it.

For some reason I’m not able to assign a physical entity to the boolean difference. I also didn’t find how to marking the surface of the sphere.

"""http://jsdokken.com/converted_files/tutorial_pygmsh.html
"""

import pygmsh
import meshio


size_box = 1e-3
radius_bubble = 0.1e-3

geom = pygmsh.occ.Geometry()
model3D = geom.__enter__()


box = model3D.add_box(x0=(0-size_box/2, 0-size_box/2, 0-size_box/2), extents=(size_box, size_box, size_box))
bubble = model3D.add_ball(center=(0, 0, 0), radius=radius_bubble)
difference = model3D.boolean_difference(box, bubble)

model3D.synchronize()

model3D.add_physical(difference, "salt")
# mark the surface

geom.generate_mesh()


pygmsh.write("mesh.msh")

def create_mesh(mesh, cell_type, prune_z=False):
    cells = mesh.get_cells_type(cell_type)
    cell_data = mesh.get_cell_data("gmsh:physical", cell_type)
    points = mesh.points[:, :2] if prune_z else mesh.points
    out_mesh = meshio.Mesh(points=points, cells={cell_type: cells}, cell_data={"name_to_read":[cell_data]})
    return out_mesh


mesh_from_file = meshio.read("mesh.msh")
meshio.write("mesh.xdmf", mesh_from_file)

boundary_mesh = create_mesh(mesh_from_file, "triangle")
meshio.write("mesh_bubble_surfaces.xdmf", boundary_mesh)

volume_mesh = create_mesh(mesh_from_file, "tetra")
meshio.write("mesh_bubble_volumes.xdmf", volume_mesh)

This is the error

Traceback (most recent call last):
  File "d:/Projets/MIT/FESTIM/mesh_bubble_pygmsh.py", line 40, in <module>
    boundary_mesh = create_mesh(mesh_from_file, "triangle")
  File "d:/Projets/MIT/FESTIM/mesh_bubble_pygmsh.py", line 31, in create_mesh
    cell_data = mesh.get_cell_data("gmsh:physical", cell_type)
  File "C:\Users\user\AppData\Local\Programs\Python\Python38\lib\site-packages\meshio\_mesh.py", line 178, in get_cell_data
    return np.concatenate(
  File "<__array_function__ internals>", line 180, in concatenate
ValueError: need at least one array to concatenate

I would strongly suggest using the gmsh python interface, as it gives you more freedom to do marking of boundaries.

See for instance: Test problem 2: Flow past a cylinder (DFG 2D-3 benchmark) — FEniCSx tutorial

1 Like

Thanks for the answer!
Do you even recommend it in 3D?

Yes, I do. It has alot more functionality than pygmsh (as pygmsh is simply a wrapper around the GMSH python API).
You can of course mix the two of them, you can start creating your geometry with pygmsh, and then use the gmsh python API to fetch boundaries so you can mark them.