]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added GridTools::get_subdomain_association() for CellId collections. 11366/head
authorMarc Fehling <mafehling.git@gmail.com>
Sun, 13 Dec 2020 06:33:52 +0000 (23:33 -0700)
committerMarc Fehling <mafehling.git@gmail.com>
Mon, 11 Jan 2021 04:13:06 +0000 (21:13 -0700)
doc/news/changes/minor/20210107Fehling [new file with mode: 0644]
include/deal.II/grid/grid_tools.h
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/grid/grid_tools_subdomain_association_01.cc [new file with mode: 0644]
tests/grid/grid_tools_subdomain_association_01.with_p4est=true.mpirun=4.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20210107Fehling b/doc/news/changes/minor/20210107Fehling
new file mode 100644 (file)
index 0000000..03b702a
--- /dev/null
@@ -0,0 +1,5 @@
+New: Function GridTools::get_subdomain_association() determines
+the owning process of any valid cell on a Triangulation that is
+represented by a globally unique CellId.
+<br>
+(Marc Fehling, 2021/01/07)
index ae933ff775500fdebd3aedb9d061b0363ccba72b..93698a6c0f9c57899537f5c0997a0a3b00de5b11 100644 (file)
@@ -2100,6 +2100,18 @@ namespace GridTools
   void
   partition_multigrid_levels(Triangulation<dim, spacedim> &triangulation);
 
