]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added 'set_fe' functionality to DoFHandlers.
authorMarc Fehling <marc.fehling@gmx.net>
Sun, 3 Feb 2019 00:24:43 +0000 (01:24 +0100)
committerMarc Fehling <marc.fehling@gmx.net>
Thu, 28 Feb 2019 12:57:39 +0000 (13:57 +0100)
12 files changed:
doc/news/changes/minor/20190210MarcFehling [new file with mode: 0644]
include/deal.II/distributed/solution_transfer.h
include/deal.II/dofs/dof_handler.h
include/deal.II/hp/dof_handler.h
source/dofs/dof_handler.cc
source/hp/dof_handler.cc
tests/mpi/hp_active_fe_indices_transfer_01.cc
tests/mpi/hp_active_fe_indices_transfer_02.cc
tests/mpi/hp_cell_weights_01.cc
tests/mpi/hp_cell_weights_02.cc
tests/mpi/hp_cell_weights_03.cc
tests/mpi/p4est_save_06.cc

diff --git a/doc/news/changes/minor/20190210MarcFehling b/doc/news/changes/minor/20190210MarcFehling
new file mode 100644 (file)
index 0000000..9a3e1e9
--- /dev/null
@@ -0,0 +1,6 @@
+New: Member functions DoFHandler::set_fe() and hp::DoFHandler::set_fe()
+that will register a finite element object without enumerating all
+degrees of freedom. In the hp case, the active_fe_indices will be
+initialized and communicated amongst all processors, additionally.
+<br>
+(2019/02/10, Marc Fehling)
index 558ebeed29663ba89e1376c38f201846e648cb93..9d61878b2770cdbe0f354059b8ca47e8ceee0087 100644 (file)
@@ -203,7 +203,7 @@ namespace parallel
      * hp::DoFHandler<dim,spacedim> hp_dof_handler(triangulation);
      * // We need to introduce our dof_handler to the fe_collection
      * // before setting all active_fe_indices.
-     * hp_dof_handler.distribute_dofs(fe_collection);
+     * hp_dof_handler.set_fe(fe_collection);
      * hp_dof_handler.deserialize_active_fe_indices();
      * hp_dof_handler.distribute_dofs(fe_collection);
      *
index 638e98b535eef00536755d17bd87aaafbacd4ab2..1da5931739ac8b6f47d3119ed3d3890c2e77b2d9 100644 (file)
@@ -472,6 +472,25 @@ public:
   initialize(const Triangulation<dim, spacedim> &tria,
              const FiniteElement<dim, spacedim> &fe);
 
+  /**
+   * Assign a FiniteElement @p fe to this object.
+   *
+   * @note This function makes a copy of the finite element given as
+   * argument, and stores it as a member variable. Consequently, it is
+   * possible to write code such as
+   * @code
+   *   dof_handler.set_fe(FE_Q<dim>(2));
+   * @endcode
+   * You can then access the finite element later on by calling
+   * DoFHandler::get_fe(). However, it is often more convenient to
+   * keep a named finite element object as a member variable in your
+   * main class and refer to it directly whenever you need to access
+   * properties of the finite element (such as
+   * FiniteElementData::dofs_per_cell). This is what all tutorial programs do.
+   */
+  virtual void
+  set_fe(const FiniteElement<dim, spacedim> &fe);
+
   /**
    * Go through the triangulation and "distribute" the degrees of
    * freedom needed for the given finite element. "Distributing"
@@ -498,17 +517,8 @@ public:
    * step-2 tutorial program.
    *
    * @note This function makes a copy of the finite element given as
-   * argument, and stores it as a member variable. Consequently, it is
-   * possible to write code such as
-   * @code
-   *   dof_handler.distribute_dofs (FE_Q<dim>(2));
-   * @endcode
-   * You can then access the finite element later on by calling
-   * DoFHandler::get_fe().  However, it is often more convenient to
-   * keep a named finite element object as a member variable in your
-   * main class and refer to it directly whenever you need to access
-   * properties of the finite element (such as FiniteElement::dofs_per_cell).
-   * This is what all tutorial programs do.
+   * argument, and stores it as a member variable, similarly to the above
+   * function set_fe().
    */
   virtual void
   distribute_dofs(const FiniteElement<dim, spacedim> &fe);
