[PATCH] ldbl-128: Fix y0 and y1 for -Inf input

[PATCH] ldbl-128: Fix y0 and y1 for -Inf input

[Gabriel F. T. Gomes]

The code that this patch changes is not executed, because the wrappers
(in math/w_j0_compat.c and math/w_j1_compat.c) call __kernel_standard
and return.

However, in the float128 branch (not on master), a new wrapper that
does not use __kernel_standard is in use, so this code path gets
executed.  This patch is tested in that branch.

---8<---
The Bessel functions of the second type (Yn) are not defined for
negative x and should return NAN in these cases.  However, current
code checks for infinity return zero, regardless of the sign.  This
error is not exposed for float, double, and long double, because the
wrappers for these functions, which use __kernel_standard
functionality, return the correct value.  For float128 wrappers, which
must not rely on __kernel_standard functionality, this error is
exposed.  This patch fixes it.

2017-01-06  Gabriel F. T. Gomes  <gftg@linux.vnet.ibm.com>

	* sysdeps/ieee754/ldbl-128/e_j0l.c (__ieee754_y0l): Return NAN
	when x is -Inf.
	* sysdeps/ieee754/ldbl-128/e_j1l.c (__ieee754_y1l): Likewise.
---
 sysdeps/ieee754/ldbl-128/e_j0l.c | 10 +++++++++-
 sysdeps/ieee754/ldbl-128/e_j1l.c | 10 +++++++++-
 2 files changed, 18 insertions(+), 2 deletions(-)

diff --git a/sysdeps/ieee754/ldbl-128/e_j0l.c b/sysdeps/ieee754/ldbl-128/e_j0l.c
index d711007..3bf3569 100644
--- a/sysdeps/ieee754/ldbl-128/e_j0l.c
+++ b/sysdeps/ieee754/ldbl-128/e_j0l.c
@@ -833,7 +833,15 @@ _Float128
       if (x != x)
 	return x + x;
       else
-	return 0;
+	{
+	  /* x = +/-Inf */
+	  if (!signbit (x))
+	    /* y0 (Inf) = 0.  */
+	    return 0;
+	  else
+	    /* y0 (-Inf) = NaN.  */
+	    return NAN;
+	}
     }
   if (x <= 0)
     {
diff --git a/sysdeps/ieee754/ldbl-128/e_j1l.c b/sysdeps/ieee754/ldbl-128/e_j1l.c
index 9e78230..99b6534 100644
--- a/sysdeps/ieee754/ldbl-128/e_j1l.c
+++ b/sysdeps/ieee754/ldbl-128/e_j1l.c
@@ -851,7 +851,15 @@ __ieee754_y1l (_Float128 x)
       if (x != x)
 	return x + x;
       else
-	return 0;
+	{
+	  /* x = +/-Inf */
+	  if (!signbit (x))
+	    /* y0 (Inf) = 0.  */
+	    return 0;
+	  else
+	    /* y0 (-Inf) = NaN.  */
+	    return NAN;
+	}
     }
   if (x <= 0)
     {
-- 
2.4.11

Re: [PATCH] ldbl-128: Fix y0 and y1 for -Inf input

[Joseph Myers]

On Wed, 8 Feb 2017, Gabriel F. T. Gomes wrote:

> The code that this patch changes is not executed, because the wrappers
> (in math/w_j0_compat.c and math/w_j1_compat.c) call __kernel_standard
> and return.

It's executed with -lieee (meaning this is user-visible and should have a 
bug filed in Bugzilla accordingly).

> The Bessel functions of the second type (Yn) are not defined for
> negative x and should return NAN in these cases.  However, current

Not just NaN, NaN with the "invalid" exception raised.  (The logical 
design is that the main function implementations should get both return 
value and exceptions right, with the wrappers only needing to deal with 
setting errno and returning the return value from the main implementation 
unchanged.)

I think the appropriate fix is for these functions to return 1 / (x + x * 
x) for all non-finite input, like those for other floating-point types do.  
That achieves the desired results: propagation of NaN payloads for NaN 
inputs, with sNaN converted to qNaN, +Inf returning 0 quietly and -Inf 
returning a default NaN with "invalid" raised.

-- 
Joseph S. Myers
joseph@codesourcery.com

