]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix confilct with Simplex::GridGenerator 10769/head
authorPeter Munch <peterrmuench@gmail.com>
Thu, 30 Jul 2020 17:57:41 +0000 (19:57 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Fri, 31 Jul 2020 12:23:49 +0000 (14:23 +0200)
include/deal.II/base/exceptions.h
include/deal.II/simplex/grid_generator.h
source/simplex/CMakeLists.txt
source/simplex/grid_generator.cc [new file with mode: 0644]
source/simplex/grid_generator.inst.in [new file with mode: 0644]

index b244dd3b70e293850f8c5b474e726abc021d3a17..171411eb079852f2017542e695b28f07931e2e6f 100644 (file)
@@ -1091,6 +1091,14 @@ namespace StandardExceptions
     "You are attempting to use functionality that is only available "
     "if deal.II was configured to use MPI.");
 
+  /**
+   * This function requires simplex support.
+   */
+  DeclExceptionMsg(
+    ExcNeedsSimplexSupport,
+    "You are attempting to use functionality that is only available "
+    "if deal.II was configured with DEAL_II_WITH_SIMPLEX_SUPPORT enabled.");
+
   /**
    * This function requires support for the FunctionParser library.
    */
index de0cf2d7db896f08e77f2984482927b41a51b053..2e68f9e1bd485428658d8d04e0254c25ce85bd03 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
-namespace Simplex
+#ifdef DEAL_II_WITH_SIMPLEX_SUPPORT
+
+namespace GridGenerator
 {
   /**
-   * This namespace provides a collection of functions to generate simplex
-   * triangulations for some basic geometries.
+   * Create a coordinate-parallel brick from the two diagonally opposite
+   * corner points @p p1 and @p p2. The number of vertices in coordinate
+   * direction @p i is given by <tt>repetitions[i]+1</tt>.
+   *
+   * @note This function connects internally 4/8 vertices to
+   *   quadrilateral/hexahedral cells and subdivides these into 2/5
+   * triangular/tetrahedral cells.
+   *
+   * @note Currently, this function only works for `dim==spacedim`.
    *
-   *  @ingroup simplex
+   * @ingroup simplex
    */
-  namespace GridGenerator
-  {
-    /**
-     * Create a coordinate-parallel brick from the two diagonally opposite
-     * corner points @p p1 and @p p2. The number of vertices in coordinate
-     * direction @p i is given by <tt>repetitions[i]+1</tt>.
-     *
-     * @note This function connects internally 4/8 vertices to quadrilateral/
-     *   hexahedral cells and subdivides these into 2/5 triangular/
-     *   tetrahedral cells.
-     *
-     * @note Currently, this function only works for `dim==spacedim`.
-     */
-    template <int dim, int spacedim>
-    void
-    subdivided_hyper_rectangle(Triangulation<dim, spacedim> &   tria,
-                               const std::vector<unsigned int> &repetitions,
-                               const Point<dim> &               p1,
-                               const Point<dim> &               p2,
-                               const bool colorize = false)
-    {
-      AssertDimension(dim, spacedim);
-
-      AssertThrow(colorize == false, ExcNotImplemented());
-
-      std::vector<Point<spacedim>> vertices;
-      std::vector<CellData<dim>>   cells;
-
-      if (dim == 2)
-        {
-          // determine cell sizes
-          const Point<dim> dx((p2[0] - p1[0]) / repetitions[0],
-                              (p2[1] - p1[1]) / repetitions[1]);
-
-          // create vertices
-          for (unsigned int j = 0; j <= repetitions[1]; ++j)
-            for (unsigned int i = 0; i <= repetitions[0]; ++i)
-              vertices.push_back(
-                Point<spacedim>(p1[0] + dx[0] * i, p1[1] + dx[1] * j));
-
-          // create cells
-          for (unsigned int j = 0; j < repetitions[1]; ++j)
-            for (unsigned int i = 0; i < repetitions[0]; ++i)
-              {
-                // create reference QUAD cell
-                std::array<unsigned int, 4> quad{
-                  (j + 0) * (repetitions[0] + 1) + i + 0, //
-                  (j + 0) * (repetitions[0] + 1) + i + 1, //
-                  (j + 1) * (repetitions[0] + 1) + i + 0, //
-                  (j + 1) * (repetitions[0] + 1) + i + 1  //
-                };                                        //
-
-                // TRI cell 0
-                {
-                  CellData<dim> tri;
-                  tri.vertices = {quad[0], quad[1], quad[2]};
-                  cells.push_back(tri);
-                }
-
-                // TRI cell 1
-                {
-                  CellData<dim> tri;
-                  tri.vertices = {quad[3], quad[2], quad[1]};
-                  cells.push_back(tri);
-                }
-              }
-        }
-      else if (dim == 3)
-        {
-          // determine cell sizes
-          const Point<dim> dx((p2[0] - p1[0]) / repetitions[0],
-                              (p2[1] - p1[1]) / repetitions[1],
-                              (p2[2] - p1[2]) / repetitions[2]);
-
-          // create vertices
-          for (unsigned int k = 0; k <= repetitions[2]; ++k)
-            for (unsigned int j = 0; j <= repetitions[1]; ++j)
-              for (unsigned int i = 0; i <= repetitions[0]; ++i)
-                vertices.push_back(Point<spacedim>(p1[0] + dx[0] * i,
-                                                   p1[1] + dx[1] * j,
-                                                   p1[2] + dx[2] * k));
-
-          // create cells
-          for (unsigned int k = 0; k < repetitions[2]; ++k)
-            for (unsigned int j = 0; j < repetitions[1]; ++j)
-              for (unsigned int i = 0; i < repetitions[0]; ++i)
-                {
-                  // create reference HEX cell
-                  std::array<unsigned int, 8> quad{
-                    (k + 0) * (repetitions[0] + 1) * (repetitions[1] + 1) +
-                      (j + 0) * (repetitions[0] + 1) + i + 0,
-                    (k + 0) * (repetitions[0] + 1) * (repetitions[1] + 1) +
-                      (j + 0) * (repetitions[0] + 1) + i + 1,
-                    (k + 0) * (repetitions[0] + 1) * (repetitions[1] + 1) +
-                      (j + 1) * (repetitions[0] + 1) + i + 0,
-                    (k + 0) * (repetitions[0] + 1) * (repetitions[1] + 1) +
-                      (j + 1) * (repetitions[0] + 1) + i + 1,
-                    (k + 1) * (repetitions[0] + 1) * (repetitions[1] + 1) +
-                      (j + 0) * (repetitions[0] + 1) + i + 0,
-                    (k + 1) * (repetitions[0] + 1) * (repetitions[1] + 1) +
-                      (j + 0) * (repetitions[0] + 1) + i + 1,
-                    (k + 1) * (repetitions[0] + 1) * (repetitions[1] + 1) +
-                      (j + 1) * (repetitions[0] + 1) + i + 0,
-                    (k + 1) * (repetitions[0] + 1) * (repetitions[1] + 1) +
-                      (j + 1) * (repetitions[0] + 1) + i + 1};
+  template <int dim, int spacedim>
+  void
+  subdivided_hyper_rectangle_with_simplices(
+    Triangulation<dim, spacedim> &   tria,
+    const std::vector<unsigned int> &repetitions,
+    const Point<dim> &               p1,
+    const Point<dim> &               p2,
+    const bool                       colorize = false);
 
-                  // TET cell 0
-                  {
-                    CellData<dim> cell;
-                    if (((i % 2) + (j % 2) + (k % 2)) % 2 == 0)
-                      cell.vertices = {quad[0], quad[1], quad[2], quad[4]};
-                    else
-                      cell.vertices = {quad[0], quad[1], quad[3], quad[5]};
-
-                    cells.push_back(cell);
-                  }
-
-                  // TET cell 1
-                  {
-                    CellData<dim> cell;
-                    if (((i % 2) + (j % 2) + (k % 2)) % 2 == 0)
-                      cell.vertices = {quad[2], quad[1], quad[3], quad[7]};
-                    else
-                      cell.vertices = {quad[0], quad[3], quad[2], quad[6]};
-                    cells.push_back(cell);
-                  }
-
-                  // TET cell 2
-                  {
-                    CellData<dim> cell;
-                    if (((i % 2) + (j % 2) + (k % 2)) % 2 == 0)
-                      cell.vertices = {quad[1], quad[4], quad[5], quad[7]};
-                    else
-                      cell.vertices = {quad[0], quad[4], quad[5], quad[6]};
-                    cells.push_back(cell);
-                  }
-
-                  // TET cell 3
-                  {
-                    CellData<dim> cell;
-                    if (((i % 2) + (j % 2) + (k % 2)) % 2 == 0)
-                      cell.vertices = {quad[2], quad[4], quad[7], quad[6]};
-                    else
-                      cell.vertices = {quad[3], quad[5], quad[7], quad[6]};
-                    cells.push_back(cell);
-                  }
-
-                  // TET cell 4
-                  {
-                    CellData<dim> cell;
-                    if (((i % 2) + (j % 2) + (k % 2)) % 2 == 0)
-                      cell.vertices = {quad[1], quad[2], quad[4], quad[7]};
-                    else
-                      cell.vertices = {quad[0], quad[3], quad[6], quad[5]};
-                    cells.push_back(cell);
-                  }
-                }
-        }
-      else
-        {
-          AssertThrow(colorize == false, ExcNotImplemented());
-        }
-
-      // actually create triangulation
-      tria.create_triangulation(vertices, cells, SubCellData());
-    }
-
-    /**
-     * Initialize the given triangulation with a hypercube (square in 2D and
-     * cube in 3D) consisting of @p repetitions cells in each direction.
-     * The hypercube volume is the tensor product interval
-     * $[left,right]^{\text{dim}}$ in the present number of dimensions, where
-     * the limits are given as arguments. They default to zero and unity, then
-     * producing the unit hypercube.
-     *
-     * @note This function connects internally 4/8 vertices to quadrilateral/
-     *   hexahedral cells and subdivides these into 2/5 triangular/
-     *   tetrahedral cells.
-     */
-    template <int dim, int spacedim>
-    void
-    subdivided_hyper_cube(Triangulation<dim, spacedim> &tria,
-                          const unsigned int            repetitions,
-                          const double                  p1       = 0.0,
-                          const double                  p2       = 1.0,
-                          const bool                    colorize = false)
-    {
-      if (dim == 2)
-        {
-          subdivided_hyper_rectangle(
-            tria, {repetitions, repetitions}, {p1, p1}, {p2, p2}, colorize);
-        }
-      else if (dim == 3)
-        {
-          subdivided_hyper_rectangle(tria,
-                                     {repetitions, repetitions, repetitions},
-                                     {p1, p1, p1},
-                                     {p2, p2, p2},
-                                     colorize);
-        }
-      else
-        {
-          AssertThrow(false, ExcNotImplemented())
-        }
-    }
-  } // namespace GridGenerator
-} // namespace Simplex
+  /**
+   * Initialize the given triangulation with a hypercube (square in 2D and
+   * cube in 3D) consisting of @p repetitions cells in each direction.
+   * The hypercube volume is the tensor product interval
+   * $[left,right]^{\text{dim}}$ in the present number of dimensions, where
+   * the limits are given as arguments. They default to zero and unity, then
+   * producing the unit hypercube.
+   *
+   * @note This function connects internally 4/8 vertices to
+   * quadrilateral/hexahedral cells and subdivides these into 2/5
+   * triangular/tetrahedral cells.
+   *
+   * @ingroup simplex
+   */
+  template <int dim, int spacedim>
+  void
+  subdivided_hyper_cube_with_simplices(Triangulation<dim, spacedim> &tria,
+                                       const unsigned int repetitions,
+                                       const double       p1       = 0.0,
+                                       const double       p2       = 1.0,
+                                       const bool         colorize = false);
 
