]> https://gitweb.dealii.org/ - dealii.git/commitdiff
copy boundary ids to remaining boundary faces 14116/head
authorJohannes Heinz <johannes.heinz@tuwien.ac.at>
Mon, 4 Jul 2022 07:31:06 +0000 (09:31 +0200)
committerJohannes Heinz <johannes.heinz@tuwien.ac.at>
Mon, 11 Jul 2022 14:21:59 +0000 (16:21 +0200)
Co-authored-by: Peter Munch <peterrmuench@gmail.com>
doc/news/changes/minor/20220707Heinz [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/grid_tools_07.cc
tests/grid/merge_triangulations_10.cc [new file with mode: 0644]
tests/grid/merge_triangulations_10.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20220707Heinz b/doc/news/changes/minor/20220707Heinz
new file mode 100644 (file)
index 0000000..26dfdb4
--- /dev/null
@@ -0,0 +1,3 @@
+Improved: The function GridGenerator::merge_triangulations() can now copy boundary ids to remaining boundary faces.
+<br>
+(Johannes Heinz, 2022/07/07)
index d6de5006ab08ace235f3c8f14883e818c08c36b3..73865b5f29ca23252e89116148d42a5beb5ee87e 100644 (file)
@@ -1884,10 +1884,8 @@ namespace GridGenerator
    *
    * @note The function copies the material ids of the cells of the two input
    * triangulations into the output triangulation. If @p copy_manifold_ids is
-   * set to @p true, manifold ids will be copied. Boundary indicators are never
-   * copied. In other words, if the two coarse meshes have anything but the
-   * default boundary indicators, then you will have to set boundary indicators
-   * again by hand in the output triangulation.
+   * set to @p true, manifold ids will be copied. If @p copy_boundary_ids is
+   * set to @p true, boundary_ids are copied to all remaining faces at the boundary.
    *
    * @note This function does not attach any manifolds to @p result, nor does
    * it set any manifold ids. In particular, manifolds attached to the two
@@ -1903,7 +1901,8 @@ namespace GridGenerator
                        const Triangulation<dim, spacedim> &triangulation_2,
                        Triangulation<dim, spacedim> &      result,
                        const double duplicated_vertex_tolerance = 1.0e-12,
-                       const bool   copy_manifold_ids           = false);
+                       const bool   copy_manifold_ids           = false,
+                       const bool   copy_boundary_ids           = false);
 
   /**
    * Same as above but allows to merge more than two triangulations at once.
@@ -1916,6 +1915,7 @@ namespace GridGenerator
    *   GridGenerator::merge_triangulations({&tria_1, &tria_2, &tria_3},
    *                                       merged_triangulation,
    *                                       1.0e-10,
+   *                                       false,
    *                                       false);
    * @endcode
    */
@@ -1925,7 +1925,8 @@ namespace GridGenerator
     const std::vector<const Triangulation<dim, spacedim> *> &triangulations,
     Triangulation<dim, spacedim> &                           result,
     const double duplicated_vertex_tolerance = 1.0e-12,
