Node ordering correction is automatic, user does not need to enter th…

…e node ordering himself. Just needs to specify the name of elements to check.

docs/geos-mesh.rst

@@ -117,41 +117,35 @@ Cells with negative volumes will typically be an issue for ``geos`` and should b
 """"""""""""""""""""""""""
 
 It sometimes happens that an exported mesh does not abide by the ``vtk`` orderings.
-The ``fix_elements_orderings`` module can rearrange the nodes of given types of elements.
+The ``fix_elements_orderings`` module can rearrange the nodes of given types of elements to respect VTK convention.
 This can be convenient if you cannot regenerate the mesh.
 
 .. code-block::
 
       $ python src/geos/mesh/doctor/mesh_doctor.py fix_elements_orderings --help
-      usage: mesh_doctor.py fix_elements_orderings [-h] [--Tetrahedron 2,0,3,1] [--Pyramid 3,4,0,2,1] [--Wedge 3,5,4,0,2,1]
-                                                   [--Hexahedron 1,6,5,4,7,0,2,3] [--Prism5 8,2,0,7,6,9,5,1,4,3] [--Prism6 11,2,8,10,5,0,9,7,6,1,4,3]
-                                                   --volume_to_reorder all,positive,negative --output OUTPUT [--data-mode binary, ascii]
+      usage: mesh_doctor.py fix_elements_orderings [-h] [--cell_names Hexahedron, Tetrahedron, Pyramid, Wedge, Prism5, Prism6]
+                                                   [--volume_to_reorder all, positive, negative]
+                                                   --output OUTPUT [--data-mode binary, ascii]
 
       options:
-      -h, --help                  show this help message and exit
-      --Tetrahedron 2,0,3,1       [list of integers]: node permutation for "Tetrahedron".
-      --Pyramid 3,4,0,2,1         [list of integers]: node permutation for "Pyramid".
-      --Wedge 3,5,4,0,2,1         [list of integers]: node permutation for "Wedge".
-      --Hexahedron 1,6,5,4,7,0,2,3
-                                  [list of integers]: node permutation for "Hexahedron".
-      --Prism5 8,2,0,7,6,9,5,1,4,3
-                                  [list of integers]: node permutation for "Prism5".
-      --Prism6 11,2,8,10,5,0,9,7,6,1,4,3
-                                  [list of integers]: node permutation for "Prism6".
-      --volume_to_reorder all,positive,negative
-                                  [str]: Select which element volume is invalid and needs reordering. 'all' will allow reordering of nodes for every element, regarding of their volume.
-                                  'positive' or 'negative' will only reorder the element with the corresponding volume.
-      --output OUTPUT             [string]: The vtk output file destination.
-      --data-mode binary, ascii   [string]: For ".vtu" output format, the data mode can be binary or ascii. Defaults to binary.
+      -h, --help              show this help message and exit
+      --cell_names Hexahedron, Tetrahedron, Pyramid, Wedge, Prism5, Prism6
+                              [list of str]: Cell names that can be reordered in your grid. You can use multiple names.Defaults to all cell names being used.
+      --volume_to_reorder all, positive, negative
+                              [str]: Select which element volume is invalid and needs reordering. 'all' will allow reordering of nodes for every element, regarding of their volume.
+                              'positive' or 'negative' will only reorder the element with the corresponding volume. Defaults to 'negative'.
+      --output OUTPUT         [string]: The vtk output file destination.
+      --data-mode binary, ascii
+                              [string]: For ".vtu" output format, the data mode can be binary or ascii. Defaults to binary.
 
 For example, assume that you have a mesh that contains 3 different cell types like Hexahedrons, Pyramids and Tetrahedrons.
 After checking ``element_volumes``, you found that all your Hexahedrons and half of your Tetrahedrons have a negative volume.
-To correct that, assuming you know the correct node ordering to correct these volumes, you can use this command:
+To correct that, you can use this command:
 
 .. code-block::
 
