Add more helper functions to number.h

diff --git a/include/deal.II/base/numbers.h b/include/deal.II/base/numbers.h
index 542853c45a4b5fe8d71e0dd74a1628db6953d672..a87894e737118d3dba77f58675d31c46e7b2980a 100644 (file)
--- a/include/deal.II/base/numbers.h
+++ b/include/deal.II/base/numbers.h
@@ -192,6 +192,29 @@ namespace numbers
    */
   bool is_finite (const std::complex<long double> &x);
 
+  /**
+   * Return whether two numbers are equal to one another. For intricate data
+   * types (e.g. some automatically differentiable numbers), this function
+   * returns only whether the scalar values stored by the input values are
+   * equal.
+   *
+   * @note This function expects that @p value_2 is castable to the type
+   * of @p value_1.
+   */
+  template <typename Number1,typename Number2>
+  bool
+  values_are_equal (const Number1 &value_1, const Number2 &value_2);
+
+  /**
+   * Return whether or not a value is equal to zero. For intricate data
+   * types (e.g. some automatically differentiable numbers), this function
+   * returns only whether the scalar value stored by the input value is
+   * equal to zero.
+   */
+  template <typename Number>
+  bool
+  value_is_zero (const Number &value);
+
   /**
    * A structure that, together with its partial specializations
    * NumberTraits<std::complex<number> >, provides traits and member functions
@@ -542,6 +565,85 @@ namespace internal
   };
 }
 
+namespace numbers
+{
+
+#ifdef DEAL_II_ADOLC_WITH_ADVANCED_BRANCHING
+
+  /**
+   * Return whether two numbers are equal to one another. For intricate data
+   * types (e.g. some automatically differentiable numbers), this function
+   * returns only whether the scalar values stored by the input values are
+   * equal.
+   *
+   * @note When Adol-C is compiled with the "advanced branching" feature, then
+   * this specialization is only intended for use in assertions and
+   * other code paths that do not affect the end result of a computation.
+   */
+  // Defined in differentiation/ad/adolc_number_types.cc
+  bool
+  values_are_equal (const adouble &value_1,
+                    const adouble &value_2);
+
+
+  /**
+   * Return whether two numbers are equal to one another. For intricate data
+   * types (e.g. some automatically differentiable numbers), this function
+   * returns only whether the scalar values stored by the input values are
+   * equal.
+   *
+   * @note When Adol-C is compiled with the "advanced branching" feature, then
+   * this specialization is only intended for use in assertions and
+   * other code paths that do not affect the end result of a computation.
+   */
+  template <typename Number>
+  bool
+  values_are_equal (const adouble &value_1,
+                    const Number  &value_2)
+  {
+    // Use the specialized definition for two Adol-C taped types
+    return values_are_equal(value_1, internal::NumberType<adouble>::value(value_2));
+  }
+
+
+  /**
+   * Return whether two numbers are equal to one another. For intricate data
+   * types (e.g. some automatically differentiable numbers), this function
+   * returns only whether the scalar values stored by the input values are
+   * equal.
+   *
+   * @note When Adol-C is compiled with the "advanced branching" feature, then
+   * this specialization is only intended for use in assertions and
+   * other code paths that do not affect the end result of a computation.
+   */
+  template <typename Number>
+  bool
+  values_are_equal (const Number  &value_1,
+                    const adouble &value_2)
+  {
+    // Use the above definition
+    return values_are_equal(value_2, value_1);
+  }
+
+#endif
+
+
+  template <typename Number1,typename Number2>
+  bool
+  values_are_equal (const Number1 &value_1, const Number2 &value_2)
+  {
+    return (value_1 == internal::NumberType<Number1>::value(value_2));
+  }
+
+
+  template <typename Number>
+  bool
+  value_is_zero (const Number &value)
+  {
+    return values_are_equal(value, 0.0);
+  }
+}
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif

diff --git a/source/differentiation/ad/adolc_number_types.cc b/source/differentiation/ad/adolc_number_types.cc
index 073e421082a72c2259d0fa3c87e669090040d445..431ad64219091aac719a66f608d8f63858ea0c9e 100644 (file)
--- a/source/differentiation/ad/adolc_number_types.cc
+++ b/source/differentiation/ad/adolc_number_types.cc
@@ -19,9 +19,55 @@
 #ifdef DEAL_II_WITH_ADOLC
 
 #include <deal.II/differentiation/ad/adolc_number_types.h>
+#include <deal.II/differentiation/ad/ad_number_traits.h>
+
+#include <utility>
 
 DEAL_II_NAMESPACE_OPEN
 
+
+#ifdef DEAL_II_ADOLC_WITH_ADVANCED_BRANCHING
+
+namespace numbers
+{
+
+  namespace internal
+  {
+    namespace
+    {
+      // Apply some comparator and extract the boolean result of the operation,
+      // instead of the "adub" return type tpyically returned by Adol-C for
+      // such a comparison. This is implemented as a general function so that
+      // the list of implemented comparative operations can be easily extended.
+      bool
+      adouble_boolean_comparator (const adouble &value_1,
+                                  const adouble &value_2,
+                                  const std::function<adouble(const adouble &, const adouble &)> &comparator)
+      {
+        typedef typename Differentiation::AD::NumberTraits<double,Differentiation::AD::NumberTypes::adolc_taped>::ad_type ad_type;
+        static_assert(std::is_same<adouble,ad_type>::value, "The type of the AD number is not that which was expected.");
+        const ad_type result = comparator(value_1,value_2);
+        return !(Differentiation::AD::ADNumberTraits<ad_type>::get_scalar_value(result) == 0.0);
+      }
+    }
+  }
+
+  bool
+  values_are_equal (const adouble &value_1,
+                    const adouble &value_2)
+  {
+    return internal::adouble_boolean_comparator(
+             value_1, value_2,
+             [](const adouble &a,const adouble &b) -> adouble
+    {
+      return dealii::internal::NumberType<adouble>::value(a == b);
+    });
+  }
+}
+
+#endif
+
+
 /*---------------------- Explicit Instantiations ----------------------*/
 
 #include "adolc_number_types.inst"