Fix powf overflow handling in non-nearest rounding mode

Message ID 75658fad-e3b2-aa0d-039b-2318a9dd5e1d@arm.com
State New
Headers show
Series
  • Fix powf overflow handling in non-nearest rounding mode
Related show

Commit Message

Szabolcs Nagy Dec. 10, 2018, 3:04 p.m.
Follows the fix in https://github.com/ARM-software/optimized-routines

Comments

Corinna Vinschen Dec. 11, 2018, 11:56 a.m. | #1
On Dec 10 15:04, Szabolcs Nagy wrote:
> Follows the fix in https://github.com/ARM-software/optimized-routines


> From 8ff59d2f4e0209f03c30462e9319baf4c6ba79f7 Mon Sep 17 00:00:00 2001

> From: Szabolcs Nagy <szabolcs.nagy@arm.com>

> Date: Mon, 10 Dec 2018 14:40:01 +0000

> Subject: [PATCH] Fix powf overflow handling in non-nearest rounding mode

> 

> The threshold value at which powf overflows depends on the rounding mode

> and the current check did not take this into account. So when the result

> was rounded away from zero it could become infinity without setting

> errno to ERANGE.

> 

> Example: pow(0x1.7ac7cp+5, 23) is 0x1.fffffep+127 + 0.1633ulp

> 

> If the result goes above 0x1.fffffep+127 + 0.5ulp then errno is set,

> which is fine in nearest rounding mode, but

> 

>   powf(0x1.7ac7cp+5, 23) is inf in upward rounding mode

>   powf(-0x1.7ac7cp+5, 23) is -inf in downward rounding mode

> 

> and the previous implementation did not set errno in these cases.

> 

> The fix tries to avoid affecting the common code path or calling a

> function that may introduce a stack frame, so float arithmetics is used

> to check the rounding mode and the threshold is selected accordingly.

> ---

>  newlib/libm/common/sf_pow.c | 10 ++++++++++

>  1 file changed, 10 insertions(+)


Pushed.


Thanks,
Corinna

-- 
Corinna Vinschen
Cygwin Maintainer
Red Hat

Patch

From 8ff59d2f4e0209f03c30462e9319baf4c6ba79f7 Mon Sep 17 00:00:00 2001
From: Szabolcs Nagy <szabolcs.nagy@arm.com>
Date: Mon, 10 Dec 2018 14:40:01 +0000
Subject: [PATCH] Fix powf overflow handling in non-nearest rounding mode

The threshold value at which powf overflows depends on the rounding mode
and the current check did not take this into account. So when the result
was rounded away from zero it could become infinity without setting
errno to ERANGE.

Example: pow(0x1.7ac7cp+5, 23) is 0x1.fffffep+127 + 0.1633ulp

If the result goes above 0x1.fffffep+127 + 0.5ulp then errno is set,
which is fine in nearest rounding mode, but

  powf(0x1.7ac7cp+5, 23) is inf in upward rounding mode
  powf(-0x1.7ac7cp+5, 23) is -inf in downward rounding mode

and the previous implementation did not set errno in these cases.

The fix tries to avoid affecting the common code path or calling a
function that may introduce a stack frame, so float arithmetics is used
to check the rounding mode and the threshold is selected accordingly.
---
 newlib/libm/common/sf_pow.c | 10 ++++++++++
 1 file changed, 10 insertions(+)

diff --git a/newlib/libm/common/sf_pow.c b/newlib/libm/common/sf_pow.c
index f7b22b612..2946c611b 100644
--- a/newlib/libm/common/sf_pow.c
+++ b/newlib/libm/common/sf_pow.c
@@ -220,7 +220,17 @@  powf (float x, float y)
     {
       /* |y*log(x)| >= 126.  */
       if (ylogx > 0x1.fffffffd1d571p+6 * POWF_SCALE)
+	/* |x^y| > 0x1.ffffffp127.  */
 	return __math_oflowf (sign_bias);
+      if (WANT_ROUNDING && WANT_ERRNO
+	  && ylogx > 0x1.fffffffa3aae2p+6 * POWF_SCALE)
+	/* |x^y| > 0x1.fffffep127, check if we round away from 0.  */
+	if ((!sign_bias
+	     && eval_as_float (1.0f + opt_barrier_float (0x1p-25f)) != 1.0f)
+	    || (sign_bias
+		&& eval_as_float (-1.0f - opt_barrier_float (0x1p-25f))
+		     != -1.0f))
+	  return __math_oflowf (sign_bias);
       if (ylogx <= -150.0 * POWF_SCALE)
 	return __math_uflowf (sign_bias);
 #if WANT_ERRNO_UFLOW
-- 
2.14.1