]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Rename FPArrayComparator -> FloatingPointComparator 14296/head
authorPeter Munch <peterrmuench@gmail.com>
Wed, 21 Sep 2022 07:26:29 +0000 (09:26 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Wed, 21 Sep 2022 21:08:10 +0000 (23:08 +0200)
include/deal.II/base/floating_point_comparator.h [new file with mode: 0644]
include/deal.II/matrix_free/dof_info.h
include/deal.II/matrix_free/dof_info.templates.h
include/deal.II/matrix_free/mapping_info.h
include/deal.II/matrix_free/mapping_info.templates.h
source/matrix_free/mapping_info.cc

diff --git a/include/deal.II/base/floating_point_comparator.h b/include/deal.II/base/floating_point_comparator.h
new file mode 100644 (file)
index 0000000..f6d6e82
--- /dev/null
@@ -0,0 +1,203 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 1998 - 2022 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.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_base_floating_point_copmerator_h
+#define dealii_base_floating_point_copmerator_h
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/vectorization.h>
+
+
+DEAL_II_NAMESPACE_OPEN
+
+/**
+ * A class that is used to compare floating point arrays (e.g. std::vector,
+ * Tensor<1,dim>, etc.). The most common use case of this comparator
+ * is for detecting arrays with the same content for the sake of
+ * compression. The idea of this class is to consider two arrays as
+ * equal if they are the same within a given tolerance. We use this
+ * comparator class within a std::map<> of the given arrays. Note that this
+ * comparison operator does not satisfy all the mathematical properties one
+ * usually wants to have (consider e.g. the numbers a=0, b=0.1, c=0.2 with
+ * tolerance 0.15; the operator gives a<c, but neither a<b? nor b<c? is
+ * satisfied). This is not a problem in the use cases for this class, but be
+ * careful when using it in other contexts.
+ */
+template <typename VectorizedArrayType>
+struct FloatingPointComparator
+{
+  using Number = typename dealii::internal::VectorizedArrayTrait<
+    VectorizedArrayType>::value_type;
+  static constexpr std::size_t width =
+    dealii::internal::VectorizedArrayTrait<VectorizedArrayType>::width;
+
+  FloatingPointComparator(const Number scaling);
+
+  /**
+   * Compare two vectors of numbers (not necessarily of the same length),
+   * where vectors of different lengths are first sorted by their length and
+   * then by the entries.
+   */
+  bool
+  operator()(const std::vector<Number> &v1,
+             const std::vector<Number> &v2) const;
+
+  /**
+   * Compare two vectorized arrays (stored as tensors to avoid alignment
+   * issues).
+   */
+  bool
+  operator()(const Tensor<1, width, Number> &t1,
+             const Tensor<1, width, Number> &t2) const;
+
+  /**
+   * Compare two rank-1 tensors of vectorized arrays (stored as tensors to
+   * avoid alignment issues).
+   */
+  template <int dim>
+  bool
+  operator()(const Tensor<1, dim, Tensor<1, width, Number>> &t1,
+             const Tensor<1, dim, Tensor<1, width, Number>> &t2) const;
+
+  /**
+   * Compare two rank-2 tensors of vectorized arrays (stored as tensors to
+   * avoid alignment issues).
+   */
+  template <int dim>
+  bool
+  operator()(const Tensor<2, dim, Tensor<1, width, Number>> &t1,
+             const Tensor<2, dim, Tensor<1, width, Number>> &t2) const;
+
+  /**
+   * Compare two arrays of tensors.
+   */
+  template <int dim>
+  bool
+  operator()(const std::array<Tensor<2, dim, Number>, dim + 1> &t1,
+             const std::array<Tensor<2, dim, Number>, dim + 1> &t2) const;
+
+  Number tolerance;
+};
+
+
+/* ------------------------------------------------------------------ */
+
+
+template <typename VectorizedArrayType>
+FloatingPointComparator<VectorizedArrayType>::FloatingPointComparator(
+  const Number scaling)
+  : tolerance(scaling * std::numeric_limits<double>::epsilon() * 1024.)
+{}
+
+
+
+template <typename VectorizedArrayType>
+bool
+FloatingPointComparator<VectorizedArrayType>::operator()(
+  const std::vector<Number> &v1,
+  const std::vector<Number> &v2) const
+{
+  const unsigned int s1 = v1.size(), s2 = v2.size();
+  if (s1 < s2)
+    return true;
+  else if (s1 > s2)
+    return false;
+  else
+    for (unsigned int i = 0; i < s1; ++i)
+      if (v1[i] < v2[i] - tolerance)
+        return true;
+      else if (v1[i] > v2[i] + tolerance)
+        return false;
+  return false;
+}
+
+
+
+template <typename VectorizedArrayType>
+bool
+FloatingPointComparator<VectorizedArrayType>::operator()(
+  const Tensor<1, width, Number> &t1,
+  const Tensor<1, width, Number> &t2) const
+{
+  for (unsigned int k = 0; k < width; ++k)
+    if (t1[k] < t2[k] - tolerance)
+      return true;
+    else if (t1[k] > t2[k] + tolerance)
+      return false;
+  return false;
+}
+
+
+
+template <typename VectorizedArrayType>
+template <int dim>
+bool
+FloatingPointComparator<VectorizedArrayType>::operator()(
+  const Tensor<1, dim, Tensor<1, width, Number>> &t1,
+  const Tensor<1, dim, Tensor<1, width, Number>> &t2) const
+{
+  for (unsigned int d = 0; d < dim; ++d)
+    for (unsigned int k = 0; k < width; ++k)
+      if (t1[d][k] < t2[d][k] - tolerance)
+        return true;
+      else if (t1[d][k] > t2[d][k] + tolerance)
+        return false;
+  return false;
+}
+
+
+
+template <typename VectorizedArrayType>
+template <int dim>
+bool
+FloatingPointComparator<VectorizedArrayType>::operator()(
+  const Tensor<2, dim, Tensor<1, width, Number>> &t1,
+  const Tensor<2, dim, Tensor<1, width, Number>> &t2) const
+{
+  for (unsigned int d = 0; d < dim; ++d)
+    for (unsigned int e = 0; e < dim; ++e)
+      for (unsigned int k = 0; k < width; ++k)
+        if (t1[d][e][k] < t2[d][e][k] - tolerance)
+          return true;
+        else if (t1[d][e][k] > t2[d][e][k] + tolerance)
+          return false;
+  return false;
+}
+
+
+
+template <typename VectorizedArrayType>
+template <int dim>
+bool
+FloatingPointComparator<VectorizedArrayType>::operator()(
+  const std::array<Tensor<2, dim, Number>, dim + 1> &t1,
+  const std::array<Tensor<2, dim, Number>, dim + 1> &t2) const
+{
+  for (unsigned int i = 0; i < t1.size(); ++i)
+    for (unsigned int d = 0; d < dim; ++d)
+      for (unsigned int e = 0; e < dim; ++e)
+        if (t1[i][d][e] < t2[i][d][e] - tolerance)
+          return true;
+        else if (t1[i][d][e] > t2[i][d][e] + tolerance)
+          return false;
+  return false;
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index 1ea1adec85991e955ce5e2f17d509545b7a2d675..f60288cff093af614f179ebd0741496738e1a651 100644 (file)
@@ -21,6 +21,7 @@
 #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>
@@ -44,9 +45,6 @@ namespace internal
     template <int dim>
     class HangingNodes;
 
-    template <typename>
-    struct FPArrayComparator;
-
     struct TaskInfo;
   } // namespace MatrixFreeFunctions
 } // namespace internal
