[libc-commits] [libc] [libc][math] Short-circuit certain powf range errors (PR #223719)

[via libc-commits]

llvmorg-github-actions[bot] wrote:


<!--LLVM PR SUMMARY COMMENT-->

@llvm/pr-subscribers-libc

Author: sriramshastry

<details>
<summary>Changes</summary>

### Summary

Return overflow and underflow results directly when float bounds prove
that the exact `powf` result is outside the representable range. This
avoids the log2 and exp2 approximations for those inputs.

Rounding-mode results, `errno`, and floating-point exceptions are
preserved.

### Correctness

- Compared 16,214,688 records against main across all four rounding
  modes.
- Result bits, `errno`, and floating-point exception flags matched.
- LLVM libc smoke tests passed: 20/20 across default, no-FMA, no-round,
  and combined configurations.

### Performance

LLVM has no checked-in `powf_perf` target. Measurements used a minimal
driver that directly instantiates LLVM's unmodified
`BINARY_INPUT_SINGLE_OUTPUT_PERF` macro from `PerfTest.h`.

Three alternating CPU 0 comparisons produced these median savings:

- Denormal range: 3.53 ns (2.01%).
- Normal range: 53.13 ns (78.64%).
- Close-exponent range: 2.95 ns (15.89%).

CPU frequency scaling was enabled, so results are reported in ns/op and
are not converted to cycles.

---
Full diff: https://github.com/llvm/llvm-project/pull/223719.diff


1 Files Affected:

- (modified) libc/src/__support/math/powf.h (+36) 


``````````diff
diff --git a/libc/src/__support/math/powf.h b/libc/src/__support/math/powf.h
index 944b477ac44612..b95d3d2bd10d7c 100644
--- a/libc/src/__support/math/powf.h
+++ b/libc/src/__support/math/powf.h
@@ -37,6 +37,7 @@
 #include "src/__support/FPUtil/double_double.h"
 #include "src/__support/FPUtil/multiply_add.h"
 #include "src/__support/FPUtil/nearest_integer.h"
+#include "src/__support/FPUtil/rounding_mode.h"
 #include "src/__support/FPUtil/sqrt.h" // Speedup for powf(x, 1/2) = sqrtf(x)
 #include "src/__support/FPUtil/triple_double.h"
 #include "src/__support/common.h"
@@ -693,6 +694,26 @@ LIBC_INLINE float powf(float x, float y) {
       return FloatBits::quiet_nan().get_val();
     }
 
+    if (y_abs > 0x4f17'0000U && y_abs < 0x7f80'0000U && x_u != 0 &&
+        x_u < 0x7f80'0000U && x_u != 0x3f80'0000U &&
+        x_u != 0x4000'0000U && x_u != 0x4120'0000U) {
+      const bool overflow = (x_u > 0x3f80'0000U) != (y_u > 0x8000'0000U);
+      if (overflow) {
+        fputil::raise_except_if_required(FE_OVERFLOW | FE_INEXACT);
+        static volatile float rounding_test = 0x1.0p-25f;
+        if (1.0f - rounding_test != 1.0f)
+          return FloatBits::max_normal().get_val();
+        fputil::set_errno_if_required(ERANGE);
+        return FloatBits::inf().get_val();
+      }
+
+      fputil::raise_except_if_required(FE_UNDERFLOW | FE_INEXACT);
+      if (fputil::fenv_is_round_up())
+        return FloatBits::min_subnormal().get_val();
+      fputil::set_errno_if_required(ERANGE);
+      return 0.0f;
+    }
+
     // Exceptional exponents.
     if (y == 0.0f)
       return 1.0f;
@@ -781,6 +802,21 @@ LIBC_INLINE float powf(float x, float y) {
     }
   }
 
+  if (LIBC_UNLIKELY(y_u >= 0x3f80'0000U && y_u < 0x7f80'0000U && x_u > 0 &&
+                    x_u < 0x3f00'0000U)) {
+    // If x is in [2^(e - 127), 2^(e - 126)), then
+    // log2(x) < -(126 - e).  Use this bound to identify results below
+    // 2^-150 without evaluating log2(x).
+    int log2_bound = 126 - static_cast<int>(x_u >> FloatBits::FRACTION_LEN);
+    if (LIBC_UNLIKELY(static_cast<double>(y) * log2_bound >= 150.0)) {
+      fputil::raise_except_if_required(FE_UNDERFLOW | FE_INEXACT);
+      if (fputil::fenv_is_round_up())
+        return FloatBits::min_subnormal().get_val();
+      fputil::set_errno_if_required(ERANGE);
+      return 0.0f;
+    }
+  }
+
   int ex = -FloatBits::EXP_BIAS;
   uint64_t sign = 0;
 

``````````

</details>


https://github.com/llvm/llvm-project/pull/223719