]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement GridGenerator::replicate_triangulation(). 9388/head
authorDavid Wells <drwells@email.unc.edu>
Tue, 14 Jan 2020 04:21:45 +0000 (23:21 -0500)
committerDavid Wells <drwells@email.unc.edu>
Sun, 26 Jan 2020 18:42:45 +0000 (13:42 -0500)
doc/doxygen/images/replicated_tria_2d.png [new file with mode: 0644]
doc/doxygen/images/replicated_tria_3d.png [new file with mode: 0644]
doc/news/changes/minor/20200113DavidWellsVictorZheng [new file with mode: 0644]
include/deal.II/grid/grid_generator.h
source/grid/grid_generator.cc
source/grid/grid_generator.inst.in
tests/grid/replicate_triangulation_01.cc [new file with mode: 0644]
tests/grid/replicate_triangulation_01.output [new file with mode: 0644]
tests/grid/replicate_triangulation_02.cc [new file with mode: 0644]
tests/grid/replicate_triangulation_02.output [new file with mode: 0644]

diff --git a/doc/doxygen/images/replicated_tria_2d.png b/doc/doxygen/images/replicated_tria_2d.png
new file mode 100644 (file)
index 0000000..a80aeab
Binary files /dev/null and b/doc/doxygen/images/replicated_tria_2d.png differ
diff --git a/doc/doxygen/images/replicated_tria_3d.png b/doc/doxygen/images/replicated_tria_3d.png
new file mode 100644 (file)
index 0000000..e185620
Binary files /dev/null and b/doc/doxygen/images/replicated_tria_3d.png differ
diff --git a/doc/news/changes/minor/20200113DavidWellsVictorZheng b/doc/news/changes/minor/20200113DavidWellsVictorZheng
new file mode 100644 (file)
index 0000000..4ec5392
--- /dev/null
@@ -0,0 +1,5 @@
+New: Added a new function GridGenerator::replicate_triangulation() for creating
+a new triangulation by copying a given triangulation repeatedly along the
+coordinate axes.
+<br>
+(David Wells, Victor Zheng, 2020/01/13)
index 8126aae9f56f04f2fd10b815ae4bc94ee1d81437..21b0cd408719292f4d79984ed571a0f91fce609d 100644 (file)
@@ -1551,6 +1551,59 @@ namespace GridGenerator
     const double                  duplicated_vertex_tolerance = 1.0e-12,
     const bool                    copy_manifold_ids           = false);
 