@@ -111,7 +109,7 @@ namespace internal
       std::pair<std::vector<Number>, types::global_dof_index> next_constraint;
       std::map<std::vector<Number>,
                types::global_dof_index,
-               FPArrayComparator<VectorizedArray<Number>>>
+               FloatingPointComparator<VectorizedArray<Number>>>
         constraints;
     };
 
index 8c3693d66ae9406f93056ca9101ad8b51f658f36..47e964a92ed88acefc10226291ba15a919cfb599 100644 (file)
@@ -19,6 +19,7 @@
 
 #include <deal.II/base/config.h>
 
+#include <deal.II/base/floating_point_comparator.h>
 #include <deal.II/base/memory_consumption.h>
 #include <deal.II/base/parallel.h>
 
@@ -28,7 +29,6 @@
 
 #include <deal.II/matrix_free/dof_info.h>
 #include <deal.II/matrix_free/hanging_nodes_internal.h>
-#include <deal.II/matrix_free/mapping_info.h>
 #include <deal.II/matrix_free/task_info.h>
 
 #include <limits>
@@ -41,7 +41,7 @@ namespace internal
   {
     template <typename Number>
     ConstraintValues<Number>::ConstraintValues()
-      : constraints(FPArrayComparator<VectorizedArray<Number>>(1.))
+      : constraints(FloatingPointComparator<VectorizedArray<Number>>(1.))
     {}
 
 
index 8efa107154452a5002530e73f4df7f15e6f36140..6113903f85e782406784b1dd1afeb3a0cda689b0 100644 (file)
@@ -298,74 +298,6 @@ namespace internal
     };
 
 
