]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Bugfix: DoFHandler connected to Triangulation::Signals only in hp-mode. 11798/head
authorMarc Fehling <mafehling.git@gmail.com>
Wed, 24 Feb 2021 18:22:33 +0000 (11:22 -0700)
committerMarc Fehling <mafehling.git@gmail.com>
Tue, 2 Mar 2021 05:40:59 +0000 (22:40 -0700)
include/deal.II/dofs/dof_handler.h
source/distributed/solution_transfer.cc
source/dofs/dof_handler.cc
tests/dofs/clear_signal.cc [new file with mode: 0644]
tests/dofs/clear_signal.with_mpi=true.with_p4est=true.mpirun=1.output [new file with mode: 0644]
tests/fullydistributed_grids/solution_transfer_01.cc

index 99601c45c7e136182d4552b151d7e58f66af22d1..d3dd609226142f5e531e2289a7690e38be9dd30b 100644 (file)
@@ -1666,6 +1666,13 @@ private:
    */
   std::vector<boost::signals2::connection> tria_listeners;
 
+  /**
+   * A list of connections with which this object connects to the
+   * triangulation. They get triggered specifially when data needs to be
+   * transferred due to refinement or repartitioning. Only active in hp-mode.
+   */
+  std::vector<boost::signals2::connection> tria_listeners_for_transfer;
+
   /**
    * Free all memory used for non-multigrid data structures.
    */
@@ -1700,14 +1707,13 @@ private:
                 const types::global_dof_index global_index) const;
 
   /**
-   * Setup DoFHandler policy.
+   * Set up DoFHandler policy.
    */
   void
   setup_policy();
 
   /**
-   * Setup connections to refinement signals of the underlying triangulation.
-   * Necessary for the hp-mode.
+   * Set up connections to signals of the underlying triangulation.
    */
   void
   connect_to_triangulation_signals();
index 6f8ac3e08fe4a7dda9c834bd3a468561ba21137f..50104891cdffcfb17f1b467041c90506e85ace7a 100644 (file)
@@ -137,6 +137,10 @@ namespace parallel
       prepare_for_coarsening_and_refinement(
         const std::vector<const VectorType *> &all_in)
     {
+      for (unsigned int i = 0; i < all_in.size(); ++i)
+        Assert(all_in[i]->size() == dof_handler->n_dofs(),
+               ExcDimensionMismatch(all_in[i]->size(), dof_handler->n_dofs()));
+
       input_vectors = all_in;
       register_data_attach();
     }
