From 08f43c0acb891f71e02f47c86e5bbaf4df1e162d Mon Sep 17 00:00:00 2001
From: Daniel Arndt <arndtd@ornl.gov>
Date: Sun, 2 Jul 2023 09:14:49 -0400
Subject: [PATCH] Remove deprecated
 Triangulation::create_triangulation_compatibility

---
 examples/step-1/doc/intro.dox                |  4 +-
 include/deal.II/grid/persistent_tria.h       | 13 ----
 include/deal.II/grid/tria.h                  | 15 ----
 source/grid/persistent_tria.cc               | 12 ---
 source/grid/tria.cc                          | 72 ------------------
 tests/bits/neighboring_cells_at_two_faces.cc | 10 +--
 tests/bits/q_points.cc                       | 48 +++++++++++-
 tests/data_out/data_out_curved_cells.cc      |  8 +-
 tests/grid/grid_invert_02.cc                 |  5 +-
 tests/grid/line_coarsening_3d.cc             | 14 ++--
 tests/grid/mesh_3d.h                         |  4 +-
 tests/grid/subcelldata.cc                    | 48 +++++++++++-
 tests/lac/constraints.cc                     | 79 ++++++++++++++++----
 13 files changed, 176 insertions(+), 156 deletions(-)

diff --git a/examples/step-1/doc/intro.dox b/examples/step-1/doc/intro.dox
index ead99932b8..c0a2b20bae 100644
--- a/examples/step-1/doc/intro.dox
+++ b/examples/step-1/doc/intro.dox
@@ -17,8 +17,8 @@ different levels:
   every single (member) function in deal.II. You get there if, for
   example, you click on the "Main page" or "Classes" tab at the top of
   this page. This is the place where you would look up what the second
-  argument of Triangulation::create_triangulation_compatibility means,
-  to give just one slightly obscure example. You need this level of
+  argument of Triangulation::create_triangulation means,
+  to give just one example. You need this level of
   documentation for when you know what you want to do, but forgot how
   exactly the function was named, what its arguments are, or what it
   returns. Note that you also get into the manual whenever you read
diff --git a/include/deal.II/grid/persistent_tria.h b/include/deal.II/grid/persistent_tria.h
index 8c99d49514..7ee795b52e 100644
--- a/include/deal.II/grid/persistent_tria.h
+++ b/include/deal.II/grid/persistent_tria.h
@@ -211,19 +211,6 @@ public:
     const TriangulationDescription::Description<dim, spacedim>
       &construction_data) override;
 
-  /**
-   * An overload of the respective function of the base class.
-   *
-   * Throw an error, since this function is not useful in the context of this
-   * class.
-   */
-  DEAL_II_DEPRECATED
-  virtual void
-  create_triangulation_compatibility(
-    const std::vector<Point<spacedim>> &vertices,
-    const std::vector<CellData<dim>> &  cells,
-    const SubCellData &                 subcelldata) override;
-
   /**
    * Write all refine and coarsen flags to the ostream @p out.
    */
diff --git a/include/deal.II/grid/tria.h b/include/deal.II/grid/tria.h
index f3b01b3976..726e3bcad7 100644
--- a/include/deal.II/grid/tria.h
+++ b/include/deal.II/grid/tria.h
@@ -1893,21 +1893,6 @@ public:
     const TriangulationDescription::Description<dim, spacedim>
       &construction_data);
 