-
-    /**
-     * A class that is used to compare floating point arrays (e.g. std::vectors,
-     * Tensor<1,dim>, etc.). The idea of this class is to consider two arrays as
-     * equal if they are the same within a given tolerance. We use this
-     * comparator class within a std::map<> of the given arrays. Note that this
-     * comparison operator does not satisfy all the mathematical properties one
-     * usually wants to have (consider e.g. the numbers a=0, b=0.1, c=0.2 with
-     * tolerance 0.15; the operator gives a<c, but neither a<b? nor b<c? is
-     * satisfied). This is not a problem in the use cases for this class, but be
-     * careful when using it in other contexts.
-     */
-    template <typename VectorizedArrayType>
-    struct FPArrayComparator
-    {
-      using Number = typename dealii::internal::VectorizedArrayTrait<
-        VectorizedArrayType>::value_type;
-      static constexpr std::size_t width =
-        dealii::internal::VectorizedArrayTrait<VectorizedArrayType>::width;
-
-      FPArrayComparator(const Number scaling);
-
-      /**
-       * Compare two vectors of numbers (not necessarily of the same length)
-       */
-      bool
-      operator()(const std::vector<Number> &v1,
-                 const std::vector<Number> &v2) const;
-
-      /**
-       * Compare two vectorized arrays (stored as tensors to avoid alignment
-       * issues).
-       */
-      bool
-      operator()(const Tensor<1, width, Number> &t1,
-                 const Tensor<1, width, Number> &t2) const;
-
-      /**
-       * Compare two rank-1 tensors of vectorized arrays (stored as tensors to
-       * avoid alignment issues).
-       */
-      template <int dim>
-      bool
-      operator()(const Tensor<1, dim, Tensor<1, width, Number>> &t1,
-                 const Tensor<1, dim, Tensor<1, width, Number>> &t2) const;
-
-      /**
-       * Compare two rank-2 tensors of vectorized arrays (stored as tensors to
-       * avoid alignment issues).
-       */
-      template <int dim>
-      bool
-      operator()(const Tensor<2, dim, Tensor<1, width, Number>> &t1,
-                 const Tensor<2, dim, Tensor<1, width, Number>> &t2) const;
-
-      /**
-       * Compare two arrays of tensors.
-       */
-      template <int dim>
-      bool
-      operator()(const std::array<Tensor<2, dim, Number>, dim + 1> &t1,
-                 const std::array<Tensor<2, dim, Number>, dim + 1> &t2) const;
-
-      Number tolerance;
-    };
-
-
-
     /* ------------------- inline functions ----------------------------- */
 
     template <int dim, typename Number, typename VectorizedArrayType>
index 57e777008f62ff2415d2e7a961096eb4ea58d2c4..75faac421c014603442c9c69034dbf8ab819ea23 100644 (file)
@@ -18,6 +18,7 @@
 
 #include <deal.II/base/config.h>
 
+#include <deal.II/base/floating_point_comparator.h>
 #include <deal.II/base/memory_consumption.h>
 #include <deal.II/base/multithread_info.h>
 #include <deal.II/base/thread_management.h>
