]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added get/set functions for future FE indices. 14208/head
authorMarc Fehling <mafehling.git@gmail.com>
Sun, 21 Aug 2022 03:03:55 +0000 (21:03 -0600)
committerMarc Fehling <mafehling.git@gmail.com>
Sun, 21 Aug 2022 03:41:35 +0000 (21:41 -0600)
doc/news/changes/minor/20220820Fehling [new file with mode: 0644]
include/deal.II/dofs/dof_handler.h
source/dofs/dof_handler.cc
tests/mpi/hp_getset_future_fe_indices_01.cc [new file with mode: 0644]
tests/mpi/hp_getset_future_fe_indices_01.with_p4est=true.mpirun=2.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20220820Fehling b/doc/news/changes/minor/20220820Fehling
new file mode 100644 (file)
index 0000000..d181b8d
--- /dev/null
@@ -0,0 +1,6 @@
+New: Functions DoFHandler::get_future_fe_indices() and
+DoFHandler::set_future_fe_indices() that store and recover future FE
+indices on a DoFHandler. They can also be used to transfer future FE
+indices between different DoFHandler objects.
+<br>
+(Marc Fehling, 2022/08/20)
index 2c9417b7be6b5ecb4b557cc3cc33dad4fd74aa18..a38e256b49617b421cf3d5d28106e9a257c9d592 100644 (file)
@@ -586,6 +586,22 @@ public:
   void
   get_active_fe_indices(std::vector<unsigned int> &active_fe_indices) const;
 
+  /**
+   * Go through the triangulation and set the future FE indices of all
+   * locally owned cells to the values given in @p future_fe_indices.
+   * Cells corresponding to invalid_active_fe_index will be skipped.
+   */
+  void
+  set_future_fe_indices(const std::vector<unsigned int> &future_fe_indices);
+
+  /**
+   * Go through the triangulation and return a vector of future FE indices of
+   * all locally owned cells. If no future FE index has been set on a cell,
+   * its value will be invalid_active_fe_index.
+   */
+  std::vector<unsigned int>
+  get_future_fe_indices() const;
+
   /**
    * Assign a Triangulation to the DoFHandler.
    *
index bd2b7bc949c5c66242ee8d1f62e849cc3a170a91..9f15f7a2b36a0367c7909cfe0a09a87795146978 100644 (file)
@@ -2661,6 +2661,47 @@ DoFHandler<dim, spacedim>::get_active_fe_indices(
 
 
 
+template <int dim, int spacedim>
+void
+DoFHandler<dim, spacedim>::set_future_fe_indices(
+  const std::vector<unsigned int> &future_fe_indices)
+{
+  Assert(future_fe_indices.size() == this->get_triangulation().n_active_cells(),
+         ExcDimensionMismatch(future_fe_indices.size(),
+                              this->get_triangulation().n_active_cells()));
+
+  this->create_active_fe_table();
+  // we could set the values directly, since they are stored as
+  // protected data of this object, but for simplicity we use the
+  // cell-wise access. this way we also have to pass some debug-mode
+  // tests which we would have to duplicate ourselves otherwise
+  for (const auto &cell : this->active_cell_iterators())
+    if (cell->is_locally_owned() &&
+        future_fe_indices[cell->active_cell_index()] != invalid_active_fe_index)
+      cell->set_future_fe_index(future_fe_indices[cell->active_cell_index()]);
+}
+
+
+
+template <int dim, int spacedim>
+std::vector<unsigned int>
+DoFHandler<dim, spacedim>::get_future_fe_indices() const
+{
+  std::vector<unsigned int> future_fe_indices(
+    this->get_triangulation().n_active_cells(), invalid_active_fe_index);
+
+  // we could try to extract the values directly, since they are
+  // stored as protected data of this object, but for simplicity we
+  // use the cell-wise access.
+  for (const auto &cell : this->active_cell_iterators())
+    if (cell->is_locally_owned() && cell->future_fe_index_set())
+      future_fe_indices[cell->active_cell_index()] = cell->future_fe_index();
+
+  return future_fe_indices;
+}
+
+
+
 template <int dim, int spacedim>
 void
 DoFHandler<dim, spacedim>::connect_to_triangulation_signals()
diff --git a/tests/mpi/hp_getset_future_fe_indices_01.cc b/tests/mpi/hp_getset_future_fe_indices_01.cc
new file mode 100644 (file)
index 0000000..6ac3b61
--- /dev/null
@@ -0,0 +1,95 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2022 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// Transfer future FE indices between DoFHandler objects.
+
+
+#include <deal.II/base/mpi.h>
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <algorithm>
+
+#include "../tests.h"
+
+#include "../test_grids.h"
+
+
+template <int dim, int spacedim>
+void
+test()
+{
+  MPI_Comm mpi_comm = MPI_COMM_WORLD;
+
+  const unsigned int n_procs = Utilities::MPI::n_mpi_processes(mpi_comm);
+
+  parallel::distributed::Triangulation<dim, spacedim> tria(mpi_comm);
+  TestGrids::hyper_line(tria, n_procs);
+  tria.refine_global(1);
+
+  DoFHandler<dim, spacedim> dofh_src(tria);
+  DoFHandler<dim, spacedim> dofh_dst(tria);
+
+  // set future fe index on first cell of every process
+  for (const auto &cell : dofh_src.active_cell_iterators())
+    if (cell->is_locally_owned())
+      {
+        cell->set_future_fe_index(1);
+        break;
+      }
+
+  // transfer future FE indices
+  dofh_dst.set_future_fe_indices(dofh_src.get_future_fe_indices());
+
+  // count copied future FE indices
+  const unsigned int n_indices_copied =
+    std::count_if(dofh_dst.active_cell_iterators().begin(),
+                  dofh_dst.active_cell_iterators().end(),
+                  [](const auto &cell) {
+                    return (cell->is_locally_owned() &&
+                            cell->future_fe_index_set());
+                  });
+  deallog << "indices copied: " << n_indices_copied << std::endl;
+
+  // verify that future FE indices match
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->is_locally_owned())
+      {
+        const auto cell_src = cell->as_dof_handler_iterator(dofh_src);
+        const auto cell_dst = cell->as_dof_handler_iterator(dofh_dst);
+
+        AssertThrow(cell_src->future_fe_index_set() ==
+                      cell_dst->future_fe_index_set(),
+                    ExcInternalError());
+        AssertThrow(cell_src->future_fe_index() == cell_dst->future_fe_index(),
+                    ExcInternalError());
+      }
+
+  deallog << "OK" << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll                    log;
+
+  test<2, 2>();
+}
diff --git a/tests/mpi/hp_getset_future_fe_indices_01.with_p4est=true.mpirun=2.output b/tests/mpi/hp_getset_future_fe_indices_01.with_p4est=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..a6f21ad
--- /dev/null
@@ -0,0 +1,7 @@
+
+DEAL:0::indices copied: 1
+DEAL:0::OK
+
+DEAL:1::indices copied: 1
+DEAL:1::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.