]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Address issue #228 (https://code.google.com/p/dealii/issues/detail?id=228): Teach...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 15 Jul 2014 02:24:32 +0000 (02:24 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 15 Jul 2014 02:24:32 +0000 (02:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@33179 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/grid/grid_reordering.h
deal.II/source/grid/grid_generator.cc
deal.II/source/grid/grid_reordering.cc

index 41acf13069897292bebe5ddad02e8ebb6c4785e4..211a2a0595229f21ca8b0bee7ba463e27c168a13 100644 (file)
@@ -172,6 +172,17 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
+  <li> New: The GridReordering::reorder_cells() function used a
+  numbering format for the vertices in a cell that was last used in
+  deal.II version 5.2. This format is still used internally, but
+  the function now also understands the numbering that has been
+  used in deal.II ever since. The choice is made by an additional
+  argument to the function that defaults to the old-style
+  format for backward compatibility.
+  <br>
+  (Wolfgang Bangerth, 2014/07/14)
+  </li>
+
   <li> New: There are now functions GridOut::write_vtk() and
    GridOut::write_vtu() that can
   write a mesh in VTK/VTU format.
index 652a59ce98103184025debba465a4267ea6bdecd..e14d6097ff851c86c4500afc09eaccf1397e16eb 100644 (file)
@@ -34,9 +34,12 @@ DEAL_II_NAMESPACE_OPEN
  * is mainly used when reading in grids from files and converting them to
  * deal.II triangulations.
  *
- * Note: In contrast to the rest of the deal.II library this class
+ * @note In contrast to the rest of the deal.II library, by default this class
  * uses the old deal.II numbering scheme, which was used up to deal.II
- * version 5.2. That is, the vertex and face ordering in 2d is assumed
+ * version 5.2 (but the main function of this class takes a flag that
+ * specifies whether it should do an implicit conversion from the new
+ * to the old format before doing its work, and then back again after
+ * reordering). In this old format, the vertex and face ordering in 2d is assumed
  * to be
  * @verbatim
  *          2
@@ -73,7 +76,6 @@ DEAL_II_NAMESPACE_OPEN
  *     |/       /       |       |/
  *     *-------*        *-------*
  * @endverbatim
- *
  * After calling the GridReordering::reorder_cells() function the
  * CellData is still in this old numbering scheme. Hence, for creating
  * a Triangulation based on the resulting CellData the
@@ -495,9 +497,9 @@ DEAL_II_NAMESPACE_OPEN
  *
  * Prior to the implementation of the algorithms developed by Michael Anderson
  * and described above, we used a branch-and-cut algorithm initially
- * implemented in 2000 by Wolfgang Bangerth. Although it is no longer used
+ * implemented in 2000 by Wolfgang Bangerth. Although it is no longer used,
  * here is how it works, and why it doesn't always work for large meshes since
- * its run-time can exponential in bad cases.
+ * its run-time can be exponential in bad cases.
  *
  * The first observation is that although there are counterexamples,
  * problems are usually local. For example, in the second example
@@ -662,8 +664,17 @@ public:
    *  If a consistent reordering is not
    *  possible in dim=3, the original
    *  connectivity data is restored.
+   *
+   * @param original_cells An object that contains the data that describes
+   *   the mesh.
+   * @param use_new_style_ordering If true, then use the standard
+   *   ordering of 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.
    */