@@ -347,12 +348,12 @@ namespace internal
       struct CompressedCellData
       {
         CompressedCellData(const double expected_size)
-          : data(FPArrayComparator<VectorizedArrayType>(expected_size))
+          : data(FloatingPointComparator<VectorizedArrayType>(expected_size))
         {}
 
         std::map<Tensor<2, dim, Tensor<1, VectorizedArrayType::size(), Number>>,
                  unsigned int,
-                 FPArrayComparator<VectorizedArrayType>>
+                 FloatingPointComparator<VectorizedArrayType>>
           data;
       };
 
@@ -1077,11 +1078,11 @@ namespace internal
         // first quadrature point on the cell - we use a relatively coarse
         // tolerance to account for some inaccuracies in the manifold
         // evaluation
-        const FPArrayComparator<VectorizedArray<double>> comparator(
+        const FloatingPointComparator<VectorizedArray<double>> comparator(
           1e4 * jacobian_size);
         std::map<std::array<Tensor<2, dim>, dim + 1>,
                  unsigned int,
-                 FPArrayComparator<VectorizedArray<double>>>
+                 FloatingPointComparator<VectorizedArray<double>>>
           compressed_jacobians(comparator);
 
         unsigned int n_data_buckets = 0;
@@ -1499,12 +1500,13 @@ namespace internal
       template <int dim, typename Number, typename VectorizedArrayType>
       struct CompressedFaceData
       {
-        // Constructor. As a scaling factor for the FPArrayComparator, we
+        // Constructor. As a scaling factor for the FloatingPointComparator, we
         // select the inverse of the Jacobian (not the Jacobian as in the
         // CompressedCellData) and add another factor of 512 to account for
         // some roundoff effects.
         CompressedFaceData(const Number jacobian_size)
-          : data(FPArrayComparator<VectorizedArrayType>(512. / jacobian_size))
+          : data(FloatingPointComparator<VectorizedArrayType>(512. /
+                                                              jacobian_size))
           , jacobian_size(jacobian_size)
         {}
 
@@ -1518,7 +1520,7 @@ namespace internal
                         2 * dim * dim + dim + 1,
                         Tensor<1, VectorizedArrayType::size(), Number>>,
                  unsigned int,
-                 FPArrayComparator<VectorizedArrayType>>
+                 FloatingPointComparator<VectorizedArrayType>>
           data;
 
         // Store the scaling factor
@@ -3347,110 +3349,6 @@ namespace internal
         }
     }
 
