From 74c69dcda6927121b517fe0da99384f5afed34c5 Mon Sep 17 00:00:00 2001
From: Martin Kronbichler <martin.kronbichler@uni-a.de>
Date: Mon, 26 Sep 2022 12:04:27 +0200
Subject: [PATCH] Restructure includes of FloatingPointComparator

---
 include/deal.II/matrix_free/constraint_info.h | 104 ++++++++++++++++++
 include/deal.II/matrix_free/dof_info.h        |  37 +------
 .../deal.II/matrix_free/dof_info.templates.h  |  68 +-----------
 .../matrix_free/mapping_info.templates.h      |   4 +-
 .../matrix_free/matrix_free.templates.h       |   1 +
 source/matrix_free/dof_info.cc                |  10 --
 6 files changed, 111 insertions(+), 113 deletions(-)

diff --git a/include/deal.II/matrix_free/constraint_info.h b/include/deal.II/matrix_free/constraint_info.h
index 91c735c36f..b11e18aeb4 100644
--- a/include/deal.II/matrix_free/constraint_info.h
+++ b/include/deal.II/matrix_free/constraint_info.h
@@ -20,6 +20,8 @@
 
 #include <deal.II/base/config.h>
 
+#include <deal.II/base/floating_point_comparator.h>
+
 #include <deal.II/matrix_free/dof_info.h>
 #include <deal.II/matrix_free/evaluation_template_factory.h>
 #include <deal.II/matrix_free/hanging_nodes_internal.h>