-  static void reorder_cells (std::vector<CellData<dim> > &original_cells);
+  static void reorder_cells (std::vector<CellData<dim> > &original_cells,
+                             const bool use_new_style_ordering = false);
 
   /**
    * Grids generated by grid
@@ -704,15 +715,18 @@ public:
 // declaration of explicit specializations
 template<>
 void
-GridReordering<2>::reorder_cells (std::vector<CellData<2> > &original_cells);
+GridReordering<2>::reorder_cells (std::vector<CellData<2> > &original_cells,
+                                  const bool);
 
 template<>
 void
-GridReordering<2,3>::reorder_cells (std::vector<CellData<2> > &original_cells);
+GridReordering<2,3>::reorder_cells (std::vector<CellData<2> > &original_cells,
+                                    const bool);
 
 template<>
 void
-GridReordering<3>::reorder_cells (std::vector<CellData<3> > &original_cells);
+GridReordering<3>::reorder_cells (std::vector<CellData<3> > &original_cells,
+                                  const bool);
 
 template<>
 void
index 5233e55acfc4718a37378c8f0aac8837f3f922ab..39fa5409cf4cd4b1e0f3208ba5718547336df063 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
-// anonymous namespace for internal helper functions
-namespace
-{
-  /**
-   * A set of three 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_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 (unsigned int cell=0; cell<cells.size(); ++cell)
-      std::swap(cells[cell].vertices[2],cells[cell].vertices[3]);
-  }
-
-
-  void
-  reorder_compatibility (std::vector<CellData<3> > &cells,
-                         SubCellData               &subcelldata)
-  {
-    unsigned int tmp[GeometryInfo<3>::vertices_per_cell];
-    for (unsigned int cell=0; cell<cells.size(); ++cell)
-      {
-        for (unsigned int i=0; i<GeometryInfo<3>::vertices_per_cell; ++i)
-          tmp[i] = cells[cell].vertices[i];
-        for (unsigned int i=0; i<GeometryInfo<3>::vertices_per_cell; ++i)
-          cells[cell].vertices[i] = tmp[GeometryInfo<3>::ucd_to_deal[i]];
-      }
-
-    // now points in boundary quads
-    std::vector<CellData<2> >::iterator boundary_quad
-      = subcelldata.boundary_quads.begin();
-    std::vector<CellData<2> >::iterator end_quad
-      = subcelldata.boundary_quads.end();
-    for (unsigned int quad_no=0; boundary_quad!=end_quad; ++boundary_quad, ++quad_no)
-      std::swap(boundary_quad->vertices[2], boundary_quad->vertices[3]);
-  }
-}
 
 namespace GridGenerator
 {
@@ -3480,12 +3428,12 @@ namespace GridGenerator
     GridTools::delete_duplicated_vertices (vertices, cells,
                                           subcell_data,
                                           considered_vertices);
-    // reorder_cells expects data in the "classic" format, transform
-    // the cell data before and after the call to that function 
-    reorder_compatibility(cells, subcell_data);
-    GridReordering<dim, spacedim>::reorder_cells(cells);
+
+    // reorder the cells to ensure that they satisfy the convention for
+    // edge and face directions
+    GridReordering<dim, spacedim>::reorder_cells(cells, true);
     result.clear ();
-    result.create_triangulation_compatibility (vertices, cells, subcell_data);
+    result.create_triangulation (vertices, cells, subcell_data);
   }
 
 
index 9b478ad08c953b7e6110b514a028ffdebc7b0968..07dcf6c947bbe1e9edf8bc4eb02e5e2201d0627c 100644 (file)
@@ -32,7 +32,8 @@ DEAL_II_NAMESPACE_OPEN
 
 template<>
 void
-GridReordering<1>::reorder_cells (std::vector<CellData<1> > &)
+GridReordering<1>::reorder_cells (std::vector<CellData<1> > &,
+                                  const bool)
 {
   // there should not be much to do
   // in 1d...
@@ -49,7 +50,8 @@ GridReordering<1>::invert_all_cells_of_negative_grid(const std::vector<Point<1>
 
 template<>
 void
-GridReordering<1,2>::reorder_cells (std::vector<CellData<1> > &)
+GridReordering<1,2>::reorder_cells (std::vector<CellData<1> > &,
+                                    const bool)
 {
   // there should not be much to do
   // in 1d...
@@ -67,7 +69,8 @@ GridReordering<1,2>::invert_all_cells_of_negative_grid(const std::vector<Point<2
 
 template<>
 void
-GridReordering<1,3>::reorder_cells (std::vector<CellData<1> > &)
+GridReordering<1,3>::reorder_cells (std::vector<CellData<1> > &,
+                                    const bool)
 {
   // there should not be much to do
   // in 1d...
@@ -614,28 +617,121 @@ namespace internal
 } // namespace internal
 
 
+// anonymous namespace for internal helper functions
+namespace
+{
+  /**
+   * A set of three 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_new_to_old_style (const std::vector<CellData<1> > &)
+  {
+    // nothing to do here: the format
+    // hasn't changed for 1d
+  }
+
+
+  void
+  reorder_new_to_old_style (std::vector<CellData<2> > &cells)
+  {
+    for (unsigned int cell=0; cell<cells.size(); ++cell)
+      std::swap(cells[cell].vertices[2], cells[cell].vertices[3]);
+  }
+
+
+  void
+  reorder_new_to_old_style (std::vector<CellData<3> > &cells)
+  {
+    unsigned int tmp[GeometryInfo<3>::vertices_per_cell];
+    for (unsigned int cell=0; cell<cells.size(); ++cell)
+      {
+        for (unsigned int i=0; i<GeometryInfo<3>::vertices_per_cell; ++i)
+          tmp[i] = cells[cell].vertices[i];
+        for (unsigned int i=0; i<GeometryInfo<3>::vertices_per_cell; ++i)
+          cells[cell].vertices[i] = tmp[GeometryInfo<3>::ucd_to_deal[i]];
+      }
+  }
+
+
+  /**
+   * And now also in the opposite direction.
+   */
+  void
+  reorder_old_to_new_style (const std::vector<CellData<1> > &)
+  {
+    // nothing to do here: the format
+    // hasn't changed for 1d
+  }
+
+
+  void
+  reorder_old_to_new_style (std::vector<CellData<2> > &cells)
+  {
+    // just invert the permutation:
+    reorder_new_to_old_style(cells);
+  }
+
+
+  void
+  reorder_old_to_new_style (std::vector<CellData<3> > &cells)
+  {
+    // undo the ordering above
+    unsigned int tmp[GeometryInfo<3>::vertices_per_cell];
+    for (unsigned int cell=0; cell<cells.size(); ++cell)
+      {
+        for (unsigned int i=0; i<GeometryInfo<3>::vertices_per_cell; ++i)
+          tmp[i] = cells[cell].vertices[i];
+        for (unsigned int i=0; i<GeometryInfo<3>::vertices_per_cell; ++i)
+          cells[cell].vertices[GeometryInfo<3>::ucd_to_deal[i]] = tmp[i];
+      }
+  }
+}
 
 
 template<>
 void