-
-
-    /* ------------------------------------------------------------------ */
-
-    template <typename VectorizedArrayType>
-    FPArrayComparator<VectorizedArrayType>::FPArrayComparator(
-      const Number scaling)
-      : tolerance(scaling * std::numeric_limits<double>::epsilon() * 1024.)
-    {}
-
-
-
-    template <typename VectorizedArrayType>
-    bool
-    FPArrayComparator<VectorizedArrayType>::operator()(
-      const std::vector<Number> &v1,
-      const std::vector<Number> &v2) const
-    {
-      const unsigned int s1 = v1.size(), s2 = v2.size();
-      if (s1 < s2)
-        return true;
-      else if (s1 > s2)
-        return false;
-      else
-        for (unsigned int i = 0; i < s1; ++i)
-          if (v1[i] < v2[i] - tolerance)
-            return true;
-          else if (v1[i] > v2[i] + tolerance)
-            return false;
-      return false;
-    }
-
-
-
-    template <typename VectorizedArrayType>
-    bool
-    FPArrayComparator<VectorizedArrayType>::operator()(
-      const Tensor<1, width, Number> &t1,
-      const Tensor<1, width, Number> &t2) const
-    {
-      for (unsigned int k = 0; k < width; ++k)
-        if (t1[k] < t2[k] - tolerance)
-          return true;
-        else if (t1[k] > t2[k] + tolerance)
-          return false;
-      return false;
-    }
-
-
-
-    template <typename VectorizedArrayType>
-    template <int dim>
-    bool
-    FPArrayComparator<VectorizedArrayType>::operator()(
-      const Tensor<1, dim, Tensor<1, width, Number>> &t1,
-      const Tensor<1, dim, Tensor<1, width, Number>> &t2) const
-    {
-      for (unsigned int d = 0; d < dim; ++d)
-        for (unsigned int k = 0; k < width; ++k)
-          if (t1[d][k] < t2[d][k] - tolerance)
-            return true;
-          else if (t1[d][k] > t2[d][k] + tolerance)
-            return false;
-      return false;
-    }
-
-
-
-    template <typename VectorizedArrayType>
-    template <int dim>
-    bool
-    FPArrayComparator<VectorizedArrayType>::operator()(
-      const Tensor<2, dim, Tensor<1, width, Number>> &t1,
-      const Tensor<2, dim, Tensor<1, width, Number>> &t2) const
-    {
-      for (unsigned int d = 0; d < dim; ++d)
-        for (unsigned int e = 0; e < dim; ++e)
-          for (unsigned int k = 0; k < width; ++k)
-            if (t1[d][e][k] < t2[d][e][k] - tolerance)
-              return true;
-            else if (t1[d][e][k] > t2[d][e][k] + tolerance)
-              return false;
-      return false;
-    }
-
-
-
-    template <typename VectorizedArrayType>
-    template <int dim>
-    bool
-    FPArrayComparator<VectorizedArrayType>::operator()(
-      const std::array<Tensor<2, dim, Number>, dim + 1> &t1,
-      const std::array<Tensor<2, dim, Number>, dim + 1> &t2) const
-    {
-      for (unsigned int i = 0; i < t1.size(); ++i)
-        for (unsigned int d = 0; d < dim; ++d)
-          for (unsigned int e = 0; e < dim; ++e)
-            if (t1[i][d][e] < t2[i][d][e] - tolerance)
-              return true;
-            else if (t1[i][d][e] > t2[i][d][e] + tolerance)
-              return false;
-      return false;
-    }
-
   } // namespace MatrixFreeFunctions
 } // end of namespace internal
 
index 666353c281d4b7c7c1267a3242a559da3d2d5ac1..36f78a1b4131efb7e6f60c4b0fb0e3b6708b2810 100644 (file)
@@ -29,38 +29,4 @@ DEAL_II_NAMESPACE_OPEN
 #endif
 #include "mapping_info.inst"
 
-#if SPLIT_INSTANTIATIONS_INDEX == 0
-
-template struct internal::MatrixFreeFunctions::FPArrayComparator<double>;
-template struct internal::MatrixFreeFunctions::FPArrayComparator<float>;
-
-template struct internal::MatrixFreeFunctions::FPArrayComparator<
-  VectorizedArray<double, 1>>;
-template struct internal::MatrixFreeFunctions::FPArrayComparator<
-  VectorizedArray<float, 1>>;
-
-#  if (DEAL_II_VECTORIZATION_WIDTH_IN_BITS >= 128 && defined(__SSE2__)) || \
-    (DEAL_II_VECTORIZATION_WIDTH_IN_BITS >= 128 && defined(__ALTIVEC__))
-template struct internal::MatrixFreeFunctions::FPArrayComparator<
-  VectorizedArray<double, 2>>;
-template struct internal::MatrixFreeFunctions::FPArrayComparator<
-  VectorizedArray<float, 4>>;
-#  endif
-
-#  if DEAL_II_VECTORIZATION_WIDTH_IN_BITS >= 256 && defined(__AVX__)
-template struct internal::MatrixFreeFunctions::FPArrayComparator<
-  VectorizedArray<double, 4>>;
-template struct internal::MatrixFreeFunctions::FPArrayComparator<
-  VectorizedArray<float, 8>>;
-#  endif
-
-#  if DEAL_II_VECTORIZATION_WIDTH_IN_BITS >= 512 && defined(__AVX512F__)
-template struct internal::MatrixFreeFunctions::FPArrayComparator<
-  VectorizedArray<double, 8>>;
-template struct internal::MatrixFreeFunctions::FPArrayComparator<
-  VectorizedArray<float, 16>>;
-#  endif
-
-#endif
-
 DEAL_II_NAMESPACE_CLOSE

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.