]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Move GridReordering<dim>::invert_all_negative_measure_cells().
authorDavid Wells <drwells@email.unc.edu>
Wed, 21 Apr 2021 15:40:49 +0000 (11:40 -0400)
committerDavid Wells <drwells@email.unc.edu>
Wed, 21 Apr 2021 17:46:29 +0000 (13:46 -0400)
include/deal.II/grid/grid_reordering.h
include/deal.II/grid/grid_tools.h
source/grid/grid_generator.cc
source/grid/grid_in.cc
source/grid/grid_reordering.cc
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in

index 9a9386a8893051aec2609fbf1b4fd37a75daed18..6d7d76ac78c6e8b1165aa562995a4ebf1feb8d07 100644 (file)
@@ -81,7 +81,10 @@ public:
    * vertices within a cell. If false (the default), then use the "old-style"
    * ordering of vertices within cells used by deal.II before version 5.2 and
    * as explained in the documentation of this class.
+   *
+   * @deprecated Use GridTools::invert_all_negative_measure_cells() instead.
    */
+  DEAL_II_DEPRECATED_EARLY
   static void
   invert_all_cells_of_negative_grid(
     const std::vector<Point<spacedim>> &all_vertices,
index a9a1666399f0782d3a39f23091c23fa80b1c6810..723ed0f80235ce69b2e05b31440feb5ea0d19c44 100644 (file)
@@ -436,6 +436,30 @@ namespace GridTools
                              std::vector<unsigned int> &   considered_vertices,
                              const double                  tol = 1e-12);
 