index 29a3746ec89af241505dbfeff6d4461bdaec3cf2..3b73789941c6437b1f7867e36820bfc7a4d95d03 100644 (file)
@@ -417,6 +417,19 @@ namespace hp
     initialize(const Triangulation<dim, spacedim> &   tria,
                const hp::FECollection<dim, spacedim> &fe);
 
+    /**
+     * Assign a hp::FECollection @p fe to this object.
+     *
+     * In case a parallel::Triangulation is assigned to this object,
+     * the active_fe_indices will be exchanged between processors so that
+     * each one knows the indices on its own cells and all ghost cells.
+     *
+     * @note In accordance with dealii::DoFHandler::set_fe(),
+     * this function also makes a copy of the object given as argument.
+     */
+    virtual void
+    set_fe(const hp::FECollection<dim, spacedim> &fe);
+
     /**
      * Go through the triangulation and "distribute" the degrees of
      * freedom needed for the given finite element. "Distributing"
index 03603ee5a05eed3309d3eb0ea441a4d6f6b42df4..f5b0a0fd206422e1245457416219cb264304bf4c 100644 (file)
@@ -1197,6 +1197,18 @@ DoFHandler<dim, spacedim>::memory_consumption() const
 
 
 
+template <int dim, int spacedim>
+void
+DoFHandler<dim, spacedim>::set_fe(const FiniteElement<dim, spacedim> &ff)
+{
+  // Only recreate the FECollection if we don't already store
+  // the exact same FiniteElement object.
+  if (fe_collection.size() == 0 || fe_collection[0] != ff)
+    fe_collection = hp::FECollection<dim, spacedim>(ff);
+}
+
+
+
 template <int dim, int spacedim>
 void
 DoFHandler<dim, spacedim>::distribute_dofs(
@@ -1210,10 +1222,8 @@ DoFHandler<dim, spacedim>::distribute_dofs(
   Assert(tria->n_levels() > 0,
          ExcMessage("The Triangulation you are using is empty!"));
 
-  // Only recreate the FECollection if we don't already store
-  // the exact same FiniteElement object.
-  if (fe_collection.size() == 0 || fe_collection[0] != ff)
-    fe_collection = hp::FECollection<dim, spacedim>(ff);
+  // first, assign the finite_element
+  set_fe(ff);
 
   // delete all levels and set them
   // up newly. note that we still
index 01900459d885d3857aa5a833cd36dc57b223a6a9..35e4a176664916ae943219f016dfb8a57d913a7c 100644 (file)
@@ -1390,6 +1390,33 @@ namespace hp
 
 
 
+  template <int dim, int spacedim>
+  void
+  DoFHandler<dim, spacedim>::set_fe(const hp::FECollection<dim, spacedim> &ff)
+  {
+    // don't create a new object if the one we have is already appropriate
+    if (fe_collection != ff)
+      fe_collection = hp::FECollection<dim, spacedim>(ff);
+
+    // ensure that the active_fe_indices vectors are initialized correctly
+    create_active_fe_table();
+
+    // make sure every processor knows the active_fe_indices
+    // on both its own cells and all ghost cells
+    dealii::internal::hp::DoFHandlerImplementation::Implementation::
+      communicate_active_fe_indices(*this);
+
+    // make sure that the fe collection is large enough to
+    // cover all fe indices presently in use on the mesh
+    for (const auto &cell : active_cell_iterators())
+      if (cell->is_locally_owned())
+        Assert(cell->active_fe_index() < fe_collection.size(),
+               ExcInvalidFEIndex(cell->active_fe_index(),
+                                 fe_collection.size()));
+  }
+
+
+
   template <int dim, int spacedim>
   void
   DoFHandler<dim, spacedim>::distribute_dofs(
@@ -1404,14 +1431,8 @@ namespace hp
            ExcMessage("The Triangulation you are using is empty!"));
     Assert(ff.size() > 0, ExcMessage("The hp::FECollection given is empty!"));
 
-    // don't create a new object if the one we have is already appropriate
-    if (fe_collection != ff)
-      fe_collection = hp::FECollection<dim, spacedim>(ff);
-
-    // at the beginning, make sure every processor knows the
-    // active_fe_indices on both its own cells and all ghost cells
-    dealii::internal::hp::DoFHandlerImplementation::Implementation::
-      communicate_active_fe_indices(*this);
+    // first, assign the fe_collection
+    set_fe(ff);
 
     // If an underlying shared::Tria allows artificial cells,
     // then save the current set of subdomain ids, and set
@@ -1425,51 +1446,29 @@ namespace hp
         {
           saved_subdomain_ids.resize(shared_tria->n_active_cells());
 
-          typename parallel::shared::Triangulation<dim, spacedim>::
-            active_cell_iterator cell = get_triangulation().begin_active(),
-                                 endc = get_triangulation().end();
-
           const std::vector<types::subdomain_id> &true_subdomain_ids =
             shared_tria->get_true_subdomain_ids_of_cells();
 
-          for (unsigned int index = 0; cell != endc; ++cell, ++index)
+          for (const auto &cell : shared_tria->active_cell_iterators())
             {
+              const unsigned int index   = cell->active_cell_index();
               saved_subdomain_ids[index] = cell->subdomain_id();
               cell->set_subdomain_id(true_subdomain_ids[index]);
             }
         }
 
-    // ensure that the active_fe_indices vectors are
-    // initialized correctly.
-    create_active_fe_table();
-
-    // up front make sure that the fe collection is large enough to
-    // cover all fe indices presently in use on the mesh
-    for (active_cell_iterator cell = begin_active(); cell != end(); ++cell)
-      if (cell->is_locally_owned())
-        Assert(cell->active_fe_index() < fe_collection.size(),
-               ExcInvalidFEIndex(cell->active_fe_index(),
-                                 fe_collection.size()));
-
-
-    // then allocate space for all the other tables
+    // then allocate space for all tables
     dealii::internal::hp::DoFHandlerImplementation::Implementation::
       reserve_space(*this);
 
-
     // now undo the subdomain modification
     if (const parallel::shared::Triangulation<dim, spacedim> *shared_tria =
           (dynamic_cast<const parallel::shared::Triangulation<dim, spacedim> *>(
             &get_triangulation())))
       if (shared_tria->with_artificial_cells())
-        {
-          typename parallel::shared::Triangulation<dim, spacedim>::
-            active_cell_iterator cell = get_triangulation().begin_active(),
-                                 endc = get_triangulation().end();
-
-          for (unsigned int index = 0; cell != endc; ++cell, ++index)
-            cell->set_subdomain_id(saved_subdomain_ids[index]);
-        }
+        for (const auto &cell : shared_tria->active_cell_iterators())
+          cell->set_subdomain_id(
+            saved_subdomain_ids[cell->active_cell_index()]);
 
 
     // Clear user flags because we will need them. But first we save
index 2d5cb00ab231e3b7d73238eded0614ada3c31e76..5fd0cc5108dc98a80b79d640be91655259a94943 100644 (file)
@@ -50,9 +50,8 @@ test()
   for (unsigned int i = 0; i < max_degree; ++i)
     fe_collection.push_back(FE_Q<dim>(max_degree - i));
 
-  // this distribute_dofs() call is necessary
   // we need to introduce dof_handler to its fe_collection first
-  dh.distribute_dofs(fe_collection);
+  dh.set_fe(fe_collection);
 
   unsigned int i = 0;
   for (auto &cell : dh.active_cell_iterators())
index 65393457401ee65135039e296a1cc130b386aa93..e40dc29a2b67db6da863894bb204f3333e4c0b82 100644 (file)
@@ -48,10 +48,9 @@ test()
     GridGenerator::subdivided_hyper_cube(tria, 2);
     tria.refine_global(1);
 
-    // this distribute_dofs() call is necessary
-    // we need to introduce dof_handler to its fe_collection first
     hp::DoFHandler<dim> dh(tria);
-    dh.distribute_dofs(fe_collection);
+    // we need to introduce dof_handler to its fe_collection first
+    dh.set_fe(fe_collection);
 
     unsigned int i = 0;
     for (auto &cell : dh.active_cell_iterators())
@@ -85,10 +84,9 @@ test()
     GridGenerator::subdivided_hyper_cube(tria, 2);
     // triangulation has to be initialized with correct coarse cells
 
-    // this distribute_dofs() call is necessary
     // we need to introduce dof_handler to its fe_collection first
     hp::DoFHandler<dim> dh(tria);
-    dh.distribute_dofs(fe_collection);
+    dh.set_fe(fe_collection);
 
     // ----- transfer -----
     tria.load("file");
index 6d15d5841b328828b8c3f240ab77dc2924aefc00..6edf11841294cadeb359d1e96e6b72df7a33d97b 100644 (file)
@@ -57,17 +57,12 @@ test()
   fe_collection.push_back(FE_Q<dim>(5));
 
   hp::DoFHandler<dim> dh(tria);
+  dh.set_fe(fe_collection);
   // default: active_fe_index = 0
   for (auto &cell : dh.active_cell_iterators())
     if (cell->is_locally_owned())
       if (cell->id().to_string() == "0_2:00")
         cell->set_active_fe_index(1);
-  dh.distribute_dofs(fe_collection);
-
-
-  parallel::CellWeights<dim> cell_weights(dh);
-  cell_weights.register_ndofs_weighting(100000);
-
 
   deallog << "Number of cells before repartitioning: "
           << tria.n_locally_owned_active_cells() << std::endl;
@@ -80,8 +75,10 @@ test()
   }
 
 
+  parallel::CellWeights<dim> cell_weights(dh);
+  cell_weights.register_ndofs_weighting(100000);
+
   tria.repartition();
-  dh.distribute_dofs(fe_collection);
 
 
   deallog << "Number of cells after repartitioning: "
index 5d13dff7746911f7e66fd77969cb0f838f9db522..e116925d70dec14b84aeb3a57c47def2bea98e3f 100644 (file)
@@ -63,17 +63,12 @@ test()
   fe_collection.push_back(FE_Q<dim>(7));
 
   hp::DoFHandler<dim> dh(tria);
+  dh.set_fe(fe_collection);
   // default: active_fe_index = 0
   for (auto &cell : dh.active_cell_iterators())
     if (cell->is_locally_owned())
       if (cell->id().to_string() == "0_3:000")
         cell->set_active_fe_index(1);
-  dh.distribute_dofs(fe_collection);
-
-
-  parallel::CellWeights<dim> cell_weights(dh);
-  cell_weights.register_ndofs_weighting(100000);
-
 
   deallog << "Number of cells before repartitioning: "
           << tria.n_locally_owned_active_cells() << std::endl;
@@ -86,8 +81,10 @@ test()
   }
 
 
+  parallel::CellWeights<dim> cell_weights(dh);
+  cell_weights.register_ndofs_weighting(100000);
+
   tria.repartition();
-  dh.distribute_dofs(fe_collection);
 
 
   deallog << "Number of cells after repartitioning: "
index c2eccc249f9831335405c616192c3071abfcebfd..afe9d5be68d917a87ba3a7f0fc2d19b7af5a9961 100644 (file)
@@ -63,17 +63,12 @@ test()
   fe_collection.push_back(FE_Q<dim>(5));
 
   hp::DoFHandler<dim> dh(tria);
+  dh.set_fe(fe_collection);
   // default: active_fe_index = 0
   for (auto &cell : dh.active_cell_iterators())
     if (cell->is_locally_owned())
       if (cell->id().to_string() == "0_2:00")
         cell->set_active_fe_index(1);
-  dh.distribute_dofs(fe_collection);
-
-
-  parallel::CellWeights<dim> cell_weights(dh);
-  cell_weights.register_ndofs_weighting(100000);
-
 
   deallog << "Number of cells before repartitioning: "
           << tria.n_locally_owned_active_cells() << std::endl;
@@ -86,9 +81,11 @@ test()
   }
 
 
+  parallel::CellWeights<dim> cell_weights(dh);
+  cell_weights.register_ndofs_weighting(100000);
+
   // we didn't mark any cells, but we want to repartition our domain
   tria.execute_coarsening_and_refinement();
-  dh.distribute_dofs(fe_collection);
 
 
   deallog << "Number of cells after repartitioning: "
index 5af46a86964b7a871d5d88c2a2b0043f83c1472a..453ed583896d3190f28726cb1f185edea8230d2a 100644 (file)
@@ -134,7 +134,7 @@ test()
     for (unsigned int i = 0; i < max_degree; ++i)
       fe_collection.push_back(FE_Q<dim>(max_degree - i));
 
-    dh.distribute_dofs(fe_collection);
+    dh.set_fe(fe_collection);
     dh.deserialize_active_fe_indices();
     dh.distribute_dofs(fe_collection);
 

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.