Re: [PATCH] Implement "P0631R4 Math Constants" for C++20

[Jonathan Wakely <jwakely at redhat dot com>]

On 31/07/19 19:07 +0200, Marc Glisse wrote:
> On Wed, 31 Jul 2019, Jonathan Wakely wrote:
>
> > The values of the constants are taken from Glibc where the equivalent
> > constant exists, or by rounding the actual constant to the same number
> > of digits as the Glibc constants have.
>
> How does it behave with __float128? I think that with -std=gnu++2a, it 
> counts as a floating point type, but the constant first becomes a long 
> double and thus loses precision.

Indeed it does:

#if !defined(__STRICT_ANSI__) && defined(_GLIBCXX_USE_FLOAT128)
 template<>
   struct __is_floating_point_helper<__float128>
   : public true_type { };
#endif

We could add explicit specializations for __float128 and use the Q
suffix. That's not valid when __STRICT_ANSI__ is defined (unless 
-fext-numeric-literals is also used) but as the snippet above shows,
__float128 isn't a floating-point type when __STRICT_ANSI__ is
defined.

So something like the attached patch.  The glibc <math.h> says the
values I used have enough digits for IEEE quad-precision:

/* The above constants are not adequate for computation using `long double's.
  Therefore we provide as an extension constants with similar names as a
  GNU extension.  Provide enough digits for the 128-bit IEEE quad.  */

I don't know if we need more accuracy for IBM double double.


diff --git a/libstdc++-v3/include/std/numbers b/libstdc++-v3/include/std/numbers
index b8e38dd8080..314b130fbd5 100644
--- a/libstdc++-v3/include/std/numbers
+++ b/libstdc++-v3/include/std/numbers
@@ -133,6 +133,72 @@ namespace numbers
   inline constexpr double egamma = egamma_v<double>;
   inline constexpr double phi = phi_v<double>;
 
+#if !defined(__STRICT_ANSI__) && defined(_GLIBCXX_USE_FLOAT128)
+  template<>
+    inline constexpr __float128 e_v<__float128>
+      = 2.718281828459045235360287471352662498Q;
+
+  /// log_2 e
+  template<>
+    inline constexpr __float128 log2e_v<__float128>
+      = 1.442695040888963407359924681001892137Q;
+
+  /// log_10 e
+  template<>
+    inline constexpr __float128 log10e_v<__float128>
+      = 0.434294481903251827651128918916605082Q;
+
+  /// pi
+  template<>
+    inline constexpr __float128 pi_v<__float128>
+      = 3.141592653589793238462643383279502884Q;
+
+  /// 1/pi
+  template<>
+    inline constexpr __float128 inv_pi_v<__float128>
+      = 0.318309886183790671537767526745028724Q;
+
+  /// 1/sqrt(pi)
+  template<>
+    inline constexpr __float128 inv_sqrtpi_v<__float128>
+      = 0.564189583547756286948079451560772586Q;
+
+  /// log_e 2
+  template<>
+    inline constexpr __float128 ln2_v<__float128>
+      = 0.693147180559945309417232121458176568Q;
+
+  /// log_e 10
+  template<>
+    inline constexpr __float128 ln10_v<__float128>
+      = 2.302585092994045684017991454684364208Q;
+
+  /// sqrt(2)
+  template<>
+    inline constexpr __float128 sqrt2_v<__float128>
+      = 1.414213562373095048801688724209698079Q;
+
+  /// sqrt(3)
+  template<>
+    inline constexpr __float128 sqrt3_v<__float128>
+      = 1.732050807568877293527446341505872367Q;
+
+  /// 1/sqrt(3)
+  template<>
+    inline constexpr __float128 inv_sqrt3_v<__float128>
+      = 0.577350269189625764509148780501957456Q;
+
+  /// The Euler-Mascheroni constant
+  template<>
+    inline constexpr __float128 egamma_v<__float128>
+      = 0.577215664901532860606512090082402431Q;
+
+  /// The golden ratio, (1+sqrt(5))/2
+  template<>
+    inline constexpr __float128 phi_v<__float128>
+      = 1.618033988749894848204586834365638118Q;
+#endif // USE_FLOAT128
+
 } // namespace numbers
 /// @}
 _GLIBCXX_END_NAMESPACE_VERSION