-GridReordering<2>::reorder_cells (std::vector<CellData<2> > &original_cells)
+GridReordering<2>::reorder_cells (std::vector<CellData<2> > &cells,
+                                  const bool use_new_style_ordering)
 {
+  // if necessary, convert to old-style format
+  if (use_new_style_ordering)
+    reorder_new_to_old_style(cells);
+
   // check if grids are already
   // consistent. if so, do
   // nothing. if not, then do the
   // reordering
-  if (internal::GridReordering2d::is_consistent (original_cells))
-    return;
+  if (!internal::GridReordering2d::is_consistent (cells))
+    internal::GridReordering2d::GridReordering().reorient(cells);
 
-  internal::GridReordering2d::GridReordering().reorient(original_cells);
+
+  // and convert back if necessary
+  if (use_new_style_ordering)
+    reorder_old_to_new_style(cells);
 }
 
 
 template<>
 void
-GridReordering<2,3>::reorder_cells (std::vector<CellData<2> > &original_cells)
+GridReordering<2,3>::reorder_cells (std::vector<CellData<2> > &cells,
+                                    const bool use_new_style_ordering)
 {
-  GridReordering<2>::reorder_cells(original_cells);
+  // if necessary, convert to old-style format
+  if (use_new_style_ordering)
+    reorder_new_to_old_style(cells);
+
+  GridReordering<2>::reorder_cells(cells);
+
+
+  // and convert back if necessary
+  if (use_new_style_ordering)
+    reorder_old_to_new_style(cells);
 }
 
 
@@ -1551,26 +1647,34 @@ namespace internal
 
 template<>
 void
-GridReordering<3>::reorder_cells (std::vector<CellData<3> > &incubes)
+GridReordering<3>::reorder_cells (std::vector<CellData<3> > &cells,
+                                  const bool use_new_style_ordering)
 {
+  Assert (cells.size() != 0,
+          ExcMessage("List of elements to orient must have at least one cell"));
 
-  Assert (incubes.size() != 0,
-          ExcMessage("List of elements to orient was of zero length"));
+  // if necessary, convert to old-style format
+  if (use_new_style_ordering)
+    reorder_new_to_old_style(cells);
 
   // create a backup to use if GridReordering
   // was not successful
-  std::vector<CellData<3> > backup=incubes;
+  std::vector<CellData<3> > backup=cells;
 
   // This does the real work
-  bool success=
-    internal::GridReordering3d::Orienter::orient_mesh (incubes);
+  const bool success=
+    internal::GridReordering3d::Orienter::orient_mesh (cells);
 
   // if reordering was not successful use
-  // original connectivity, otherwiese do
+  // original connectivity, otherwise do
   // nothing (i.e. use the reordered
   // connectivity)
   if (!success)
-    incubes=backup;
+    cells=backup;
+
+  // and convert back if necessary
+  if (use_new_style_ordering)
+    reorder_old_to_new_style(cells);
 }
 
 

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.