BZ 5857: Fix rint/nearbyint

BZ 5857: Fix rint/nearbyint

[Andreas Jaeger]

[-- Attachment #1: Type: text/plain, Size: 3274 bytes --]


rint and nearbyint gave in some cases wrong values.

Ok to commit?

A patch for ldbl-* might be needed as well - anybody interested to look
at it?

Andreas

2008-03-09  Andreas Jaeger  <aj@suse.de>

	[BZ #5857]
	* sysdeps/ieee754/dbl-64/s_rint.c (__rint): Handle j0==18.
	* sysdeps/ieee754/dbl-64/s_nearbyint.c (__nearbyint): Likewise.
  	Patch by Mark Elliott <mark.h.elliott@lmco.com>.

	* math/libm-test.inc (nearbyint_test): Add new test cases from #5857.
	(rint_test): Likewise.

============================================================
Index: sysdeps/ieee754/dbl-64/s_rint.c
--- sysdeps/ieee754/dbl-64/s_rint.c	13 Jul 1999 23:53:26 -0000	1.1
+++ sysdeps/ieee754/dbl-64/s_rint.c	9 Mar 2008 14:20:59 -0000
@@ -67,8 +67,11 @@ TWO52[2]={
 		if(((i0&i)|i1)==0) return x; /* x is integral */
 		i>>=1;
 		if(((i0&i)|i1)!=0) {
-		    if(j0==19) i1 = 0x40000000; else
-		    i0 = (i0&(~i))|((0x20000)>>j0);
+		  if(j0==19) i1 = 0x40000000; else if (j0<18)
+		    i0 = (i0&(~i))|((0x20000)>>j0); else {
+		    i0 &= ~i;
+		    i1 = 0x80000000;
+		  }
 		}
 	    }
 	} else if (j0>51) {
============================================================
Index: sysdeps/ieee754/dbl-64/s_nearbyint.c
--- sysdeps/ieee754/dbl-64/s_nearbyint.c	13 Jul 1999 23:53:03 -0000	1.1
+++ sysdeps/ieee754/dbl-64/s_nearbyint.c	9 Mar 2008 14:21:00 -0000
@@ -71,8 +71,11 @@ TWO52[2]={
 		if(((i0&i)|i1)==0) return x; /* x is integral */
 		i>>=1;
 		if(((i0&i)|i1)!=0) {
-		    if(j0==19) i1 = 0x40000000; else
-		    i0 = (i0&(~i))|((0x20000)>>j0);
+		  if(j0==19) i1 = 0x40000000; else if (j0<18)
+		    i0 = (i0&(~i))|((0x20000)>>j0); else {
+		    i0 &= ~i;
+		    i1 = 0x80000000;
+		  }
 		}
 	    }
 	} else if (j0>51) {
============================================================
Index: math/libm-test.inc
--- math/libm-test.inc	28 Oct 2007 01:14:25 -0000	1.77
+++ math/libm-test.inc	9 Mar 2008 14:21:00 -0000
@@ -1,4 +1,4 @@
-/* Copyright (C) 1997-2006, 2007 Free Software Foundation, Inc.
+/* Copyright (C) 1997-2006, 2007, 2008 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Andreas Jaeger <aj@suse.de>, 1997.
 
@@ -4481,6 +4481,11 @@ nearbyint_test (void)
   TEST_f_f (nearbyint, -0.5, minus_zero);
   TEST_f_f (nearbyint, -1.5, -2.0);
 
+  TEST_f_f (nearbyint, 262144.75, 262145.0);
+  TEST_f_f (nearbyint, 262142.75, 262143.0);
+  TEST_f_f (nearbyint, 524286.75, 524287.0);
+  TEST_f_f (nearbyint, 524288.75, 524289.0);
+
   END (nearbyint);
 }
 
@@ -4852,6 +4857,10 @@ rint_test (void)
   TEST_f_f (rint, -0.1, -0.0);
   TEST_f_f (rint, -0.25, -0.0);
   TEST_f_f (rint, -0.625, -1.0);
+  TEST_f_f (rint, 262144.75, 262145.0);
+  TEST_f_f (rint, 262142.75, 262143.0);
+  TEST_f_f (rint, 524286.75, 524287.0);
+  TEST_f_f (rint, 524288.75, 524289.0);
 #ifdef TEST_LDOUBLE
   /* The result can only be represented in long double.  */
   TEST_f_f (rint, 4503599627370495.5L, 4503599627370496.0L);

-- 
 Andreas Jaeger, Director Platform / openSUSE, aj@suse.de
  SUSE LINUX Products GmbH, GF: Markus Rex, HRB 16746 (AG Nürnberg)
   Maxfeldstr. 5, 90409 Nürnberg, Germany
    GPG fingerprint = 93A3 365E CE47 B889 DF7F  FED1 389A 563C C272 A126

[-- Attachment #2: Type: application/pgp-signature, Size: 193 bytes --]