-  /**
-   * For backward compatibility, only. This function takes the cell data in
-   * the ordering as requested by deal.II versions up to 5.2, converts it to
-   * the new (lexicographic) ordering and calls create_triangulation().
-   *
-   * @note This function internally calls create_triangulation and therefore
-   * can throw the same exception as the other function.
-   */
-  DEAL_II_DEPRECATED
-  virtual void
-  create_triangulation_compatibility(
-    const std::vector<Point<spacedim>> &vertices,
-    const std::vector<CellData<dim>> &  cells,
-    const SubCellData &                 subcelldata);
-
   /**
    * Revert or flip the direction_flags of a dim<spacedim triangulation, see
    * @ref GlossDirectionFlag.
diff --git a/source/grid/persistent_tria.cc b/source/grid/persistent_tria.cc
index 21189ab35f..a5913ee39e 100644
--- a/source/grid/persistent_tria.cc
+++ b/source/grid/persistent_tria.cc
@@ -157,18 +157,6 @@ PersistentTriangulation<dim, spacedim>::create_triangulation(
 
 
 
-template <int dim, int spacedim>
-void
-PersistentTriangulation<dim, spacedim>::create_triangulation_compatibility(
-  const std::vector<Point<spacedim>> &,
-  const std::vector<CellData<dim>> &,
-  const SubCellData &)
-{
-  Assert(false, ExcImpossibleInDim(dim));
-}
-
-
-
 template <int dim, int spacedim>
 void
 PersistentTriangulation<dim, spacedim>::write_flags(std::ostream &out) const
diff --git a/source/grid/tria.cc b/source/grid/tria.cc
index c568252229..73570ab6c3 100644
--- a/source/grid/tria.cc
+++ b/source/grid/tria.cc
@@ -414,59 +414,6 @@ namespace
 
 
 
-  /**
-   * A set of three functions that
-   * reorder the data given to
-   * create_triangulation_compatibility
-   * from the "classic" to the
-   * "current" format of vertex
-   * numbering of cells and
-   * faces. These functions do the
-   * reordering of their arguments
-   * in-place.
-   */
-  void
-  reorder_compatibility(const std::vector<CellData<1>> &, const SubCellData &)
-  {
-    // nothing to do here: the format
-    // hasn't changed for 1d
-  }
-
-
-  void
-  reorder_compatibility(std::vector<CellData<2>> &cells, const SubCellData &)
-  {
-    for (auto &cell : cells)
-      if (cell.vertices.size() == GeometryInfo<2>::vertices_per_cell)
-        std::swap(cell.vertices[2], cell.vertices[3]);
-  }
-
-
-  void
-  reorder_compatibility(std::vector<CellData<3>> &cells,
-                        SubCellData &             subcelldata)
-  {
-    unsigned int tmp[GeometryInfo<3>::vertices_per_cell];
-    static constexpr std::array<unsigned int,
-                                GeometryInfo<3>::vertices_per_cell>
-      local_vertex_numbering{{0, 1, 5, 4, 2, 3, 7, 6}};
-    for (auto &cell : cells)
-      if (cell.vertices.size() == GeometryInfo<3>::vertices_per_cell)
-        {
-          for (const unsigned int i : GeometryInfo<3>::vertex_indices())
-            tmp[i] = cell.vertices[i];
-          for (const unsigned int i : GeometryInfo<3>::vertex_indices())
-            cell.vertices[local_vertex_numbering[i]] = tmp[i];
-        }
-
-    // now points in boundary quads
-    for (auto &boundary_quad : subcelldata.boundary_quads)
-      if (boundary_quad.vertices.size() == GeometryInfo<2>::vertices_per_cell)
-        std::swap(boundary_quad.vertices[2], boundary_quad.vertices[3]);
-  }
-
-
-
   /**
    * Return the index of the vertex
    * in the middle of this object,
@@ -11548,25 +11495,6 @@ void Triangulation<dim, spacedim>::copy_triangulation(
 
 
 
-template <int dim, int spacedim>
-DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
-void Triangulation<dim, spacedim>::create_triangulation_compatibility(
-  const std::vector<Point<spacedim>> &v,
-  const std::vector<CellData<dim>> &  cells,
-  const SubCellData &                 subcelldata)
-{
-  std::vector<CellData<dim>> reordered_cells(cells);             // NOLINT
-  SubCellData                reordered_subcelldata(subcelldata); // NOLINT
-
-  // in-place reordering of data
-  reorder_compatibility(reordered_cells, reordered_subcelldata);
-
-  // now create triangulation from
-  // reordered data
-  create_triangulation(v, reordered_cells, reordered_subcelldata);
-}
-
-
 template <int dim, int spacedim>
 DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
 void Triangulation<dim, spacedim>::reset_policy()
diff --git a/tests/bits/neighboring_cells_at_two_faces.cc b/tests/bits/neighboring_cells_at_two_faces.cc
index 4ad6151b44..b675f4aa30 100644
--- a/tests/bits/neighboring_cells_at_two_faces.cc
+++ b/tests/bits/neighboring_cells_at_two_faces.cc
@@ -60,19 +60,19 @@ create_grid(Triangulation<2> &tria)
   std::vector<CellData<2>> cells(2);
   cells[0].vertices[0] = 0;
   cells[0].vertices[1] = 1;
-  cells[0].vertices[2] = 2;
-  cells[0].vertices[3] = 3;
+  cells[0].vertices[2] = 3;
+  cells[0].vertices[3] = 2;
   cells[1].vertices[0] = 0;
   cells[1].vertices[1] = 4;
-  cells[1].vertices[2] = 2;
-  cells[1].vertices[3] = 1;
+  cells[1].vertices[2] = 1;
+  cells[1].vertices[3] = 2;
 
   // generate a triangulation
   // out of this
   GridTools::consistently_order_cells(cells);
   try
     {
-      tria.create_triangulation_compatibility(vertices, cells, SubCellData());
+      tria.create_triangulation(vertices, cells, SubCellData());
     }
   catch (Triangulation<2>::DistortedCellList &dcv)
     {
diff --git a/tests/bits/q_points.cc b/tests/bits/q_points.cc
index d5b04dd97a..377856ea09 100644
--- a/tests/bits/q_points.cc
+++ b/tests/bits/q_points.cc
@@ -36,7 +36,48 @@
 
 #include "../tests.h"
 
-
+namespace
+{
+  /**
+   * This file uses a different ordering for the vertices in a hex
+   * cell than we usually do in deal.II. The different convention used
+   * here originates in what we believed the ordering to be in UCD
+   * format, until it was discovered in 2022 that UCD will interpret
+   * this ordering to correspond to inverted cells -- as a
+   * consequence, the UCD ordering was fixed, but the current file is
+   * stuck on the old ordering.
+   */
+  constexpr std::array<unsigned int, 8> local_vertex_numbering{
+    {0, 1, 5, 4, 2, 3, 7, 6}};
+
+  /**
+   * Following is a set of functions that reorder the data from the
+   * "current" to the "classic" format of vertex numbering of cells
+   * and faces. These functions do the reordering of their arguments
+   * in-place.
+   */
+  void
+  reorder_old_to_new_style(std::vector<CellData<2>> &cells)
+  {
+    for (auto &cell : cells)
+      std::swap(cell.vertices[2], cell.vertices[3]);
+  }
+
+
+  void
+  reorder_old_to_new_style(std::vector<CellData<3>> &cells)
+  {
+    // undo the ordering above
+    unsigned int tmp[GeometryInfo<3>::vertices_per_cell];
+    for (auto &cell : cells)
+      {
+        for (const unsigned int i : GeometryInfo<3>::vertex_indices())
+          tmp[i] = cell.vertices[i];
+        for (const unsigned int i : GeometryInfo<3>::vertex_indices())
+          cell.vertices[local_vertex_numbering[i]] = tmp[i];
+      }
+  }
+} // namespace
 
 void
 create_two_cubes(Triangulation<3> &coarse_grid)