-      $ python src/geos/mesh/doctor/mesh_doctor.py -i /path/to/your/mesh/file fix_elements_orderings --Hexahedron 1,6,5,4,7,0,2,3
-      --Tetrahedron 2,0,3,1 --volume_to_reorder negative --output /new/path/for/your/new/mesh/reordered/file --data-mode binary
+      $ python src/geos/mesh/doctor/mesh_doctor.py -i /path/to/your/mesh/file fix_elements_orderings --cell_names Hexahedron,Tetrahedron
+      --volume_to_reorder negative --output /new/path/for/your/new/mesh/reordered/file --data-mode binary
 
 ``generate_cube``
 """""""""""""""""

geos-mesh/src/geos/mesh/doctor/checks/fix_elements_orderings.py

@@ -1,18 +1,18 @@
-from dataclasses import dataclass
 import numpy as np
 import logging
+from dataclasses import dataclass
+from itertools import permutations
 from vtk import vtkCellSizeFilter
-from vtkmodules.vtkCommonCore import vtkIdList
 from vtkmodules.util.numpy_support import vtk_to_numpy
-from vtkmodules.vtkCommonDataModel import ( vtkDataSet, VTK_HEXAHEDRON, VTK_TETRA, VTK_PYRAMID, VTK_WEDGE,
+from vtkmodules.vtkCommonDataModel import ( vtkDataSet, vtkCell, VTK_HEXAHEDRON, VTK_TETRA, VTK_PYRAMID, VTK_WEDGE,
                                             VTK_PENTAGONAL_PRISM, VTK_HEXAGONAL_PRISM )
-from .vtk_utils import VtkOutput, to_vtk_id_list, write_mesh, read_mesh
+from geos.mesh.doctor.checks.vtk_utils import VtkOutput, vtk_iter, to_vtk_id_list, write_mesh, read_mesh
 
 
 @dataclass( frozen=True )
 class Options:
     vtk_output: VtkOutput
-    cell_name_to_ordering: dict[ str, list[ int ] ]
+    cell_names_to_reorder: tuple[ str ]
     volume_to_reorder: str
 
 
@@ -103,6 +103,90 @@ def is_cell_to_reorder( cell_volume: str, options: Options ) -> bool:
     return False
 
 
+def are_face_nodes_counterclockwise( face: any ) -> bool:
+    """Checks if the nodes of a face are ordered counterclockwise when looking at the plan created by the face.
+
+    Args:
+        face (any): Face of a vtkCell.
+
+    Returns:
+        bool: True if counterclockwise, False instead.
+    """
+    face_points = face.GetPoints()
+    number_points = face_points.GetNumberOfPoints()
+    if number_points < 3:
+        err_msg = f"The face has less than 3 nodes which is invalid."
+        logging.error( err_msg )
+        raise ValueError( err_msg )
+    # first calculate the normal vector of the face
+    a = np.array( face_points.GetPoint( 0 ) )
+    b = np.array( face_points.GetPoint( 1 ) )
+    c = np.array( face_points.GetPoint( 2 ) )
+    AB = b - a
+    AC = c - a
+    normal = np.cross(AB, AC)
+
+    # then calculate the vector cross products sum of all points of the face with a random point P
+    # from discussion https://math.stackexchange.com/questions/2152623/determine-the-order-of-a-3d-polygon
+    P = np.array( [ 0.0, 0.0, 0.0 ] )  # P position does not change the value of sum
+    all_points = [ np.array( face_points.GetPoint( v ) ) for v in range( number_points ) ]
+    vector_sum = np.array( [ 0.0, 0.0, 0.0 ] )
+    for i in range( number_points ):
+        PAi = all_points[ i % number_points ] - P
+        PAiplus1 = all_points[ ( i + 1 ) % number_points ] - P
+        vector_sum += np.cross( PAi, PAiplus1 )
+    vector_sum = vector_sum / 2  # needs to be half
+
+    # if dot product is positive, the nodes are ordered counterclockwise 
+    if np.dot( vector_sum, normal ) > 0.0:
+        return True
+    return False
+
+
+def valid_cell_point_ids_ordering( cell: vtkCell ) -> list[ int ]:
+    """Returns the valid order of point ids of a cell that respect counterclockwise convention of each face.
+
+    Args:
+        cell (vtkCell): A cell of vtk grid.
+
+    Raises:
+        ValueError: "No node ids permutation made the face nodes to be ordered counterclockwise."
+
+    Returns:
+        list[ int ]: [ pt_id0, pt_id2, pt_id1, ..., pt_idN ]
+    """
+    initial_ids_order: list[ int ] = list( vtk_iter( cell.GetPointIds() ) )
+    valid_points_id: list[ int ] = initial_ids_order.copy()
+
+    for face_id in range( cell.GetNumberOfFaces() ):
+        face = cell.GetFace( face_id )
+        if are_face_nodes_counterclockwise( face ):
+            continue  # the face nodes respect the convention so continue to next face
+
+        reordered = False
+        initial_ids: list[ int ] = [ face.GetPointId( i ) for i in range( face.GetNumberOfPoints() ) ]
+        initial_coords: list[ float ] = [ face.GetPoints().GetPoint( v ) for v in range( face.GetNumberOfPoints() ) ]
+        for permutation_ids in permutations( initial_ids ):
+            for i, node_id in enumerate( permutation_ids ):
+                face.GetPointIds().SetId( i, node_id )
+                face.GetPoints().SetPoint( i, initial_coords[ initial_ids.index( node_id ) ] )
+
+            if are_face_nodes_counterclockwise( face ):  # the correct permutation was found
+                for j, initial_id in enumerate(initial_ids):
+                    if initial_id != permutation_ids[ j ]:
+                        valid_points_id[ initial_ids_order.index( initial_id ) ] = permutation_ids[ j ]
+                reordered = True
+                break
+
+        if not reordered:
+            err_msg = ( f"For face with nodes id '{initial_ids}', that corresponds to coordinates '{initial_coords}'" +
+                        ", no node ids permutation made the face nodes to be ordered counterclockwise." )
+            logging.error( err_msg )
+            raise ValueError( err_msg )
+
+    return valid_points_id
+
+
 def reorder_nodes_to_new_mesh( old_mesh: vtkDataSet, options: Options ) -> tuple:
     """From an old mesh, creates a new one where certain cell nodes are reordered to obtain a correct volume.
 
