Merge pull request #13912 from luca-heltai/cgal-fix-orientation-hex

Fix orientation issue with hex cells.

source/cgal/surface_mesh.cc

@@ -114,17 +114,30 @@ namespace CGALWrappers
     else
       // in 3d, we build a surface mesh containing all the faces of the 3d cell.
       // Simplices, Tetrahedrons, and Pyramids have their faces numbered in the
-      // same way as CGAL does (all faces are numbered clockwise). Hexahedrons,
-      // instead, have their faces numbered lexicographically, and one cannot
+      // same way as CGAL does (all faces of a bounding polyhedron are numbered
+      // counter-clockwise, so that their normal points outwards). Hexahedrons
+      // in deal.II have their faces numbered lexicographically, and one cannot
       // deduce the direction of the normals by just looking at the vertices.
+      //
       // In order for CGAL to be able to produce the right orientation, we need
-      // to revers the order of the vertices for faces with even index.
+      // to reverse the order of the vertices for faces with even index.
+      // However, in order to allow for all kinds of meshes in 3d, including
+      // Moebius-loops, a deal.II face might even be rotated looking from one
+      // cell, whereas it is according to the standard when looking at it from
+      // the neighboring cell sharing that particular face. Therefore, when
+      // building a cgal face we must also take into account the fact that a
+      // face may have a non-standard orientation.
       for (const auto &f : cell->face_indices())
-        add_facet(cell->face(f),
-                  deal2cgal,
-                  mesh,
-                  cell->reference_cell() != ReferenceCells::Hexahedron ||
-                    (f % 2 == 0));
+        {
+          // Check for standard orientation of faces
+          bool face_is_clockwise_oriented =
+            cell->reference_cell() != ReferenceCells::Hexahedron ||
+            (f % 2 == 0);
+          // Make sure that we revert the orientation if required
+          if (cell->face_orientation(f) == false)
+            face_is_clockwise_oriented = !face_is_clockwise_oriented;
+          add_facet(cell->face(f), deal2cgal, mesh, face_is_clockwise_oriented);
+        }
   }
 
 

tests/cgal/cgal_surface_mesh_06.cc

@@ -0,0 +1,50 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2022 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Create a Surface_mesh from each cell of an hyper ball, to check for
+// orientation issues that may arise in Hexaedrons with non-standard face
+// orientations.
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/fe/mapping_q1.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/cgal/surface_mesh.h>
+#include <deal.II/cgal/utilities.h>
+
+#include "../tests.h"
+
+using K         = CGAL::Exact_predicates_inexact_constructions_kernel;
+using CGALPoint = CGAL::Point_3<K>;
+
+int
+main()
+{
+  initlog();
+  Triangulation<3> tria;
+  GridGenerator::hyper_ball(tria);
+
+  for (const auto &cell : tria.active_cell_iterators())
+    {
+      CGAL::Surface_mesh<CGALPoint> mesh;
+      CGALWrappers::dealii_cell_to_cgal_surface_mesh(
+        cell, StaticMappingQ1<3>::mapping, mesh);
+      deallog << "cell: " << cell << " is valid: " << (int)mesh.is_valid()
+              << std::endl;
+    }
+}

tests/cgal/cgal_surface_mesh_06.output

@@ -0,0 +1,8 @@
+
+DEAL::cell: 0.0 is valid: 1
+DEAL::cell: 0.1 is valid: 1
+DEAL::cell: 0.2 is valid: 1
+DEAL::cell: 0.3 is valid: 1
+DEAL::cell: 0.4 is valid: 1
+DEAL::cell: 0.5 is valid: 1
+DEAL::cell: 0.6 is valid: 1

tests/cgal/cgal_surface_mesh_07.cc

@@ -0,0 +1,55 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2022 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Create a Surface_mesh from each cell of an hyper ball, to check for
+// orientation issues that may arise in Hexaedrons with non-standard face
+// orientations. Make sure this works also when tranforming the grid to a
+// simplex one.
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/fe/mapping_q1.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/cgal/surface_mesh.h>
+#include <deal.II/cgal/utilities.h>
+
+#include "../tests.h"
+
+using K         = CGAL::Exact_predicates_inexact_constructions_kernel;
+using CGALPoint = CGAL::Point_3<K>;
+
+int
+main()
+{
+  initlog();
+  Triangulation<3> tria;
+  {
+    Triangulation<3> tria_tmp;
+    GridGenerator::hyper_ball(tria_tmp);
+    GridGenerator::convert_hypercube_to_simplex_mesh(tria_tmp, tria);
+  }
+
+  for (const auto &cell : tria.active_cell_iterators())
+    {
+      CGAL::Surface_mesh<CGALPoint> mesh;
+      CGALWrappers::dealii_cell_to_cgal_surface_mesh(
+        cell, StaticMappingQ1<3>::mapping, mesh);
+      deallog << "cell: " << cell << " is valid: " << (int)mesh.is_valid()
+              << std::endl;
+    }
+}

