]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Move active_fe_index_type to types. 14210/head
authorMarc Fehling <mafehling.git@gmail.com>
Sun, 21 Aug 2022 23:56:55 +0000 (17:56 -0600)
committerMarc Fehling <mafehling.git@gmail.com>
Sun, 21 Aug 2022 23:56:55 +0000 (17:56 -0600)
include/deal.II/base/types.h
include/deal.II/dofs/dof_handler.h
source/dofs/dof_handler.cc
source/hp/refinement.cc

index c025f887d475f319517fabbab12078de41b28b98..0b8efbc4812a029d08d8648c60d025609b77fd7c 100644 (file)
@@ -53,6 +53,11 @@ namespace types
    */
 #define DEAL_II_VERTEX_INDEX_MPI_TYPE MPI_UINT64_T
 
+  /**
+   * The type in which we store the active and future FE indices.
+   */
+  using fe_index = unsigned short int;
+
   /**
    * The type used to denote the global index of degrees of freedom. This
    * type is then also used for querying the global *number* of degrees
index 2c9417b7be6b5ecb4b557cc3cc33dad4fd74aa18..f3665f3c3cd898264a28b544a568b9165c5859c3 100644 (file)
@@ -25,6 +25,7 @@
 #include <deal.II/base/index_set.h>
 #include <deal.II/base/iterator_range.h>
 #include <deal.II/base/smartpointer.h>
+#include <deal.II/base/types.h>
 
 #include <deal.II/dofs/block_info.h>
 #include <deal.II/dofs/dof_accessor.h>
@@ -523,8 +524,10 @@ public:
 
   /**
    * The type in which we store the active FE index.
+   *
+   * @deprecated Use types::fe_index instead.
    */
-  using active_fe_index_type = unsigned short int;
+  using active_fe_index_type DEAL_II_DEPRECATED = types::fe_index;
 
   /**
    * The type in which we store the offsets in the CRS data structures.
@@ -535,8 +538,8 @@ public:
    * Invalid active FE index which will be used as a default value to determine
    * whether a future FE index has been set or not.
    */
-  static const active_fe_index_type invalid_active_fe_index =
-    static_cast<active_fe_index_type>(-1);
+  static const types::fe_index invalid_active_fe_index =
+    static_cast<types::fe_index>(-1);
 
   /**
    * Standard constructor, not initializing any data. After constructing an
@@ -1489,7 +1492,7 @@ private:
    * of the appropriate position of a cell in the vectors is done via
    * hp_object_fe_ptr (CRS scheme).
    */
-  mutable std::array<std::vector<active_fe_index_type>, dim + 1>
+  mutable std::array<std::vector<types::fe_index>, dim + 1>
     hp_object_fe_indices;
 
   /**
@@ -1501,15 +1504,13 @@ private:
    * Active FE index of an active cell (identified by level and level index).
    * This vector is only used in hp-mode.
    */
-  mutable std::vector<std::vector<active_fe_index_type>>
-    hp_cell_active_fe_indices;
+  mutable std::vector<std::vector<types::fe_index>> hp_cell_active_fe_indices;
 
   /**
    * Future FE index of an active cell (identified by level and level index).
    * This vector is only used in hp-mode.
    */