[RFC] Add math tests for libieee

[Tulio Magno Quites Machado Filho]

Joseph Myers <joseph@codesourcery.com> writes:

> On Wed, 8 Feb 2017, Gabriel F. T. Gomes wrote:
>
>> The code that this patch changes is not executed, because the wrappers
>> (in math/w_j0_compat.c and math/w_j1_compat.c) call __kernel_standard
>> and return.
>
> It's executed with -lieee (meaning this is user-visible and should have a 
> bug filed in Bugzilla accordingly).

In this case, should we be testing libieee too?

--- 8< ---

Reuse the math tests to test IEEE error handling rules.
Integrating this patch will cause 148 errors on each math type test for
the following ppc64, ppc64le and x86_64.

2017-02-13  Tulio Magno Quites Machado Filho  <tuliom@linux.vnet.ibm.com>

	* math/Makefile (libm-tests): Add ieee tests.
	(LDLIBS-test-float-ieee): New.
	(LDLIBS-test-double-ieee): Likewise.
	(LDLIBS-test-ldouble-ieee): Likewise.
	* math/libm-test-support.c (check_float_internal): Avoid testing
	errno under IEEE error handling rules.
	* math/test-double-ieee.c: New file.
	* math/test-float-ieee.c: Likewise.
	* math/test-ldouble-ieee.c: Likewise.
---
 math/Makefile            |  7 ++++++-
 math/libm-test-support.c |  4 +++-
 math/test-double-ieee.c  | 18 ++++++++++++++++++
 math/test-float-ieee.c   | 18 ++++++++++++++++++
 math/test-ldouble-ieee.c | 18 ++++++++++++++++++
 5 files changed, 63 insertions(+), 2 deletions(-)
 create mode 100644 math/test-double-ieee.c
 create mode 100644 math/test-float-ieee.c
 create mode 100644 math/test-ldouble-ieee.c

diff --git a/math/Makefile b/math/Makefile
index 2735212..c1d7676 100644
--- a/math/Makefile
+++ b/math/Makefile
@@ -193,7 +193,8 @@ endif
 ifneq (no,$(PERL))
 libm-vec-tests = $(addprefix test-,$(libmvec-tests))
 
-libm-tests = $(foreach t,$(types),test-$(t) test-$(t)-finite test-i$(t)) \
+libm-tests = $(foreach t,$(types),test-$(t) test-$(t)-finite test-i$(t) \
+				  test-$(t)-ieee) \
 	     $(libm-vec-tests)
 
 libm-tests.o = $(addsuffix .o,$(libm-tests))
@@ -334,6 +335,10 @@ CFLAGS-test-flt-eval-method.c = -fexcess-precision=standard
 
 CFLAGS-test-fe-snans-always-signal.c = -fsignaling-nans
 
+LDLIBS-test-float-ieee = -lieee
+LDLIBS-test-double-ieee = -lieee
+LDLIBS-test-ldouble-ieee = -lieee
+
 # The -lieee module sets the _LIB_VERSION_ switch to IEEE mode
 # for error handling in the -lm functions.
 install-lib += libieee.a