tests/cgal/cgal_surface_mesh_07.output

@@ -0,0 +1,169 @@
+
+DEAL::cell: 0.0 is valid: 1
+DEAL::cell: 0.1 is valid: 1
+DEAL::cell: 0.2 is valid: 1
+DEAL::cell: 0.3 is valid: 1
+DEAL::cell: 0.4 is valid: 1
+DEAL::cell: 0.5 is valid: 1
+DEAL::cell: 0.6 is valid: 1
+DEAL::cell: 0.7 is valid: 1
+DEAL::cell: 0.8 is valid: 1
+DEAL::cell: 0.9 is valid: 1
+DEAL::cell: 0.10 is valid: 1
+DEAL::cell: 0.11 is valid: 1
+DEAL::cell: 0.12 is valid: 1
+DEAL::cell: 0.13 is valid: 1
+DEAL::cell: 0.14 is valid: 1
+DEAL::cell: 0.15 is valid: 1
+DEAL::cell: 0.16 is valid: 1
+DEAL::cell: 0.17 is valid: 1
+DEAL::cell: 0.18 is valid: 1
+DEAL::cell: 0.19 is valid: 1
+DEAL::cell: 0.20 is valid: 1
+DEAL::cell: 0.21 is valid: 1
+DEAL::cell: 0.22 is valid: 1
+DEAL::cell: 0.23 is valid: 1
+DEAL::cell: 0.24 is valid: 1
+DEAL::cell: 0.25 is valid: 1
+DEAL::cell: 0.26 is valid: 1
+DEAL::cell: 0.27 is valid: 1
+DEAL::cell: 0.28 is valid: 1
+DEAL::cell: 0.29 is valid: 1
+DEAL::cell: 0.30 is valid: 1
+DEAL::cell: 0.31 is valid: 1
+DEAL::cell: 0.32 is valid: 1
+DEAL::cell: 0.33 is valid: 1
+DEAL::cell: 0.34 is valid: 1
+DEAL::cell: 0.35 is valid: 1
+DEAL::cell: 0.36 is valid: 1
+DEAL::cell: 0.37 is valid: 1
+DEAL::cell: 0.38 is valid: 1
+DEAL::cell: 0.39 is valid: 1
+DEAL::cell: 0.40 is valid: 1
+DEAL::cell: 0.41 is valid: 1
+DEAL::cell: 0.42 is valid: 1
+DEAL::cell: 0.43 is valid: 1
+DEAL::cell: 0.44 is valid: 1
+DEAL::cell: 0.45 is valid: 1
+DEAL::cell: 0.46 is valid: 1
+DEAL::cell: 0.47 is valid: 1
+DEAL::cell: 0.48 is valid: 1
+DEAL::cell: 0.49 is valid: 1
+DEAL::cell: 0.50 is valid: 1
+DEAL::cell: 0.51 is valid: 1
+DEAL::cell: 0.52 is valid: 1
+DEAL::cell: 0.53 is valid: 1
+DEAL::cell: 0.54 is valid: 1
+DEAL::cell: 0.55 is valid: 1
+DEAL::cell: 0.56 is valid: 1
+DEAL::cell: 0.57 is valid: 1
+DEAL::cell: 0.58 is valid: 1
+DEAL::cell: 0.59 is valid: 1
+DEAL::cell: 0.60 is valid: 1
+DEAL::cell: 0.61 is valid: 1
+DEAL::cell: 0.62 is valid: 1
+DEAL::cell: 0.63 is valid: 1
+DEAL::cell: 0.64 is valid: 1
+DEAL::cell: 0.65 is valid: 1
+DEAL::cell: 0.66 is valid: 1
+DEAL::cell: 0.67 is valid: 1
+DEAL::cell: 0.68 is valid: 1
+DEAL::cell: 0.69 is valid: 1
+DEAL::cell: 0.70 is valid: 1
+DEAL::cell: 0.71 is valid: 1
+DEAL::cell: 0.72 is valid: 1
+DEAL::cell: 0.73 is valid: 1
+DEAL::cell: 0.74 is valid: 1
+DEAL::cell: 0.75 is valid: 1
+DEAL::cell: 0.76 is valid: 1
+DEAL::cell: 0.77 is valid: 1
+DEAL::cell: 0.78 is valid: 1
+DEAL::cell: 0.79 is valid: 1
+DEAL::cell: 0.80 is valid: 1
+DEAL::cell: 0.81 is valid: 1
+DEAL::cell: 0.82 is valid: 1
+DEAL::cell: 0.83 is valid: 1
+DEAL::cell: 0.84 is valid: 1
+DEAL::cell: 0.85 is valid: 1
+DEAL::cell: 0.86 is valid: 1
+DEAL::cell: 0.87 is valid: 1
+DEAL::cell: 0.88 is valid: 1
+DEAL::cell: 0.89 is valid: 1
+DEAL::cell: 0.90 is valid: 1
+DEAL::cell: 0.91 is valid: 1
+DEAL::cell: 0.92 is valid: 1
+DEAL::cell: 0.93 is valid: 1
+DEAL::cell: 0.94 is valid: 1
+DEAL::cell: 0.95 is valid: 1
+DEAL::cell: 0.96 is valid: 1
+DEAL::cell: 0.97 is valid: 1
+DEAL::cell: 0.98 is valid: 1
+DEAL::cell: 0.99 is valid: 1
+DEAL::cell: 0.100 is valid: 1
+DEAL::cell: 0.101 is valid: 1
+DEAL::cell: 0.102 is valid: 1
+DEAL::cell: 0.103 is valid: 1
+DEAL::cell: 0.104 is valid: 1
+DEAL::cell: 0.105 is valid: 1
+DEAL::cell: 0.106 is valid: 1
+DEAL::cell: 0.107 is valid: 1
+DEAL::cell: 0.108 is valid: 1
+DEAL::cell: 0.109 is valid: 1
+DEAL::cell: 0.110 is valid: 1
+DEAL::cell: 0.111 is valid: 1
+DEAL::cell: 0.112 is valid: 1
+DEAL::cell: 0.113 is valid: 1
+DEAL::cell: 0.114 is valid: 1
+DEAL::cell: 0.115 is valid: 1
+DEAL::cell: 0.116 is valid: 1
+DEAL::cell: 0.117 is valid: 1
+DEAL::cell: 0.118 is valid: 1
+DEAL::cell: 0.119 is valid: 1
+DEAL::cell: 0.120 is valid: 1
+DEAL::cell: 0.121 is valid: 1
+DEAL::cell: 0.122 is valid: 1
+DEAL::cell: 0.123 is valid: 1
+DEAL::cell: 0.124 is valid: 1
+DEAL::cell: 0.125 is valid: 1
+DEAL::cell: 0.126 is valid: 1
+DEAL::cell: 0.127 is valid: 1
+DEAL::cell: 0.128 is valid: 1
+DEAL::cell: 0.129 is valid: 1
+DEAL::cell: 0.130 is valid: 1
+DEAL::cell: 0.131 is valid: 1
+DEAL::cell: 0.132 is valid: 1
+DEAL::cell: 0.133 is valid: 1
+DEAL::cell: 0.134 is valid: 1
+DEAL::cell: 0.135 is valid: 1
+DEAL::cell: 0.136 is valid: 1
+DEAL::cell: 0.137 is valid: 1
+DEAL::cell: 0.138 is valid: 1
+DEAL::cell: 0.139 is valid: 1
+DEAL::cell: 0.140 is valid: 1
+DEAL::cell: 0.141 is valid: 1
+DEAL::cell: 0.142 is valid: 1
+DEAL::cell: 0.143 is valid: 1
+DEAL::cell: 0.144 is valid: 1
+DEAL::cell: 0.145 is valid: 1
+DEAL::cell: 0.146 is valid: 1
+DEAL::cell: 0.147 is valid: 1
+DEAL::cell: 0.148 is valid: 1
+DEAL::cell: 0.149 is valid: 1
+DEAL::cell: 0.150 is valid: 1
+DEAL::cell: 0.151 is valid: 1
+DEAL::cell: 0.152 is valid: 1
+DEAL::cell: 0.153 is valid: 1
+DEAL::cell: 0.154 is valid: 1
+DEAL::cell: 0.155 is valid: 1
+DEAL::cell: 0.156 is valid: 1
+DEAL::cell: 0.157 is valid: 1
+DEAL::cell: 0.158 is valid: 1
+DEAL::cell: 0.159 is valid: 1
+DEAL::cell: 0.160 is valid: 1
+DEAL::cell: 0.161 is valid: 1
+DEAL::cell: 0.162 is valid: 1
+DEAL::cell: 0.163 is valid: 1
+DEAL::cell: 0.164 is valid: 1
+DEAL::cell: 0.165 is valid: 1
+DEAL::cell: 0.166 is valid: 1
+DEAL::cell: 0.167 is valid: 1