add hanging node tests before writing mesh to file

Author: prashjha

tools/python_utils/geom_util.py

@@ -71,15 +71,15 @@ def get_ref_triangle_points(center, radius, add_center = False):
 
   Reference triangle:
               
-      + v2
-      | \
-      |      \
-      |            \  
-      |    o           + v1
-      |            /
-      |      /
-      | /
-      + v3
+  #    + v2
+  #    | \
+  #    |      \
+  #    |            \  
+  #    |    o           + v1
+  #    |            /
+  #    |      /
+  #    | /
+  #    + v3
   
   Radius: o-v1 
   Axis vector: o-v1
@@ -378,4 +378,39 @@ def gmsh_translate(xc):
         new_coord = old_coord + translation_vector
         # Handle paramCoord properly - use empty list if paramCoord is None or empty
         param = paramCoord[3*i:3*(i+1)] if paramCoord is not None and len(paramCoord) > 0 else []
-        gmsh.model.mesh.setNode(nodeTag, new_coord, param)
+        gmsh.model.mesh.setNode(nodeTag, new_coord, param)
+
+def check_hanging_nodes():
+    """Check for hanging/orphan nodes in the GMSH mesh that are not connected to any elements.
+    
+    Returns:
+    -------
+    list
+        List of node tags that are hanging (not connected to any elements)
+    bool
+        True if hanging nodes were found, False otherwise
+    """
+    # Get all nodes in the mesh
+    all_node_tags, _, _ = gmsh.model.mesh.getNodes()
+    all_node_tags = set(all_node_tags)
+    
+    # Get all nodes that are part of elements
+    connected_nodes = set()
+    
+    # Iterate through all entities in the model
+    for entity in gmsh.model.getEntities():
+        # Get elements for this entity
+        elem_types, elem_tags, elem_node_tags = gmsh.model.mesh.getElements(entity[0], entity[1])
+        
+        # Add all nodes from elements to connected set
+        for nodes in elem_node_tags:
+            connected_nodes.update(nodes)
+    
+    # Find hanging nodes (nodes that exist but aren't in any elements)
+    hanging_nodes = list(all_node_tags - connected_nodes)
+    
+    has_hanging = len(hanging_nodes) > 0
+    if has_hanging:
+        print(f"Warning: Found {len(hanging_nodes)} hanging nodes with tags: {hanging_nodes}")
+    
+    return hanging_nodes, has_hanging

tools/python_utils/gmsh_particles.py

@@ -8,6 +8,7 @@
 from geom_util import *
 
 