@@ -33,6 +35,41 @@ namespace internal
 {
   namespace MatrixFreeFunctions
   {
+    /**
+     * A struct that takes entries describing a constraint and puts them into
+     * a sorted list where duplicates are filtered out
+     */
+    template <typename Number>
+    struct ConstraintValues
+    {
+      ConstraintValues();
+
+      /**
+       * This function inserts some constrained entries to the collection of
+       * all values. It stores the (reordered) numbering of the dofs
+       * (according to the ordering that matches with the function) in
+       * new_indices, and returns the storage position the double array for
+       * access later on.
+       */
+      template <typename number2>
+      unsigned short
+      insert_entries(
+        const std::vector<std::pair<types::global_dof_index, number2>>
+          &entries);
+
+      std::vector<std::pair<types::global_dof_index, double>>
+                                           constraint_entries;
+      std::vector<types::global_dof_index> constraint_indices;
+
+      std::pair<std::vector<Number>, types::global_dof_index> next_constraint;
+      std::map<std::vector<Number>,
+               types::global_dof_index,
+               FloatingPointComparator<Number>>
+        constraints;
+    };
+
+
+
     /**
      * A helper class to apply constraints in matrix-free loops in
      * user code. It combines constraint related functionalties from
@@ -131,6 +168,73 @@ namespace internal
 
 
 
+    // ------------------------- inline functions --------------------------
+
+    template <typename Number>
+    ConstraintValues<Number>::ConstraintValues()
+      : constraints(FloatingPointComparator<Number>(
+          1. * std::numeric_limits<double>::epsilon() * 1024.))
+    {}
+
+
+
+    template <typename Number>
+    template <typename number2>
+    unsigned short
+    ConstraintValues<Number>::insert_entries(
+      const std::vector<std::pair<types::global_dof_index, number2>> &entries)
+    {
+      next_constraint.first.resize(entries.size());
+      if (entries.size() > 0)
+        {
+          constraint_indices.resize(entries.size());
+          // Use assign so that values for nonmatching Number / number2 are
+          // converted:
+          constraint_entries.assign(entries.begin(), entries.end());
+          std::sort(constraint_entries.begin(),
+                    constraint_entries.end(),
+                    [](const std::pair<types::global_dof_index, double> &p1,
+                       const std::pair<types::global_dof_index, double> &p2) {
+                      return p1.second < p2.second;
+                    });
+          for (types::global_dof_index j = 0; j < constraint_entries.size();
+               j++)
+            {
+              // copy the indices of the constraint entries after sorting.
+              constraint_indices[j] = constraint_entries[j].first;
+
+              // one_constraint takes the weights of the constraint
+              next_constraint.first[j] = constraint_entries[j].second;
+            }
+        }
+
+      // check whether or not constraint is already in pool. the initial
+      // implementation computed a hash value based on the truncated array (to
+      // given accuracy around 1e-13) in order to easily detect different
+      // arrays and then made a fine-grained check when the hash values were
+      // equal. this was quite lengthy and now we use a std::map with a
+      // user-defined comparator to compare floating point arrays to a
+      // tolerance 1e-13.
+      types::global_dof_index insert_position = numbers::invalid_dof_index;
+      const auto position = constraints.find(next_constraint.first);
+      if (position != constraints.end())
+        insert_position = position->second;
+      else
+        {
+          next_constraint.second = constraints.size();
+          constraints.insert(next_constraint);
+          insert_position = next_constraint.second;
+        }
+
+      // we want to store the result as a short variable, so we have to make
+      // sure that the result does not exceed the limits when casting.
+      Assert(insert_position < (1 << (8 * sizeof(unsigned short))),
+             ExcInternalError());
+      return static_cast<unsigned short>(insert_position);
+    }
+
+
+
     template <int dim, typename Number>
     inline void
     ConstraintInfo<dim, Number>::reinit(
diff --git a/include/deal.II/matrix_free/dof_info.h b/include/deal.II/matrix_free/dof_info.h
index f60288cff0..1e45a1bac8 100644
--- a/include/deal.II/matrix_free/dof_info.h
+++ b/include/deal.II/matrix_free/dof_info.h
@@ -21,7 +21,6 @@
 #include <deal.II/base/config.h>
 
 #include <deal.II/base/exceptions.h>
-#include <deal.II/base/floating_point_comparator.h>
 #include <deal.II/base/vectorization.h>
 
 #include <deal.II/matrix_free/face_info.h>
@@ -46,6 +45,9 @@ namespace internal
     class HangingNodes;
 
     struct TaskInfo;
+
+    template <typename Number>
+    struct ConstraintValues;
   } // namespace MatrixFreeFunctions
 } // namespace internal
 
@@ -80,39 +82,6 @@ namespace internal
      */
     using compressed_constraint_kind = std::uint8_t;
 
-    /**
-     * A struct that takes entries describing a constraint and puts them into
-     * a sorted list where duplicates are filtered out
-     */
-    template <typename Number>
-    struct ConstraintValues
-    {
-      ConstraintValues();
-
-      /**
-       * This function inserts some constrained entries to the collection of
-       * all values. It stores the (reordered) numbering of the dofs
-       * (according to the ordering that matches with the function) in
-       * new_indices, and returns the storage position the double array for
-       * access later on.
-       */
-      template <typename number2>
-      unsigned short
-      insert_entries(
-        const std::vector<std::pair<types::global_dof_index, number2>>
-          &entries);
-
-      std::vector<std::pair<types::global_dof_index, double>>
-                                           constraint_entries;
-      std::vector<types::global_dof_index> constraint_indices;
-
-      std::pair<std::vector<Number>, types::global_dof_index> next_constraint;
-      std::map<std::vector<Number>,
-               types::global_dof_index,
-               FloatingPointComparator<VectorizedArray<Number>>>
-        constraints;
-    };
-
     /**
      * The class that stores the indices of the degrees of freedom for all the
      * cells. Essentially, this is a smart number cache in the style of a
diff --git a/include/deal.II/matrix_free/dof_info.templates.h b/include/deal.II/matrix_free/dof_info.templates.h
index 7acb48b8b9..d3f84e0cdb 100644
--- a/include/deal.II/matrix_free/dof_info.templates.h
+++ b/include/deal.II/matrix_free/dof_info.templates.h
@@ -27,6 +27,7 @@
 #include <deal.II/lac/dynamic_sparsity_pattern.h>
 #include <deal.II/lac/sparsity_pattern.h>
 
+#include <deal.II/matrix_free/constraint_info.h>
 #include <deal.II/matrix_free/dof_info.h>
 #include <deal.II/matrix_free/hanging_nodes_internal.h>
 #include <deal.II/matrix_free/task_info.h>
@@ -39,73 +40,6 @@ namespace internal
 {
   namespace MatrixFreeFunctions
   {
-    template <typename Number>
-    ConstraintValues<Number>::ConstraintValues()
-      : constraints(FloatingPointComparator<VectorizedArray<Number>>(
-          1. * std::numeric_limits<double>::epsilon() * 1024.))
-    {}
-
-
-
-    template <typename Number>
-    template <typename number2>
-    unsigned short
-    ConstraintValues<Number>::insert_entries(
-      const std::vector<std::pair<types::global_dof_index, number2>> &entries)
-    {
-      next_constraint.first.resize(entries.size());
-      if (entries.size() > 0)
-        {
-          constraint_indices.resize(entries.size());
-          // Use assign so that values for nonmatching Number / number2 are
-          // converted:
-          constraint_entries.assign(entries.begin(), entries.end());
-          std::sort(constraint_entries.begin(),
-                    constraint_entries.end(),
-                    [](const std::pair<types::global_dof_index, double> &p1,
-                       const std::pair<types::global_dof_index, double> &p2) {
-                      return p1.second < p2.second;
-                    });
-          for (types::global_dof_index j = 0; j < constraint_entries.size();
-               j++)
-            {
-              // copy the indices of the constraint entries after sorting.
-              constraint_indices[j] = constraint_entries[j].first;
-
-              // one_constraint takes the weights of the constraint
-              next_constraint.first[j] = constraint_entries[j].second;
-            }
-        }
-
-      // check whether or not constraint is already in pool. the initial
-      // implementation computed a hash value based on the truncated array (to
-      // given accuracy around 1e-13) in order to easily detect different
-      // arrays and then made a fine-grained check when the hash values were
-      // equal. this was quite lengthy and now we use a std::map with a
-      // user-defined comparator to compare floating point arrays to a
-      // tolerance 1e-13.
-      types::global_dof_index insert_position = numbers::invalid_dof_index;
-      const auto position = constraints.find(next_constraint.first);
-      if (position != constraints.end())
-        insert_position = position->second;
-      else
-        {
-          next_constraint.second = constraints.size();
-          constraints.insert(next_constraint);
-          insert_position = next_constraint.second;
-        }
-
-      // we want to store the result as a short variable, so we have to make
-      // sure that the result does not exceed the limits when casting.
-      Assert(insert_position < (1 << (8 * sizeof(unsigned short))),
-             ExcInternalError());
-      return static_cast<unsigned short>(insert_position);
-    }
-
-
-
-    // ----------------- actual DoFInfo functions -----------------------------
-
     template <typename number>
     void
     DoFInfo::read_dof_indices(
diff --git a/include/deal.II/matrix_free/mapping_info.templates.h b/include/deal.II/matrix_free/mapping_info.templates.h
index bb584128e9..62ffd3e98f 100644
--- a/include/deal.II/matrix_free/mapping_info.templates.h
+++ b/include/deal.II/matrix_free/mapping_info.templates.h
@@ -1081,8 +1081,8 @@ namespace internal
         // evaluation
         std::map<std::array<Tensor<2, dim>, dim + 1>,
                  unsigned int,
-                 FloatingPointComparator<VectorizedArray<double>>>
-          compressed_jacobians(FloatingPointComparator<VectorizedArray<double>>(
+                 FloatingPointComparator<double>>
+          compressed_jacobians(FloatingPointComparator<double>(
             1e4 * jacobian_size * std::numeric_limits<double>::epsilon() *
             1024.));
 
diff --git a/include/deal.II/matrix_free/matrix_free.templates.h b/include/deal.II/matrix_free/matrix_free.templates.h
index 1b43fc4d03..f5c1baf106 100644
--- a/include/deal.II/matrix_free/matrix_free.templates.h
+++ b/include/deal.II/matrix_free/matrix_free.templates.h
@@ -40,6 +40,7 @@
 
 #include <deal.II/lac/dynamic_sparsity_pattern.h>
 
+#include <deal.II/matrix_free/constraint_info.h>
 #include <deal.II/matrix_free/face_info.h>
 #include <deal.II/matrix_free/face_setup_internal.h>
 #include <deal.II/matrix_free/hanging_nodes_internal.h>
diff --git a/source/matrix_free/dof_info.cc b/source/matrix_free/dof_info.cc
index 6fe7abaebd..7ba077ef9e 100644
--- a/source/matrix_free/dof_info.cc
+++ b/source/matrix_free/dof_info.cc
@@ -1543,16 +1543,6 @@ namespace internal
 {
   namespace MatrixFreeFunctions
   {
-    template struct ConstraintValues<double>;
-
-    template unsigned short
-    ConstraintValues<double>::insert_entries(
-      const std::vector<std::pair<types::global_dof_index, float>> &);
-    template unsigned short
-    ConstraintValues<double>::insert_entries(
-      const std::vector<std::pair<types::global_dof_index, double>> &);
-
-
     template void
     DoFInfo::read_dof_indices<double>(
       const std::vector<types::global_dof_index> &,
-- 
2.39.5