]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix get_boundary_ids() and get_manifold_ids() for parallel triangulations 9812/head
authorPeter Munch <peterrmuench@gmail.com>
Fri, 3 Apr 2020 22:30:24 +0000 (00:30 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Wed, 8 Apr 2020 11:02:26 +0000 (13:02 +0200)
include/deal.II/distributed/tria_base.h
include/deal.II/grid/tria.h
source/distributed/tria_base.cc
source/grid/tria.cc
tests/distributed_grids/get_boundary_ids_01.cc [new file with mode: 0644]
tests/distributed_grids/get_boundary_ids_01.mpirun=2.output [new file with mode: 0644]

index 458c72699c046fd95207806e8b9f6f2fdcc5c3a4..0203f6cc7732b4d32bf5d2935a3fa95cbb4d3bc6 100644 (file)
@@ -211,6 +211,24 @@ namespace parallel
     virtual std::map<unsigned int, std::set<dealii::types::subdomain_id>>
     compute_vertices_with_ghost_neighbors() const;
 
+    /**
+     * @copydoc dealii::Triangulation::get_boundary_ids()
+     *
+     * @note This function involves a global communication gathering all current
+     *   IDs from all processes.
+     */
+    virtual std::vector<types::boundary_id>
+    get_boundary_ids() const override;
+
+    /**
+     * @copydoc dealii::Triangulation::get_manifold_ids()
+     *
+     * @note This function involves a global communication gathering all current
+     *   IDs from all processes.
+     */
+    virtual std::vector<types::manifold_id>
+    get_manifold_ids() const override;
+
   protected:
     /**
      * MPI communicator to be used for the triangulation. We create a unique
index 5866fd21884bd1cdf4a6b05bbc033729e7336206..0c40ea2ef4c8eec1fe41316bef705a7387e17191 100644 (file)
@@ -1708,30 +1708,32 @@ public:
 
   /**
    * Return a vector containing all boundary indicators assigned to boundary
-   * faces of this Triangulation object. Note, that each boundary indicator is
-   * reported only once. The size of the return vector will represent the
-   * number of different indicators (which is greater or equal one).
+   * faces of active cells of this Triangulation object. Note, that each
+   * boundary indicator is reported only once. The size of the return vector
+   * will represent the number of different indicators (which is greater or
+   * equal one).
    *
    * @ingroup boundary
    *
    * @see
    * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators"
    */
-  std::vector<types::boundary_id>
+  virtual std::vector<types::boundary_id>
   get_boundary_ids() const;
 
   /**
    * Return a vector containing all manifold indicators assigned to the
-   * objects of this Triangulation. Note, that each manifold indicator is
-   * reported only once. The size of the return vector will represent the
-   * number of different indicators (which is greater or equal one).
+   * objects of the active cells of this Triangulation. Note, that each
+   * manifold indicator is reported only once. The size of the return vector
+   * will represent the number of different indicators (which is greater or
+   * equal one).
    *
    * @ingroup manifold
    *
    * @see
    * @ref GlossManifoldIndicator "Glossary entry on manifold indicators"
    */
-  std::vector<types::manifold_id>
+  virtual std::vector<types::manifold_id>
   get_manifold_ids() const;
 
   /**
index 7117d6f2ac8c8cf1b4e788a4e579da32b37216dd..beecfd94862b1acad6d0fb038b02c60721337721 100644 (file)
@@ -403,6 +403,28 @@ namespace parallel
 
 
 
+  template <int dim, int spacedim>
+  std::vector<types::boundary_id>
+  TriangulationBase<dim, spacedim>::get_boundary_ids() const
+  {
+    return Utilities::MPI::compute_set_union(
+      dealii::Triangulation<dim, spacedim>::get_boundary_ids(),
+      this->mpi_communicator);
+  }
+
+
+
+  template <int dim, int spacedim>
+  std::vector<types::manifold_id>
+  TriangulationBase<dim, spacedim>::get_manifold_ids() const
+  {
+    return Utilities::MPI::compute_set_union(
+      dealii::Triangulation<dim, spacedim>::get_manifold_ids(),
+      this->mpi_communicator);
+  }
+
+
+
   template <int dim, int spacedim>
   DistributedTriangulationBase<dim, spacedim>::DistributedTriangulationBase(
     MPI_Comm mpi_communicator,
index 98f8d49917499cda190d4f2bc03d723b9e49d8de..3bc8fa488a120748b6e389993ee3f8782a7beb79 100644 (file)
@@ -10397,11 +10397,11 @@ Triangulation<dim, spacedim>::get_boundary_ids() const
   else
     {
       std::set<types::boundary_id> b_ids;
-      active_cell_iterator         cell = begin_active();
-      for (; cell != end(); ++cell)
-        for (const unsigned int face : GeometryInfo<dim>::face_indices())
-          if (cell->at_boundary(face))
-            b_ids.insert(cell->face(face)->boundary_id());
+      for (auto cell : active_cell_iterators())
+        if (cell->is_locally_owned())
+          for (const unsigned int face : GeometryInfo<dim>::face_indices())
+            if (cell->at_boundary(face))
+              b_ids.insert(cell->face(face)->boundary_id());
       std::vector<types::boundary_id> boundary_ids(b_ids.begin(), b_ids.end());
       return boundary_ids;
     }
@@ -10414,15 +10414,15 @@ std::vector<types::manifold_id>
 Triangulation<dim, spacedim>::get_manifold_ids() const
 {
   std::set<types::manifold_id> m_ids;
-  active_cell_iterator         cell = begin_active();
-  for (; cell != end(); ++cell)
-    {
-      m_ids.insert(cell->manifold_id());
-      if (dim > 1)
-        for (const unsigned int face : GeometryInfo<dim>::face_indices())
-          if (cell->at_boundary(face))
-            m_ids.insert(cell->face(face)->manifold_id());
-    }
+  for (auto cell : active_cell_iterators())
+    if (cell->is_locally_owned())
+      {
+        m_ids.insert(cell->manifold_id());
+        if (dim > 1)
+          for (const unsigned int face : GeometryInfo<dim>::face_indices())
+            if (cell->at_boundary(face))
+              m_ids.insert(cell->face(face)->manifold_id());
+      }
   std::vector<types::manifold_id> manifold_indicators(m_ids.begin(),
                                                       m_ids.end());
   return manifold_indicators;
diff --git a/tests/distributed_grids/get_boundary_ids_01.cc b/tests/distributed_grids/get_boundary_ids_01.cc
new file mode 100644 (file)
index 0000000..3020c84
--- /dev/null
@@ -0,0 +1,58 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+// create a tria mesh and copy it
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include "../tests.h"
+
+template <int dim>
+void
+test()
+{
+  parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD);
+  GridGenerator::subdivided_hyper_cube(tria, 4);
+
+  for (auto cell : tria.active_cell_iterators())
+    if (cell->is_locally_owned())
+      for (const unsigned int face : GeometryInfo<dim>::face_indices())
+        {
+          if (std::fabs(cell->face(face)->center()(0) - 0.0) < 1e-12)
+            cell->face(face)->set_all_boundary_ids(1);
+          if (std::fabs(cell->face(face)->center()(0) - 1.0) < 1e-12)
+            cell->face(face)->set_all_boundary_ids(2);
+          if (std::fabs(cell->face(face)->center()(1) - 0.0) < 1e-12)
+            cell->face(face)->set_all_boundary_ids(3);
+          if (std::fabs(cell->face(face)->center()(1) - 1.0) < 1e-12)
+            cell->face(face)->set_all_boundary_ids(4);
+        }
+
+  for (auto i : tria.get_boundary_ids())
+    deallog << i << " ";
+  deallog << std::endl;
+}
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll                    all;
+
+  test<2>();
+}
diff --git a/tests/distributed_grids/get_boundary_ids_01.mpirun=2.output b/tests/distributed_grids/get_boundary_ids_01.mpirun=2.output
new file mode 100644 (file)
index 0000000..8bdee82
--- /dev/null
@@ -0,0 +1,5 @@
+
+DEAL:0::1 2 3 4 
+
+DEAL:1::1 2 3 4 
+

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.