diff --git a/math/libm-test-support.c b/math/libm-test-support.c
index d387f81..dabf6bc 100644
--- a/math/libm-test-support.c
+++ b/math/libm-test-support.c
@@ -619,7 +619,9 @@ check_float_internal (const char *test_name, FLOAT computed, FLOAT expected,
   int errno_value = errno;
 
   test_exceptions (test_name, exceptions);
-  test_errno (test_name, errno_value, exceptions);
+  /* IEEE error handling (via libieee) doesn't support errno.  */
+  if (_LIB_VERSION != _IEEE_)
+    test_errno (test_name, errno_value, exceptions);
   if (exceptions & IGNORE_RESULT)
     goto out;
   if (issignaling (computed) && issignaling (expected))
diff --git a/math/test-double-ieee.c b/math/test-double-ieee.c
new file mode 100644
index 0000000..f638c45
--- /dev/null
+++ b/math/test-double-ieee.c
@@ -0,0 +1,18 @@
+/* Copyright (C) 2017 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+
+   The GNU C Library is free software; you can redistribute it and/or
+   modify it under the terms of the GNU Lesser General Public
+   License as published by the Free Software Foundation; either
+   version 2.1 of the License, or (at your option) any later version.
+
+   The GNU C Library is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   Lesser General Public License for more details.
+
+   You should have received a copy of the GNU Lesser General Public
+   License along with the GNU C Library; if not, see
+   <http://www.gnu.org/licenses/>.  */
+
+#include "test-double.c"
diff --git a/math/test-float-ieee.c b/math/test-float-ieee.c
new file mode 100644
index 0000000..88b68e1
--- /dev/null
+++ b/math/test-float-ieee.c
@@ -0,0 +1,18 @@
+/* Copyright (C) 2017 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+
+   The GNU C Library is free software; you can redistribute it and/or
+   modify it under the terms of the GNU Lesser General Public
+   License as published by the Free Software Foundation; either
+   version 2.1 of the License, or (at your option) any later version.
+
+   The GNU C Library is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   Lesser General Public License for more details.
+
+   You should have received a copy of the GNU Lesser General Public
+   License along with the GNU C Library; if not, see
+   <http://www.gnu.org/licenses/>.  */
+
+#include "test-float.c"
diff --git a/math/test-ldouble-ieee.c b/math/test-ldouble-ieee.c
new file mode 100644
index 0000000..e8c3cc6
--- /dev/null
+++ b/math/test-ldouble-ieee.c
@@ -0,0 +1,18 @@
+/* Copyright (C) 2017 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+
+   The GNU C Library is free software; you can redistribute it and/or
+   modify it under the terms of the GNU Lesser General Public
+   License as published by the Free Software Foundation; either
+   version 2.1 of the License, or (at your option) any later version.
+
+   The GNU C Library is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   Lesser General Public License for more details.
+
+   You should have received a copy of the GNU Lesser General Public
+   License along with the GNU C Library; if not, see
+   <http://www.gnu.org/licenses/>.  */
+
+#include "test-ldouble.c"
-- 
2.1.0

Re: [RFC] Add math tests for libieee

[Tulio Magno Quites Machado Filho]

Tulio Magno Quites Machado Filho <tuliom@linux.vnet.ibm.com> writes:

> diff --git a/math/libm-test-support.c b/math/libm-test-support.c
> index d387f81..dabf6bc 100644
> --- a/math/libm-test-support.c
> +++ b/math/libm-test-support.c
> @@ -619,7 +619,9 @@ check_float_internal (const char *test_name, FLOAT computed, FLOAT expected,
>    int errno_value = errno;
>
>    test_exceptions (test_name, exceptions);
> -  test_errno (test_name, errno_value, exceptions);
> +  /* IEEE error handling (via libieee) doesn't support errno.  */
> +  if (_LIB_VERSION != _IEEE_)
> +    test_errno (test_name, errno_value, exceptions);

Please, ignore this part.
I've just noticed the recent inclusion of TEST_ERRNO.
I can change this later, if it really makes sense to have these tests.

-- 
Tulio Magno

Re: [RFC] Add math tests for libieee

[Joseph Myers]

On Mon, 13 Feb 2017, Tulio Magno Quites Machado Filho wrote:

> Joseph Myers <joseph@codesourcery.com> writes:
> 
> > On Wed, 8 Feb 2017, Gabriel F. T. Gomes wrote:
> >
> >> The code that this patch changes is not executed, because the wrappers
> >> (in math/w_j0_compat.c and math/w_j1_compat.c) call __kernel_standard
> >> and return.
> >
> > It's executed with -lieee (meaning this is user-visible and should have a 
> > bug filed in Bugzilla accordingly).
> 
> In this case, should we be testing libieee too?

Since -lieee is an obsolescent feature I'd be wary of adding such tests.

In a world where _LIB_VERSION and matherr are compat symbols and libieee 
is removed, but without new symbol versions for existing libm functions in 
existing configurations of glibc, there is the issue of how to find issues 
hidden by the existing wrappers when testing on architectures that still 
use them in shared libm.  But the obvious approach for that would be 
testing static libm which in such a context would use the new wrappers for 
all types rather than the old ones.  (Right now testing static libm 
wouldn't achieve much, it's only in a case where static libm and shared 
libm use different wrappers for the same type that it becomes relevant.)

-- 
Joseph S. Myers
joseph@codesourcery.com