[libc][math] Fix signed zeros for acosf, acoshf, and atanf in FE_DOWNWARD mode.

[lntue]

Summary

Fix signed zeros for acosf, acoshf, and atanf in FE_DOWNWARD mode.

[sivachandra]

A nit marked inline.

libc/src/math/generic/acosf.cpp

79

You should add comments for each branch taken. For example:

if (x_abs == 0x3f80'0000U) { // |x| == 1
    return x_sign ? /* x == -1 */ ... :
                             /* x == 1*/ ...;
}

[lntue]

Address comment.

libc/src/math/generic/acosf.cpp

if (LIBC_UNLIKELY(x_abs <= 0x3f00'0000U)) {
     //                 [|1, D...|], [0, 0.5]);
     double xd = static_cast<double>(x);
     double xsq = xd * xd;
     double x3 = xd * xsq;
     double r = asin_eval(xsq);
     return static_cast<float>(fputil::multiply_add(-x3, r, M_MATH_PI_2 - xd));
   }
 
-  // |x| > 1, return NaNs.
-  if (LIBC_UNLIKELY(x_abs > 0x3f80'0000U)) {
+  // |x| >= 1, return 0, 2pi, or NaNs.
+  if (LIBC_UNLIKELY(x_abs >= 0x3f80'0000U)) {
+    if (x_abs == 0x3f80'0000U)
+      return x_sign ? /* x == -1.0f */ fputil::round_result_slightly_down(
+                          0x1.921fb6p+1f)
+                    : /* x == 1.0f */ 0.0f;
+
     if (x_abs <= 0x7f80'0000U) {
       fputil::set_errno_if_required(EDOM);
       fputil::raise_except_if_required(FE_INVALID);
     }
     return x + FPBits::build_quiet_nan(0);
   }
 
-  // When 0.5 < |x| <= 1, we perform range reduction as follow:
+  // When 0.5 < |x| < 1, we perform range reduction as follow:
   //
   // Assume further that 0.5 < x <= 1, and let:
   //   y = acos(x)
   // We use the double angle formula:
   //   x = cos(y) = 1 - 2 sin^2(y/2)
   // So:
   //   sin(y/2) = sqrt( (1 - x)/2 )
   // And hence:
   //   y = 2 * asin( sqrt( (1 - x)/2 ) )
   // Let u = (1 - x)/2, then
   //   acos(x) = 2 * asin( sqrt(u) )
   // Moreover, since 0.5 < x <= 1,
   //   0 <= u < 1/4, and 0 <= sqrt(u) < 0.5,
   // And hence we can reuse the same polynomial approximation of asin(x) when
   // |x| <= 0.5:
   //   acos(x) ~ 2 * sqrt(u) * P(u).
   //
-  // When -1 <= x <= -0.5, we use the identity:
+  // When -1 < x <= -0.5, we use the identity:
   //   acos(x) = pi - acos(-x)
   // which is reduced to the postive case.
 
   xbits.set_sign(false);
   double xd = static_cast<double>(xbits.get_val());
   double u = fputil::multiply_add(-0.5, xd, 0.5);
   double cv = 2 * fputil::sqrt(u);
 
   double r3 = asin_eval(u);
   double r = fputil::multiply_add(cv * u, r3, cv);
   return static_cast<float>(x_sign ? M_MATH_PI - r : r);
 }
 
 } // namespace __llvm_libc

libc/src/math/generic/acoshf.cpp

 
 namespace __llvm_libc {
 
 LLVM_LIBC_FUNCTION(float, acoshf, (float x)) {
   using FPBits_t = typename fputil::FPBits<float>;
   FPBits_t xbits(x);
   uint32_t x_u = xbits.uintval();
 
-  if (LIBC_UNLIKELY(x < 1.0f)) {
+  if (LIBC_UNLIKELY(x <= 1.0f)) {
+    if (x == 1.0f)
+      return 0.0f;
     // x < 1.
     fputil::set_errno_if_required(EDOM);
     fputil::raise_except_if_required(FE_INVALID);
     return FPBits_t::build_quiet_nan(0);
   }
 
   if (LIBC_UNLIKELY(x_u >= 0x4f8ffb03)) {
     // Check for exceptional values.

libc/src/math/generic/atanf.cpp

 #include "src/__support/FPUtil/rounding_mode.h"
 #include "src/__support/macros/optimization.h" // LIBC_UNLIKELY
 #include "src/math/generic/inv_trigf_utils.h"
 
 namespace __llvm_libc {
 
 LLVM_LIBC_FUNCTION(float, atanf, (float x)) {
   using FPBits = typename fputil::FPBits<float>;
+
+  // x == 0.0
+  if (LIBC_UNLIKELY(x == 0.0f))
+    return x;
+
   FPBits xbits(x);
   bool sign = xbits.get_sign();
   xbits.set_sign(false);
 
   if (LIBC_UNLIKELY(xbits.is_inf_or_nan())) {
     if (xbits.is_inf())
       return static_cast<float>(opt_barrier(sign ? -M_MATH_PI_2 : M_MATH_PI_2));
     else

libc/test/src/math/acosf_test.cpp

 
 namespace mpfr = __llvm_libc::testing::mpfr;
 
 DECLARE_SPECIAL_CONSTANTS(float)
 
 TEST(LlvmLibcAcosfTest, SpecialNumbers) {
   libc_errno = 0;
 
-  EXPECT_FP_EQ(aNaN, __llvm_libc::acosf(aNaN));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::acosf(aNaN));
   EXPECT_MATH_ERRNO(0);
 
-  EXPECT_FP_EQ(aNaN, __llvm_libc::acosf(inf));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::acosf(inf));
   EXPECT_MATH_ERRNO(EDOM);
 
-  EXPECT_FP_EQ(aNaN, __llvm_libc::acosf(neg_inf));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::acosf(neg_inf));
+  EXPECT_MATH_ERRNO(EDOM);
+
+  EXPECT_FP_EQ_ALL_ROUNDING(zero, __llvm_libc::acosf(1.0f));
+  EXPECT_MATH_ERRNO(0);
+
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::acosf(2.0f));
+  EXPECT_MATH_ERRNO(EDOM);
+
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::acosf(-2.0f));
   EXPECT_MATH_ERRNO(EDOM);
 }
 
 TEST(LlvmLibcAcosfTest, InFloatRange) {
   constexpr uint32_t COUNT = 100'000;
   constexpr uint32_t STEP = UINT32_MAX / COUNT;
   for (uint32_t i = 0, v = 0; i <= COUNT; ++i, v += STEP) {
     float x = float(FPBits(v));

libc/test/src/math/acoshf_test.cpp

 
 namespace mpfr = __llvm_libc::testing::mpfr;
 
 DECLARE_SPECIAL_CONSTANTS(float)
 
 TEST(LlvmLibcAcoshfTest, SpecialNumbers) {
   libc_errno = 0;
 
-  EXPECT_FP_EQ(aNaN, __llvm_libc::acoshf(aNaN));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::acoshf(aNaN));
   EXPECT_MATH_ERRNO(0);
 
-  EXPECT_FP_EQ(aNaN, __llvm_libc::acoshf(0.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::acoshf(0.0f));
   EXPECT_MATH_ERRNO(EDOM);
 
-  EXPECT_FP_EQ(0.0f, __llvm_libc::acoshf(1.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(0.0f, __llvm_libc::acoshf(1.0f));
   EXPECT_MATH_ERRNO(0);
 
-  EXPECT_FP_EQ(inf, __llvm_libc::acoshf(inf));
+  EXPECT_FP_EQ_ALL_ROUNDING(inf, __llvm_libc::acoshf(inf));
   EXPECT_MATH_ERRNO(0);
 
-  EXPECT_FP_EQ(aNaN, __llvm_libc::acoshf(neg_inf));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::acoshf(neg_inf));
   EXPECT_MATH_ERRNO(EDOM);
 }
 
 TEST(LlvmLibcAcoshfTest, InFloatRange) {
   constexpr uint32_t COUNT = 100'000;
   constexpr uint32_t STEP = UINT32_MAX / COUNT;
   for (uint32_t i = 0, v = 0; i <= COUNT; ++i, v += STEP) {
     float x = float(FPBits_t(v));

libc/test/src/math/asinf_test.cpp

 
 namespace mpfr = __llvm_libc::testing::mpfr;
 
 DECLARE_SPECIAL_CONSTANTS(float)
 
 TEST(LlvmLibcAsinfTest, SpecialNumbers) {
   libc_errno = 0;
 
-  EXPECT_FP_EQ(aNaN, __llvm_libc::asinf(aNaN));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::asinf(aNaN));
   EXPECT_MATH_ERRNO(0);
 
-  EXPECT_FP_EQ(0.0f, __llvm_libc::asinf(0.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(0.0f, __llvm_libc::asinf(0.0f));
   EXPECT_MATH_ERRNO(0);
 
-  EXPECT_FP_EQ(-0.0f, __llvm_libc::asinf(-0.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(-0.0f, __llvm_libc::asinf(-0.0f));
   EXPECT_MATH_ERRNO(0);
 
-  EXPECT_FP_EQ(aNaN, __llvm_libc::asinf(inf));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::asinf(inf));
   EXPECT_MATH_ERRNO(EDOM);
 
-  EXPECT_FP_EQ(aNaN, __llvm_libc::asinf(neg_inf));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::asinf(neg_inf));
   EXPECT_MATH_ERRNO(EDOM);
 }
 
 TEST(LlvmLibcAsinfTest, InFloatRange) {
   constexpr uint32_t COUNT = 100'000;
   constexpr uint32_t STEP = UINT32_MAX / COUNT;
   for (uint32_t i = 0, v = 0; i <= COUNT; ++i, v += STEP) {
     float x = float(FPBits(v));

libc/test/src/math/asinhf_test.cpp

 
 namespace mpfr = __llvm_libc::testing::mpfr;
 
 DECLARE_SPECIAL_CONSTANTS(float)
 
 TEST(LlvmLibcAsinhfTest, SpecialNumbers) {
   libc_errno = 0;
 
-  EXPECT_FP_EQ(aNaN, __llvm_libc::asinhf(aNaN));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::asinhf(aNaN));
   EXPECT_MATH_ERRNO(0);
 
-  EXPECT_FP_EQ(0.0f, __llvm_libc::asinhf(0.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(0.0f, __llvm_libc::asinhf(0.0f));
   EXPECT_MATH_ERRNO(0);
 
-  EXPECT_FP_EQ(-0.0f, __llvm_libc::asinhf(-0.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(-0.0f, __llvm_libc::asinhf(-0.0f));
   EXPECT_MATH_ERRNO(0);
 
-  EXPECT_FP_EQ(inf, __llvm_libc::asinhf(inf));
+  EXPECT_FP_EQ_ALL_ROUNDING(inf, __llvm_libc::asinhf(inf));
   EXPECT_MATH_ERRNO(0);
 
-  EXPECT_FP_EQ(neg_inf, __llvm_libc::asinhf(neg_inf));
+  EXPECT_FP_EQ_ALL_ROUNDING(neg_inf, __llvm_libc::asinhf(neg_inf));
   EXPECT_MATH_ERRNO(0);
 }
 
 TEST(LlvmLibcAsinhfTest, InFloatRange) {
   constexpr uint32_t COUNT = 1'001;
   constexpr uint32_t STEP = UINT32_MAX / COUNT;
   for (uint32_t i = 0, v = 0; i <= COUNT; ++i, v += STEP) {
     float x = float(FPBits_t(v));

libc/test/src/math/atanf_test.cpp

 
 namespace mpfr = __llvm_libc::testing::mpfr;
 
 DECLARE_SPECIAL_CONSTANTS(float)
 
 TEST(LlvmLibcAtanfTest, SpecialNumbers) {
   libc_errno = 0;
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(aNaN, __llvm_libc::atanf(aNaN));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::atanf(aNaN));
   EXPECT_FP_EXCEPTION(0);
   EXPECT_MATH_ERRNO(0);
 
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(0.0f, __llvm_libc::atanf(0.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(0.0f, __llvm_libc::atanf(0.0f));
   EXPECT_FP_EXCEPTION(0);
   EXPECT_MATH_ERRNO(0);
 
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(-0.0f, __llvm_libc::atanf(-0.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(-0.0f, __llvm_libc::atanf(-0.0f));
   EXPECT_FP_EXCEPTION(0);
   EXPECT_MATH_ERRNO(0);
 }
 
 TEST(LlvmLibcAtanfTest, InFloatRange) {
   constexpr uint32_t COUNT = 100'000;
   const uint32_t STEP = FPBits(inf).uintval() / COUNT;
   for (uint32_t i = 0, v = 0; i <= COUNT; ++i, v += STEP) {

libc/test/src/math/atanhf_test.cpp

 
 namespace mpfr = __llvm_libc::testing::mpfr;
 
 DECLARE_SPECIAL_CONSTANTS(float)
 
 TEST(LlvmLibcAtanhfTest, SpecialNumbers) {
   libc_errno = 0;
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(aNaN, __llvm_libc::atanhf(aNaN));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::atanhf(aNaN));
   EXPECT_FP_EXCEPTION(0);
   EXPECT_MATH_ERRNO(0);
 
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(0.0f, __llvm_libc::atanhf(0.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(0.0f, __llvm_libc::atanhf(0.0f));
   EXPECT_FP_EXCEPTION(0);
   EXPECT_MATH_ERRNO(0);
 
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(-0.0f, __llvm_libc::atanhf(-0.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(-0.0f, __llvm_libc::atanhf(-0.0f));
   EXPECT_FP_EXCEPTION(0);
   EXPECT_MATH_ERRNO(0);
 
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(inf, __llvm_libc::atanhf(1.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(inf, __llvm_libc::atanhf(1.0f));
   EXPECT_FP_EXCEPTION(FE_DIVBYZERO);
   EXPECT_MATH_ERRNO(ERANGE);
 
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(neg_inf, __llvm_libc::atanhf(-1.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(neg_inf, __llvm_libc::atanhf(-1.0f));
   EXPECT_FP_EXCEPTION(FE_DIVBYZERO);
   EXPECT_MATH_ERRNO(ERANGE);
 
   auto bt = FPBits(1.0f);
   bt.bits += 1;
 
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(aNaN, __llvm_libc::atanhf(bt.get_val()));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::atanhf(bt.get_val()));
   EXPECT_FP_EXCEPTION(FE_INVALID);
   EXPECT_MATH_ERRNO(EDOM);
 
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
   bt.set_sign(true);
-  EXPECT_FP_EQ(aNaN, __llvm_libc::atanhf(bt.get_val()));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::atanhf(bt.get_val()));
   EXPECT_FP_EXCEPTION(FE_INVALID);
   EXPECT_MATH_ERRNO(EDOM);
 
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(aNaN, __llvm_libc::atanhf(2.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::atanhf(2.0f));
   EXPECT_FP_EXCEPTION(FE_INVALID);
   EXPECT_MATH_ERRNO(EDOM);
 
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(aNaN, __llvm_libc::atanhf(-2.0f));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::atanhf(-2.0f));
   EXPECT_FP_EXCEPTION(FE_INVALID);
   EXPECT_MATH_ERRNO(EDOM);
 
   __llvm_libc::fputil::clear_except(FE_ALL_EXCEPT);
-  EXPECT_FP_EQ(aNaN, __llvm_libc::atanhf(inf));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::atanhf(inf));
   EXPECT_FP_EXCEPTION(FE_INVALID);
   EXPECT_MATH_ERRNO(EDOM);
 
   bt.set_sign(true);
-  EXPECT_FP_EQ(aNaN, __llvm_libc::atanhf(neg_inf));
+  EXPECT_FP_EQ_ALL_ROUNDING(aNaN, __llvm_libc::atanhf(neg_inf));
   EXPECT_FP_EXCEPTION(FE_INVALID);
   EXPECT_MATH_ERRNO(EDOM);
 }
 
 TEST(LlvmLibcAtanhfTest, InFloatRange) {
   constexpr uint32_t COUNT = 100'000;
   const uint32_t STEP = FPBits(1.0f).uintval() / COUNT;
   for (uint32_t i = 0, v = 0; i <= COUNT; ++i, v += STEP) {