@@ -69,10 +110,9 @@ create_two_cubes(Triangulation<3> &coarse_grid)
 
   // finally generate a triangulation
   // out of this
+  reorder_old_to_new_style(cells);
   GridTools::consistently_order_cells(cells);
-  coarse_grid.create_triangulation_compatibility(vertices,
-                                                 cells,
-                                                 SubCellData());
+  coarse_grid.create_triangulation(vertices, cells, SubCellData());
 }
 
 
diff --git a/tests/data_out/data_out_curved_cells.cc b/tests/data_out/data_out_curved_cells.cc
index cea603b299..de957fb092 100644
--- a/tests/data_out/data_out_curved_cells.cc
+++ b/tests/data_out/data_out_curved_cells.cc
@@ -124,14 +124,12 @@ curved_grid(std::ostream &out)
 
         cells[index].vertices[0] = p_iplus;
         cells[index].vertices[1] = p;
-        cells[index].vertices[2] = p_jplus;
-        cells[index].vertices[3] = p_iplus_jplus;
+        cells[index].vertices[2] = p_iplus_jplus;
+        cells[index].vertices[3] = p_jplus;
       }
   // create triangulation
   Triangulation<2> triangulation;
-  triangulation.create_triangulation_compatibility(vertices,
-                                                   cells,
-                                                   SubCellData());
+  triangulation.create_triangulation(vertices, cells, SubCellData());
   // now provide everything that is
   // needed for solving a Laplace
   // equation.
diff --git a/tests/grid/grid_invert_02.cc b/tests/grid/grid_invert_02.cc
index 5cf5a6bc5b..0a3ce5743f 100644
--- a/tests/grid/grid_invert_02.cc
+++ b/tests/grid/grid_invert_02.cc
@@ -58,8 +58,11 @@ test()
 
   deallog << "We inverted " << n_cells_inverted << " cells." << std::endl;
 
+  std::swap(cells[0].vertices[2], cells[0].vertices[3]);
+  std::swap(cells[1].vertices[2], cells[1].vertices[3]);
+
   Triangulation<dim> tria;
-  tria.create_triangulation_compatibility(vertices, cells, subcelldata);
+  tria.create_triangulation(vertices, cells, subcelldata);
 
   std::ostream &logfile = deallog.get_file_stream();
   logfile << "---------------------------------------------" << std::endl