@@ -230,6 +234,9 @@ namespace parallel
     {
       Assert(input_vectors.size() == all_out.size(),
              ExcDimensionMismatch(input_vectors.size(), all_out.size()));
+      for (unsigned int i = 0; i < all_out.size(); ++i)
+        Assert(all_out[i]->size() == dof_handler->n_dofs(),
+               ExcDimensionMismatch(all_out[i]->size(), dof_handler->n_dofs()));
 
       // TODO: casting away constness is bad
       auto *tria = (dynamic_cast<parallel::DistributedTriangulationBase<
index 16564dcd148fe3332f02c29d7f6bfd047dfa1b02..5e440212febb765bccc1acffcf709b20f099a83a 100644 (file)
@@ -1999,12 +1999,15 @@ DoFHandler<dim, spacedim>::DoFHandler(const Triangulation<dim, spacedim> &tria)
 template <int dim, int spacedim>
 DoFHandler<dim, spacedim>::~DoFHandler()
 {
-  // unsubscribe all attachments to refinement signals of the underlying
-  // triangulation
+  // unsubscribe all attachments to signals of the underlying triangulation
   for (auto &connection : this->tria_listeners)
     connection.disconnect();
   this->tria_listeners.clear();
 
+  for (auto &connection : this->tria_listeners_for_transfer)
+    connection.disconnect();
+  this->tria_listeners_for_transfer.clear();
+
   // release allocated memory
   // virtual functions called in constructors and destructors never use the
   // override in a derived class
@@ -2052,6 +2055,10 @@ DoFHandler<dim, spacedim>::reinit(const Triangulation<dim, spacedim> &tria)
     connection.disconnect();
   this->tria_listeners.clear();
 
+  for (auto &connection : this->tria_listeners_for_transfer)
+    connection.disconnect();
+  this->tria_listeners_for_transfer.clear();
+
   // release allocated memory and policy
   DoFHandler<dim, spacedim>::clear();
   this->policy.reset();
@@ -2431,10 +2438,11 @@ DoFHandler<dim, spacedim>::set_fe(const hp::FECollection<dim, spacedim> &ff)
         {
           hp_capability_enabled = false;
 
-          // unsubscribe connections to signals
-          for (auto &connection : this->tria_listeners)
+          // unsubscribe connections to signals that are only relevant for
+          // hp-mode, since we only have a single element here
+          for (auto &connection : this->tria_listeners_for_transfer)
             connection.disconnect();
-          this->tria_listeners.clear();
+          this->tria_listeners_for_transfer.clear();
 
           // release active and future finite element tables
           this->hp_cell_active_fe_indices.clear();
@@ -2510,10 +2518,11 @@ DoFHandler<dim, spacedim>::distribute_dofs(
         {
           hp_capability_enabled = false;
 
-          // unsubscribe connections to signals
-          for (auto &connection : this->tria_listeners)
+          // unsubscribe connections to signals that are only relevant for
+          // hp-mode, since we only have a single element here
+          for (auto &connection : this->tria_listeners_for_transfer)
             connection.disconnect();
-          this->tria_listeners.clear();
+          this->tria_listeners_for_transfer.clear();
 
           // release active and future finite element tables
           this->hp_cell_active_fe_indices.clear();
@@ -3025,6 +3034,9 @@ template <int dim, int spacedim>
 void
 DoFHandler<dim, spacedim>::connect_to_triangulation_signals()
 {
+  // make sure this is called during initialization in hp-mode
+  Assert(hp_capability_enabled, ExcOnlyAvailableWithHP());
+
   // connect functions to signals of the underlying triangulation
   this->tria_listeners.push_back(this->tria->signals.create.connect(
     [this]() { this->reinit(*(this->tria)); }));
@@ -3038,31 +3050,31 @@ DoFHandler<dim, spacedim>::connect_to_triangulation_signals()
         &this->get_triangulation()))
     {
       // repartitioning signals
-      this->tria_listeners.push_back(
+      this->tria_listeners_for_transfer.push_back(
         this->tria->signals.pre_distributed_repartition.connect([this]() {
           internal::hp::DoFHandlerImplementation::Implementation::
             ensure_absence_of_future_fe_indices<dim, spacedim>(*this);
         }));
-      this->tria_listeners.push_back(
+      this->tria_listeners_for_transfer.push_back(
         this->tria->signals.pre_distributed_repartition.connect(
           [this]() { this->pre_distributed_transfer_action(); }));
-      this->tria_listeners.push_back(
+      this->tria_listeners_for_transfer.push_back(
         this->tria->signals.post_distributed_repartition.connect(
           [this] { this->post_distributed_transfer_action(); }));
 
       // refinement signals
-      this->tria_listeners.push_back(
+      this->tria_listeners_for_transfer.push_back(
         this->tria->signals.pre_distributed_refinement.connect(
           [this]() { this->pre_distributed_transfer_action(); }));
-      this->tria_listeners.push_back(
+      this->tria_listeners_for_transfer.push_back(
         this->tria->signals.post_distributed_refinement.connect(
           [this]() { this->post_distributed_transfer_action(); }));
 
       // serialization signals
-      this->tria_listeners.push_back(
+      this->tria_listeners_for_transfer.push_back(
         this->tria->signals.post_distributed_save.connect(
           [this]() { this->active_fe_index_transfer.reset(); }));
-      this->tria_listeners.push_back(
+      this->tria_listeners_for_transfer.push_back(
         this->tria->signals.post_distributed_load.connect(
           [this]() { this->update_active_fe_table(); }));
     }
@@ -3071,25 +3083,27 @@ DoFHandler<dim, spacedim>::connect_to_triangulation_signals()
              &this->get_triangulation()) != nullptr)
     {
       // partitioning signals
-      this->tria_listeners.push_back(
+      this->tria_listeners_for_transfer.push_back(
         this->tria->signals.pre_partition.connect([this]() {
           internal::hp::DoFHandlerImplementation::Implementation::
             ensure_absence_of_future_fe_indices(*this);
         }));
 
       // refinement signals
-      this->tria_listeners.push_back(this->tria->signals.pre_refinement.connect(
-        [this] { this->pre_transfer_action(); }));
-      this->tria_listeners.push_back(
+      this->tria_listeners_for_transfer.push_back(
+        this->tria->signals.pre_refinement.connect(
+          [this] { this->pre_transfer_action(); }));
+      this->tria_listeners_for_transfer.push_back(
         this->tria->signals.post_refinement.connect(
           [this] { this->post_transfer_action(); }));
     }
   else
     {
       // refinement signals
-      this->tria_listeners.push_back(this->tria->signals.pre_refinement.connect(
-        [this] { this->pre_transfer_action(); }));
-      this->tria_listeners.push_back(
+      this->tria_listeners_for_transfer.push_back(
+        this->tria->signals.pre_refinement.connect(
+          [this] { this->pre_transfer_action(); }));
+      this->tria_listeners_for_transfer.push_back(
         this->tria->signals.post_refinement.connect(
           [this] { this->post_transfer_action(); }));
     }
diff --git a/tests/dofs/clear_signal.cc b/tests/dofs/clear_signal.cc
new file mode 100644 (file)
index 0000000..1ed76cb
--- /dev/null
@@ -0,0 +1,118 @@
+// ---------------------------------------------------------------------
+//
+// 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 that DoFHandler::clear() will be called
+// whenever Triangulation::Signals::clear is triggered
+
+
+#include <deal.II/distributed/fully_distributed_tria.h>
+#include <deal.II/distributed/shared_tria.h>
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include "../tests.h"
+
+
+void
+test_serial()
+{
+  constexpr const unsigned dim = 2;
+
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria);
+
+  DoFHandler<dim> dh(tria);
+  dh.distribute_dofs(FE_Q<dim>(1));
+
+  tria.clear();
+  deallog << "ndofs:" << dh.n_dofs() << std::endl;
+}
+
+
+void
+test_parallel_shared()
+{
+  constexpr const unsigned dim = 2;
+
+  parallel::shared::Triangulation<dim> tria(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(tria);
+
+  DoFHandler<dim> dh(tria);
+  dh.distribute_dofs(FE_Q<dim>(1));
+
+  tria.clear();
+  deallog << "ndofs:" << dh.n_dofs() << std::endl;
+}
+
+
+void
+test_parallel_distributed()
+{
+  constexpr const unsigned dim = 2;
+
+  parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(tria);
+
+  DoFHandler<dim> dh(tria);
+  dh.distribute_dofs(FE_Q<dim>(1));
+
+  tria.clear();
+  deallog << "ndofs:" << dh.n_dofs() << std::endl;
+}
+
+
+void
+test_parallel_fullydistributed()
+{
+  constexpr const unsigned dim = 2;
+
+  parallel::distributed::Triangulation<dim> tria_base(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(tria_base);
+
+  const auto description =
+    TriangulationDescription::Utilities::create_description_from_triangulation(
+      tria_base, MPI_COMM_WORLD);
+
+  parallel::fullydistributed::Triangulation<dim> tria(MPI_COMM_WORLD);
+  tria.create_triangulation(description);
+
+  DoFHandler<dim> dh(tria);
+  dh.distribute_dofs(FE_Q<dim>(1));
+
+  tria.clear();
+  deallog << "ndofs:" << dh.n_dofs() << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+  initlog();
+
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  test_serial();
+  test_parallel_shared();
+  test_parallel_distributed();
+  test_parallel_fullydistributed();
+}
diff --git a/tests/dofs/clear_signal.with_mpi=true.with_p4est=true.mpirun=1.output b/tests/dofs/clear_signal.with_mpi=true.with_p4est=true.mpirun=1.output
new file mode 100644 (file)
index 0000000..6c2ac32
--- /dev/null
@@ -0,0 +1,5 @@
+
+DEAL::ndofs:0
+DEAL::ndofs:0
+DEAL::ndofs:0
+DEAL::ndofs:0
index 1750f1a1c12307927e553a1a21f6304d23f6708a..0aedae0de23ff7c0516a09b7a418f8be1763aab4 100644 (file)
@@ -89,6 +89,7 @@ test(TriangulationType &triangulation)
 
   {
     triangulation.load(filename);
+    dof_handler.distribute_dofs(FE_Q<dim>(2));
 
     parallel::distributed::SolutionTransfer<dim, VectorType> solution_transfer(
       dof_handler);

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.