+  /**
+   * \brief Replicate a given triangulation in multiple coordinate axes
+   *
+   * @param input The triangulation which will be replicated along the
+   * coordinate axes.
+   *
+   * @param extents A vector with <tt>dim</tt> entries specifying how many
+   * copies of a triangulation should be present along each coordinate axis.
+   *
+   * @param result The triangulation to be created. It needs to be empty upon
+   * calling this function.
+   *
+   * This function creates a new Triangulation equal to a
+   * <tt>dim</tt>-dimensional array of copies of @p input. Copies of @p input
+   * are created by translating @p input along the coordinate axes. Boundary
+   * ids of faces (but not lines in 3D) and all manifold ids are copied but
+   * Manifold objects are not since most Manifold objects do not work
+   * correctly when a Triangulation has been translated.
+   *
+   * To see how this works, consider the following code:
+   * @code
+   * Triangulation<2> input;
+   * GridGenerator::hyper_cube_with_cylindrical_hole(input);
+   * Triangulation<2> output;
+   * GridGenerator::replicate_triangulation(input, {3, 2}, output);
+   * @endcode
+   * results in
+   *
+   * @image html replicated_tria_2d.png
+   *
+   * And, similarly, in 3D:
+   * @code
+   * Triangulation<3> input;
+   * GridGenerator::hyper_cross(1, 1, 1, 2, 1, 2);
+   * Triangulation<3> output;
+   * GridGenerator::replicate_triangulation(input, {3, 2, 1}, output);
+   * @endcode
+   * results in
+   *
+   * @image html replicated_tria_3d.png
+   *
+   * @note This function determines the spacing of the copies of @p input
+   * based on the BoundingBox of @p input. If the boundary faces of @p input
+   * are not aligned with the coordinate axes then the copies might not share
+   * common faces; i.e., this function is intended for simple geometries with
+   * boundary faces aligned along the coordinate axes.
+   */
+  template <int dim, int spacedim = dim>
+  void
+  replicate_triangulation(const Triangulation<dim, spacedim> &input,
+                          const std::vector<unsigned int> &   extents,
+                          Triangulation<dim, spacedim> &      result);
+
   /**
    * Given the two triangulations specified as the first two arguments, create
    * the triangulation that contains the finest cells of both triangulation
index a09b1c086ca849cd0baa391d6ef62f814a315a82..9e19f3f72405b1b069adeae8b86738466286d456 100644 (file)
@@ -5845,6 +5845,232 @@ namespace GridGenerator
 
 
 
+  namespace
+  {
+    /**
+     * Merging or replicating triangulations usually results in duplicated
+     * boundary objects - in particular, faces that used to be boundary faces
+     * will now be internal faces.
+     *
+     * This function modifies @p subcell_data by detecting duplicated objects,
+     * marking said duplicates as internal faces, and then deleting all but
+     * one of the duplicates.
+     *
+     * This function relies on some implementation details of
+     * create_triangulation to uniquify objects in SubCellData - in
+     * particular, quadrilaterals are only identified by their lines and not
+     * by their orientation or volume, so rotating and flipping a
+     * quadrilateral doesn't effect the way said quadrilateral is read by that
+     * function.
+     *
+     * @warning even though this function is implemented for structdim 1 and
+     * structdim 2, it will produce <em>wrong</em> results when called for
+     * boundary lines in 3D in most cases since a boundary line can be shared
+     * by an arbitrary number of cells in 3D.
+     */
+    template <int structdim>
+    void
+    delete_duplicated_objects(std::vector<CellData<structdim>> &subcell_data)
+    {
+      static_assert(structdim == 1 || structdim == 2,
+                    "This function is only implemented for lines and "
+                    "quadrilaterals.");
+      // start by making sure that all objects representing the same vertices
+      // are numbered in the same way by canonicalizing the numberings. This
+      // makes it possible to detect duplicates.
+      for (CellData<structdim> &cell_data : subcell_data)
+        {
+          if (structdim == 1)
+            std::sort(std::begin(cell_data.vertices),
+                      std::end(cell_data.vertices));
+          else if (structdim == 2)
+            {
+              // rotate the vertex numbers so that the lowest one is first
+              std::array<unsigned int, 4> renumbering;
+              std::copy(std::begin(cell_data.vertices),
+                        std::end(cell_data.vertices),
+                        renumbering.begin());
+
+              // convert to old style vertex numbering. This makes the
+              // permutations easy since the valid configurations are
+              //
+              // 3  2   2  1   1  0   0  3
+              // 0  1   3  0   2  3   1  2
+              // (0123) (3012) (2310) (1230)
+              //
+              // rather than the lexical ordering which is harder to permute
+              // by rotation.
+              std::swap(renumbering[2], renumbering[3]);
+              std::rotate(renumbering.begin(),
+                          std::min_element(renumbering.begin(),
+                                           renumbering.end()),
+                          renumbering.end());
+              // convert to new style
+              std::swap(renumbering[2], renumbering[3]);
+              // deal with cases where we might have
+              //
+              // 3 2   1 2
+              // 0 1   0 3
+              //
+              // by forcing the second vertex (in lexical ordering) to be
+              // smaller than the third
+              if (renumbering[1] > renumbering[2])
+                std::swap(renumbering[1], renumbering[2]);
+              std::copy(renumbering.begin(),
+                        renumbering.end(),
+                        std::begin(cell_data.vertices));
+            }
+        }
+
+      // Now that all cell objects have been canonicalized they can be sorted:
+      auto compare = [](const CellData<structdim> &a,
+                        const CellData<structdim> &b) {
+        return std::lexicographical_compare(std::begin(a.vertices),
+                                            std::end(a.vertices),
+                                            std::begin(b.vertices),
+                                            std::end(b.vertices));
+      };
+      std::sort(subcell_data.begin(), subcell_data.end(), compare);
+
+      // Finally, determine which objects are duplicates. Duplicates are
+      // assumed to be interior objects, so delete all but one and change the
+      // boundary id:
+      auto left = subcell_data.begin();
+      while (left != subcell_data.end())
+        {
+          const auto right =
+            std::upper_bound(left, subcell_data.end(), *left, compare);
+          // if the range has more than one item, then there are duplicates -
+          // set all boundary ids in the range to the internal boundary id
+          if (left + 1 != right)
+            for (auto it = left; it != right; ++it)
+              {
+                it->boundary_id = numbers::internal_face_boundary_id;
+                Assert(it->manifold_id == left->manifold_id,
+                       ExcMessage(
+                         "In the process of grid generation a single "
+                         "line or quadrilateral has been assigned two "
+                         "different manifold ids. This can happen when "
+                         "a Triangulation is copied, e.g., via "
+                         "GridGenerator::replicate_triangulation() and "
+                         "not all external boundary faces have the same "
+                         "manifold id. Double check that all faces "
+                         "which you expect to be merged together have "
+                         "the same manifold id."));
+              }
+          left = right;
+        }
+
+      subcell_data.erase(std::unique(subcell_data.begin(), subcell_data.end()),
+                         subcell_data.end());
+    }
+  } // namespace
+
+
+
+  template <int dim, int spacedim>
+  void
+  replicate_triangulation(const Triangulation<dim, spacedim> &input,
+                          const std::vector<unsigned int> &   extents,
+                          Triangulation<dim, spacedim> &      result)
+  {
+    AssertDimension(dim, extents.size());
+#ifdef DEBUG
+    for (const auto &extent : extents)
+      Assert(0 < extent,
+             ExcMessage("The Triangulation must be copied at least one time in "
+                        "each coordinate dimension."));
+#endif
+    const BoundingBox<spacedim> bbox(input.get_vertices());
+    const auto &                min = bbox.get_boundary_points().first;
+    const auto &                max = bbox.get_boundary_points().second;
+
+    std::array<Tensor<1, spacedim>, dim> offsets;
+    for (unsigned int d = 0; d < dim; ++d)
+      offsets[d][d] = max[d] - min[d];
+
+    Triangulation<dim, spacedim> tria_to_replicate;
+    tria_to_replicate.copy_triangulation(input);
+    for (unsigned int d = 0; d < dim; ++d)
+      {
+        std::vector<Point<spacedim>> input_vertices;
+        std::vector<CellData<dim>>   input_cell_data;
+        SubCellData                  input_subcell_data;
+        std::tie(input_vertices, input_cell_data, input_subcell_data) =
+          GridTools::get_coarse_mesh_description(tria_to_replicate);
+        std::vector<Point<spacedim>> output_vertices     = input_vertices;
+        std::vector<CellData<dim>>   output_cell_data    = input_cell_data;
+        SubCellData                  output_subcell_data = input_subcell_data;
+
+        for (unsigned int k = 1; k < extents[d]; ++k)
+          {
+            const std::size_t vertex_offset = k * input_vertices.size();
+            // vertices
+            for (const Point<spacedim> &point : input_vertices)
+              output_vertices.push_back(point + double(k) * offsets[d]);
+            // cell data
+            for (const CellData<dim> &cell_data : input_cell_data)
+              {
+                output_cell_data.push_back(cell_data);
+                for (unsigned int &vertex : output_cell_data.back().vertices)
+                  vertex += vertex_offset;
+              }
+            // subcell data
+            for (const CellData<1> &boundary_line :
+                 input_subcell_data.boundary_lines)
+              {
+                output_subcell_data.boundary_lines.push_back(boundary_line);
+                for (unsigned int &vertex :
+                     output_subcell_data.boundary_lines.back().vertices)
+                  vertex += vertex_offset;
+              }
+            for (const CellData<2> &boundary_quad :
+                 input_subcell_data.boundary_quads)
+              {
+                output_subcell_data.boundary_quads.push_back(boundary_quad);
+                for (unsigned int &vertex :
+                     output_subcell_data.boundary_quads.back().vertices)
+                  vertex += vertex_offset;
+              }
+          }
+        // check all vertices: since the grid is coarse, most will be on the
+        // boundary anyway
+        std::vector<unsigned int> boundary_vertices;
+        GridTools::delete_duplicated_vertices(
+          output_vertices,
+          output_cell_data,
+          output_subcell_data,
+          boundary_vertices,
+          1e-6 * input.begin_active()->diameter());
+        // delete_duplicated_vertices also deletes any unused vertices
+        // deal with any reordering issues created by delete_duplicated_vertices
+        GridReordering<dim>::reorder_cells(output_cell_data, true);
+        // clean up the boundary ids of the boundary objects: note that we
+        // have to do this after delete_duplicated_vertices so that boundary
+        // objects are actually duplicated at this point
+        if (dim == 2)
+          delete_duplicated_objects(output_subcell_data.boundary_lines);
+        else if (dim == 3)
+          {
+            delete_duplicated_objects(output_subcell_data.boundary_quads);
+            for (CellData<1> &boundary_line :
+                 output_subcell_data.boundary_lines)
+              // set boundary lines to the default value - let
+              // create_triangulation figure out the rest.
+              boundary_line.boundary_id = numbers::internal_face_boundary_id;
+          }
+
+        tria_to_replicate.clear();
+        tria_to_replicate.create_triangulation(output_vertices,
+                                               output_cell_data,
+                                               output_subcell_data);
+      }
+
+    result.copy_triangulation(tria_to_replicate);
+  }
+
+
+
   template <int dim, int spacedim>
   void
   create_union_triangulation(
index 76522efd9f1698d0be5cba6a740a8ce828480002..0e6c15512d60ff0079cce8a3caa84a95dca3a0d1 100644 (file)
@@ -97,6 +97,12 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
         const double duplicated_vertex_tolerance,
         const bool   copy_manifold_ids);
 
