This is the mail archive of the
libc-alpha@sources.redhat.com
mailing list for the glibc project.
sinhl patch
- To: libc-alpha at sources dot redhat dot com
- Subject: sinhl patch
- From: Stephen L Moshier <moshier at mediaone dot net>
- Date: Tue, 29 May 2001 08:07:14 -0400 (EDT)
- Reply-To: moshier at moshier dot ne dot mediaone dot net
This patch fixes a domain boundary bug in the 80-bit hyperbolic sine.
The program sets sinh(x) = x when x < 2^-28, but that threshold
should be lower, about 2^-32. Here is a test case, for x = 2^-29:
long double x = 1.86264514923095703125e-9L;
long double y = 1.862645149230957032327058e-9L;
extern long double sinhl (long double);
int
main()
{
long double z;
z = sinhl(x);
printf ("sinhl(2^-29)=%.20Le,\n should be %.20Le\n", z, y);
exit(0);
}
* sysdeps/ieee754/ldbl_96/e_sinhl.c (__ieee754_sinhl): sinhl(x) = x
when x < 2^-32.
*** e_sinhl.c 1999/07/14 00:11:45 1.1
--- e_sinhl.c 2001/05/27 21:39:27
*************** static long double one = 1.0, shuge = 1.
*** 66,72 ****
if (jx & 0x8000) h = -h;
/* |x| in [0,25], return sign(x)*0.5*(E+E/(E+1))) */
if (ix < 0x4003 || (ix == 0x4003 && i0 <= 0xc8000000)) { /* |x|<25 */
! if (ix<0x3fe3) /* |x|<2**-28 */
if(shuge+x>one) return x;/* sinh(tiny) = tiny with inexact */
t = __expm1l(fabsl(x));
if(ix<0x3fff) return h*(2.0*t-t*t/(t+one));
--- 66,72 ----
if (jx & 0x8000) h = -h;
/* |x| in [0,25], return sign(x)*0.5*(E+E/(E+1))) */
if (ix < 0x4003 || (ix == 0x4003 && i0 <= 0xc8000000)) { /* |x|<25 */
! if (ix<0x3fdf) /* |x|<2**-32 */
if(shuge+x>one) return x;/* sinh(tiny) = tiny with inexact */
t = __expm1l(fabsl(x));
if(ix<0x3fff) return h*(2.0*t-t*t/(t+one));