Skip to content

Commit 9f2215a

Browse files
authored
[libc][math] Fix signed zeros for erff. (#97742)
The inexact exception flag was raised for the exact cases of signed zeros. This was reported by Paul Zimmermann using the CORE-MATH test suites.
1 parent d4216b5 commit 9f2215a

File tree

1 file changed

+15
-9
lines changed

1 file changed

+15
-9
lines changed

libc/src/math/generic/erff.cpp

Lines changed: 15 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -127,15 +127,6 @@ LLVM_LIBC_FUNCTION(float, erff, (float x)) {
127127
uint32_t x_u = xbits.uintval();
128128
uint32_t x_abs = x_u & 0x7fff'ffffU;
129129

130-
// Exceptional values
131-
if (LIBC_UNLIKELY(x_abs == 0x3f65'9229U)) // |x| = 0x1.cb2452p-1f
132-
return x < 0.0f ? fputil::round_result_slightly_down(-0x1.972ea8p-1f)
133-
: fputil::round_result_slightly_up(0x1.972ea8p-1f);
134-
if (LIBC_UNLIKELY(x_abs == 0x4004'1e6aU)) // |x| = 0x1.083cd4p+1f
135-
return x < 0.0f ? fputil::round_result_slightly_down(-0x1.fe3462p-1f)
136-
: fputil::round_result_slightly_up(0x1.fe3462p-1f);
137-
138-
// if (LIBC_UNLIKELY(x_abs > 0x407a'd444U)) {
139130
if (LIBC_UNLIKELY(x_abs >= 0x4080'0000U)) {
140131
const float ONE[2] = {1.0f, -1.0f};
141132
const float SMALL[2] = {-0x1.0p-25f, 0x1.0p-25f};
@@ -149,6 +140,21 @@ LLVM_LIBC_FUNCTION(float, erff, (float x)) {
149140
return ONE[sign] + SMALL[sign];
150141
}
151142

143+
// Exceptional mask = common 0 bits of 2 exceptional values.
144+
constexpr uint32_t EXCEPT_MASK = 0x809a'6184U;
145+
146+
if (LIBC_UNLIKELY((x_abs & EXCEPT_MASK) == 0)) {
147+
// Exceptional values
148+
if (LIBC_UNLIKELY(x_abs == 0x3f65'9229U)) // |x| = 0x1.cb2452p-1f
149+
return x < 0.0f ? fputil::round_result_slightly_down(-0x1.972ea8p-1f)
150+
: fputil::round_result_slightly_up(0x1.972ea8p-1f);
151+
if (LIBC_UNLIKELY(x_abs == 0x4004'1e6aU)) // |x| = 0x1.083cd4p+1f
152+
return x < 0.0f ? fputil::round_result_slightly_down(-0x1.fe3462p-1f)
153+
: fputil::round_result_slightly_up(0x1.fe3462p-1f);
154+
if (x_abs == 0U)
155+
return x;
156+
}
157+
152158
// Polynomial approximation:
153159
// erf(x) ~ x * (c0 + c1 * x^2 + c2 * x^4 + ... + c7 * x^14)
154160
double xd = static_cast<double>(x);

0 commit comments

Comments
 (0)