@@ -118,7 +202,7 @@ def reorder_nodes_to_new_mesh( old_mesh: vtkDataSet, options: Options ) -> tuple
     names_with_totals: dict[ str, int ] = { n: v for n, v in zip( unique_cell_names, total_per_cell_types ) }
     # sorted dict allow for sorted output of reordering_stats
     names_with_totals_sorted: dict[ str, int ] = dict( sorted( names_with_totals.items() ) )
-    useful_VTK_TYPEs: list[ int ] = [ NAME_TO_VTK_TYPE[ name ] for name in options.cell_name_to_ordering.keys() ]
+    useful_VTK_TYPEs: list[ int ] = [ NAME_TO_VTK_TYPE[ name ] for name in options.cell_names_to_reorder ]
     all_cells_volume: np.array = compute_mesh_cells_volume( old_mesh )
     # a new mesh with the same data is created from the old mesh
     new_mesh: vtkDataSet = old_mesh.NewInstance()
@@ -132,7 +216,7 @@ def reorder_nodes_to_new_mesh( old_mesh: vtkDataSet, options: Options ) -> tuple
         "Types non reordered": list( names_with_totals_sorted.keys() ),
         "Number of cells non reordered": list( names_with_totals_sorted.values() )
     }
-    counter_cells_reordered: dict[ str, int ] = { name: 0 for name in options.cell_name_to_ordering.keys() }
+    counter_cells_reordered: dict[ str, int ] = { name: 0 for name in options.cell_names_to_reorder }
     # sorted dict allow for sorted output of reordering_stats
     ounter_cells_reordered_sorted: dict[ str, int ] = dict( sorted( counter_cells_reordered.items() ) )
     # Reordering operations