-  mutable std::vector<std::vector<active_fe_index_type>>
-    hp_cell_future_fe_indices;
+  mutable std::vector<std::vector<types::fe_index>> hp_cell_future_fe_indices;
 
   /**
    * An array to store the indices for level degrees of freedom located at
index bd2b7bc949c5c66242ee8d1f62e849cc3a170a91..1e2a06f6b8686b9a038497bd71bd2efa8e3be398 100644 (file)
@@ -47,8 +47,7 @@ const unsigned int DoFHandler<dim, spacedim>::default_fe_index;
 
 
 template <int dim, int spacedim>
-const typename DoFHandler<dim, spacedim>::active_fe_index_type
-  DoFHandler<dim, spacedim>::invalid_active_fe_index;
+const types::fe_index DoFHandler<dim, spacedim>::invalid_active_fe_index;
 
 
 namespace internal
@@ -1322,9 +1321,6 @@ namespace internal
             dof_handler.hp_capability_enabled == true,
             (typename DoFHandler<dim, spacedim>::ExcOnlyAvailableWithHP()));
 
-          using active_fe_index_type =
-            typename dealii::DoFHandler<dim, spacedim>::active_fe_index_type;
-
           if (const dealii::parallel::shared::Triangulation<dim, spacedim> *tr =
                 dynamic_cast<
                   const dealii::parallel::shared::Triangulation<dim, spacedim>
@@ -1340,7 +1336,7 @@ namespace internal
               // on the other cells to zero. then we add all of these vectors
               // up, and because every vector entry has exactly one processor
               // that owns it, the sum is correct
-              std::vector<active_fe_index_type> active_fe_indices(
+              std::vector<types::fe_index> active_fe_indices(
                 tr->n_active_cells(), 0u);
               for (const auto &cell : dof_handler.active_cell_iterators())
                 if (cell->is_locally_owned())
@@ -1376,17 +1372,15 @@ namespace internal
               // to have functions that can pack and unpack the data we want to
               // transport -- namely, the single unsigned int active_fe_index
               // objects
-              auto pack =
-                [](const typename dealii::DoFHandler<dim, spacedim>::
-                     active_cell_iterator &cell) -> active_fe_index_type {
+              auto pack = [](const typename dealii::DoFHandler<dim, spacedim>::
+                               active_cell_iterator &cell) -> types::fe_index {
                 return cell->active_fe_index();
               };
 
-              auto unpack =
-                [&dof_handler](
-                  const typename dealii::DoFHandler<dim, spacedim>::
-                    active_cell_iterator &   cell,
-                  const active_fe_index_type active_fe_index) -> void {
+              auto unpack = [&dof_handler](
+                              const typename dealii::DoFHandler<dim, spacedim>::
+                                active_cell_iterator &cell,
+                              const types::fe_index   active_fe_index) -> void {
                 // we would like to say
                 //   cell->set_active_fe_index(active_fe_index);
                 // but this is not allowed on cells that are not
@@ -1397,7 +1391,7 @@ namespace internal
               };
 
               GridTools::exchange_cell_data_to_ghosts<
-                active_fe_index_type,
+                types::fe_index,
                 dealii::DoFHandler<dim, spacedim>>(dof_handler, pack, unpack);
             }
           else
@@ -1434,15 +1428,12 @@ namespace internal
             dof_handler.hp_capability_enabled == true,
             (typename DoFHandler<dim, spacedim>::ExcOnlyAvailableWithHP()));
 
-          using active_fe_index_type =
-            typename dealii::DoFHandler<dim, spacedim>::active_fe_index_type;
-
           if (const dealii::parallel::shared::Triangulation<dim, spacedim> *tr =
                 dynamic_cast<
                   const dealii::parallel::shared::Triangulation<dim, spacedim>
                     *>(&dof_handler.get_triangulation()))
             {
-              std::vector<active_fe_index_type> future_fe_indices(
+              std::vector<types::fe_index> future_fe_indices(
                 tr->n_active_cells(), 0u);
               for (const auto &cell : dof_handler.active_cell_iterators() |
                                         IteratorFilters::LocallyOwnedCell())
@@ -1467,26 +1458,24 @@ namespace internal
                            DistributedTriangulationBase<dim, spacedim> *>(
                          &dof_handler.get_triangulation()))
             {
-              auto pack =
-                [&dof_handler](
-                  const typename dealii::DoFHandler<dim, spacedim>::
-                    active_cell_iterator &cell) -> active_fe_index_type {
+              auto pack = [&dof_handler](
+                            const typename dealii::DoFHandler<dim, spacedim>::
+                              active_cell_iterator &cell) -> types::fe_index {
                 return dof_handler
                   .hp_cell_future_fe_indices[cell->level()][cell->index()];
               };
 
-              auto unpack =
-                [&dof_handler](
-                  const typename dealii::DoFHandler<dim, spacedim>::
-                    active_cell_iterator &   cell,
-                  const active_fe_index_type future_fe_index) -> void {
+              auto unpack = [&dof_handler](
+                              const typename dealii::DoFHandler<dim, spacedim>::
+                                active_cell_iterator &cell,
+                              const types::fe_index   future_fe_index) -> void {
                 dof_handler
                   .hp_cell_future_fe_indices[cell->level()][cell->index()] =
                   future_fe_index;
               };
 
               GridTools::exchange_cell_data_to_ghosts<
-                active_fe_index_type,
+                types::fe_index,
                 dealii::DoFHandler<dim, spacedim>>(dof_handler, pack, unpack);
             }
           else
index cd3b3a70a6ab0cc26bfabe7630613efea25146ae..a2b8c9fa96bfd5e4a04ddd8cc2e5980c32ca64c5 100644 (file)
@@ -811,8 +811,7 @@ namespace hp
 
       // there can be as many levels in the hierarchy as active FE indices are
       // possible
-      using level_type =
-        typename DoFHandler<dim, spacedim>::active_fe_index_type;
+      using level_type         = types::fe_index;
       const auto invalid_level = static_cast<level_type>(-1);
 
       // map from FE index to level in hierarchy

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.