+  /**
+   * Grids generated by grid generators may have an orientation of cells which
+   * is the inverse of the orientation required by deal.II.
+   *
+   * In 2d and 3d this function checks whether all cells have negative or
+   * positive measure/volume. In the former case, all cells are inverted. It
+   * does nothing in 1d.
+   *
+   * The inversion of cells might also work when only a subset of all cells
+   * have negative volume. However, grids consisting of a mixture of negative
+   * and positively oriented cells are very likely to be broken. Therefore, an
+   * exception is thrown, in case cells are not uniformly oriented.
+   *
+   * @note This function should be called before GridTools::consistently_order_cells().
+   *
+   * @param all_vertices The vertices of the mesh.
+   * @param cells The array of CellData objects that describe the mesh's topology.
+   */
+  template <int dim, int spacedim>
+  void
+  invert_all_negative_measure_cells(
+    const std::vector<Point<spacedim>> &all_vertices,
+    std::vector<CellData<dim>> &        cells);
+
   /*@}*/
   /**
    * @name Rotating, stretching and otherwise transforming meshes
index 2330327063491f5e562034c53b53714705869c71..f7e99e478c6c17c8bb0868069e845f0e2c5d1741 100644 (file)
@@ -1860,9 +1860,7 @@ namespace GridGenerator
     cells[n_cells - 1].vertices[GeometryInfo<3>::ucd_to_deal[7]] =
       (1 + n_rotations) % 4;
 
-    GridReordering<dim>::invert_all_cells_of_negative_grid(vertices,
-                                                           cells,
-                                                           true);
+    GridTools::invert_all_negative_measure_cells(vertices, cells);
     tria.create_triangulation(vertices, cells, SubCellData());
   }
 
index 4696f6a679fce54769b1f22186670f427b25d3da..0ee5e9b82e7c5c55d8eec420f7b98995d9704fdd 100644 (file)
@@ -560,7 +560,7 @@ GridIn<dim, spacedim>::read_vtk(std::istream &in)
 
 
       // TODO: the functions below (GridTools::delete_unused_vertices(),
-      // GridTools::invert_all_cells_of_negative_grid(),
+      // GridTools::invert_all_negative_measure_cells(),
       // GridReordering::reorder_cells()) need to be
       // revisited for simplex/mixed meshes
 
@@ -569,8 +569,7 @@ GridIn<dim, spacedim>::read_vtk(std::istream &in)
           GridTools::delete_unused_vertices(vertices, cells, subcelldata);
 
           if (dim == spacedim)
-            GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid(
-              vertices, cells, true);
+            GridTools::invert_all_negative_measure_cells(vertices, cells);
 
           GridReordering<dim, spacedim>::reorder_cells(cells, true);
           tria->create_triangulation(vertices, cells, subcelldata);
@@ -861,9 +860,7 @@ GridIn<dim, spacedim>::read_unv(std::istream &in)
   GridTools::delete_unused_vertices(vertices, cells, subcelldata);
 
   if (dim == spacedim)
-    GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid(vertices,
-                                                                     cells,
-                                                                     true);
+    GridTools::invert_all_negative_measure_cells(vertices, cells);
 
   GridReordering<dim, spacedim>::reorder_cells(cells, true);
 
@@ -1097,9 +1094,7 @@ GridIn<dim, spacedim>::read_ucd(std::istream &in,
   GridTools::delete_unused_vertices(vertices, cells, subcelldata);
   // ... and cells
   if (dim == spacedim)
-    GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid(vertices,
-                                                                     cells,
-                                                                     true);
+    GridTools::invert_all_negative_measure_cells(vertices, cells);
   GridReordering<dim, spacedim>::reorder_cells(cells, true);
   tria->create_triangulation(vertices, cells, subcelldata);
 }
@@ -1352,9 +1347,7 @@ GridIn<dim, spacedim>::read_dbmesh(std::istream &in)
   // do some clean-up on vertices...
   GridTools::delete_unused_vertices(vertices, cells, subcelldata);
   // ...and cells
-  GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid(vertices,
-                                                                   cells,
-                                                                   true);
+  GridTools::invert_all_negative_measure_cells(vertices, cells);
   GridReordering<dim, spacedim>::reorder_cells(cells, true);
   tria->create_triangulation(vertices, cells, subcelldata);
 }
@@ -1422,9 +1415,7 @@ GridIn<dim, spacedim>::read_xda(std::istream &in)
   // do some clean-up on vertices...
   GridTools::delete_unused_vertices(vertices, cells, subcelldata);
   // ... and cells
-  GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid(vertices,
-                                                                   cells,
-                                                                   true);
+  GridTools::invert_all_negative_measure_cells(vertices, cells);
   GridReordering<dim>::reorder_cells(cells, true);
   tria->create_triangulation(vertices, cells, subcelldata);
 }
@@ -2104,7 +2095,7 @@ GridIn<dim, spacedim>::read_msh(std::istream &in)
   AssertThrow(cells.size() > 0, ExcGmshNoCellInformation());
 
   // TODO: the functions below (GridTools::delete_unused_vertices(),
-  // GridTools::invert_all_cells_of_negative_grid(),
+  // GridTools::invert_all_negative_measure_cells(),
   // GridReordering::reorder_cells()) need to be revisited
   // for simplex/mixed meshes
 
@@ -2114,8 +2105,7 @@ GridIn<dim, spacedim>::read_msh(std::istream &in)
       GridTools::delete_unused_vertices(vertices, cells, subcelldata);
       // ... and cells
       if (dim == spacedim)
-        GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid(
-          vertices, cells, true);
+        GridTools::invert_all_negative_measure_cells(vertices, cells);
       GridReordering<dim, spacedim>::reorder_cells(cells, true);
     }
   tria->create_triangulation(vertices, cells, subcelldata);
@@ -2808,9 +2798,7 @@ GridIn<2>::read_tecplot(std::istream &in)
   AssertThrow(in, ExcIO());
 
   // do some cleanup on cells
-  GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid(vertices,
-                                                                   cells,
-                                                                   true);
+  GridTools::invert_all_negative_measure_cells(vertices, cells);
   GridReordering<dim, spacedim>::reorder_cells(cells, true);
   tria->create_triangulation(vertices, cells, subcelldata);
 }
@@ -2965,9 +2953,7 @@ GridIn<dim, spacedim>::read_assimp(const std::string &filename,
 
   GridTools::delete_unused_vertices(vertices, cells, subcelldata);
   if (dim == spacedim)
-    GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid(vertices,
-                                                                     cells,
-                                                                     true);
+    GridTools::invert_all_negative_measure_cells(vertices, cells);
   GridReordering<dim, spacedim>::reorder_cells(cells, true);
   tria->create_triangulation(vertices, cells, subcelldata);
 
index 922af99583743a0052d10b48bb42a029948424ec..76f532e1d53ea625f879559135aaee7c5f53d17b 100644 (file)
@@ -1154,66 +1154,13 @@ GridReordering<2>::invert_all_cells_of_negative_grid(
   std::vector<CellData<2>> &   cells,
   const bool                   use_new_style_ordering)
 {
-  unsigned int vertices_lex[GeometryInfo<2>::vertices_per_cell];
-  unsigned int n_negative_cells = 0;
-
-  // GridTools::cell_measure requires the vertices to be in lexicographic
-  // ordering
-  auto copy_vertices_to_temp = [&](const CellData<2> &cell) {
-    if (use_new_style_ordering)
-      {
-        std::copy(cell.vertices.begin(),
-                  cell.vertices.end(),
-                  std::begin(vertices_lex));
-      }
-    else
-      {
-        for (const unsigned int i : GeometryInfo<2>::vertex_indices())
-          vertices_lex[GeometryInfo<2>::ucd_to_deal[i]] = cell.vertices[i];
-      }
-  };
+  if (!use_new_style_ordering)
+    reorder_old_to_new_style(cells);
 
-  for (auto &cell : cells)
-    {
-      copy_vertices_to_temp(cell);
-      if (GridTools::cell_measure<2>(all_vertices, vertices_lex) < 0)
-        {
-          ++n_negative_cells;
-          if (use_new_style_ordering)
-            std::swap(cell.vertices[1], cell.vertices[2]);
-          else
-            std::swap(cell.vertices[1], cell.vertices[3]);
-
-          // Check whether the resulting cell is now ok.
-          // If not, then the grid is seriously broken and
-          // we just give up.
-          copy_vertices_to_temp(cell);
-          AssertThrow(GridTools::cell_measure<2>(all_vertices, vertices_lex) >
-                        0,
-                      ExcInternalError());
-        }
-    }
+  GridTools::invert_all_negative_measure_cells(all_vertices, cells);
 
-  // We assume that all cells of a grid have
-  // either positive or negative volumes but
-  // not both mixed. Although above reordering
-  // might work also on single cells, grids
-  // with both kind of cells are very likely to
-  // be broken. Check for this here.
-  AssertThrow(
-    n_negative_cells == 0 || n_negative_cells == cells.size(),
-    ExcMessage(
-      std::string("This class assumes that either all cells have positive "
-                  "volume, or that all cells have been specified in an "
-                  "inverted vertex order so that their volume is negative. "
-                  "(In the latter case, this class automatically inverts "
-                  "every cell.) However, the mesh you have specified "
-                  "appears to have both cells with positive and cells with "
-                  "negative volume. You need to check your mesh which "
-                  "cells these are and how they got there.\n"
-                  "As a hint, of the total ") +
-      std::to_string(cells.size()) + " cells in the mesh, " +
-      std::to_string(n_negative_cells) + " appear to have a negative volume."));
+  if (!use_new_style_ordering)
+    reorder_new_to_old_style(cells);
 }
 
 
@@ -1237,82 +1184,13 @@ GridReordering<3>::invert_all_cells_of_negative_grid(
   std::vector<CellData<3>> &   cells,
   const bool                   use_new_style_ordering)
 {
-  unsigned int vertices_lex[GeometryInfo<3>::vertices_per_cell];
-  unsigned int n_negative_cells = 0;
-
-  // GridTools::cell_measure requires the vertices to be in lexicographic
-  // ordering
-  auto copy_vertices_to_temp = [&](const CellData<3> &cell) {
-    if (use_new_style_ordering)
-      {
-        std::copy(cell.vertices.begin(),
-                  cell.vertices.end(),
-                  std::begin(vertices_lex));
-      }
-    else
-      {
-        for (const unsigned int i : GeometryInfo<3>::vertex_indices())
-          vertices_lex[GeometryInfo<3>::ucd_to_deal[i]] = cell.vertices[i];
-      }
-  };
+  if (!use_new_style_ordering)
+    reorder_old_to_new_style(cells);
 
-  for (auto &cell : cells)
-    {
-      copy_vertices_to_temp(cell);
-      if (GridTools::cell_measure<3>(all_vertices, vertices_lex) < 0)
-        {
-          ++n_negative_cells;
-          // reorder vertices: swap front and back face
-          if (use_new_style_ordering)
-            {
-              std::swap(cell.vertices[0], cell.vertices[2]);
-              std::swap(cell.vertices[1], cell.vertices[3]);
-              std::swap(cell.vertices[4], cell.vertices[6]);
-              std::swap(cell.vertices[5], cell.vertices[7]);
-            }
-          else
-            {
-              for (unsigned int i = 0; i < 4; ++i)
-                std::swap(cell.vertices[i], cell.vertices[i + 4]);
-            }
-          copy_vertices_to_temp(cell);
-
-          // Check whether the resulting cell is now ok.
-          // If not, then the grid is seriously broken and
-          // we just give up.
-          AssertThrow(GridTools::cell_measure<3>(all_vertices, vertices_lex) >
-                        0,
-                      ExcInternalError());
-        }
-    }
+  GridTools::invert_all_negative_measure_cells(all_vertices, cells);
 
-  // We assume that all cells of a
-  // grid have either positive or
-  // negative volumes but not both
-  // mixed. Although above reordering
-  // might work also on single cells,
-  // grids with both kind of cells
-  // are very likely to be
-  // broken. Check for this here.
-  AssertThrow(n_negative_cells == 0 || n_negative_cells == cells.size(),
-              ExcMessage("While sorting the cells that will be passed for "
-                         "creating a Triangulation object, deal.II found that "
-                         "some but not all cells have a negative volume. (If "
-                         "all cells had a negative volume, they would simply "
-                         "all have been inverted.) This usually happens in "
-                         "hand-generated meshes if one accidentally uses an "
-                         "incorrect convention for ordering the vertices in "
-                         "one or more cells; in that case, you may want to "
-                         "double check that you specified the vertex indices "
-                         "in their correct order. If you are reading a mesh "
-                         "that was created by a mesh generator, then this "
-                         "exception indicates that some of the cells created "
-                         "are so badly distorted that their volume becomes "
-                         "negative; this commonly occurs at complex geometric "
-                         "features, and you may see if the problem can be "
-                         "fixed by playing with the parameters that control "
-                         "mesh properties in your mesh generator, such as "
-                         "the number of cells, the mesh density, etc."));
+  if (!use_new_style_ordering)
+    reorder_new_to_old_style(cells);
 }
 
 
index 0b7195776d2644f8a746e71a74ba1a6998022a56..8b125a912f692eab35aaa2de99848942722a59bb 100644 (file)
@@ -858,6 +858,77 @@ namespace GridTools
 
 
 
+  template <int dim, int spacedim>
+  void
+  invert_all_negative_measure_cells(
+    const std::vector<Point<spacedim>> &all_vertices,
+    std::vector<CellData<dim>> &        cells)
+  {
+    if (dim == 1)
+      return;
+    if (dim == 2 && spacedim == 3)
+      Assert(false, ExcNotImplemented());
+
+    std::size_t n_negative_cells = 0;
+    for (auto &cell : cells)
+      {
+        Assert(cell.vertices.size() ==
+                 ReferenceCells::get_hypercube<dim>().n_vertices(),
+               ExcNotImplemented());
+        const ArrayView<const unsigned int> vertices(cell.vertices);
+        if (GridTools::cell_measure<dim>(all_vertices, vertices) < 0)
+          {
+            ++n_negative_cells;
+
+            // TODO: this only works for quads and hexes
+            if (dim == 2)
+              {
+                // flip the cell across the y = x line in 2D
+                std::swap(cell.vertices[1], cell.vertices[2]);
+              }
+            else if (dim == 3)
+              {
+                // swap the front and back faces in 3D
+                std::swap(cell.vertices[0], cell.vertices[2]);
+                std::swap(cell.vertices[1], cell.vertices[3]);
+                std::swap(cell.vertices[4], cell.vertices[6]);
+                std::swap(cell.vertices[5], cell.vertices[7]);
+              }
+
+            // Check whether the resulting cell is now ok.
+            // If not, then the grid is seriously broken and
+            // we just give up.
+            AssertThrow(GridTools::cell_measure<dim>(all_vertices, vertices) >
+                          0,
+                        ExcInternalError());
+          }
+      }
+
+    // We assume that all cells of a grid have
+    // either positive or negative volumes but
+    // not both mixed. Although above reordering
+    // might work also on single cells, grids
+    // with both kind of cells are very likely to
+    // be broken. Check for this here.
+    AssertThrow(n_negative_cells == 0 || n_negative_cells == cells.size(),
+                ExcMessage(
+                  std::string(
+                    "This function assumes that either all cells have positive "
+                    "volume, or that all cells have been specified in an "
+                    "inverted vertex order so that their volume is negative. "
+                    "(In the latter case, this class automatically inverts "
+                    "every cell.) However, the mesh you have specified "
+                    "appears to have both cells with positive and cells with "
+                    "negative volume. You need to check your mesh which "
+                    "cells these are and how they got there.\n"
+                    "As a hint, of the total ") +
+                  std::to_string(cells.size()) + " cells in the mesh, " +
+                  std::to_string(n_negative_cells) +
+                  " appear to have a negative volume."));
+  }
+
+
+
   // define some transformations
   namespace internal
   {
index b972191dc9b61df7c4653b2171d44618998f025d..145b5922a9f8c32162521240529e025a48c905c9 100644 (file)
@@ -245,6 +245,11 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
                                  std::vector<unsigned int> &,
                                  double);
 
+      template void
+      invert_all_negative_measure_cells(
+        const std::vector<Point<deal_II_space_dimension>> &,
+        std::vector<CellData<deal_II_dimension>> &);
+
       template void
       shift<deal_II_dimension>(
         const Tensor<1, deal_II_space_dimension> &,

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.