@@ -141,13 +225,8 @@ def reorder_nodes_to_new_mesh( old_mesh: vtkDataSet, options: Options ) -> tuple
         if vtk_type in useful_VTK_TYPEs:
             if is_cell_to_reorder( float( all_cells_volume[ cell_id ] ), options ):
                 cell_name: str = VTK_TYPE_TO_NAME[ vtk_type ]
-                support_point_ids = vtkIdList()
-                cells.GetCellAtId( cell_id, support_point_ids )
-                new_support_point_ids: list[ int ] = list()
-                node_ordering: list[ int ] = options.cell_name_to_ordering[ cell_name ]
-                for i in range( len( node_ordering ) ):
-                    new_support_point_ids.append( support_point_ids.GetId( node_ordering[ i ] ) )
-                cells.ReplaceCellAtId( cell_id, to_vtk_id_list( new_support_point_ids ) )
+                point_ids_ordering: list[ int ] = valid_cell_point_ids_ordering( new_mesh.GetCell( cell_id ) )
+                cells.ReplaceCellAtId( cell_id, to_vtk_id_list( point_ids_ordering ) )
                 ounter_cells_reordered_sorted[ cell_name ] += 1
     # Calculation of stats
     for cell_name_reordered, amount in ounter_cells_reordered_sorted.items():

geos-mesh/src/geos/mesh/doctor/checks/vtk_utils.py

@@ -103,7 +103,7 @@ def read_mesh( vtk_input_file: str ) -> vtkUnstructuredGrid:
     :return: A unstructured grid.
     """
     if not is_filepath_valid( vtk_input_file ):
-        err_msg: str = f"Invalid file path. Could not read \"{vtk_input_file}\". Dying..."
+        err_msg: str = f"Invalid file path. Could not read \"{vtk_input_file}\"."
         logging.error( err_msg )
         raise ValueError( err_msg )
     file_extension = os.path.splitext( vtk_input_file )[ -1 ]
@@ -119,7 +119,7 @@ def read_mesh( vtk_input_file: str ) -> vtkUnstructuredGrid:
         if output_mesh:
             return output_mesh
     # No reader did work. Dying.
-    err_msg = f"Could not find the appropriate VTK reader for file \"{vtk_input_file}\". Dying..."
+    err_msg = f"Could not find the appropriate VTK reader for file \"{vtk_input_file}\"."
     logging.error( err_msg )
     raise ValueError( err_msg )
 
@@ -160,7 +160,7 @@ def write_mesh( mesh: vtkUnstructuredGrid, vtk_output: VtkOutput ) -> int:
         success_code = __write_vtu( mesh, vtk_output.output, vtk_output.is_data_mode_binary )
     else:
         # No writer found did work. Dying.
-        err_msg = f"Could not find the appropriate VTK writer for extension \"{file_extension}\". Dying..."
+        err_msg = f"Could not find the appropriate VTK writer for extension \"{file_extension}\"."
         logging.error( err_msg )
         raise ValueError( err_msg )
     return 0 if success_code else 2  # the Write member function return 1 in case of success, 0 otherwise.

geos-mesh/src/geos/mesh/doctor/parsing/fix_elements_orderings_parsing.py

@@ -1,41 +1,32 @@
 import logging
 import random
-from geos.mesh.doctor.checks.fix_elements_orderings import Options, Result
+from geos.mesh.doctor.checks.fix_elements_orderings import Options, Result, NAME_TO_VTK_TYPE
 from . import vtk_output_parsing, FIX_ELEMENTS_ORDERINGS
 
-__CELL_NAME_WITH_NUMBER_NODES = {
-    "Tetrahedron": 4,
-    "Pyramid": 5,
-    "Wedge": 6,
-    "Hexahedron": 8,
-    "Prism5": 10,
-    "Prism6": 12
-}
+__CELL_NAMES = "cell_names"
+__CELL_NAMES_CHOICES = list( NAME_TO_VTK_TYPE.keys() )
 
 __VOLUME_TO_REORDER = "volume_to_reorder"
-__VOLUME_TO_REORDER_DEFAULT = "all"
+__VOLUME_TO_REORDER_DEFAULT = "negative"
 __VOLUME_TO_REORDER_CHOICES = [ "all", "positive", "negative" ]
 
 
 def fill_subparser( subparsers ) -> None:
     p = subparsers.add_parser( FIX_ELEMENTS_ORDERINGS, help="Reorders the support nodes for the given cell types." )
-    for element_name, size in __CELL_NAME_WITH_NUMBER_NODES.items():
-        tmp = list( range( size ) )
-        random.Random( 4 ).shuffle( tmp )
-        p.add_argument( '--' + element_name,
-                        type=str,
-                        metavar=",".join( map( str, tmp ) ),
-                        default=None,
-                        required=False,
-                        help=f"[list of integers]: node permutation for \"{element_name}\"." )
+    p.add_argument( '--' + __CELL_NAMES,
+                    type=str,
+                    metavar=", ".join( map( str, __CELL_NAMES_CHOICES ) ),
+                    default=", ".join( map( str, __CELL_NAMES_CHOICES ) ),
+                    help=f"[list of str]: Cell names that can be reordered in your grid. You can use multiple names." +
+                    "Defaults to all cell names being used."  )
     p.add_argument( '--' + __VOLUME_TO_REORDER,
                     type=str,
                     default=__VOLUME_TO_REORDER_DEFAULT,
-                    metavar=",".join( map( str, __VOLUME_TO_REORDER_CHOICES ) ),
-                    required=True,
+                    metavar=", ".join( __VOLUME_TO_REORDER_CHOICES ),
                     help="[str]: Select which element volume is invalid and needs reordering." +
                     " 'all' will allow reordering of nodes for every element, regarding of their volume." +
-                    " 'positive' or 'negative' will only reorder the element with the corresponding volume." )
+                    " 'positive' or 'negative' will only reorder the element with the corresponding volume." +
+                    " Defaults to 'negative'."  )
     vtk_output_parsing.fill_vtk_output_subparser( p )
 
 
@@ -45,22 +36,19 @@ def convert( parsed_options ) -> Options:
     :param options_str: Parsed cli options.
     :return: Options instance.
     """
