]> https://gitweb.dealii.org/ - dealii.git/commitdiff
MF: allow to get active_fe_index and category for specific dof index 18350/head
authorNils Much <nils.much@tum.de>
Thu, 10 Apr 2025 16:52:24 +0000 (18:52 +0200)
committerNils Much <nils.much@tum.de>
Thu, 10 Apr 2025 16:52:24 +0000 (18:52 +0200)
include/deal.II/matrix_free/fe_evaluation.h
include/deal.II/matrix_free/matrix_free.h
include/deal.II/matrix_free/tools.h

index 496bae9dd7d683740724113480e485fc4f5327d1..5764fc3d996d436c811667f0c7ed75b6cb20e9aa 100644 (file)
@@ -6160,9 +6160,7 @@ inline FEEvaluation<dim,
                  dof_no,
                  quad_no,
                  first_selected_component,
-                 (matrix_free.first_hp_dof_handler_index == dof_no) ?
-                   matrix_free.get_cell_active_fe_index(range) :
-                   0)
+                 matrix_free.get_cell_active_fe_index(range, dof_no))
 {}
 
 
index aef158d9e58fe3b0aeac0e0a82114de000904b15..a7b6fa85ad08853ef0a3a80e81231e6c1a2cd15c 100644 (file)
@@ -20,6 +20,7 @@
 
 #include <deal.II/base/aligned_vector.h>
 #include <deal.II/base/exceptions.h>
+#include <deal.II/base/numbers.h>
 #include <deal.II/base/quadrature.h>
 #include <deal.II/base/template_constraints.h>
 #include <deal.II/base/thread_local_storage.h>
@@ -1607,14 +1608,17 @@ public:
    */
   unsigned int
   get_cell_active_fe_index(
-    const std::pair<unsigned int, unsigned int> range) const;
+    const std::pair<unsigned int, unsigned int> range,
+    const unsigned int dof_handler_index = numbers::invalid_unsigned_int) const;
 
   /**
    * In the hp-adaptive case, return the active FE index of a face range.
    */
   unsigned int
-  get_face_active_fe_index(const std::pair<unsigned int, unsigned int> range,
-                           const bool is_interior_face = true) const;
+  get_face_active_fe_index(
+    const std::pair<unsigned int, unsigned int> range,
+    const bool                                  is_interior_face = true,
+    const unsigned int dof_handler_index = numbers::invalid_unsigned_int) const;
 
   /** @} */
 
@@ -2008,7 +2012,8 @@ public:
    */
   unsigned int
   get_cell_range_category(
-    const std::pair<unsigned int, unsigned int> cell_batch_range) const;
+    const std::pair<unsigned int, unsigned int> cell_batch_range,
+    const unsigned int dof_handler_index = numbers::invalid_unsigned_int) const;
 
   /**
    * Return the category of the cells on the two sides of the current batch
@@ -2016,7 +2021,8 @@ public:
    */
   std::pair<unsigned int, unsigned int>
   get_face_range_category(
-    const std::pair<unsigned int, unsigned int> face_batch_range) const;
+    const std::pair<unsigned int, unsigned int> face_batch_range,
+    const unsigned int dof_handler_index = numbers::invalid_unsigned_int) const;
 
   /**
    * Return the category the current batch of cells was assigned to. Categories
@@ -2030,14 +2036,18 @@ public:
    * enabled, it is guaranteed that all cells have the same category.
    */
   unsigned int
-  get_cell_category(const unsigned int cell_batch_index) const;
+  get_cell_category(
+    const unsigned int cell_batch_index,
+    const unsigned int dof_handler_index = numbers::invalid_unsigned_int) const;
 
   /**
    * Return the category of the cells on the two sides of the current batch of
    * faces.
    */
   std::pair<unsigned int, unsigned int>
-  get_face_category(const unsigned int face_batch_index) const;
+  get_face_category(
+    const unsigned int face_batch_index,
+    const unsigned int dof_handler_index = numbers::invalid_unsigned_int) const;
 
   /**
    * Queries whether or not the indexation has been set.
@@ -2353,7 +2363,6 @@ private:
    * Stores the index of the first DoFHandler that is in hp-mode. If no
    * DoFHandler is in hp-mode, the value is 0.
    */
-public:
   unsigned int first_hp_dof_handler_index;
 };
 