+
 def circle_mesh_symmetric(xc = [0., 0., 0.], r = 1., h = 0.1, filename = 'mesh', vtk_out = False, symmetric_mesh = True):
   """
   Create a circular mesh, either symmetric or full.
@@ -78,6 +79,9 @@ def circle_mesh_symmetric(xc = [0., 0., 0.], r = 1., h = 0.1, filename = 'mesh',
     # generate mesh
     gmsh.model.mesh.generate(3)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+
   # write
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -158,6 +162,9 @@ def ellipse_mesh_symmetric(xc = [0., 0., 0.], rx = 1., ry = 0.5, h = 0.1, filena
     # generate mesh
     gmsh.model.mesh.generate(3)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+
   # Write output files
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -243,6 +250,9 @@ def sphere_mesh_symmetric(xc = [0., 0., 0.], r = 1., h = 0.1, filename = 'mesh',
     # Generate 3D mesh
     gmsh.model.mesh.generate(3)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+
   # Write output files
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -328,6 +338,9 @@ def rectangle_mesh_symmetric(xc = [0., 0., 0.], Lx = 1., Ly = 1., h = 0.1, filen
     # generate mesh
     gmsh.model.mesh.generate(3)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+
   # write
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -409,6 +422,9 @@ def hexagon_mesh_symmetric(xc = [0., 0., 0.], r = 1., h = 0.1, filename = 'mesh'
     # generate mesh
     gmsh.model.mesh.generate(3)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+
   # Write output files
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -497,6 +513,9 @@ def drum2d_mesh_symmetric(xc = [0., 0., 0.], r = 1., width = 1., h = 0.1, filena
     # generate mesh
     gmsh.model.mesh.generate(3)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+
   # Write output files
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -565,6 +584,9 @@ def polygon_mesh_symmetric(points, theta, xc=[0., 0., 0.], h=0.1, filename='mesh
 
   # Translate to final position if needed
   gmsh_translate(xc)
+
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
 
   # Write output files
   gmsh.write(filename + '.msh')
@@ -636,10 +658,18 @@ def cylindrical2d_wall_mesh(center=[0., 0., 0.], outer_radius=1.0, inner_radius=
   # Create surface with hole
   s1 = gmsh.model.geo.addPlaneSurface([outer_loop, inner_loop])
 
-  # Synchronize and generate mesh
+  # Synchronize geometry before adding physical groups
   gmsh.model.geo.synchronize()
+  
+  # Add physical group for the surface
+  gmsh.model.addPhysicalGroup(2, [s1], 1)
+  
+  # Generate mesh
   gmsh.model.mesh.generate(2)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+  
   # Write output files
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -729,6 +759,9 @@ def triangle_mesh_symmetric(xc=[0., 0., 0.], r=1., h=0.1, filename='mesh', vtk_o
     # Generate mesh
     gmsh.model.mesh.generate(3)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+  
   # Write output files
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -807,6 +840,9 @@ def cuboid_mesh_symmetric(xc = [0., 0., 0.], Lx = 1., Ly = 1., Lz = 1., h = 0.1,
     # Generate 3D mesh
     gmsh.model.mesh.generate(3)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+
   # Write output files
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -891,6 +927,9 @@ def ellipsoid_mesh_symmetric(xc = [0., 0., 0.], rx = 1., ry = 0.5, rz = 0.3, h =
     # Generate 3D mesh
     gmsh.model.mesh.generate(3)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+  
   # Write output files
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -979,6 +1018,9 @@ def cylinder_mesh_symmetric(xc = [0., 0., 0.], r = 1., h = 1., mesh_size = 0.1,
     # Generate 3D mesh
     gmsh.model.mesh.generate(3)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+  
   # Write output files
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -1071,12 +1113,17 @@ def annulus_circle_mesh_symmetric(xc = [0., 0., 0.], r_outer = 1., r_inner = 0.5
     # surface
     s1 = gmsh.model.geo.addPlaneSurface([outer_loop, inner_loop])
 
-    # designate surface as physical surface
-    gmsh.model.addPhysicalGroup(2, [s1], 1)
-    
-    # Synchronize and generate mesh
+    # Synchronize geometry before adding physical groups
     gmsh.model.geo.synchronize()
+
+    # Add physical group for the surface
+    gmsh.model.addPhysicalGroup(2, [s1], 1)
+
+    # Generate mesh
     gmsh.model.mesh.generate(2)
+
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
 
   # Write output files
   gmsh.write(filename + '.msh')
@@ -1194,6 +1241,9 @@ def annulus_rectangle_mesh(xc=[0., 0., 0.], Lx=1., Ly=1., hole_Lx=0.3, hole_Ly=0
     gmsh.model.geo.synchronize()
     gmsh.model.mesh.generate(2)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+
   # Write mesh files
   gmsh.write(filename + '.msh')
   if vtk_out:
@@ -1257,6 +1307,9 @@ def open_rectangle_mesh(xc = [0., 0., 0.], Lx = 1., Ly = 1., hole_Lx = 0.5, hole
   gmsh.model.mesh.removeDuplicateNodes()
   gmsh_translate(xc)
 
+  # Check for hanging nodes before writing
+  check_hanging_nodes()
+
   gmsh.write(filename + '.msh')
   if vtk_out:
     gmsh.write(filename + '.vtk')
@@ -1340,6 +1393,9 @@ def open_pipe_mesh(xc=[0., 0., 0.], axis=[0., 0., 1.], length=2., outer_radius=1
 
     # Generate 3D mesh
     gmsh.model.mesh.generate(3)
+
+    # Check for hanging nodes before writing
+    check_hanging_nodes()
 
     # Write mesh files
     gmsh.write(filename + '.msh')
@@ -1353,7 +1409,9 @@ def open_pipe_mesh(xc=[0., 0., 0.], axis=[0., 0., 1.], length=2., outer_radius=1
 if __name__ == "__main__":
   inp_dir = './'
 
-  test_meshes = ['circle', 'ellipse', 'sphere', 'cuboid', 'ellipsoid', 'rectangle', 'hexagon', 'drum2d', 'triangle', 'polygon', 'cylindrical2d_wall', 'cylinder', 'annulus_circle', 'annulus_rectangle', 'open_rectangle', 'open_pipe']
+  test_meshes = ['circle', 'ellipse', 'sphere', 'cuboid', 'ellipsoid', 'rectangle', \
+                 'hexagon', 'drum2d', 'triangle', 'polygon', 'cylindrical2d_wall', \
+                 'cylinder', 'annulus_circle', 'annulus_rectangle', 'open_rectangle', 'open_pipe']
   test_meshes = ['open_pipe']
   symm_flag = 0 # 0 - both, 1 - symmetric, 2 - non-symmetric