diff --git a/tests/grid/line_coarsening_3d.cc b/tests/grid/line_coarsening_3d.cc
index ed7994d831..20a12d6bc6 100644
--- a/tests/grid/line_coarsening_3d.cc
+++ b/tests/grid/line_coarsening_3d.cc
@@ -55,16 +55,16 @@ create_star_structured_cylinder(Triangulation<3> & coarse_grid,
     {
       cells[c].vertices[0] = 0;
       cells[c].vertices[1] = 1 + 2 * c;
-      cells[c].vertices[2] = 2 + 2 * c;
-      cells[c].vertices[3] = (3 + 2 * c) % (2 * n_cells);
-      cells[c].vertices[4] = 0 + 2 * n_cells + 1;
-      cells[c].vertices[5] = 1 + 2 * c + 2 * n_cells + 1;
-      cells[c].vertices[6] = 2 + 2 * c + 2 * n_cells + 1;
-      cells[c].vertices[7] = (3 + 2 * c) % (2 * n_cells) + 2 * n_cells + 1;
+      cells[c].vertices[2] = 0 + 2 * n_cells + 1;
+      cells[c].vertices[3] = 1 + 2 * c + 2 * n_cells + 1;
+      cells[c].vertices[4] = (3 + 2 * c) % (2 * n_cells);
+      cells[c].vertices[5] = 2 + 2 * c;
+      cells[c].vertices[6] = (3 + 2 * c) % (2 * n_cells) + 2 * n_cells + 1;
+      cells[c].vertices[7] = 2 + 2 * c + 2 * n_cells + 1;
     }
   // finally generate a triangulation
   // out of this
-  coarse_grid.create_triangulation_compatibility(points, cells, SubCellData());
+  coarse_grid.create_triangulation(points, cells, SubCellData());
 }
 
 
diff --git a/tests/grid/mesh_3d.h b/tests/grid/mesh_3d.h
index a7d5d3b317..f4564b3c8a 100644
--- a/tests/grid/mesh_3d.h
+++ b/tests/grid/mesh_3d.h
@@ -139,9 +139,7 @@ create_two_cubes_rotation(Triangulation<3> & coarse_grid,
     }
   // finally generate a triangulation
   // out of this
-  coarse_grid.create_triangulation_compatibility(vertices,
-                                                 cells,
-                                                 SubCellData());
+  coarse_grid.create_triangulation(vertices, cells, SubCellData());
 }
 
 
diff --git a/tests/grid/subcelldata.cc b/tests/grid/subcelldata.cc
index 23dce0cc9f..cc26ed861e 100644
--- a/tests/grid/subcelldata.cc
+++ b/tests/grid/subcelldata.cc
@@ -23,6 +23,51 @@
 #include "../tests.h"
 
 
+namespace
+{
+  /**
+   * This file uses a different ordering for the vertices in a hex
+   * cell than we usually do in deal.II. The different convention used
+   * here originates in what we believed the ordering to be in UCD
+   * format, until it was discovered in 2022 that UCD will interpret
+   * this ordering to correspond to inverted cells -- as a
+   * consequence, the UCD ordering was fixed, but the current file is
+   * stuck on the old ordering.
+   */
+  constexpr std::array<unsigned int, 8> local_vertex_numbering{
+    {0, 1, 5, 4, 2, 3, 7, 6}};
+
+  /**
+   * And now also in the opposite direction.
+   */
+  void
+  reorder_old_to_new_style(std::vector<CellData<3>> &cells)
+  {
+    // undo the ordering above
+    unsigned int tmp[GeometryInfo<3>::vertices_per_cell];
+    for (auto &cell : cells)
+      {
+        for (const unsigned int i : GeometryInfo<3>::vertex_indices())
+          tmp[i] = cell.vertices[i];
+        for (const unsigned int i : GeometryInfo<3>::vertex_indices())
+          cell.vertices[local_vertex_numbering[i]] = tmp[i];
+      }
+  }
+
+  /**
+   * Following is a set of functions that reorder the data from the
+   * "current" to the "classic" format of vertex numbering of cells
+   * and faces. These functions do the reordering of their arguments
+   * in-place.
+   */
+  void
+  reorder_old_to_new_style(std::vector<CellData<2>> &cells)
+  {
+    for (auto &cell : cells)
+      std::swap(cell.vertices[2], cell.vertices[3]);
+  }
+} // namespace
+
 
 static unsigned subcells[6][4] = {{0, 1, 2, 3},
                                   {4, 5, 6, 7},
@@ -93,8 +138,9 @@ test()
         }
     }
 
+  reorder_old_to_new_style(cells);
   Triangulation<dim> tria;