+      template void
+      replicate_triangulation(
+        const Triangulation<deal_II_dimension, deal_II_space_dimension> &input,
+        const std::vector<unsigned int> &                          extents,
+        Triangulation<deal_II_dimension, deal_II_space_dimension> &result);
+
       template void
       create_union_triangulation(
         const Triangulation<deal_II_dimension, deal_II_space_dimension>
diff --git a/tests/grid/replicate_triangulation_01.cc b/tests/grid/replicate_triangulation_01.cc
new file mode 100644 (file)
index 0000000..f89acb7
--- /dev/null
@@ -0,0 +1,76 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 - 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 at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// test GridGenerator::replicate_triangulation
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+
+#include "../tests.h"
+
+void
+test()
+{
+  {
+    Triangulation<1, 1> tr1;
+    GridGenerator::hyper_cube(tr1);
+    std::vector<unsigned int> reps1 = {2};
+    Triangulation<1>          res1;
+    GridGenerator::replicate_triangulation(tr1, reps1, res1);
+    deallog << "1d triangulation has the following cell centers:" << std::endl;
+    for (const auto &cell : res1.active_cell_iterators())
+      deallog << '(' << cell->center() << ')' << std::endl;
+
+    deallog << "Number of global vertices: " << res1.n_vertices() << std::endl;
+  }
+
+  {
+    Triangulation<2, 2> tr2;
+    GridGenerator::hyper_cube_with_cylindrical_hole(tr2, 0.3, 0.5);
+    const std::vector<unsigned int> reps2 = {{2, 2}};
+    Triangulation<2, 2>             res2;
+    GridGenerator::replicate_triangulation(tr2, reps2, res2);
+    deallog << "2d triangulation has the following cell centers:" << std::endl;
+    for (const auto &cell : res2.active_cell_iterators())
+      deallog << '(' << cell->center() << ')' << std::endl;
+
+    deallog << "Number of global vertices: " << res2.n_vertices() << std::endl;
+  }
+
+  {
+    Triangulation<3, 3>       tr3;
+    std::vector<unsigned int> reps3 = {2, 2, 2};
+    GridGenerator::subdivided_hyper_rectangle(tr3,
+                                              reps3,
+                                              Point<3>(-0.7, -0.5, -1.5),
+                                              Point<3>(0.5, 0.5, 2.));
+    Triangulation<3, 3> res3;
+    GridGenerator::replicate_triangulation(tr3, reps3, res3);
+    deallog << "3d triangulation has the following cell centers:" << std::endl;
+    for (const auto &cell : res3.active_cell_iterators())
+      deallog << '(' << cell->center() << ')' << std::endl;
+
+    deallog << "Number of global vertices: " << res3.n_vertices() << std::endl;
+  }
+}
+
+int
+main()
+{
+  initlog();
+  test();
+}
diff --git a/tests/grid/replicate_triangulation_01.output b/tests/grid/replicate_triangulation_01.output
new file mode 100644 (file)
index 0000000..1a154bf
--- /dev/null
@@ -0,0 +1,105 @@
+
+DEAL::1d triangulation has the following cell centers:
+DEAL::(0.500000)
+DEAL::(1.50000)
+DEAL::Number of global vertices: 3
+DEAL::2d triangulation has the following cell centers:
+DEAL::(0.378033 0.178033)
+DEAL::(0.178033 0.378033)
+DEAL::(-0.178033 0.378033)
+DEAL::(-0.378033 0.178033)
+DEAL::(-0.378033 -0.178033)
+DEAL::(-0.178033 -0.378033)
+DEAL::(0.178033 -0.378033)
+DEAL::(0.378033 -0.178033)
+DEAL::(1.37803 0.178033)
+DEAL::(1.17803 0.378033)
+DEAL::(0.821967 0.378033)
+DEAL::(0.621967 0.178033)
+DEAL::(0.621967 -0.178033)
+DEAL::(0.821967 -0.378033)
+DEAL::(1.17803 -0.378033)
+DEAL::(1.37803 -0.178033)
+DEAL::(0.378033 1.17803)
+DEAL::(0.178033 1.37803)
+DEAL::(-0.178033 1.37803)
+DEAL::(-0.378033 1.17803)
+DEAL::(-0.378033 0.821967)
+DEAL::(-0.178033 0.621967)
+DEAL::(0.178033 0.621967)
+DEAL::(0.378033 0.821967)
+DEAL::(1.37803 1.17803)
+DEAL::(1.17803 1.37803)
+DEAL::(0.821967 1.37803)
+DEAL::(0.621967 1.17803)
+DEAL::(0.621967 0.821967)
+DEAL::(0.821967 0.621967)
+DEAL::(1.17803 0.621967)
+DEAL::(1.37803 0.821967)
+DEAL::Number of global vertices: 53
+DEAL::3d triangulation has the following cell centers:
+DEAL::(-0.400000 -0.250000 -0.625000)
+DEAL::(0.200000 -0.250000 -0.625000)
+DEAL::(-0.400000 0.250000 -0.625000)
+DEAL::(0.200000 0.250000 -0.625000)
+DEAL::(-0.400000 -0.250000 1.12500)
+DEAL::(0.200000 -0.250000 1.12500)
+DEAL::(-0.400000 0.250000 1.12500)
+DEAL::(0.200000 0.250000 1.12500)
+DEAL::(0.800000 -0.250000 -0.625000)
+DEAL::(1.40000 -0.250000 -0.625000)
+DEAL::(0.800000 0.250000 -0.625000)
+DEAL::(1.40000 0.250000 -0.625000)
+DEAL::(0.800000 -0.250000 1.12500)
+DEAL::(1.40000 -0.250000 1.12500)
+DEAL::(0.800000 0.250000 1.12500)
+DEAL::(1.40000 0.250000 1.12500)
+DEAL::(-0.400000 0.750000 -0.625000)
+DEAL::(0.200000 0.750000 -0.625000)
+DEAL::(-0.400000 1.25000 -0.625000)
+DEAL::(0.200000 1.25000 -0.625000)
+DEAL::(-0.400000 0.750000 1.12500)
+DEAL::(0.200000 0.750000 1.12500)
+DEAL::(-0.400000 1.25000 1.12500)
+DEAL::(0.200000 1.25000 1.12500)
+DEAL::(0.800000 0.750000 -0.625000)
+DEAL::(1.40000 0.750000 -0.625000)
+DEAL::(0.800000 1.25000 -0.625000)
+DEAL::(1.40000 1.25000 -0.625000)
+DEAL::(0.800000 0.750000 1.12500)
+DEAL::(1.40000 0.750000 1.12500)
+DEAL::(0.800000 1.25000 1.12500)
+DEAL::(1.40000 1.25000 1.12500)
+DEAL::(-0.400000 -0.250000 2.87500)
+DEAL::(0.200000 -0.250000 2.87500)
+DEAL::(-0.400000 0.250000 2.87500)
+DEAL::(0.200000 0.250000 2.87500)
+DEAL::(-0.400000 -0.250000 4.62500)
+DEAL::(0.200000 -0.250000 4.62500)
+DEAL::(-0.400000 0.250000 4.62500)
+DEAL::(0.200000 0.250000 4.62500)
+DEAL::(0.800000 -0.250000 2.87500)
+DEAL::(1.40000 -0.250000 2.87500)
+DEAL::(0.800000 0.250000 2.87500)
+DEAL::(1.40000 0.250000 2.87500)
+DEAL::(0.800000 -0.250000 4.62500)
+DEAL::(1.40000 -0.250000 4.62500)
+DEAL::(0.800000 0.250000 4.62500)
+DEAL::(1.40000 0.250000 4.62500)
+DEAL::(-0.400000 0.750000 2.87500)
+DEAL::(0.200000 0.750000 2.87500)
+DEAL::(-0.400000 1.25000 2.87500)
+DEAL::(0.200000 1.25000 2.87500)
+DEAL::(-0.400000 0.750000 4.62500)
+DEAL::(0.200000 0.750000 4.62500)
+DEAL::(-0.400000 1.25000 4.62500)
+DEAL::(0.200000 1.25000 4.62500)
+DEAL::(0.800000 0.750000 2.87500)
+DEAL::(1.40000 0.750000 2.87500)
+DEAL::(0.800000 1.25000 2.87500)
+DEAL::(1.40000 1.25000 2.87500)
+DEAL::(0.800000 0.750000 4.62500)
+DEAL::(1.40000 0.750000 4.62500)
+DEAL::(0.800000 1.25000 4.62500)
+DEAL::(1.40000 1.25000 4.62500)
+DEAL::Number of global vertices: 125
diff --git a/tests/grid/replicate_triangulation_02.cc b/tests/grid/replicate_triangulation_02.cc
new file mode 100644 (file)
index 0000000..5968abf
--- /dev/null
@@ -0,0 +1,81 @@
+// ---------------------------------------------------------------------
+//
+// 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 at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// test GridGenerator::replicate_triangulation's ability to copy boundary ids
+// correctly
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+
+#include "../tests.h"
+
+void
+test()
+{
+  // no subcell data in 1D so nothing to test
+
+  {
+    Triangulation<2, 2> tr2;
+    GridGenerator::hyper_cube(tr2, 0.0, 1.0, true);
+    const std::vector<unsigned int> reps2 = {{2, 2}};
+    Triangulation<2, 2>             res2;
+    GridGenerator::replicate_triangulation(tr2, reps2, res2);
+    for (const auto &face : res2.active_face_iterators())
+      if (face->at_boundary())
+        deallog << "face: " << face->center()
+                << " boundary id: " << face->boundary_id() << std::endl;
+  }
+
+  {
+    Triangulation<3, 3> tr3;
+    GridGenerator::hyper_cube(tr3, 0.0, 1.0, true);
+    // actually set some boundary line boundary ids to something besides zero
+    // for testing purposes
+    for (auto &face : tr3.active_face_iterators())
+      {
+        if (face->boundary_id() == 4)
+          face->set_all_boundary_ids(4);
+        else if (face->boundary_id() == 5)
+          face->set_all_boundary_ids(5);
+      }
+    // set all (i.e., overwrite two set above) lines on the x = 1 face to
+    // boundary id 1
+    for (auto &face : tr3.active_face_iterators())
+      if (face->boundary_id() == 1)
+        face->set_all_boundary_ids(1);
+    std::vector<unsigned int> reps3 = {2, 1, 1};
+    Triangulation<3, 3>       res3;
+    GridGenerator::replicate_triangulation(tr3, reps3, res3);
+    for (const auto &face : res3.active_face_iterators())
+      if (face->at_boundary())
+        {
+          deallog << "face: " << face->center()
+                  << " boundary id: " << face->boundary_id() << std::endl;
+          for (unsigned int l = 0; l < GeometryInfo<2>::lines_per_cell; ++l)
+            deallog << "line: " << face->line(l)->center()
+                    << " boundary id: " << face->line(l)->boundary_id()
+                    << std::endl;
+        }
+  }
+}
+
+int
+main()
+{
+  initlog();
+  test();
+}
diff --git a/tests/grid/replicate_triangulation_02.output b/tests/grid/replicate_triangulation_02.output
new file mode 100644 (file)
index 0000000..2fc72ee
--- /dev/null
@@ -0,0 +1,59 @@
+
+DEAL::face: 0.500000 0.00000 boundary id: 2
+DEAL::face: 0.00000 0.500000 boundary id: 0
+DEAL::face: 1.50000 0.00000 boundary id: 2
+DEAL::face: 0.00000 1.50000 boundary id: 0
+DEAL::face: 2.00000 0.500000 boundary id: 1
+DEAL::face: 2.00000 1.50000 boundary id: 1
+DEAL::face: 0.500000 2.00000 boundary id: 3
+DEAL::face: 1.50000 2.00000 boundary id: 3
+DEAL::face: 0.500000 0.00000 0.500000 boundary id: 2
+DEAL::line: 0.500000 0.00000 0.00000 boundary id: 0
+DEAL::line: 0.500000 0.00000 1.00000 boundary id: 0
+DEAL::line: 0.00000 0.00000 0.500000 boundary id: 0
+DEAL::line: 1.00000 0.00000 0.500000 boundary id: 0
+DEAL::face: 0.500000 0.500000 0.00000 boundary id: 4
+DEAL::line: 0.00000 0.500000 0.00000 boundary id: 0
+DEAL::line: 1.00000 0.500000 0.00000 boundary id: 0
+DEAL::line: 0.500000 0.00000 0.00000 boundary id: 0
+DEAL::line: 0.500000 1.00000 0.00000 boundary id: 0
+DEAL::face: 0.00000 0.500000 0.500000 boundary id: 0
+DEAL::line: 0.00000 0.00000 0.500000 boundary id: 0
+DEAL::line: 0.00000 1.00000 0.500000 boundary id: 0
+DEAL::line: 0.00000 0.500000 0.00000 boundary id: 0
+DEAL::line: 0.00000 0.500000 1.00000 boundary id: 0
+DEAL::face: 1.50000 0.500000 0.00000 boundary id: 4
+DEAL::line: 1.00000 0.500000 0.00000 boundary id: 0
+DEAL::line: 2.00000 0.500000 0.00000 boundary id: 0
+DEAL::line: 1.50000 0.00000 0.00000 boundary id: 0
+DEAL::line: 1.50000 1.00000 0.00000 boundary id: 0
+DEAL::face: 1.50000 0.00000 0.500000 boundary id: 2
+DEAL::line: 1.50000 0.00000 0.00000 boundary id: 0
+DEAL::line: 1.50000 0.00000 1.00000 boundary id: 0
+DEAL::line: 1.00000 0.00000 0.500000 boundary id: 0
+DEAL::line: 2.00000 0.00000 0.500000 boundary id: 0
+DEAL::face: 0.500000 1.00000 0.500000 boundary id: 3
+DEAL::line: 0.500000 1.00000 0.00000 boundary id: 0
+DEAL::line: 0.500000 1.00000 1.00000 boundary id: 0
+DEAL::line: 0.00000 1.00000 0.500000 boundary id: 0
+DEAL::line: 1.00000 1.00000 0.500000 boundary id: 0
+DEAL::face: 1.50000 1.00000 0.500000 boundary id: 3
+DEAL::line: 1.50000 1.00000 0.00000 boundary id: 0
+DEAL::line: 1.50000 1.00000 1.00000 boundary id: 0
+DEAL::line: 1.00000 1.00000 0.500000 boundary id: 0
+DEAL::line: 2.00000 1.00000 0.500000 boundary id: 0
+DEAL::face: 0.500000 0.500000 1.00000 boundary id: 5
+DEAL::line: 0.00000 0.500000 1.00000 boundary id: 0
+DEAL::line: 1.00000 0.500000 1.00000 boundary id: 0
+DEAL::line: 0.500000 0.00000 1.00000 boundary id: 0
+DEAL::line: 0.500000 1.00000 1.00000 boundary id: 0
+DEAL::face: 1.50000 0.500000 1.00000 boundary id: 5
+DEAL::line: 1.00000 0.500000 1.00000 boundary id: 0
+DEAL::line: 2.00000 0.500000 1.00000 boundary id: 0
+DEAL::line: 1.50000 0.00000 1.00000 boundary id: 0
+DEAL::line: 1.50000 1.00000 1.00000 boundary id: 0
+DEAL::face: 2.00000 0.500000 0.500000 boundary id: 1
+DEAL::line: 2.00000 0.00000 0.500000 boundary id: 0
+DEAL::line: 2.00000 1.00000 0.500000 boundary id: 0
+DEAL::line: 2.00000 0.500000 0.00000 boundary id: 0
+DEAL::line: 2.00000 0.500000 1.00000 boundary id: 0

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.