-    const bool   copy_manifold_ids           = false);
+    const bool   copy_manifold_ids           = false,
+    const bool   copy_boundary_ids           = false);
 
   /**
    * \brief Replicate a given triangulation in multiple coordinate axes
index 41fa3a00f452e69101cea19121957878b88a0979..637a714d0ee0d17079dea65d15306ba941efcc47 100644 (file)
@@ -6469,7 +6469,8 @@ namespace GridGenerator
     const std::vector<const Triangulation<dim, spacedim> *> &triangulations,
     Triangulation<dim, spacedim> &                           result,
     const double duplicated_vertex_tolerance,
-    const bool   copy_manifold_ids)
+    const bool   copy_manifold_ids,
+    const bool   copy_boundary_ids)
   {
     std::vector<Point<spacedim>> vertices;
     std::vector<CellData<dim>>   cells;
@@ -6545,6 +6546,22 @@ namespace GridGenerator
     result.clear();
     result.create_triangulation(vertices, cells, subcell_data);
 
+    if (copy_boundary_ids)
+      {
+        auto result_cell = result.begin();
+        for (const auto &tria : triangulations)
+          {
+            for (const auto &cell : tria->cell_iterators())
+              {
+                for (auto const &f : cell->face_indices())
+                  if (result_cell->face(f)->at_boundary())
+                    result_cell->face(f)->set_boundary_id(
+                      cell->face(f)->boundary_id());
+                result_cell++;
+              }
+          }
+      }
+
     Assert(duplicated_vertex_tolerance > 0.0 ||
              n_accumulated_vertices == result.n_vertices(),
            ExcInternalError());
@@ -6558,7 +6575,8 @@ namespace GridGenerator
                        const Triangulation<dim, spacedim> &triangulation_2,
                        Triangulation<dim, spacedim> &      result,
                        const double duplicated_vertex_tolerance,
-                       const bool   copy_manifold_ids)
+                       const bool   copy_manifold_ids,
+                       const bool   copy_boundary_ids)
   {
     // if either Triangulation is empty then merging is just a copy.
     if (triangulation_1.n_cells() == 0)
@@ -6574,7 +6592,8 @@ namespace GridGenerator
     merge_triangulations({&triangulation_1, &triangulation_2},
                          result,
                          duplicated_vertex_tolerance,
-                         copy_manifold_ids);
+                         copy_manifold_ids,
+                         copy_boundary_ids);
   }
 
 
index 05d1fff1c24cf1f39894b5854ebce97b3d725031..c0b6cce0bdf6e4f3c8d915c7def24ae9f0cca99e 100644 (file)
@@ -84,6 +84,7 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
           const Triangulation<deal_II_dimension, deal_II_space_dimension> *> &,
         Triangulation<deal_II_dimension, deal_II_space_dimension> &,
         const double,
+        const bool,
         const bool);
 
       template void
@@ -94,7 +95,8 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
           &triangulation_2,
         Triangulation<deal_II_dimension, deal_II_space_dimension> &result,
         const double duplicated_vertex_tolerance,
-        const bool   copy_manifold_ids);
+        const bool   copy_manifold_ids,
+        const bool   copy_boundary_ids);
 
       template void
       replicate_triangulation(
index a40c3c36679f29c00d020407c2ad47c0d2c6a9ea..141a680ed5d3a9c84557f6fd69e2eec07dc798cc 100644 (file)
@@ -65,21 +65,7 @@ generate_grid(Triangulation<2> &triangulation)
 
   // tolerance 0. to ensure vertices are not merged
   GridGenerator::merge_triangulations(
-    triangulation1, triangulation2, triangulation, 0., false);
-
-  // make sure boundary ids are kept
-  std::vector<types::boundary_id> boundary_ids;
-  for (const auto &face : triangulation1.active_face_iterators())
-    if (face->at_boundary())
-      boundary_ids.emplace_back(face->boundary_id());
-  for (const auto &face : triangulation2.active_face_iterators())
-    if (face->at_boundary())
-      boundary_ids.emplace_back(face->boundary_id());
-
-  unsigned int i = 0;
-  for (const auto &face : triangulation.active_face_iterators())
-    if (face->at_boundary())
-      face->set_boundary_id(boundary_ids[i++]);
+    triangulation1, triangulation2, triangulation, 0., false, true);
 }
 
 std::vector<unsigned int>
diff --git a/tests/grid/merge_triangulations_10.cc b/tests/grid/merge_triangulations_10.cc
new file mode 100644 (file)
index 0000000..d701127
--- /dev/null
@@ -0,0 +1,73 @@
+// ---------------------------------------------------------------------
+//
+// 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Test for flag copy_boundary_ids in merge_triangulations()
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include "../tests.h"
+
+template <int dim>
+void
+print_boundary_ids(const dealii::Triangulation<dim> &tria)
+{
+  for (const auto &face : tria.active_face_iterators())
+    if (face->at_boundary())
+      deallog << "Face at " << face->center()
+              << " has id: " << face->boundary_id() << "\n";
+  deallog << "-----------------" << std::endl;
+}
+
+
+int
+main()
+{
+  using namespace dealii;
+
+  initlog();
+
+  const unsigned int dim = 2;
+
+  Triangulation<dim> tria1;
+  Triangulation<dim> tria2;
+  Triangulation<dim> tria3;
+  Triangulation<dim> tria;
+  GridGenerator::hyper_rectangle(tria1, {-1., -1.}, {0., 1.}, true);
+  GridGenerator::hyper_rectangle(tria2, {0., -1.}, {1., 1.}, true);
+  for (auto &face : tria2.active_face_iterators())
+    if (face->at_boundary())
+      face->set_boundary_id(face->boundary_id() + dim * 2);
+
+  deallog << "Initial IDs of Triangulation 1" << std::endl;
+  print_boundary_ids(tria1);
+  deallog << "Initial IDs of Triangulation 2" << std::endl;
+  print_boundary_ids(tria2);
+
+  deallog << "Merge two rectangles, one face is eliminated" << std::endl;
+  GridGenerator::merge_triangulations(tria1, tria2, tria, 1e-6, false, true);
+  print_boundary_ids(tria);
+
+  deallog
+    << "Merge two rectangles, due to hanging node there is no internal face"
+    << std::endl;
+  tria.clear();
+  tria2.refine_global(1);
+  GridGenerator::flatten_triangulation(tria2, tria3);
+  GridGenerator::merge_triangulations(tria1, tria3, tria, 1e-6, false, true);
+  print_boundary_ids(tria);
+
+  deallog << std::endl;
+}
diff --git a/tests/grid/merge_triangulations_10.output b/tests/grid/merge_triangulations_10.output
new file mode 100644 (file)
index 0000000..7c77652
--- /dev/null
@@ -0,0 +1,36 @@
+
+DEAL::Initial IDs of Triangulation 1
+DEAL::Face at -0.500000 -1.00000 has id: 2
+Face at -1.00000 0.00000 has id: 0
+Face at 0.00000 0.00000 has id: 1
+Face at -0.500000 1.00000 has id: 3
+-----------------
+DEAL::Initial IDs of Triangulation 2
+DEAL::Face at 0.500000 -1.00000 has id: 6
+Face at 0.00000 0.00000 has id: 4
+Face at 1.00000 0.00000 has id: 5
+Face at 0.500000 1.00000 has id: 7
+-----------------
+DEAL::Merge two rectangles, one face is eliminated
+DEAL::Face at -0.500000 -1.00000 has id: 2
+Face at -1.00000 0.00000 has id: 0
+Face at 0.500000 -1.00000 has id: 6
+Face at -0.500000 1.00000 has id: 3
+Face at 0.500000 1.00000 has id: 7
+Face at 1.00000 0.00000 has id: 5
+-----------------
+DEAL::Merge two rectangles, due to hanging node there is no internal face
+DEAL::Face at -0.500000 -1.00000 has id: 2
+Face at -1.00000 0.00000 has id: 0
+Face at 0.00000 0.00000 has id: 1
+Face at 0.250000 -1.00000 has id: 6
+Face at 0.00000 -0.500000 has id: 4
+Face at -0.500000 1.00000 has id: 3
+Face at 0.250000 1.00000 has id: 7
+Face at 1.00000 -0.500000 has id: 5
+Face at 0.750000 -1.00000 has id: 6
+Face at 0.00000 0.500000 has id: 4
+Face at 1.00000 0.500000 has id: 5
+Face at 0.750000 1.00000 has id: 7
+-----------------
+DEAL::

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.