-  tria.create_triangulation_compatibility(vertices, cells, subcelldata);
+  tria.create_triangulation(vertices, cells, subcelldata);
 
   GridOutFlags::Ucd ucd_flags(true, true);
   GridOut           grid_out;
diff --git a/tests/lac/constraints.cc b/tests/lac/constraints.cc
index f342a66bdf..c9324b5b20 100644
--- a/tests/lac/constraints.cc
+++ b/tests/lac/constraints.cc
@@ -37,6 +37,49 @@
 #include "../tests.h"
 
 
+namespace
+{
+  /**
+   * This file uses a different ordering for the vertices in a hex
+   * cell than we usually do in deal.II. The different convention used
+   * here originates in what we believed the ordering to be in UCD
+   * format, until it was discovered in 2022 that UCD will interpret
+   * this ordering to correspond to inverted cells -- as a
+   * consequence, the UCD ordering was fixed, but the current file is
+   * stuck on the old ordering.
+   */
+  constexpr std::array<unsigned int, 8> local_vertex_numbering{
+    {0, 1, 5, 4, 2, 3, 7, 6}};
+
+  /**
+   * Following is a set of functions that reorder the data from the
+   * "current" to the "classic" format of vertex numbering of cells
+   * and faces. These functions do the reordering of their arguments
+   * in-place.
+   */
+  void
+  reorder_old_to_new_style(std::vector<CellData<2>> &cells)
+  {
+    for (auto &cell : cells)
+      std::swap(cell.vertices[2], cell.vertices[3]);
+  }
+
+
+  void
+  reorder_old_to_new_style(std::vector<CellData<3>> &cells)
+  {
+    // undo the ordering above
+    unsigned int tmp[GeometryInfo<3>::vertices_per_cell];
+    for (auto &cell : cells)
+      {
+        for (const unsigned int i : GeometryInfo<3>::vertex_indices())
+          tmp[i] = cell.vertices[i];
+        for (const unsigned int i : GeometryInfo<3>::vertex_indices())
+          cell.vertices[local_vertex_numbering[i]] = tmp[i];
+      }
+  }
+} // namespace
+
 
 void
 make_tria(Triangulation<3> &tria, int step)
@@ -72,10 +115,11 @@ make_tria(Triangulation<3> &tria, int step)
           cells[0].material_id = 0;
           cells[1].material_id = 0;
 
-          tria.create_triangulation_compatibility(
-            std::vector<Point<3>>(&vertices[0], &vertices[12]),
-            cells,
-            SubCellData()); // no boundary information
+          reorder_old_to_new_style(cells);
+          tria.create_triangulation(std::vector<Point<3>>(&vertices[0],
+                                                          &vertices[12]),
+                                    cells,
+                                    SubCellData()); // no boundary information
 
           if (step == 0)
             tria.last_active()->set_refine_flag();
@@ -115,10 +159,11 @@ make_tria(Triangulation<3> &tria, int step)
           cells[0].material_id = 0;
           cells[1].material_id = 0;
 
-          tria.create_triangulation_compatibility(
-            std::vector<Point<3>>(&vertices[0], &vertices[12]),
-            cells,
-            SubCellData()); // no boundary information
+          reorder_old_to_new_style(cells);
+          tria.create_triangulation(std::vector<Point<3>>(&vertices[0],
+                                                          &vertices[12]),
+                                    cells,
+                                    SubCellData()); // no boundary information
 
           if (step == 2)
             tria.last_active()->set_refine_flag();
@@ -158,10 +203,11 @@ make_tria(Triangulation<3> &tria, int step)
           cells[0].material_id = 0;
           cells[1].material_id = 0;
 
-          tria.create_triangulation_compatibility(
-            std::vector<Point<3>>(&vertices[0], &vertices[12]),
-            cells,
-            SubCellData()); // no boundary information
+          reorder_old_to_new_style(cells);
+          tria.create_triangulation(std::vector<Point<3>>(&vertices[0],
+                                                          &vertices[12]),
+                                    cells,
+                                    SubCellData()); // no boundary information
 
           if (step == 4)
             tria.last_active()->set_refine_flag();
@@ -215,10 +261,11 @@ make_tria(Triangulation<3> &tria, int step)
           cells[2].material_id = 0;
           cells[3].material_id = 0;
 
-          tria.create_triangulation_compatibility(
-            std::vector<Point<3>>(&vertices[0], &vertices[18]),
-            cells,
-            SubCellData()); // no boundary information
+          reorder_old_to_new_style(cells);
+          tria.create_triangulation(std::vector<Point<3>>(&vertices[0],
+                                                          &vertices[18]),
+                                    cells,
+                                    SubCellData()); // no boundary information
 
           switch (step)
             {
-- 
2.39.5