+  /**
+   * This function allows to ask for the owning subdomain of cells identified by
+   * CellId objects that do not have to exist on the current process.
+   *
+   * @note This function has not been implemented yet for
+   *   parallel::fullydistributed::Triangulation.
+   */
+  template <int dim, int spacedim>
+  std::vector<types::subdomain_id>
+  get_subdomain_association(const Triangulation<dim, spacedim> &triangulation,
+                            const std::vector<CellId> &         cell_ids);
+
   /**
    * For each active cell, return in the output array to which subdomain (as
    * given by the <tt>cell->subdomain_id()</tt> function) it belongs. The
@@ -2135,7 +2147,6 @@ namespace GridTools
     const Triangulation<dim, spacedim> &triangulation,
     const types::subdomain_id           subdomain);
 
-
   /**
    * For a triangulation, return a mask that represents which of its vertices
    * are "owned" by the current process in the same way as we talk about
index a4fdbb3c0cdd4fdcf1e5ba2f9ded82d6421adb72..b81947d6237a5ad1784eeba90de47f8f09c8db69 100644 (file)
@@ -18,6 +18,8 @@
 #include <deal.II/base/quadrature_lib.h>
 #include <deal.II/base/thread_management.h>
 
+#include <deal.II/distributed/fully_distributed_tria.h>
+#include <deal.II/distributed/p4est_wrappers.h>
 #include <deal.II/distributed/shared_tria.h>
 #include <deal.II/distributed/tria.h>
 
@@ -3035,6 +3037,101 @@ namespace GridTools
   }
 
 
+
+  template <int dim, int spacedim>
+  std::vector<types::subdomain_id>
+  get_subdomain_association(const Triangulation<dim, spacedim> &triangulation,
+                            const std::vector<CellId> &         cell_ids)
+  {
+    std::vector<types::subdomain_id> subdomain_ids;
+    subdomain_ids.reserve(cell_ids.size());
+
+    if (dynamic_cast<
+          const parallel::fullydistributed::Triangulation<dim, spacedim> *>(
+          &triangulation) != nullptr)
+      {
+        Assert(false, ExcNotImplemented());
+      }
+    else if (const parallel::distributed::Triangulation<dim, spacedim>
+               *parallel_tria = dynamic_cast<
+                 const parallel::distributed::Triangulation<dim, spacedim> *>(
+                 &triangulation))
+      {
+#ifndef DEAL_II_WITH_P4EST
+        Assert(
+          false,
+          ExcMessage(
+            "You are attempting to use a functionality that is only available "
+            "if deal.II was configured to use p4est, but cmake did not find a "
+            "valid p4est library."));
+#else
+        // for parallel distributed triangulations, we will ask the p4est oracle
+        // about the global partitioning of active cells since this information
+        // is stored on every process
+        for (const auto &cell_id : cell_ids)
+          {
+            // find descendent from coarse quadrant
+            typename dealii::internal::p4est::types<dim>::quadrant p4est_cell,
+              p4est_children[GeometryInfo<dim>::max_children_per_cell];
+
+            dealii::internal::p4est::init_coarse_quadrant<dim>(p4est_cell);
+            for (const auto &child_index : cell_id.get_child_indices())
+              {
+                dealii::internal::p4est::init_quadrant_children<dim>(
+                  p4est_cell, p4est_children);
+                p4est_cell =
+                  p4est_children[static_cast<unsigned int>(child_index)];
+              }
+
+            // find owning process, i.e., the subdomain id
+            const int owner =
+              dealii::internal::p4est::functions<dim>::comm_find_owner(
+                const_cast<typename dealii::internal::p4est::types<dim>::forest
+                             *>(parallel_tria->get_p4est()),
+                cell_id.get_coarse_cell_id(),
+                &p4est_cell,
+                Utilities::MPI::this_mpi_process(
+                  parallel_tria->get_communicator()));
+
+            Assert(owner >= 0, ExcMessage("p4est should know the owner."));
+
+            subdomain_ids.push_back(owner);
+          }
+#endif
+      }
+    else if (const parallel::shared::Triangulation<dim, spacedim> *shared_tria =
+               dynamic_cast<const parallel::shared::Triangulation<dim, spacedim>
+                              *>(&triangulation))
+      {
+        // for parallel shared triangulations, we need to access true subdomain
+        // ids which are also valid for artificial cells
+        const std::vector<types::subdomain_id> &true_subdomain_ids_of_cells =
+          shared_tria->get_true_subdomain_ids_of_cells();
+
+        for (const auto &cell_id : cell_ids)
+          {
+            const unsigned int active_cell_index =
+              shared_tria->create_cell_iterator(cell_id)->active_cell_index();
+            subdomain_ids.push_back(
+              true_subdomain_ids_of_cells[active_cell_index]);
+          }
+      }
+    else
+      {
+        // the most general type of triangulation is the serial one. here, all
+        // subdomain information is directly available
+        for (const auto &cell_id : cell_ids)
+          {
+            subdomain_ids.push_back(
+              triangulation.create_cell_iterator(cell_id)->subdomain_id());
+          }
+      }
+
+    return subdomain_ids;
+  }
+
+
+
   template <int dim, int spacedim>
   void
   get_subdomain_association(const Triangulation<dim, spacedim> &triangulation,
index 701a0444f6c5e1aa672a1b995963afed87e7aecc..4c42a33e840cad67131d9c3fa2409436c1ae3b3e 100644 (file)
@@ -316,6 +316,13 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
       partition_multigrid_levels(
         Triangulation<deal_II_dimension, deal_II_space_dimension> &);
 
+#  if deal_II_dimension >= 2
+      template std::vector<types::subdomain_id>
+      get_subdomain_association(
+        const Triangulation<deal_II_dimension, deal_II_space_dimension> &,
+        const std::vector<CellId> &);
+#  endif
+
       template void
       get_subdomain_association(
         const Triangulation<deal_II_dimension, deal_II_space_dimension> &,
diff --git a/tests/grid/grid_tools_subdomain_association_01.cc b/tests/grid/grid_tools_subdomain_association_01.cc
new file mode 100644 (file)
index 0000000..3d69e3d
--- /dev/null
@@ -0,0 +1,141 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+// check GridTools::get_subdomain_association
+
+
+#include <deal.II/base/mpi.h>
+
+#include <deal.II/distributed/shared_tria.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/distributed/tria_base.h>
+
+#include <deal.II/grid/cell_id.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/tria.h>
+
+#include <memory>
+#include <vector>
+
+#include "../tests.h"
+
+#include "../test_grids.h"
+
+
+enum Type
+{
+  Shared,
+  Distributed
+};
+
+
+template <int dim>
+std::unique_ptr<parallel::Triangulation<dim>>
+create_triangulation(Type type)
+{
+  if (type == Type::Shared)
+    return std::make_unique<parallel::shared::Triangulation<dim>>(
+      MPI_COMM_WORLD, Triangulation<dim>::none, true);
+  else if (type == Type::Distributed)
+    return std::make_unique<parallel::distributed::Triangulation<dim>>(
+      MPI_COMM_WORLD);
+  else
+    return nullptr;
+}
+
+
+template <int dim>
+void
+test(parallel::Triangulation<dim> &tria)
+{
+  // ----- setup -----
+  tria.clear();
+  TestGrids::hyper_line(tria, 4);
+
+  // refine the last cell
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->id().to_string() == "0_0:")
+      cell->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+
+  // ----- test -----
+  // gather global cellid objects
+  std::vector<CellId> global_cell_ids;
+  global_cell_ids.reserve(tria.n_global_active_cells());
+  {
+    std::vector<CellId> local_cell_ids;
+    local_cell_ids.reserve(tria.n_active_cells());
+    for (const auto &cell : tria.active_cell_iterators())
+      if (cell->is_locally_owned())
+        local_cell_ids.push_back(cell->id());
+
+    std::vector<std::vector<CellId>> cell_ids_per_processor =
+      Utilities::MPI::all_gather(MPI_COMM_WORLD, local_cell_ids);
+
+    for (const auto &cell_ids : cell_ids_per_processor)
+      global_cell_ids.insert(global_cell_ids.end(),
+                             cell_ids.cbegin(),
+                             cell_ids.cend());
+  }
+
+  // determine subdomain of every cellid
+  std::vector<types::subdomain_id> subdomain_ids =
+    GridTools::get_subdomain_association(tria, global_cell_ids);
+
+  // ----- verify results -----
+  AssertDimension(tria.n_global_active_cells(), global_cell_ids.size());
+  AssertDimension(tria.n_global_active_cells(), subdomain_ids.size());
+
+  // check if processor owns represented cells that have been assigned to them
+  for (unsigned int i = 0; i < tria.n_global_active_cells(); ++i)
+    if (tria.locally_owned_subdomain() == subdomain_ids[i])
+      {
+        const auto cell = tria.create_cell_iterator(global_cell_ids[i]);
+        Assert(cell->is_locally_owned(), ExcInternalError());
+      }
+
+  // check if all processors have the same list of subdomain ids
+  {
+    std::vector<std::vector<types::subdomain_id>> subdomain_ids_per_processor =
+      Utilities::MPI::gather(MPI_COMM_WORLD, subdomain_ids);
+
+    if (tria.locally_owned_subdomain() == 0)
+      for (unsigned int i = 1; i < subdomain_ids_per_processor.size(); ++i)
+        Assert(subdomain_ids_per_processor[0] == subdomain_ids_per_processor[i],
+               ExcInternalError());
+  }
+
+  deallog << "OK" << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll                    log;
+
+
+  for (const auto &type : {Type::Shared, Type::Distributed})
+    {
+      deallog.push("2d");
+      test(*create_triangulation<2>(type));
+      deallog.pop();
+      deallog.push("3d");
+      test(*create_triangulation<3>(type));
+      deallog.pop();
+    }
+}
diff --git a/tests/grid/grid_tools_subdomain_association_01.with_p4est=true.mpirun=4.output b/tests/grid/grid_tools_subdomain_association_01.with_p4est=true.mpirun=4.output
new file mode 100644 (file)
index 0000000..e7df6f3
--- /dev/null
@@ -0,0 +1,23 @@
+
+DEAL:0:2d::OK
+DEAL:0:3d::OK
+DEAL:0:2d::OK
+DEAL:0:3d::OK
+
+DEAL:1:2d::OK
+DEAL:1:3d::OK
+DEAL:1:2d::OK
+DEAL:1:3d::OK
+
+
+DEAL:2:2d::OK
+DEAL:2:3d::OK
+DEAL:2:2d::OK
+DEAL:2:3d::OK
+
+
+DEAL:3:2d::OK
+DEAL:3:3d::OK
+DEAL:3:2d::OK
+DEAL:3:3d::OK
+

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.