test_physical working in pytest but not working when individually run

[connorivy]

When I run the pytest for test_physical, I get the following output

(venv) C:\Users\civy\Music\GitHub\pygmsh\tests\built_in>pytest test_physical.py ================================== test session starts ================================== platform win32 -- Python 3.9.6, pytest-7.0.1, pluggy-1.0.0 rootdir: C:\Users\civy\Music\GitHub\pygmsh\tests\built_in collected 1 item

test_physical.py . [100%]

================================== 1 passed in 0.54s ===================================

However, when I run the test with python instead of pytest I get this error

(venv) C:\Users\civy\Music\GitHub\pygmsh\tests\built_in>python test_physical.py Info: VTU format cannot write cell_sets. Converting them to cell_data... Traceback (most recent call last): File "C:\Users\civy\Music\GitHub\pygmsh\tests\built_in\test_physical.py", line 24, in test().write("physical.vtu") File "C:\Users\civy\Music\GitHub\pygmsh\venv\lib\site-packages\meshio_mesh.py", line 241, in write write(path_or_buf, self, file_format, **kwargs) File "C:\Users\civy\Music\GitHub\pygmsh\venv\lib\site-packages\meshio_helpers.py", line 187, in write return writer(filename, mesh, **kwargs) File "C:\Users\civy\Music\GitHub\pygmsh\venv\lib\site-packages\meshio\vtu_vtu.py", line 644, in write mesh.cell_sets_to_data(key) File "C:\Users\civy\Music\GitHub\pygmsh\venv\lib\site-packages\meshio_mesh.py", line 328, in cell_sets_to_data arr[cc] = i IndexError: index 116 is out of bounds for axis 0 with size 116

I'm also getting the same error in the program that I'm building with tries to mark a polygon with a label as demonstrated in the unit test. I'm not sure if this is a me problem or a pygmsh problem.

[Jul3k]

I can confirm this issue, there is a bug in the generation of the physical labels when both pyhsicals for faces and elements are created. For example in this code, the indices for the inner elements are incorrectly ofset by the length of the physicals for the faces. When the pyhsicals for the left and right faces are commented out, the correct indices for the inner elements are created:

lx = 20e-3
ly = 4e-3

lcar = 1e-3

with pygmsh.geo.Geometry() as geom:  
    inner = geom.add_rectangle(3/8*lx, 5/8*lx, 3/8*ly, 5/8*ly, 0, mesh_size=lcar/2)
    domain = geom.add_polygon([[0,0], [0,ly], [lx,ly], [lx,0]], mesh_size=lcar, holes=[inner])  
    
    geom.add_physical(inner, label='inner')
    geom.add_physical(domain, label='outer')
        
    geom.add_physical(domain.lines[0], label='left')
    geom.add_physical(domain.lines[2], label='right')

    mesh_gmsh = geom.generate_mesh()
    

mesh_gmsh.cell_sets:

   {'left': [array([0, 1, 2, 3], dtype=uint64), None, None],
 'right': [array([4, 5, 6, 7], dtype=uint64), None, None],
 'inner': [None,
  array([ 8,  9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24,  <---- Should start at 0
         25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41,
         42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53], dtype=uint64),
  None],
 'outer': [None,
  array([ 54,  55,  56,  57,  58,  59,  60,  61,  62,  63,  64,  65,  66,
          67,  68,  69,  70,  71,  72,  73,  74,  75,  76,  77,  78,  79,
          80,  81,  82,  83,  84,  85,  86,  87,  88,  89,  90,  91,  92,
          93,  94,  95,  96,  97,  98,  99, 100, 101, 102, 103, 104, 105,
         106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118,
         119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130, 131,
         132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144,
         145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156, 157,
         158, 159, 160, 161, 162, 163, 164, 165, 166, 167, 168, 169, 170,
         171, 172, 173, 174, 175, 176, 177, 178, 179, 180, 181, 182, 183,
         184, 185, 186, 187, 188, 189, 190, 191, 192, 193, 194, 195, 196,
         197, 198, 199, 200, 201, 202, 203, 204, 205, 206, 207, 208, 209,
         210, 211, 212, 213, 214, 215, 216, 217, 218, 219, 220, 221, 222,
         223, 224, 225, 226, 227, 228, 229, 230, 231, 232, 233, 234, 235,
         236, 237, 238, 239, 240, 241, 242, 243, 244, 245, 246, 247, 248,
         249, 250, 251, 252, 253, 254, 255, 256, 257, 258, 259, 260, 261,
         262, 263, 264, 265, 266, 267, 268, 269, 270, 271, 272, 273, 274,
         275, 276, 277, 278, 279, 280, 281, 282, 283, 284, 285, 286, 287,
         288, 289, 290, 291, 292, 293, 294, 295, 296, 297, 298, 299, 300,
         301, 302, 303, 304, 305, 306, 307, 308, 309, 310, 311, 312, 313,
         314, 315, 316, 317, 318, 319, 320, 321, 322, 323, 324, 325, 326,
         327, 328, 329, 330, 331, 332, 333], dtype=uint64),
  None]}

[Jul3k]

A workaround I discovered is to write the mesh to a file and then to re-read it into the program. So it seems the issue only exists when the mesh i directly used in the program:

lx = 20e-3
ly = 4e-3

lcar = 1e-3

with pygmsh.geo.Geometry() as geom:  
    inner = geom.add_rectangle(3/8*lx, 5/8*lx, 3/8*ly, 5/8*ly, 0, mesh_size=lcar/2)
    domain = geom.add_polygon([[0,0], [0,ly], [lx,ly], [lx,0]], mesh_size=lcar, holes=[inner])  
    
    geom.add_physical(inner, label='inner')
    geom.add_physical(domain, label='outer')
        
    geom.add_physical(domain.lines[0], label='left')
    geom.add_physical(domain.lines[2], label='right')

    mesh_gmsh = geom.generate_mesh()
      
    pygmsh.write("test.msh")
    
    mesh_gmsh = meshio.read('test.msh') # Workaround