diff --git a/geos-processing/src/geos/processing/generic_processing_tools/SplitMesh.py b/geos-processing/src/geos/processing/generic_processing_tools/SplitMesh.py index f4a7a6bd..f02ba010 100644 --- a/geos-processing/src/geos/processing/generic_processing_tools/SplitMesh.py +++ b/geos-processing/src/geos/processing/generic_processing_tools/SplitMesh.py @@ -27,14 +27,14 @@ VTK_POLYHEDRON, VTK_POLYGON, ) - from vtkmodules.util.numpy_support import ( numpy_to_vtk, vtk_to_numpy ) - from geos.processing.pre_processing.CellTypeCounterEnhanced import CellTypeCounterEnhanced from geos.mesh.model.CellTypeCounts import CellTypeCounts __doc__ = """ -SplitMesh module is a vtk filter that split cells of a mesh composed of Tetrahedra, pyramids, and hexahedra. +SplitMesh module is a vtk filter that splits cells of a mesh composed of tetrahedra, pyramids, hexahedra, triangles, and quads. + +Warning: Current implementation only supports meshes composed of either polygons or polyhedra, not both together. Filter input and output types are vtkUnstructuredGrid. @@ -133,21 +133,24 @@ def RequestData( assert self.inData is not None, "Input mesh is undefined." assert output is not None, "Output mesh is undefined." + # Count the number of cells before splitting. Then we will be able to know how many new cells and points + # to allocate because each cell type is splitted in a known number of new cells and points. nbCells: int = self.inData.GetNumberOfCells() counts: CellTypeCounts = self._getCellCounts() - nbTet: int = counts.getTypeCount( VTK_TETRA ) - nbPyr: int = counts.getTypeCount( VTK_PYRAMID ) - nbHex: int = counts.getTypeCount( VTK_HEXAHEDRON ) - nbTriangles: int = counts.getTypeCount( VTK_TRIANGLE ) - nbQuad: int = counts.getTypeCount( VTK_QUAD ) + nbTet: int = counts.getTypeCount( VTK_TETRA ) # will divide into 8 tets + nbPyr: int = counts.getTypeCount( VTK_PYRAMID ) # will divide into 6 pyramids and 4 tets so 10 new cells + nbHex: int = counts.getTypeCount( VTK_HEXAHEDRON ) # will divide into 8 hexes + nbTriangles: int = counts.getTypeCount( VTK_TRIANGLE ) # will divide into 4 triangles + nbQuad: int = counts.getTypeCount( VTK_QUAD ) # will divide into 4 quads nbPolygon: int = counts.getTypeCount( VTK_POLYGON ) nbPolyhedra: int = counts.getTypeCount( VTK_POLYHEDRON ) assert counts.getTypeCount( VTK_WEDGE ) == 0, "Input mesh contains wedges that are not currently supported." + # Current implementation only supports meshes composed of either polygons or polyhedra assert nbPolyhedra * nbPolygon == 0, ( "Input mesh is composed of both polygons and polyhedra," " but it must contains only one of the two." ) nbNewPoints: int = 0 nbNewPoints = nbHex * 19 + nbTet * 6 + nbPyr * 9 if nbPolyhedra > 0 else nbTriangles * 3 + nbQuad * 5 - nbNewCells: int = nbHex * 8 + nbTet * 8 + nbPyr * 10 * nbTriangles * 4 + nbQuad * 4 + nbNewCells: int = nbHex * 8 + nbTet * 8 + nbPyr * 10 + nbTriangles * 4 + nbQuad * 4 self.points = vtkPoints() self.points.DeepCopy( self.inData.GetPoints() ) @@ -159,19 +162,20 @@ def RequestData( self.originalId.SetName( "OriginalID" ) self.originalId.Allocate( nbNewCells ) self.cellTypes = [] + # Define cell type to splitting method mapping + splitMethods = { + VTK_HEXAHEDRON: self._splitHexahedron, + VTK_TETRA: self._splitTetrahedron, + VTK_PYRAMID: self._splitPyramid, + VTK_TRIANGLE: self._splitTriangle, + VTK_QUAD: self._splitQuad, + } for c in range( nbCells ): cell: vtkCell = self.inData.GetCell( c ) cellType: int = cell.GetCellType() - if cellType == VTK_HEXAHEDRON: - self._splitHexahedron( cell, c ) - elif cellType == VTK_TETRA: - self._splitTetrahedron( cell, c ) - elif cellType == VTK_PYRAMID: - self._splitPyramid( cell, c ) - elif cellType == VTK_TRIANGLE: - self._splitTriangle( cell, c ) - elif cellType == VTK_QUAD: - self._splitQuad( cell, c ) + splitMethod = splitMethods.get( cellType ) + if splitMethod is not None: + splitMethod( cell, c ) else: raise TypeError( f"Cell type {vtkCellTypes.GetClassNameFromTypeId(cellType)} is not supported." ) # add points and cells @@ -261,7 +265,7 @@ def _splitPyramid( self: Self, cell: vtkCell, index: int ) -> None: r"""Split a pyramid. Let's suppose an input pyramid composed of nodes (0, 1, 2, 3, 4), - the cell is splitted in 8 pyramids using edge centers. + the cell is split into 6 pyramids and 4 tetrahedra using edge centers. 4 ,/|\ @@ -310,7 +314,8 @@ def _splitPyramid( self: Self, cell: vtkCell, index: int ) -> None: self.cells.InsertNextCell( 4, [ pt12, pt7, pt6, pt13 ] ) for _ in range( 10 ): self.originalId.InsertNextValue( index ) - self.cellTypes.extend( [ VTK_PYRAMID ] * 8 ) + self.cellTypes.extend( [ VTK_PYRAMID ] * 6 ) + self.cellTypes.extend( [ VTK_TETRA ] * 4 ) def _splitHexahedron( self: Self, cell: vtkCell, index: int ) -> None: r"""Split a hexahedron. @@ -378,7 +383,7 @@ def _splitTriangle( self: Self, cell: vtkCell, index: int ) -> None: r"""Split a triangle. Let's suppose an input triangle composed of nodes (0, 1, 2), - the cell is splitted in 3 triangles using edge centers. + the cell is split into 4 triangles using edge centers. 2 |\ diff --git a/geos-processing/tests/data/hexs3_tets36_pyrs18.vtu b/geos-processing/tests/data/hexs3_tets36_pyrs18.vtu new file mode 100644 index 00000000..3ee9b963 --- /dev/null +++ b/geos-processing/tests/data/hexs3_tets36_pyrs18.vtu @@ -0,0 +1,48 @@ + + + + + + mAAAAAAAAADXs10/17NdP9ezXT/Xs10/17NdP9ezXT/Xs10/17NdP9ezXT/Xs10/17NdP9ezXT/Xs10/17NdP9ezXT/Xs10/17NdP9ezXT/Xs10/17NdP9ezXT/Xs10/AAAAANezXT/Xs10/AAAAANezXT/Xs10/17NdP9ezXT/Xs10/17NdPwAAAADXs10/17NdPwAAAAAAAAAAAAAAAA== + + + mAAAAAAAAAAAAIA/AAAAQAAAAEAAAEBAAAAAQAAAQEAAAEBAAACAQAAAAEAAAABAAABAQAAAQEAAAEBAAACAQAAAAEAAAABAAABAQAAAQEAAAEBAAACAQAAAAEAAAEBAAACAPwAAAEAAAEBAAACAPwAAgD8AAABAAABAQAAAAEAAAEBAAACAQAAAgD8AAABAAABAQAAAgD8AAIA/AACAPw== + + + + + + + yAEAAAAAAAAAAAAAAAAAAAAAAAAAAIA/AAAAAAAAAAAAAAAAAACAPwAAAAAAAIA/AACAPwAAAAAAAAAAAAAAAAAAgD8AAIA/AAAAAAAAgD8AAAAAAACAPwAAgD8AAIA/AACAPwAAgD8AAABAAACAPwAAAAAAAIA/AAAAQAAAAAAAAABAAAAAQAAAAAAAAABAAACAPwAAgD8AAIA/AAAAQAAAgD8AAABAAAAAQAAAgD8AAEBAAAAAQAAAAAAAAABAAABAQAAAAAAAAEBAAABAQAAAAAAAAEBAAAAAQAAAgD8AAABAAABAQAAAgD8AAEBAAABAQAAAgD8AAABAAAAAAAAAAAAAAABAAAAAAAAAgD8AAMA/AAAAPwAAAD8AAEBAAACAPwAAAAAAAEBAAACAPwAAgD8AACBAAADAPwAAAD8AAAAAAAAAQAAAAAAAAAAAAABAQAAAAAAAAIA/AABAQAAAAAAAAAAAAAAAQAAAgD8AAAAAAABAQAAAgD8AAIA/AABAQAAAgD8AAAA/AAAgQAAAAD8AAEBAAAAAAAAAAAAAAEBAAAAAAAAAgD8AACBAAAAAPwAAAD8AAAA/AADAPwAAAD8AAMA/AAAgQAAAAD8= + + + 0 + + + 4.3588989435 + + + + + 0 + + + 4.3588989435 + + + + + + + EAgAAAAAAAAAAAAAAAAAAAEAAAAAAAAAAwAAAAAAAAACAAAAAAAAAAQAAAAAAAAABQAAAAAAAAAHAAAAAAAAAAYAAAAAAAAAAwAAAAAAAAAIAAAAAAAAAAoAAAAAAAAACQAAAAAAAAAHAAAAAAAAAAsAAAAAAAAADQAAAAAAAAAMAAAAAAAAAAoAAAAAAAAADgAAAAAAAAAQAAAAAAAAAA8AAAAAAAAADQAAAAAAAAARAAAAAAAAABMAAAAAAAAAEgAAAAAAAAABAAAAAAAAABQAAAAAAAAACAAAAAAAAAAWAAAAAAAAAAEAAAAAAAAACAAAAAAAAAADAAAAAAAAABYAAAAAAAAACwAAAAAAAAAVAAAAAAAAAAUAAAAAAAAAFgAAAAAAAAALAAAAAAAAAAUAAAAAAAAABwAAAAAAAAAWAAAAAAAAABQAAAAAAAAAFQAAAAAAAAALAAAAAAAAABYAAAAAAAAAFAAAAAAAAAALAAAAAAAAAAgAAAAAAAAAFgAAAAAAAAABAAAAAAAAAAUAAAAAAAAAFQAAAAAAAAAWAAAAAAAAAAEAAAAAAAAAFQAAAAAAAAAUAAAAAAAAABYAAAAAAAAAAQAAAAAAAAADAAAAAAAAAAcAAAAAAAAAFgAAAAAAAAABAAAAAAAAAAcAAAAAAAAABQAAAAAAAAAWAAAAAAAAAAsAAAAAAAAABwAAAAAAAAADAAAAAAAAABYAAAAAAAAACwAAAAAAAAADAAAAAAAAAAgAAAAAAAAAFgAAAAAAAAAIAAAAAAAAABcAAAAAAAAADgAAAAAAAAAZAAAAAAAAAAgAAAAAAAAADgAAAAAAAAAKAAAAAAAAABkAAAAAAAAAEQAAAAAAAAAYAAAAAAAAAAsAAAAAAAAAGQAAAAAAAAARAAAAAAAAAAsAAAAAAAAADQAAAAAAAAAZAAAAAAAAABcAAAAAAAAAGAAAAAAAAAARAAAAAAAAABkAAAAAAAAAFwAAAAAAAAARAAAAAAAAAA4AAAAAAAAAGQAAAAAAAAAIAAAAAAAAAAsAAAAAAAAAGAAAAAAAAAAZAAAAAAAAAAgAAAAAAAAAGAAAAAAAAAAXAAAAAAAAABkAAAAAAAAACAAAAAAAAAAKAAAAAAAAAA0AAAAAAAAAGQAAAAAAAAAIAAAAAAAAAA0AAAAAAAAACwAAAAAAAAAZAAAAAAAAABEAAAAAAAAADQAAAAAAAAAKAAAAAAAAABkAAAAAAAAAEQAAAAAAAAAKAAAAAAAAAA4AAAAAAAAAGQAAAAAAAAAaAAAAAAAAAAkAAAAAAAAAHAAAAAAAAAAgAAAAAAAAABoAAAAAAAAAHAAAAAAAAAAbAAAAAAAAACAAAAAAAAAAHwAAAAAAAAAMAAAAAAAAAB0AAAAAAAAAIAAAAAAAAAAfAAAAAAAAAB0AAAAAAAAAHgAAAAAAAAAgAAAAAAAAAAkAAAAAAAAADAAAAAAAAAAfAAAAAAAAACAAAAAAAAAACQAAAAAAAAAfAAAAAAAAABwAAAAAAAAAIAAAAAAAAAAaAAAAAAAAAB0AAAAAAAAADAAAAAAAAAAgAAAAAAAAABoAAAAAAAAADAAAAAAAAAAJAAAAAAAAACAAAAAAAAAAGgAAAAAAAAAbAAAAAAAAAB4AAAAAAAAAIAAAAAAAAAAaAAAAAAAAAB4AAAAAAAAAHQAAAAAAAAAgAAAAAAAAAB8AAAAAAAAAHgAAAAAAAAAbAAAAAAAAACAAAAAAAAAAHwAAAAAAAAAbAAAAAAAAABwAAAAAAAAAIAAAAAAAAAAUAAAAAAAAACEAAAAAAAAAFwAAAAAAAAAIAAAAAAAAACMAAAAAAAAAGAAAAAAAAAAiAAAAAAAAABUAAAAAAAAACwAAAAAAAAAjAAAAAAAAACEAAAAAAAAAIgAAAAAAAAAYAAAAAAAAABcAAAAAAAAAIwAAAAAAAAAUAAAAAAAAABUAAAAAAAAAIgAAAAAAAAAhAAAAAAAAACMAAAAAAAAAFAAAAAAAAAAIAAAAAAAAAAsAAAAAAAAAFQAAAAAAAAAjAAAAAAAAABgAAAAAAAAACwAAAAAAAAAIAAAAAAAAABcAAAAAAAAAIwAAAAAAAAACAAAAAAAAAAMAAAAAAAAACQAAAAAAAAAaAAAAAAAAACQAAAAAAAAADAAAAAAAAAAHAAAAAAAAAAYAAAAAAAAAHQAAAAAAAAAkAAAAAAAAAAMAAAAAAAAABwAAAAAAAAAMAAAAAAAAAAkAAAAAAAAAJAAAAAAAAAACAAAAAAAAAAYAAAAAAAAABwAAAAAAAAADAAAAAAAAACQAAAAAAAAAAgAAAAAAAAAaAAAAAAAAAB0AAAAAAAAABgAAAAAAAAAkAAAAAAAAAAwAAAAAAAAAHQAAAAAAAAAaAAAAAAAAAAkAAAAAAAAAJAAAAAAAAAAJAAAAAAAAAAoAAAAAAAAADwAAAAAAAAAcAAAAAAAAACUAAAAAAAAAEgAAAAAAAAANAAAAAAAAAAwAAAAAAAAAHwAAAAAAAAAlAAAAAAAAAAoAAAAAAAAADQAAAAAAAAASAAAAAAAAAA8AAAAAAAAAJQAAAAAAAAAJAAAAAAAAAAwAAAAAAAAADQAAAAAAAAAKAAAAAAAAACUAAAAAAAAACQAAAAAAAAAcAAAAAAAAAB8AAAAAAAAADAAAAAAAAAAlAAAAAAAAABIAAAAAAAAAHwAAAAAAAAAcAAAAAAAAAA8AAAAAAAAAJQAAAAAAAAA= + + + yAEAAAAAAAAIAAAAAAAAABAAAAAAAAAAGAAAAAAAAAAcAAAAAAAAACAAAAAAAAAAJAAAAAAAAAAoAAAAAAAAACwAAAAAAAAAMAAAAAAAAAA0AAAAAAAAADgAAAAAAAAAPAAAAAAAAABAAAAAAAAAAEQAAAAAAAAASAAAAAAAAABMAAAAAAAAAFAAAAAAAAAAVAAAAAAAAABYAAAAAAAAAFwAAAAAAAAAYAAAAAAAAABkAAAAAAAAAGgAAAAAAAAAbAAAAAAAAABwAAAAAAAAAHQAAAAAAAAAeAAAAAAAAAB8AAAAAAAAAIAAAAAAAAAAhAAAAAAAAACIAAAAAAAAAIwAAAAAAAAAkAAAAAAAAACUAAAAAAAAAJgAAAAAAAAAnAAAAAAAAACgAAAAAAAAAKQAAAAAAAAAqAAAAAAAAACtAAAAAAAAALIAAAAAAAAAtwAAAAAAAAC8AAAAAAAAAMEAAAAAAAAAxgAAAAAAAADLAAAAAAAAANAAAAAAAAAA1QAAAAAAAADaAAAAAAAAAN8AAAAAAAAA5AAAAAAAAADpAAAAAAAAAO4AAAAAAAAA8wAAAAAAAAD4AAAAAAAAAP0AAAAAAAAAAgEAAAAAAAA= + + + OQAAAAAAAAAMDAwKCgoKCgoKCgoKCgoKCgoKCgoKCgoKCgoKCgoKCgoKCgoKCgoODg4ODg4ODg4ODg4ODg4ODg4= + + + + + diff --git a/geos-processing/tests/data/quads2_tris4.vtu b/geos-processing/tests/data/quads2_tris4.vtu new file mode 100644 index 00000000..f95eb1f5 --- /dev/null +++ b/geos-processing/tests/data/quads2_tris4.vtu @@ -0,0 +1,48 @@ + + + + + + JAAAAAAAAADzBDU/8wQ1P/MENT/zBDU/8wQ1P/MENT/zBDU/8wQ1P/MENT8= + + + JAAAAAAAAAAAAIA/AAAAQAAAAEAAAEBAAAAAQAAAAEAAAEBAAAAAQAAAAEA= + + + + + + + bAAAAAAAAAAAAAAAAAAAAAAAAAAAAIA/AAAAAAAAAAAAAAAAAACAPwAAAAAAAIA/AACAPwAAAAAAAABAAACAPwAAAAAAAIA/AAAAQAAAAAAAAABAAAAAQAAAAAAAAABAAAAAAAAAAAAAAAAAAAAAQAAAAAA= + + + 0 + + + 2.8284271247 + + + + + 0 + + + 2.8284271247 + + + + + + + oAAAAAAAAAAAAAAAAAAAAAEAAAAAAAAAAwAAAAAAAAACAAAAAAAAAAMAAAAAAAAABAAAAAAAAAAGAAAAAAAAAAUAAAAAAAAAAQAAAAAAAAAHAAAAAAAAAAMAAAAAAAAABwAAAAAAAAAEAAAAAAAAAAMAAAAAAAAAAgAAAAAAAAADAAAAAAAAAAgAAAAAAAAAAwAAAAAAAAAFAAAAAAAAAAgAAAAAAAAA + + + MAAAAAAAAAAEAAAAAAAAAAgAAAAAAAAACwAAAAAAAAAOAAAAAAAAABEAAAAAAAAAFAAAAAAAAAA= + + + BgAAAAAAAAAJCQUFBQU= + + + + + diff --git a/geos-processing/tests/test_SplitMesh.py b/geos-processing/tests/test_SplitMesh.py index 0612b729..41a0bb7f 100644 --- a/geos-processing/tests/test_SplitMesh.py +++ b/geos-processing/tests/test_SplitMesh.py @@ -6,22 +6,14 @@ import numpy as np import numpy.typing as npt import pytest -from typing import ( - Iterator, ) - -from geos.mesh.utils.genericHelpers import createSingleCellMesh -from geos.processing.generic_processing_tools.SplitMesh import SplitMesh - +from typing import Iterator from vtkmodules.util.numpy_support import vtk_to_numpy - from vtkmodules.vtkCommonDataModel import ( vtkUnstructuredGrid, vtkCellArray, vtkCellData, vtkCellTypes, VTK_TRIANGLE, VTK_QUAD, VTK_TETRA, VTK_HEXAHEDRON, VTK_PYRAMID ) - -from vtkmodules.vtkCommonCore import ( - vtkPoints, - vtkIdList, - vtkDataArray, -) +from vtkmodules.vtkCommonCore import vtkPoints, vtkIdList, vtkDataArray +from geos.mesh.io.vtkIO import readUnstructuredGrid +from geos.mesh.utils.genericHelpers import createSingleCellMesh +from geos.processing.generic_processing_tools.SplitMesh import SplitMesh data_root: str = os.path.join( os.path.dirname( os.path.abspath( __file__ ) ), "data" ) @@ -98,7 +90,8 @@ ############################################################### # create multi cell mesh # ############################################################### -# TODO: add tests cases composed of multi-cell meshes of various types +multi_polygon_types_mesh_path = "quads2_tris4.vtu" +multi_polyhedron_types_mesh_path = "hexs3_tets36_pyrs18.vtu" data_filename_all = ( tetra_path, hexa_path, pyramid_path, triangle_path, quad_path ) cell_types_all = ( tetra_cell_type, hexa_cell_type, pyramid_cell_type, triangle_cell_type, quad_cell_type ) @@ -176,18 +169,34 @@ def test_single_cell_split( test_case: TestCase ) -> None: typesArray: npt.NDArray[ np.int64 ] = vtk_to_numpy( types.GetCellTypesArray() ) print( "typesArray", cellTypeName, typesArray ) - assert ( typesArray.size == 1 ) and ( typesArray[ 0 ] == test_case.cellType ), f"All cells must be {cellTypeName}" + # Pyramid splitting produces both pyramids (first 6 cells) and tetrahedra (last 4 cells) + if test_case.cellType == VTK_PYRAMID: + assert typesArray.size == 2, "Pyramid splitting should produce 2 distinct cell types" + assert VTK_PYRAMID in typesArray, "Pyramid splitting should produce pyramids" + assert VTK_TETRA in typesArray, "Pyramid splitting should produce tetrahedra" + else: + assert ( typesArray.size == 1 ) and ( typesArray[ 0 ] == test_case.cellType ), \ + f"All cells must be {cellTypeName}" for i in range( cellsOut.GetNumberOfCells() ): ptIds = vtkIdList() cellsOut.GetCellAtId( i, ptIds ) cellsOutObs: list[ int ] = [ ptIds.GetId( j ) for j in range( ptIds.GetNumberOfIds() ) ] nbPtsExp: int = len( test_case.cellsExp[ i ] ) - print( "cell type", cellTypeName, i, vtkCellTypes.GetClassNameFromTypeId( types.GetCellType( i ) ) ) + actualCellType: int = output.GetCellType( i ) + print( "cell type", cellTypeName, i, vtkCellTypes.GetClassNameFromTypeId( actualCellType ) ) print( "cellsOutObs: ", cellTypeName, i, cellsOutObs ) assert ptIds is not None, "Point ids must be defined" assert ptIds.GetNumberOfIds() == nbPtsExp, f"Cells must be defined by {nbPtsExp} points." assert cellsOutObs == test_case.cellsExp[ i ], "Cell point ids are wrong." + # For pyramids, first 6 cells should be pyramids, last 4 should be tetrahedra + if test_case.cellType == VTK_PYRAMID: + if i < 6: + assert actualCellType == VTK_PYRAMID, f"Cell {i} should be a pyramid" + else: + assert actualCellType == VTK_TETRA, f"Cell {i} should be a tetrahedron" + else: + assert actualCellType == test_case.cellType, f"Cell {i} should be {cellTypeName}" # test originalId array was created cellData: vtkCellData = output.GetCellData() @@ -202,4 +211,136 @@ def test_single_cell_split( test_case: TestCase ) -> None: nbArraySplited: int = cellData.GetNumberOfArrays() assert nbArraySplited == nbArrayInput + 1, f"Number of arrays should be {nbArrayInput + 1}" - #assert False + +def test_multi_cells_mesh_split() -> None: + """Test splitting a mesh with multiple cell types (3 hexahedra, 36 tetrahedra, 18 pyramids). + + This test verifies that the total number of cells generated matches the expected formula: + nbNewCells = nbHex * 8 + nbTet * 8 + nbPyr * 10 + """ + # Load the multi-cell mesh + input_mesh = readUnstructuredGrid( os.path.join( data_root, multi_polyhedron_types_mesh_path ) ) + assert input_mesh is not None, "Input mesh should be loaded successfully" + + # Count cells by type in input mesh + nbHex = 0 + nbTet = 0 + nbPyr = 0 + for i in range( input_mesh.GetNumberOfCells() ): + cellType = input_mesh.GetCellType( i ) + if cellType == VTK_HEXAHEDRON: + nbHex += 1 + elif cellType == VTK_TETRA: + nbTet += 1 + elif cellType == VTK_PYRAMID: + nbPyr += 1 + + print( f"Input mesh contains: {nbHex} hexahedra, {nbTet} tetrahedra, {nbPyr} pyramids" ) + assert nbHex == 3, "Expected 3 hexahedra in input mesh" + assert nbTet == 36, "Expected 36 tetrahedra in input mesh" + assert nbPyr == 18, "Expected 18 pyramids in input mesh" + + # Apply the split filter + splitMeshFilter = SplitMesh() + splitMeshFilter.SetInputDataObject( input_mesh ) + splitMeshFilter.Update() + output: vtkUnstructuredGrid = splitMeshFilter.GetOutputDataObject( 0 ) + assert output is not None, "Output mesh should be defined" + + # Calculate expected number of cells using the formula + # 1 hex -> 8 hexes, 1 tet -> 8 tets, 1 pyramid -> 6 pyramids + 4 tets = 10 cells + expectedNbCells = nbHex * 8 + nbTet * 8 + nbPyr * 10 + print( f"Expected number of cells: {expectedNbCells} (3*8 + 36*8 + 18*10)" ) + print( f"Actual number of cells: {output.GetNumberOfCells()}" ) + assert output.GetNumberOfCells() == expectedNbCells, \ + f"Expected {expectedNbCells} cells, got {output.GetNumberOfCells()}" + + # Verify cell type distribution in output + nbHexOut = 0 + nbTetOut = 0 + nbPyrOut = 0 + for i in range( output.GetNumberOfCells() ): + cellType = output.GetCellType( i ) + if cellType == VTK_HEXAHEDRON: + nbHexOut += 1 + elif cellType == VTK_TETRA: + nbTetOut += 1 + elif cellType == VTK_PYRAMID: + nbPyrOut += 1 + + print( f"Output mesh contains: {nbHexOut} hexahedra, {nbTetOut} tetrahedra, {nbPyrOut} pyramids" ) + # Expected output: 3*8=24 hexes, 36*8 + 18*4=360 tets, 18*6=108 pyramids + assert nbHexOut == 3 * 8, f"Expected {3*8} hexahedra in output, got {nbHexOut}" + assert nbTetOut == 36 * 8 + 18 * 4, f"Expected {36*8 + 18*4} tetrahedra in output, got {nbTetOut}" + assert nbPyrOut == 18 * 6, f"Expected {18*6} pyramids in output, got {nbPyrOut}" + + # Verify OriginalID array exists and has correct size + cellData: vtkCellData = output.GetCellData() + assert cellData is not None, "Cell data should be defined" + originalIdArray: vtkDataArray = cellData.GetArray( "OriginalID" ) + assert originalIdArray is not None, "OriginalID array should be defined" + assert originalIdArray.GetNumberOfTuples() == expectedNbCells, \ + f"OriginalID array should have {expectedNbCells} values" + + +def test_multi_polygon_mesh_split() -> None: + """Test splitting a mesh with multiple polygon types (2 quads, 4 triangles). + + This test verifies that the total number of cells generated matches the expected formula: + nbNewCells = nbQuad * 4 + nbTriangle * 4 + """ + # Load the multi-polygon mesh + input_mesh = readUnstructuredGrid( os.path.join( data_root, multi_polygon_types_mesh_path ) ) + assert input_mesh is not None, "Input mesh should be loaded successfully" + + # Count cells by type in input mesh + nbQuad = 0 + nbTriangle = 0 + for i in range( input_mesh.GetNumberOfCells() ): + cellType = input_mesh.GetCellType( i ) + if cellType == VTK_QUAD: + nbQuad += 1 + elif cellType == VTK_TRIANGLE: + nbTriangle += 1 + + print( f"Input mesh contains: {nbQuad} quads, {nbTriangle} triangles" ) + assert nbQuad == 2, "Expected 2 quads in input mesh" + assert nbTriangle == 4, "Expected 4 triangles in input mesh" + + # Apply the split filter + splitMeshFilter = SplitMesh() + splitMeshFilter.SetInputDataObject( input_mesh ) + splitMeshFilter.Update() + output: vtkUnstructuredGrid = splitMeshFilter.GetOutputDataObject( 0 ) + assert output is not None, "Output mesh should be defined" + + # Calculate expected number of cells using the formula + # 1 quad -> 4 quads, 1 triangle -> 4 triangles + expectedNbCells = nbQuad * 4 + nbTriangle * 4 + print( f"Expected number of cells: {expectedNbCells} (2*4 + 4*4)" ) + print( f"Actual number of cells: {output.GetNumberOfCells()}" ) + assert output.GetNumberOfCells() == expectedNbCells, \ + f"Expected {expectedNbCells} cells, got {output.GetNumberOfCells()}" + + # Verify cell type distribution in output + nbQuadOut = 0 + nbTriangleOut = 0 + for i in range( output.GetNumberOfCells() ): + cellType = output.GetCellType( i ) + if cellType == VTK_QUAD: + nbQuadOut += 1 + elif cellType == VTK_TRIANGLE: + nbTriangleOut += 1 + + print( f"Output mesh contains: {nbQuadOut} quads, {nbTriangleOut} triangles" ) + # Expected output: 2*4=8 quads, 4*4=16 triangles + assert nbQuadOut == 2 * 4, f"Expected {2*4} quads in output, got {nbQuadOut}" + assert nbTriangleOut == 4 * 4, f"Expected {4*4} triangles in output, got {nbTriangleOut}" + + # Verify OriginalID array exists and has correct size + cellData: vtkCellData = output.GetCellData() + assert cellData is not None, "Cell data should be defined" + originalIdArray: vtkDataArray = cellData.GetArray( "OriginalID" ) + assert originalIdArray is not None, "OriginalID array should be defined" + assert originalIdArray.GetNumberOfTuples() == expectedNbCells, \ + f"OriginalID array should have {expectedNbCells} values"