@@ -2685,13 +2694,18 @@ MatrixFree<dim, Number, VectorizedArrayType>::n_active_fe_indices() const
 template <int dim, typename Number, typename VectorizedArrayType>
 unsigned int
 MatrixFree<dim, Number, VectorizedArrayType>::get_cell_active_fe_index(
-  const std::pair<unsigned int, unsigned int> range) const
+  const std::pair<unsigned int, unsigned int> range,
+  const unsigned int                          dof_handler_index) const
 {
-  const auto &fe_indices =
-    dof_info[first_hp_dof_handler_index].cell_active_fe_index;
+  const unsigned int dof_no =
+    dof_handler_index == numbers::invalid_unsigned_int ?
+      first_hp_dof_handler_index :
+      dof_handler_index;
+
+  const auto &fe_indices = dof_info[dof_no].cell_active_fe_index;
 
   if (fe_indices.empty() == true ||
-      dof_handlers[first_hp_dof_handler_index]->get_fe_collection().size() == 1)
+      dof_handlers[dof_no]->get_fe_collection().size() == 1)
     return 0;
 
   const auto index = fe_indices[range.first];
@@ -2708,10 +2722,15 @@ template <int dim, typename Number, typename VectorizedArrayType>
 unsigned int
 MatrixFree<dim, Number, VectorizedArrayType>::get_face_active_fe_index(
   const std::pair<unsigned int, unsigned int> range,
-  const bool                                  is_interior_face) const
+  const bool                                  is_interior_face,
+  const unsigned int                          dof_handler_index) const
 {
-  const auto &fe_indices =
-    dof_info[first_hp_dof_handler_index].cell_active_fe_index;
+  const unsigned int dof_no =
+    dof_handler_index == numbers::invalid_unsigned_int ?
+      first_hp_dof_handler_index :
+      dof_handler_index;
+
+  const auto &fe_indices = dof_info[dof_no].cell_active_fe_index;
 
   if (fe_indices.empty() == true)
     return 0;
@@ -2923,12 +2942,13 @@ MatrixFree<dim, Number, VectorizedArrayType>::get_face_quadrature(
 template <int dim, typename Number, typename VectorizedArrayType>
 inline unsigned int
 MatrixFree<dim, Number, VectorizedArrayType>::get_cell_range_category(
-  const std::pair<unsigned int, unsigned int> range) const
+  const std::pair<unsigned int, unsigned int> range,
+  const unsigned int                          dof_handler_index) const
 {
-  auto result = get_cell_category(range.first);
+  auto result = get_cell_category(range.first, dof_handler_index);
 
   for (unsigned int i = range.first; i < range.second; ++i)
-    result = std::max(result, get_cell_category(i));
+    result = std::max(result, get_cell_category(i, dof_handler_index));
 
   return result;
 }
@@ -2938,14 +2958,17 @@ MatrixFree<dim, Number, VectorizedArrayType>::get_cell_range_category(
 template <int dim, typename Number, typename VectorizedArrayType>
 inline std::pair<unsigned int, unsigned int>
 MatrixFree<dim, Number, VectorizedArrayType>::get_face_range_category(
-  const std::pair<unsigned int, unsigned int> range) const
+  const std::pair<unsigned int, unsigned int> range,
+  const unsigned int                          dof_handler_index) const
 {
-  auto result = get_face_category(range.first);
+  auto result = get_face_category(range.first, dof_handler_index);
 
   for (unsigned int i = range.first; i < range.second; ++i)
     {
-      result.first  = std::max(result.first, get_face_category(i).first);
-      result.second = std::max(result.second, get_face_category(i).second);
+      result.first =
+        std::max(result.first, get_face_category(i, dof_handler_index).first);
+      result.second =
+        std::max(result.second, get_face_category(i, dof_handler_index).second);
     }
 
   return result;
@@ -2956,17 +2979,22 @@ MatrixFree<dim, Number, VectorizedArrayType>::get_face_range_category(
 template <int dim, typename Number, typename VectorizedArrayType>
 inline unsigned int
 MatrixFree<dim, Number, VectorizedArrayType>::get_cell_category(
-  const unsigned int cell_batch_index) const
+  const unsigned int cell_batch_index,
+  const unsigned int dof_handler_index) const
 {
   AssertIndexRange(0, dof_info.size());
-  AssertIndexRange(
-    cell_batch_index,
-    dof_info[first_hp_dof_handler_index].cell_active_fe_index.size());
-  if (dof_info[first_hp_dof_handler_index].cell_active_fe_index.empty())
+
+  const unsigned int dof_no =
+    dof_handler_index == numbers::invalid_unsigned_int ?
+      first_hp_dof_handler_index :
+      dof_handler_index;
+
+  AssertIndexRange(cell_batch_index,
+                   dof_info[dof_no].cell_active_fe_index.size());
+  if (dof_info[dof_no].cell_active_fe_index.empty())
     return 0;
   else
-    return dof_info[first_hp_dof_handler_index]
-      .cell_active_fe_index[cell_batch_index];
+    return dof_info[dof_no].cell_active_fe_index[cell_batch_index];
 }
 
 
@@ -2974,10 +3002,16 @@ MatrixFree<dim, Number, VectorizedArrayType>::get_cell_category(
 template <int dim, typename Number, typename VectorizedArrayType>
 inline std::pair<unsigned int, unsigned int>
 MatrixFree<dim, Number, VectorizedArrayType>::get_face_category(
-  const unsigned int face_batch_index) const
+  const unsigned int face_batch_index,
+  const unsigned int dof_handler_index) const
 {
+  const unsigned int dof_no =
+    dof_handler_index == numbers::invalid_unsigned_int ?
+      first_hp_dof_handler_index :
+      dof_handler_index;
+
   AssertIndexRange(face_batch_index, face_info.faces.size());
-  if (dof_info[first_hp_dof_handler_index].cell_active_fe_index.empty())
+  if (dof_info[dof_no].cell_active_fe_index.empty())
     return std::make_pair(0U, 0U);
 
   std::pair<unsigned int, unsigned int> result = std::make_pair(0U, 0U);
@@ -2986,11 +3020,11 @@ MatrixFree<dim, Number, VectorizedArrayType>::get_face_category(
        face_info.faces[face_batch_index].cells_interior[v] !=
          numbers::invalid_unsigned_int;
        ++v)
-    result.first =
-      std::max(result.first,
-               dof_info[first_hp_dof_handler_index].cell_active_fe_index
-                 [face_info.faces[face_batch_index].cells_interior[v] /
-                  VectorizedArrayType::size()]);
+    result.first = std::max(
+      result.first,
+      dof_info[dof_no].cell_active_fe_index[face_info.faces[face_batch_index]
+                                              .cells_interior[v] /
+                                            VectorizedArrayType::size()]);
   if (face_info.faces[face_batch_index].cells_exterior[0] !=
       numbers::invalid_unsigned_int)
     for (unsigned int v = 0;
@@ -2998,11 +3032,11 @@ MatrixFree<dim, Number, VectorizedArrayType>::get_face_category(
          face_info.faces[face_batch_index].cells_exterior[v] !=
            numbers::invalid_unsigned_int;
          ++v)
-      result.second =
-        std::max(result.second,
-                 dof_info[first_hp_dof_handler_index].cell_active_fe_index
-                   [face_info.faces[face_batch_index].cells_exterior[v] /
-                    VectorizedArrayType::size()]);
+      result.second = std::max(
+        result.second,
+        dof_info[dof_no].cell_active_fe_index[face_info.faces[face_batch_index]
+                                                .cells_exterior[v] /
+                                              VectorizedArrayType::size()]);
   else
     result.second = numbers::invalid_unsigned_int;
   return result;
index 2740936df3c54bd1891902150b9cd9cd380aa895..0ada9ef47266b1f1e8d827afc6037b6d6a9b65f0 100644 (file)
@@ -1352,13 +1352,13 @@ namespace MatrixFreeTools
 
       unsigned int active_fe_index = 0;
       if (!is_face)
-        active_fe_index = (matrix_free.first_hp_dof_handler_index == dof_no) ?
-                            matrix_free.get_cell_active_fe_index(range) :
-                            0;
+        active_fe_index = matrix_free.get_cell_active_fe_index(range, dof_no);
       else if (is_interior_face)
-        active_fe_index = matrix_free.get_face_range_category(range).first;
+        active_fe_index =
+          matrix_free.get_face_range_category(range, dof_no).first;
       else
-        active_fe_index = matrix_free.get_face_range_category(range).second;
+        active_fe_index =
+          matrix_free.get_face_range_category(range, dof_no).second;
 
       const auto init_data = dealii::internal::
         extract_initialization_data<is_face, dim, Number, VectorizedArrayType>(

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.