Understand surface creation after "remove overlap volume id1 id2 modify larger" command

I have the following script to create mesh from boolean operations on cylinders and bricks.

import sys

import numpy as np

sys.path.append("/opt/Coreform-Cubit-2025.3/bin")
import cubit

last_id = 0
LX = 20
LY = 20
lxs = np.arange(-0.5*LX+2.5, 0.5*LX, 5)
lys = np.arange(-0.5*LY+2.5, 0.5*LY, 5)
radius = 1.5
height = 50
for x in lxs:
    for y in lys:
        cubit.cmd(f"create cylinder radius {radius} height {height}")
        last_id += 1
        cubit.cmd(f"volume {last_id} move x {x} y {y} z 50")

cubit.cmd("create brick x 20 y 20 z 5")
last_id += 1
cubit.cmd(f"volume {last_id} move z 77.5")
cubit.cmd("unite all")
cubit.cmd("create brick x 20 y 20 z 80")
last_id = cubit.parse_cubit_list('volume', 'all')[-1]
cubit.cmd(f"volume {last_id} move z 40")
volume_ids = cubit.parse_cubit_list('volume', 'all')
cubit.cmd(f"remove overlap volume {volume_ids[0]} {volume_ids[1]} modify larger")
curves = cubit.parse_cubit_list('curve', 'all')
circle_arcs = []
for curve in curves:
    # print(cubit.get_curve_center(curve))
    if np.isclose(cubit.get_arc_length(curve), 2 * np.pi * radius):
        circle_arcs.append(curve)
cubit.cmd(f"modify curve {' '.join(map(str, circle_arcs))} blend radius 0.025")
volumes = cubit.parse_cubit_list('volume', 'all')
cubit.cmd(f"volume {volumes[0]} name 'vol1'")
cubit.cmd(f"volume {volumes[1]} name 'vol2'")
surfaces = cubit.parse_cubit_list('surface', 'all')
for surf in surfaces:
    centroid = cubit.get_surface_centroid(surf)
    area = cubit.get_surface_area(surf)
    if np.all(np.isclose(centroid[:2], [0, 0])):
        normal = cubit.get_surface_normal(surf)
        print(surf, "centroid:", centroid, "area:", area, "normal:", normal)
cubit.cmd("vol all scheme tetmesh")
cubit.cmd("mesh volume all")
filename = "mesh.bdf"
cubit.cmd(f'export nastran {filename} overwrite')

Can someone explain why at the intersection two surface ids are provided, and each have slightly different centroids?
`
55 centroid: (-3.120176742220682e-15, 4.019035726780615e-15, 75.00000000000122) area: 283.10133735992366 normal: (-0.0, -0.0, -1.0)

100 centroid: (-5.55674100982099e-16, 1.4934047829731165e-15, 75.00000000000074) area: 283.10133735992366 normal: (0.0, 0.0, 1.0)
`
I understand that the centroids are almost equal given float precision >>> np.finfo(float).eps np.float64(2.220446049250313e-16), but I did not expect that two distinct surfaces would be reported. A similar operation in GMSH produces only one distinct surface hence my confusion.

Hi @leshinka,
your script is constructing 2 unmerged volumes. This means that your geometry and also your mesh is not connected for those volumes. So every volume got their own surfaces where they overlap.

Do you need the mesh fully connected at the intersection between the two volumes?

1 Like

Yes, I need the mesh to be fully connected and to be able to label the two volumes as separate blocks so that I can specify different material properties in my FEA code.

I have figured out that the magic commands are:

cubit.cmd("imprint all")
cubit.cmd("merge all")