mirror of
https://sourceware.org/git/glibc.git
synced 2024-11-14 01:00:07 +00:00
59a63cca11
ISO C requires overflowing results from nexttoward to be the appropriate infinity independent of the rounding mode, but some implementations use a rounding-mode-dependent result (this is the same issue as was fixed for nextafter in bug 16677). This patch fixes the problem by making the nexttoward implementations discard the result from the floating-point computation that forced an overflow exception and then return the infinity previously computed with integer arithmetic. Tested for x86_64, x86, mips64 and powerpc. [BZ #19059] * math/s_nexttowardf.c (__nexttowardf): Do not return value from overflowing computation. * sysdeps/i386/fpu/s_nexttoward.c (__nexttoward): Likewise. * sysdeps/i386/fpu/s_nexttowardf.c (__nexttowardf): Likewise. * sysdeps/ieee754/ldbl-128/s_nexttoward.c (__nexttoward): Likewise. * sysdeps/ieee754/ldbl-128/s_nexttowardf.c (__nexttowardf): Likewise. * sysdeps/ieee754/ldbl-128ibm/s_nexttoward.c (__nexttoward): Likewise. * sysdeps/ieee754/ldbl-128ibm/s_nexttowardf.c (__nexttowardf): Likewise. * sysdeps/ieee754/ldbl-96/s_nexttoward.c (__nexttoward): Likewise. * sysdeps/ieee754/ldbl-96/s_nexttowardf.c (__nexttowardf): Likewise. * sysdeps/ieee754/ldbl-opt/s_nexttowardfd.c (__nldbl_nexttowardf): Likewise. * math/libm-test.inc (nexttoward_test_data): Add more tests.
73 lines
1.9 KiB
C
73 lines
1.9 KiB
C
/* Single precision version of nexttoward.c.
|
|
Conversion to IEEE single float by Jakub Jelinek, jj@ultra.linux.cz. */
|
|
/*
|
|
* ====================================================
|
|
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
|
|
*
|
|
* Developed at SunPro, a Sun Microsystems, Inc. business.
|
|
* Permission to use, copy, modify, and distribute this
|
|
* software is freely granted, provided that this notice
|
|
* is preserved.
|
|
* ====================================================
|
|
*/
|
|
|
|
/* IEEE functions
|
|
* nexttowardf(x,y)
|
|
* return the next machine floating-point number of x in the
|
|
* direction toward y.
|
|
* This is for machines which use the same binary type for double and
|
|
* long double.
|
|
* Special cases:
|
|
*/
|
|
|
|
#include <math.h>
|
|
#include <math_private.h>
|
|
#include <float.h>
|
|
|
|
float __nexttowardf(float x, long double y)
|
|
{
|
|
int32_t hx,hy,ix,iy;
|
|
u_int32_t ly;
|
|
|
|
GET_FLOAT_WORD(hx,x);
|
|
EXTRACT_WORDS(hy,ly,y);
|
|
ix = hx&0x7fffffff; /* |x| */
|
|
iy = hy&0x7fffffff; /* |y| */
|
|
|
|
if((ix>0x7f800000) || /* x is nan */
|
|
((iy>=0x7ff00000)&&((iy-0x7ff00000)|ly)!=0)) /* y is nan */
|
|
return x+y;
|
|
if((long double) x==y) return y; /* x=y, return y */
|
|
if(ix==0) { /* x == 0 */
|
|
float u;
|
|
SET_FLOAT_WORD(x,(u_int32_t)(hy&0x80000000)|1);/* return +-minsub*/
|
|
u = math_opt_barrier (x);
|
|
u = u * u;
|
|
math_force_eval (u); /* raise underflow flag */
|
|
return x;
|
|
}
|
|
if(hx>=0) { /* x > 0 */
|
|
if(x > y) /* x -= ulp */
|
|
hx -= 1;
|
|
else /* x < y, x += ulp */
|
|
hx += 1;
|
|
} else { /* x < 0 */
|
|
if(x < y) /* x -= ulp */
|
|
hx -= 1;
|
|
else /* x > y, x += ulp */
|
|
hx += 1;
|
|
}
|
|
hy = hx&0x7f800000;
|
|
if(hy>=0x7f800000) {
|
|
float u = x+x; /* overflow */
|
|
math_force_eval (u);
|
|
}
|
|
if(hy<0x00800000) {
|
|
float u = x*x; /* underflow */
|
|
math_force_eval (u); /* raise underflow flag */
|
|
}
|
|
SET_FLOAT_WORD(x,hx);
|
|
return x;
|
|
}
|
|
weak_alias (__nexttowardf, nexttowardf)
|