-    cell_name_to_ordering: dict[ str, list[ int ] ] = {}
-    for element_name, size in __CELL_NAME_WITH_NUMBER_NODES.items():
-        raw_mapping = parsed_options[ element_name ]
-        if raw_mapping:
-            nodes_ordering = tuple( map( int, raw_mapping.split( "," ) ) )
-            if not set( nodes_ordering ) == set( range( size ) ):
-                err_msg: str = f"Permutation {raw_mapping} for type {element_name} is not valid."
-                logging.error( err_msg )
-                raise ValueError( err_msg )
-            cell_name_to_ordering[ element_name ] = nodes_ordering
+    raw_mapping = parsed_options[ __CELL_NAMES ]
+    cell_names_to_reorder = tuple( raw_mapping.split( "," ) )
+    for cell_name in cell_names_to_reorder:
+        if cell_name not in __CELL_NAMES_CHOICES:
+            raise ValueError( f"Please choose names between these options for --{__CELL_NAMES_CHOICES}:" +
+                              f" {__CELL_NAMES_CHOICES}." )
     vtk_output = vtk_output_parsing.convert( parsed_options )
     volume_to_reorder: str = parsed_options[ __VOLUME_TO_REORDER ]
     if volume_to_reorder.lower() not in __VOLUME_TO_REORDER_CHOICES:
-        raise ValueError( f"Please use one of these options for --volume_to_reorder: {__VOLUME_TO_REORDER_CHOICES}." )
+        raise ValueError( f"Please use one of these options for --{__VOLUME_TO_REORDER}:" +
+                          f" {__VOLUME_TO_REORDER_CHOICES}." )
     return Options( vtk_output=vtk_output,
-                    cell_name_to_ordering=cell_name_to_ordering,
+                    cell_names_to_reorder=cell_names_to_reorder,
                     volume_to_reorder=volume_to_reorder )