+} // namespace GridGenerator
 
+#endif
 
 DEAL_II_NAMESPACE_CLOSE
 
index 604ce1d1f2d29321c69fbbc345fb14e8eae6cd03..f547691c5aae70c170468633a559cc4c9ab47cc2 100644 (file)
@@ -17,6 +17,7 @@ INCLUDE_DIRECTORIES(BEFORE ${CMAKE_CURRENT_BINARY_DIR})
 
 SET(_unity_include_src
   fe_lib.cc
+  grid_generator.cc
   polynomials.cc
   quadrature_lib.cc
   )
@@ -35,6 +36,7 @@ SETUP_SOURCE_LIST("${_unity_include_src}"
 
 SET(_inst
   fe_lib.inst.in
+  grid_generator.inst.in
 )
 
 FILE(GLOB _header
diff --git a/source/simplex/grid_generator.cc b/source/simplex/grid_generator.cc
new file mode 100644 (file)
index 0000000..e17e69e
--- /dev/null
@@ -0,0 +1,215 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 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.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/simplex/grid_generator.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+#ifdef DEAL_II_WITH_SIMPLEX_SUPPORT
+
+namespace GridGenerator
+{
+  template <int dim, int spacedim>
+  void
+  subdivided_hyper_rectangle_with_simplices(
+    Triangulation<dim, spacedim> &   tria,
+    const std::vector<unsigned int> &repetitions,
+    const Point<dim> &               p1,
+    const Point<dim> &               p2,
+    const bool                       colorize)
+  {
+    AssertDimension(dim, spacedim);
+
+    AssertThrow(colorize == false, ExcNotImplemented());
+
+    std::vector<Point<spacedim>> vertices;
+    std::vector<CellData<dim>>   cells;
+
+    if (dim == 2)
+      {
+        // determine cell sizes
+        const Point<dim> dx((p2[0] - p1[0]) / repetitions[0],
+                            (p2[1] - p1[1]) / repetitions[1]);
+
+        // create vertices
+        for (unsigned int j = 0; j <= repetitions[1]; ++j)
+          for (unsigned int i = 0; i <= repetitions[0]; ++i)
+            vertices.push_back(
+              Point<spacedim>(p1[0] + dx[0] * i, p1[1] + dx[1] * j));
+
+        // create cells
+        for (unsigned int j = 0; j < repetitions[1]; ++j)
+          for (unsigned int i = 0; i < repetitions[0]; ++i)
+            {
+              // create reference QUAD cell
+              std::array<unsigned int, 4> quad{{
+                (j + 0) * (repetitions[0] + 1) + i + 0, //
+                (j + 0) * (repetitions[0] + 1) + i + 1, //
+                (j + 1) * (repetitions[0] + 1) + i + 0, //
+                (j + 1) * (repetitions[0] + 1) + i + 1  //
+              }};                                       //
+
+              // TRI cell 0
+              {
+                CellData<dim> tri;
+                tri.vertices = {quad[0], quad[1], quad[2]};
+                cells.push_back(tri);
+              }
+
+              // TRI cell 1
+              {
+                CellData<dim> tri;
+                tri.vertices = {quad[3], quad[2], quad[1]};
+                cells.push_back(tri);
+              }
+            }
+      }
+    else if (dim == 3)
+      {
+        // determine cell sizes
+        const Point<dim> dx((p2[0] - p1[0]) / repetitions[0],
+                            (p2[1] - p1[1]) / repetitions[1],
+                            (p2[2] - p1[2]) / repetitions[2]);
+
+        // create vertices
+        for (unsigned int k = 0; k <= repetitions[2]; ++k)
+          for (unsigned int j = 0; j <= repetitions[1]; ++j)
+            for (unsigned int i = 0; i <= repetitions[0]; ++i)
+              vertices.push_back(Point<spacedim>(p1[0] + dx[0] * i,
+                                                 p1[1] + dx[1] * j,
+                                                 p1[2] + dx[2] * k));
+
+        // create cells
+        for (unsigned int k = 0; k < repetitions[2]; ++k)
+          for (unsigned int j = 0; j < repetitions[1]; ++j)
+            for (unsigned int i = 0; i < repetitions[0]; ++i)
+              {
+                // create reference HEX cell
+                std::array<unsigned int, 8> quad{
+                  {(k + 0) * (repetitions[0] + 1) * (repetitions[1] + 1) +
+                     (j + 0) * (repetitions[0] + 1) + i + 0,
+                   (k + 0) * (repetitions[0] + 1) * (repetitions[1] + 1) +
+                     (j + 0) * (repetitions[0] + 1) + i + 1,
+                   (k + 0) * (repetitions[0] + 1) * (repetitions[1] + 1) +
+                     (j + 1) * (repetitions[0] + 1) + i + 0,
+                   (k + 0) * (repetitions[0] + 1) * (repetitions[1] + 1) +
+                     (j + 1) * (repetitions[0] + 1) + i + 1,
+                   (k + 1) * (repetitions[0] + 1) * (repetitions[1] + 1) +
+                     (j + 0) * (repetitions[0] + 1) + i + 0,
+                   (k + 1) * (repetitions[0] + 1) * (repetitions[1] + 1) +
+                     (j + 0) * (repetitions[0] + 1) + i + 1,
+                   (k + 1) * (repetitions[0] + 1) * (repetitions[1] + 1) +
+                     (j + 1) * (repetitions[0] + 1) + i + 0,
+                   (k + 1) * (repetitions[0] + 1) * (repetitions[1] + 1) +
+                     (j + 1) * (repetitions[0] + 1) + i + 1}};
+
+                // TET cell 0
+                {
+                  CellData<dim> cell;
+                  if (((i % 2) + (j % 2) + (k % 2)) % 2 == 0)
+                    cell.vertices = {{quad[0], quad[1], quad[2], quad[4]}};
+                  else
+                    cell.vertices = {{quad[0], quad[1], quad[3], quad[5]}};
+
+                  cells.push_back(cell);
+                }
+
+                // TET cell 1
+                {
+                  CellData<dim> cell;
+                  if (((i % 2) + (j % 2) + (k % 2)) % 2 == 0)
+                    cell.vertices = {{quad[2], quad[1], quad[3], quad[7]}};
+                  else
+                    cell.vertices = {{quad[0], quad[3], quad[2], quad[6]}};
+                  cells.push_back(cell);
+                }
+
+                // TET cell 2
+                {
+                  CellData<dim> cell;
+                  if (((i % 2) + (j % 2) + (k % 2)) % 2 == 0)
+                    cell.vertices = {{quad[1], quad[4], quad[5], quad[7]}};
+                  else
+                    cell.vertices = {{quad[0], quad[4], quad[5], quad[6]}};
+                  cells.push_back(cell);
+                }
+
+                // TET cell 3
+                {
+                  CellData<dim> cell;
+                  if (((i % 2) + (j % 2) + (k % 2)) % 2 == 0)
+                    cell.vertices = {{quad[2], quad[4], quad[7], quad[6]}};
+                  else
+                    cell.vertices = {{quad[3], quad[5], quad[7], quad[6]}};
+                  cells.push_back(cell);
+                }
+
+                // TET cell 4
+                {
+                  CellData<dim> cell;
+                  if (((i % 2) + (j % 2) + (k % 2)) % 2 == 0)
+                    cell.vertices = {{quad[1], quad[2], quad[4], quad[7]}};
+                  else
+                    cell.vertices = {{quad[0], quad[3], quad[6], quad[5]}};
+                  cells.push_back(cell);
+                }
+              }
+      }
+    else
+      {
+        AssertThrow(colorize == false, ExcNotImplemented());
+      }
+
+    // actually create triangulation
+    tria.create_triangulation(vertices, cells, SubCellData());
+  }
+
+  template <int dim, int spacedim>
+  void
+  subdivided_hyper_cube_with_simplices(Triangulation<dim, spacedim> &tria,
+                                       const unsigned int repetitions,
+                                       const double       p1,
+                                       const double       p2,
+                                       const bool         colorize)
+  {
+    if (dim == 2)
+      {
+        subdivided_hyper_rectangle_with_simplices(
+          tria, {{repetitions, repetitions}}, {p1, p1}, {p2, p2}, colorize);
+      }
+    else if (dim == 3)
+      {
+        subdivided_hyper_rectangle_with_simplices(
+          tria,
+          {{repetitions, repetitions, repetitions}},
+          {p1, p1, p1},
+          {p2, p2, p2},
+          colorize);
+      }
+    else
+      {
+        AssertThrow(false, ExcNotImplemented())
+      }
+  }
+} // namespace GridGenerator
+
+#endif
+
+// explicit instantiations
+#include "grid_generator.inst"
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/source/simplex/grid_generator.inst.in b/source/simplex/grid_generator.inst.in
new file mode 100644 (file)
index 0000000..3c855b6
--- /dev/null
@@ -0,0 +1,37 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+#ifdef DEAL_II_WITH_SIMPLEX_SUPPORT
+
+for (deal_II_dimension : DIMENSIONS)
+  {
+    template void GridGenerator::subdivided_hyper_rectangle_with_simplices(
+      Triangulation<deal_II_dimension> & tria,
+      const std::vector<unsigned int> &repetitions,
+      const Point<deal_II_dimension> & p1,
+      const Point<deal_II_dimension> & p2,
+      const bool                       colorize);
+
+    template void GridGenerator::subdivided_hyper_cube_with_simplices(
+      Triangulation<deal_II_dimension> & tria,
+      const unsigned int repetitions,
+      const double       p1,
+      const double       p2,
+      const bool         colorize);
+  }
+
+#endif

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.