AuroraMiddleware/glibc
1
0
Fork 0
mirror of https://sourceware.org/git/glibc.git synced 2024-12-26 04:31:03 +00:00
Code Issues Packages Projects Releases Wiki Activity

Remove dbl-64/wordsize-64 (part 2)

Browse Source 
Remove the wordsize-64 implementations by merging them into the main dbl-64
directory.  The second patch just moves all wordsize-64 files and removes a
few wordsize-64 uses in comments and Implies files.

Reviewed-by: Adhemerval Zanella  <adhemerval.zanella@linaro.org>
This commit is contained in:
Wilco Dijkstra 2021-01-07 15:26:26 +00:00
parent caa884dda7
commit 9e97f239ea
40 changed files with 427 additions and 1848 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1
sysdeps/aarch64/Implies
View File

@ -1,5 +1,4 @@
wordsize-64 wordsize-64
ieee754/ldbl-128 ieee754/ldbl-128
ieee754/dbl-64/wordsize-64
ieee754/dbl-64 ieee754/dbl-64
ieee754/flt-32 ieee754/flt-32

1
sysdeps/alpha/Implies
View File

@ -1,6 +1,5 @@
wordsize-64 wordsize-64
# Alpha uses IEEE 754 single, double and quad precision floating point. # Alpha uses IEEE 754 single, double and quad precision floating point.
ieee754/ldbl-128 ieee754/ldbl-128
ieee754/dbl-64/wordsize-64
ieee754/dbl-64 ieee754/dbl-64
ieee754/flt-32 ieee754/flt-32

50
sysdeps/ieee754/dbl-64/e_acosh.c
View File

@ -1,4 +1,4 @@
/* @(#)e_acosh.c 5.1 93/09/24 */ /* Optimized for 64-bit by Ulrich Drepper <drepper@gmail.com>, 2012 */
/* /*
* ==================================================== * ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
@ -29,42 +29,40 @@
#include <libm-alias-finite.h> #include <libm-alias-finite.h>
static const double static const double
one = 1.0, one = 1.0,
ln2 = 6.93147180559945286227e-01; /* 0x3FE62E42, 0xFEFA39EF */ ln2 = 6.93147180559945286227e-01; /* 0x3FE62E42, 0xFEFA39EF */
double double
__ieee754_acosh (double x) __ieee754_acosh (double x)
{ {
double t; int64_t hx;
int32_t hx; EXTRACT_WORDS64 (hx, x);
uint32_t lx;
EXTRACT_WORDS (hx, lx, x); if (hx > INT64_C (0x4000000000000000))
if (hx < 0x3ff00000) /* x < 1 */
{ {
return (x - x) / (x - x); if (__glibc_unlikely (hx >= INT64_C (0x41b0000000000000)))
}
else if (hx >= 0x41b00000) /* x > 2**28 */
{
if (hx >= 0x7ff00000) /* x is inf of NaN */
{ {
return x + x; /* x > 2**28 */
if (hx >= INT64_C (0x7ff0000000000000))
/* x is inf of NaN */
return x + x;
else
return __ieee754_log (x) + ln2;/* acosh(huge)=log(2x) */
} }
else
return __ieee754_log (x) + ln2; /* acosh(huge)=log(2x) */ /* 2**28 > x > 2 */
} double t = x * x;
else if (((hx - 0x3ff00000) | lx) == 0)
{
return 0.0; /* acosh(1) = 0 */
}
else if (hx > 0x40000000) /* 2**28 > x > 2 */
{
t = x * x;
return __ieee754_log (2.0 * x - one / (x + sqrt (t - one))); return __ieee754_log (2.0 * x - one / (x + sqrt (t - one)));
} }
else /* 1<x<2 */ else if (__glibc_likely (hx > INT64_C (0x3ff0000000000000)))
{ {
t = x - one; /* 1<x<2 */
double t = x - one;
return __log1p (t + sqrt (2.0 * t + t * t)); return __log1p (t + sqrt (2.0 * t + t * t));
} }
else if (__glibc_likely (hx == INT64_C (0x3ff0000000000000)))
return 0.0; /* acosh(1) = 0 */
else /* x < 1 */
return (x - x) / (x - x);
} }
libm_alias_finite (__ieee754_acosh, __acosh) libm_alias_finite (__ieee754_acosh, __acosh)

75
sysdeps/ieee754/dbl-64/e_cosh.c
View File

@ -32,59 +32,54 @@
*/ */
#include <math.h> #include <math.h>
#include <math-narrow-eval.h>
#include <math_private.h> #include <math_private.h>
#include <libm-alias-finite.h> #include <libm-alias-finite.h>
static const double one = 1.0, half = 0.5, huge = 1.0e300; static const double one = 1.0, half=0.5, huge = 1.0e300;
double double
__ieee754_cosh (double x) __ieee754_cosh (double x)
{ {
double t, w; double t,w;
int32_t ix; int32_t ix;
uint32_t lx;
/* High word of |x|. */ /* High word of |x|. */
GET_HIGH_WORD (ix, x); GET_HIGH_WORD(ix,x);
ix &= 0x7fffffff; ix &= 0x7fffffff;
/* |x| in [0,22] */ /* |x| in [0,22] */
if (ix < 0x40360000) if (ix < 0x40360000) {
{ /* |x| in [0,0.5*ln2], return 1+expm1(|x|)^2/(2*exp(|x|)) */
/* |x| in [0,0.5*ln2], return 1+expm1(|x|)^2/(2*exp(|x|)) */ if(ix<0x3fd62e43) {
if (ix < 0x3fd62e43) if (ix<0x3c800000) /* cosh(tiny) = 1 */
{ return one;
if (ix < 0x3c800000) t = __expm1(fabs(x));
return one; /* cosh(tiny) = 1 */ w = one+t;
t = __expm1 (fabs (x)); return one+(t*t)/(w+w);
w = one + t; }
return one + (t * t) / (w + w);
/* |x| in [0.5*ln2,22], return (exp(|x|)+1/exp(|x|)/2; */
t = __ieee754_exp(fabs(x));
return half*t+half/t;
} }
/* |x| in [0.5*ln2,22], return (exp(|x|)+1/exp(|x|)/2; */ /* |x| in [22, log(maxdouble)] return half*exp(|x|) */
t = __ieee754_exp (fabs (x)); if (ix < 0x40862e42) return half*__ieee754_exp(fabs(x));
return half * t + half / t;
}
/* |x| in [22, log(maxdouble)] return half*exp(|x|) */ /* |x| in [log(maxdouble), overflowthresold] */
if (ix < 0x40862e42) int64_t fix;
return half * __ieee754_exp (fabs (x)); EXTRACT_WORDS64(fix, x);
fix &= UINT64_C(0x7fffffffffffffff);
if (fix <= UINT64_C(0x408633ce8fb9f87d)) {
w = __ieee754_exp(half*fabs(x));
t = half*w;
return t*w;
}
/* |x| in [log(maxdouble), overflowthresold] */ /* x is INF or NaN */
GET_LOW_WORD (lx, x); if(ix>=0x7ff00000) return x*x;
if (ix < 0x408633ce || ((ix == 0x408633ce) && (lx <= (uint32_t) 0x8fb9f87d)))
{
w = __ieee754_exp (half * fabs (x));
t = half * w;
return t * w;
}
/* x is INF or NaN */ /* |x| > overflowthresold, cosh(x) overflow */
if (ix >= 0x7ff00000) return huge*huge;
return x * x;
/* |x| > overflowthresold, cosh(x) overflow */
return math_narrow_eval (huge * huge);
} }
libm_alias_finite (__ieee754_cosh, __cosh) libm_alias_finite (__ieee754_cosh, __cosh)

212
sysdeps/ieee754/dbl-64/e_fmod.c
View File

@ -1,3 +1,4 @@
/* Rewritten for 64-bit machines by Ulrich Drepper <drepper@gmail.com>. */
/* /*
* ==================================================== * ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
@ -17,158 +18,89 @@
#include <math.h> #include <math.h>
#include <math_private.h> #include <math_private.h>
#include <stdint.h>
#include <libm-alias-finite.h> #include <libm-alias-finite.h>
static const double one = 1.0, Zero[] = { 0.0, -0.0, }; static const double one = 1.0, Zero[] = {0.0, -0.0,};
double double
__ieee754_fmod (double x, double y) __ieee754_fmod (double x, double y)
{ {
int32_t n, hx, hy, hz, ix, iy, sx, i; int32_t n,ix,iy;
uint32_t lx, ly, lz; int64_t hx,hy,hz,sx,i;
EXTRACT_WORDS (hx, lx, x); EXTRACT_WORDS64(hx,x);
EXTRACT_WORDS (hy, ly, y); EXTRACT_WORDS64(hy,y);
sx = hx & 0x80000000; /* sign of x */ sx = hx&UINT64_C(0x8000000000000000); /* sign of x */
hx ^= sx; /* |x| */ hx ^=sx; /* |x| */
hy &= 0x7fffffff; /* |y| */ hy &= UINT64_C(0x7fffffffffffffff); /* |y| */
/* purge off exception values */ /* purge off exception values */
if ((hy | ly) == 0 || (hx >= 0x7ff00000) || /* y=0,or x not finite */ if(__builtin_expect(hy==0
((hy | ((ly | -ly) >> 31)) > 0x7ff00000)) /* or y is NaN */ || hx >= UINT64_C(0x7ff0000000000000)
return (x * y) / (x * y); || hy > UINT64_C(0x7ff0000000000000), 0))
if (hx <= hy) /* y=0,or x not finite or y is NaN */
{ return (x*y)/(x*y);
if ((hx < hy) || (lx < ly)) if(__builtin_expect(hx<=hy, 0)) {
return x; /* |x|<|y| return x */ if(hx<hy) return x; /* |x|<|y| return x */
if (lx == ly) return Zero[(uint64_t)sx>>63]; /* |x|=|y| return x*0*/
return Zero[(uint32_t) sx >> 31]; /* |x|=|y| return x*0*/
}
/* determine ix = ilogb(x) */
if (__glibc_unlikely (hx < 0x00100000)) /* subnormal x */
{
if (hx == 0)
{
for (ix = -1043, i = lx; i > 0; i <<= 1)
ix -= 1;
} }
else
{
for (ix = -1022, i = (hx << 11); i > 0; i <<= 1)
ix -= 1;
}
}
else
ix = (hx >> 20) - 1023;
/* determine iy = ilogb(y) */ /* determine ix = ilogb(x) */
if (__glibc_unlikely (hy < 0x00100000)) /* subnormal y */ if(__builtin_expect(hx<UINT64_C(0x0010000000000000), 0)) {
{ /* subnormal x */
if (hy == 0) for (ix = -1022,i=(hx<<11); i>0; i<<=1) ix -=1;
{ } else ix = (hx>>52)-1023;
for (iy = -1043, i = ly; i > 0; i <<= 1)
/* determine iy = ilogb(y) */
if(__builtin_expect(hy<UINT64_C(0x0010000000000000), 0)) { /* subnormal y */
for (iy = -1022,i=(hy<<11); i>0; i<<=1) iy -=1;
} else iy = (hy>>52)-1023;
/* set up hx, hy and align y to x */
if(__builtin_expect(ix >= -1022, 1))
hx = UINT64_C(0x0010000000000000)|(UINT64_C(0x000fffffffffffff)&hx);
else { /* subnormal x, shift x to normal */
n = -1022-ix;
hx<<=n;
}
if(__builtin_expect(iy >= -1022, 1))
hy = UINT64_C(0x0010000000000000)|(UINT64_C(0x000fffffffffffff)&hy);
else { /* subnormal y, shift y to normal */
n = -1022-iy;
hy<<=n;
}
/* fix point fmod */
n = ix - iy;
while(n--) {
hz=hx-hy;
if(hz<0){hx = hx+hx;}
else {
if(hz==0) /* return sign(x)*0 */
return Zero[(uint64_t)sx>>63];
hx = hz+hz;
}
}
hz=hx-hy;
if(hz>=0) {hx=hz;}
/* convert back to floating value and restore the sign */
if(hx==0) /* return sign(x)*0 */
return Zero[(uint64_t)sx>>63];
while(hx<UINT64_C(0x0010000000000000)) { /* normalize x */
hx = hx+hx;
iy -= 1; iy -= 1;
} }
else if(__builtin_expect(iy>= -1022, 1)) { /* normalize output */
{ hx = ((hx-UINT64_C(0x0010000000000000))|((uint64_t)(iy+1023)<<52));
for (iy = -1022, i = (hy << 11); i > 0; i <<= 1) INSERT_WORDS64(x,hx|sx);
iy -= 1; } else { /* subnormal output */
n = -1022 - iy;
hx>>=n;
INSERT_WORDS64(x,hx|sx);
x *= one; /* create necessary signal */
} }
} return x; /* exact output */
else
iy = (hy >> 20) - 1023;
/* set up {hx,lx}, {hy,ly} and align y to x */
if (__glibc_likely (ix >= -1022))
hx = 0x00100000 | (0x000fffff & hx);
else /* subnormal x, shift x to normal */
{
n = -1022 - ix;
if (n <= 31)
{
hx = (hx << n) | (lx >> (32 - n));
lx <<= n;
}
else
{
hx = lx << (n - 32);
lx = 0;
}
}
if (__glibc_likely (iy >= -1022))
hy = 0x00100000 | (0x000fffff & hy);
else /* subnormal y, shift y to normal */
{
n = -1022 - iy;
if (n <= 31)
{
hy = (hy << n) | (ly >> (32 - n));
ly <<= n;
}
else
{
hy = ly << (n - 32);
ly = 0;
}
}
/* fix point fmod */
n = ix - iy;
while (n--)
{
hz = hx - hy; lz = lx - ly; if (lx < ly)
hz -= 1;
if (hz < 0)
{
hx = hx + hx + (lx >> 31); lx = lx + lx;
}
else
{
if ((hz | lz) == 0) /* return sign(x)*0 */
return Zero[(uint32_t) sx >> 31];
hx = hz + hz + (lz >> 31); lx = lz + lz;
}
}
hz = hx - hy; lz = lx - ly; if (lx < ly)
hz -= 1;
if (hz >= 0)
{
hx = hz; lx = lz;
}
/* convert back to floating value and restore the sign */
if ((hx | lx) == 0) /* return sign(x)*0 */
return Zero[(uint32_t) sx >> 31];
while (hx < 0x00100000) /* normalize x */
{
hx = hx + hx + (lx >> 31); lx = lx + lx;
iy -= 1;
}
if (__glibc_likely (iy >= -1022)) /* normalize output */
{
hx = ((hx - 0x00100000) | ((iy + 1023) << 20));
INSERT_WORDS (x, hx | sx, lx);
}
else /* subnormal output */
{
n = -1022 - iy;
if (n <= 20)
{
lx = (lx >> n) | ((uint32_t) hx << (32 - n));
hx >>= n;
}
else if (n <= 31)
{
lx = (hx << (32 - n)) | (lx >> n); hx = sx;
}
else
{
lx = hx >> (n - 32); hx = sx;
}
INSERT_WORDS (x, hx | sx, lx);
x *= one; /* create necessary signal */
}
return x; /* exact output */
} }
libm_alias_finite (__ieee754_fmod, __fmod) libm_alias_finite (__ieee754_fmod, __fmod)

42
sysdeps/ieee754/dbl-64/e_log10.c
View File

@ -44,44 +44,46 @@
*/ */
#include <math.h> #include <math.h>
#include <math_private.h>
#include <fix-int-fp-convert-zero.h> #include <fix-int-fp-convert-zero.h>
#include <math_private.h>
#include <stdint.h>
#include <libm-alias-finite.h> #include <libm-alias-finite.h>
static const double two54 = 1.80143985094819840000e+16; /* 0x43500000, 0x00000000 */ static const double two54 = 1.80143985094819840000e+16; /* 0x4350000000000000 */
static const double ivln10 = 4.34294481903251816668e-01; /* 0x3FDBCB7B, 0x1526E50E */ static const double ivln10 = 4.34294481903251816668e-01; /* 0x3FDBCB7B1526E50E */
static const double log10_2hi = 3.01029995663611771306e-01; /* 0x3FD34413, 0x509F6000 */ static const double log10_2hi = 3.01029995663611771306e-01; /* 0x3FD34413509F6000 */
static const double log10_2lo = 3.69423907715893078616e-13; /* 0x3D59FEF3, 0x11F12B36 */ static const double log10_2lo = 3.69423907715893078616e-13; /* 0x3D59FEF311F12B36 */
double double
__ieee754_log10 (double x) __ieee754_log10 (double x)
{ {
double y, z; double y, z;
int32_t i, k, hx; int64_t i, hx;
uint32_t lx; int32_t k;
EXTRACT_WORDS (hx, lx, x); EXTRACT_WORDS64 (hx, x);
k = 0; k = 0;
if (hx < 0x00100000) if (hx < INT64_C(0x0010000000000000))
{ /* x < 2**-1022 */ { /* x < 2**-1022 */
if (__glibc_unlikely (((hx & 0x7fffffff) | lx) == 0)) if (__glibc_unlikely ((hx & UINT64_C(0x7fffffffffffffff)) == 0))
return -two54 / fabs (x); /* log(+-0)=-inf */ return -two54 / fabs (x); /* log(+-0)=-inf */
if (__glibc_unlikely (hx < 0)) if (__glibc_unlikely (hx < 0))
return (x - x) / (x - x); /* log(-#) = NaN */ return (x - x) / (x - x); /* log(-#) = NaN */
k -= 54; k -= 54;
x *= two54; /* subnormal number, scale up x */ x *= two54; /* subnormal number, scale up x */
GET_HIGH_WORD (hx, x); EXTRACT_WORDS64 (hx, x);
} }
if (__glibc_unlikely (hx >= 0x7ff00000)) /* scale up resulted in a NaN number */
if (__glibc_unlikely (hx >= UINT64_C(0x7ff0000000000000)))
return x + x; return x + x;
k += (hx >> 20) - 1023; k += (hx >> 52) - 1023;
i = ((uint32_t) k & 0x80000000) >> 31; i = ((uint64_t) k & UINT64_C(0x8000000000000000)) >> 63;
hx = (hx & 0x000fffff) | ((0x3ff - i) << 20); hx = (hx & UINT64_C(0x000fffffffffffff)) | ((0x3ff - i) << 52);
y = (double) (k + i); y = (double) (k + i);
if (FIX_INT_FP_CONVERT_ZERO && y == 0.0) if (FIX_INT_FP_CONVERT_ZERO && y == 0.0)
y = 0.0; y = 0.0;
SET_HIGH_WORD (x, hx); INSERT_WORDS64 (x, hx);
z = y * log10_2lo + ivln10 * __ieee754_log (x); z = y * log10_2lo + ivln10 * __ieee754_log (x);
return z + y * log10_2hi; return z + y * log10_2hi;
} }

81
sysdeps/ieee754/dbl-64/s_frexp.c
View File

@ -1,21 +1,28 @@
/* @(#)s_frexp.c 5.1 93/09/24 */ /* Copyright (C) 2011-2021 Free Software Foundation, Inc.
/* This file is part of the GNU C Library.
* ==================================================== Contributed by Ulrich Drepper <drepper@gmail.com>, 2011.
* 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.
* ====================================================
*/
#if defined(LIBM_SCCS) && !defined(lint) The GNU C Library is free software; you can redistribute it and/or
static char rcsid[] = "$NetBSD: s_frexp.c,v 1.9 1995/05/10 20:47:24 jtc Exp $"; modify it under the terms of the GNU Lesser General Public
#endif 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
<https://www.gnu.org/licenses/>. */
#include <inttypes.h>
#include <math.h>
#include <math_private.h>
#include <libm-alias-double.h>
/* /*
* for non-zero x * for non-zero, finite x
* x = frexp(arg,&exp); * x = frexp(arg,&exp);
* return a double fp quantity x such that 0.5 <= |x| <1.0 * return a double fp quantity x such that 0.5 <= |x| <1.0
* and the corresponding binary exponent "exp". That is * and the corresponding binary exponent "exp". That is
@ -24,32 +31,36 @@ static char rcsid[] = "$NetBSD: s_frexp.c,v 1.9 1995/05/10 20:47:24 jtc Exp $";
* with *exp=0. * with *exp=0.
*/ */
#include <math.h>
#include <math_private.h>
#include <libm-alias-double.h>
static const double
two54 = 1.80143985094819840000e+16; /* 0x43500000, 0x00000000 */
double double
__frexp (double x, int *eptr) __frexp (double x, int *eptr)
{ {
int32_t hx, ix, lx; int64_t ix;
EXTRACT_WORDS (hx, lx, x); EXTRACT_WORDS64 (ix, x);
ix = 0x7fffffff & hx; int32_t ex = 0x7ff & (ix >> 52);
*eptr = 0; int e = 0;
if (ix >= 0x7ff00000 || ((ix | lx) == 0))
return x + x; /* 0,inf,nan */ if (__glibc_likely (ex != 0x7ff && x != 0.0))
if (ix < 0x00100000) /* subnormal */
{ {
x *= two54; /* Not zero and finite. */
GET_HIGH_WORD (hx, x); e = ex - 1022;
ix = hx & 0x7fffffff; if (__glibc_unlikely (ex == 0))
*eptr = -54; {
/* Subnormal. */
x *= 0x1p54;
EXTRACT_WORDS64 (ix, x);
ex = 0x7ff & (ix >> 52);
e = ex - 1022 - 54;
}
ix = (ix & INT64_C (0x800fffffffffffff)) | INT64_C (0x3fe0000000000000);
INSERT_WORDS64 (x, ix);
} }
*eptr += (ix >> 20) - 1022; else
hx = (hx & 0x800fffff) | 0x3fe00000; /* Quiet signaling NaNs. */
SET_HIGH_WORD (x, hx); x += x;
*eptr = e;
return x; return x;
} }
libm_alias_double (__frexp, frexp) libm_alias_double (__frexp, frexp)

15
sysdeps/ieee754/dbl-64/s_getpayload.c
View File

@ -1,4 +1,4 @@
/* Get NaN payload. dbl-64 version. /* Get NaN payload.
Copyright (C) 2016-2021 Free Software Foundation, Inc. Copyright (C) 2016-2021 Free Software Foundation, Inc.
This file is part of the GNU C Library. This file is part of the GNU C Library.
@ -16,22 +16,21 @@
License along with the GNU C Library; if not, see License along with the GNU C Library; if not, see
<https://www.gnu.org/licenses/>. */ <https://www.gnu.org/licenses/>. */
#include <fix-int-fp-convert-zero.h>
#include <math.h> #include <math.h>
#include <math_private.h> #include <math_private.h>
#include <libm-alias-double.h> #include <libm-alias-double.h>
#include <stdint.h> #include <stdint.h>
#include <fix-int-fp-convert-zero.h>
double double
__getpayload (const double *x) __getpayload (const double *x)
{ {
uint32_t hx, lx; uint64_t ix;
EXTRACT_WORDS (hx, lx, *x); EXTRACT_WORDS64 (ix, *x);
if ((hx & 0x7ff00000) != 0x7ff00000 if ((ix & 0x7ff0000000000000ULL) != 0x7ff0000000000000ULL
|| ((hx & 0xfffff) | lx) == 0) || (ix & 0xfffffffffffffULL) == 0)
return -1; return -1;
hx &= 0x7ffff; ix &= 0x7ffffffffffffULL;
uint64_t ix = ((uint64_t) hx << 32) | lx;
if (FIX_INT_FP_CONVERT_ZERO && ix == 0) if (FIX_INT_FP_CONVERT_ZERO && ix == 0)
return 0.0f; return 0.0f;
return (double) ix; return (double) ix;

14
sysdeps/ieee754/dbl-64/s_issignaling.c
View File

@ -23,25 +23,21 @@
int int
__issignaling (double x) __issignaling (double x)
{ {
uint64_t xi;
EXTRACT_WORDS64 (xi, x);
#if HIGH_ORDER_BIT_IS_SET_FOR_SNAN #if HIGH_ORDER_BIT_IS_SET_FOR_SNAN
uint32_t hxi;
GET_HIGH_WORD (hxi, x);
/* We only have to care about the high-order bit of x's significand, because /* We only have to care about the high-order bit of x's significand, because
having it set (sNaN) already makes the significand different from that having it set (sNaN) already makes the significand different from that
used to designate infinity. */ used to designate infinity. */
return (hxi & 0x7ff80000) == 0x7ff80000; return (xi & UINT64_C (0x7ff8000000000000)) == UINT64_C (0x7ff8000000000000);
#else #else
uint32_t hxi, lxi;
EXTRACT_WORDS (hxi, lxi, x);
/* To keep the following comparison simple, toggle the quiet/signaling bit, /* To keep the following comparison simple, toggle the quiet/signaling bit,
so that it is set for sNaNs. This is inverse to IEEE 754-2008 (as well as so that it is set for sNaNs. This is inverse to IEEE 754-2008 (as well as
common practice for IEEE 754-1985). */ common practice for IEEE 754-1985). */
hxi ^= 0x00080000; xi ^= UINT64_C (0x0008000000000000);
/* If lxi != 0, then set any suitable bit of the significand in hxi. */
hxi |= (lxi | -lxi) >> 31;
/* We have to compare for greater (instead of greater or equal), because x's /* We have to compare for greater (instead of greater or equal), because x's
significand being all-zero designates infinity not NaN. */ significand being all-zero designates infinity not NaN. */
return (hxi & 0x7fffffff) > 0x7ff80000; return (xi & UINT64_C (0x7fffffffffffffff)) > UINT64_C (0x7ff8000000000000);
#endif #endif
} }
libm_hidden_def (__issignaling) libm_hidden_def (__issignaling)

49
sysdeps/ieee754/dbl-64/s_llround.c
View File

@ -17,54 +17,43 @@
License along with the GNU C Library; if not, see License along with the GNU C Library; if not, see
<https://www.gnu.org/licenses/>. */ <https://www.gnu.org/licenses/>. */
#define lround __hidden_lround
#define __lround __hidden___lround
#include <fenv.h> #include <fenv.h>
#include <limits.h> #include <limits.h>
#include <math.h> #include <math.h>
#include <sysdep.h>
#include <math_private.h> #include <math_private.h>
#include <libm-alias-double.h> #include <libm-alias-double.h>
#include <fix-fp-int-convert-overflow.h> #include <fix-fp-int-convert-overflow.h>
long long int long long int
__llround (double x) __llround (double x)
{ {
int32_t j0; int32_t j0;
uint32_t i1, i0; int64_t i0;
long long int result; long long int result;
int sign; int sign;
EXTRACT_WORDS (i0, i1, x); EXTRACT_WORDS64 (i0, x);
j0 = ((i0 >> 20) & 0x7ff) - 0x3ff; j0 = ((i0 >> 52) & 0x7ff) - 0x3ff;
sign = (i0 & 0x80000000) != 0 ? -1 : 1; sign = i0 < 0 ? -1 : 1;
i0 &= 0xfffff; i0 &= UINT64_C(0xfffffffffffff);
i0 |= 0x100000; i0 |= UINT64_C(0x10000000000000);
if (j0 < 20) if (j0 < (int32_t) (8 * sizeof (long long int)) - 1)
{ {
if (j0 < 0) if (j0 < 0)
return j0 < -1 ? 0 : sign; return j0 < -1 ? 0 : sign;
else if (j0 >= 52)
result = i0 << (j0 - 52);
else else
{ {
i0 += 0x80000 >> j0; i0 += UINT64_C(0x8000000000000) >> j0;
result = i0 >> (20 - j0); result = i0 >> (52 - j0);
}
}
else if (j0 < (int32_t) (8 * sizeof (long long int)) - 1)
{
if (j0 >= 52)
result = (((long long int) i0 << 32) | i1) << (j0 - 52);
else
{
uint32_t j = i1 + (0x80000000 >> (j0 - 20));
if (j < i1)
++i0;
if (j0 == 20)
result = (long long int) i0;
else
result = ((long long int) i0 << (j0 - 20)) | (j >> (52 - j0));
} }
} }
else else
@ -86,3 +75,11 @@ __llround (double x)
} }
libm_alias_double (__llround, llround) libm_alias_double (__llround, llround)
/* long has the same width as long long on LP64 machines, so use an alias. */
#undef lround
#undef __lround
#ifdef _LP64
strong_alias (__llround, __lround)
libm_alias_double (__lround, lround)
#endif

57
sysdeps/ieee754/dbl-64/s_lround.c
View File

@ -1,7 +1,6 @@
/* Round double value to long int. /* Round double value to long int.
Copyright (C) 1997-2021 Free Software Foundation, Inc. Copyright (C) 1997-2021 Free Software Foundation, Inc.
This file is part of the GNU C Library. This file is part of the GNU C Library.
Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
The GNU C Library is free software; you can redistribute it and/or The GNU C Library is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public modify it under the terms of the GNU Lesser General Public
@ -25,55 +24,41 @@
#include <libm-alias-double.h> #include <libm-alias-double.h>
#include <fix-fp-int-convert-overflow.h> #include <fix-fp-int-convert-overflow.h>
/* For LP64, lround is an alias for llround. */
#ifndef _LP64
long int long int
__lround (double x) __lround (double x)
{ {
int32_t j0; int32_t j0;
uint32_t i1, i0; int64_t i0;
long int result; long int result;
int sign; int sign;
EXTRACT_WORDS (i0, i1, x); EXTRACT_WORDS64 (i0, x);
j0 = ((i0 >> 20) & 0x7ff) - 0x3ff; j0 = ((i0 >> 52) & 0x7ff) - 0x3ff;
sign = (i0 & 0x80000000) != 0 ? -1 : 1; sign = i0 < 0 ? -1 : 1;
i0 &= 0xfffff; i0 &= UINT64_C(0xfffffffffffff);
i0 |= 0x100000; i0 |= UINT64_C(0x10000000000000);
if (j0 < 20) if (j0 < (int32_t) (8 * sizeof (long int)) - 1)
{ {
if (j0 < 0) if (j0 < 0)
return j0 < -1 ? 0 : sign; return j0 < -1 ? 0 : sign;
else if (j0 >= 52)
result = i0 << (j0 - 52);
else else
{ {
i0 += 0x80000 >> j0; i0 += UINT64_C(0x8000000000000) >> j0;
result = i0 >> (20 - j0); result = i0 >> (52 - j0);
}
}
else if (j0 < (int32_t) (8 * sizeof (long int)) - 1)
{
if (j0 >= 52)
result = ((long int) i0 << (j0 - 20)) | ((long int) i1 << (j0 - 52));
else
{
uint32_t j = i1 + (0x80000000 >> (j0 - 20));
if (j < i1)
++i0;
if (j0 == 20)
result = (long int) i0;
else
{
result = ((long int) i0 << (j0 - 20)) | (j >> (52 - j0));
#ifdef FE_INVALID #ifdef FE_INVALID
if (sizeof (long int) == 4 if (sizeof (long int) == 4
&& sign == 1 && sign == 1
&& result == LONG_MIN) && result == LONG_MIN)
/* Rounding brought the value out of range. */ /* Rounding brought the value out of range. */
feraiseexcept (FE_INVALID); feraiseexcept (FE_INVALID);
#endif #endif
}
} }
} }
else else
@ -92,8 +77,8 @@ __lround (double x)
return sign == 1 ? LONG_MAX : LONG_MIN; return sign == 1 ? LONG_MAX : LONG_MIN;
} }
else if (!FIX_DBL_LONG_CONVERT_OVERFLOW else if (!FIX_DBL_LONG_CONVERT_OVERFLOW
&& sizeof (long int) == 4 && sizeof (long int) == 4
&& x <= (double) LONG_MIN - 0.5) && x <= (double) LONG_MIN - 0.5)
{ {
/* If truncation produces LONG_MIN, the cast will not raise /* If truncation produces LONG_MIN, the cast will not raise
the exception, but may raise "inexact". */ the exception, but may raise "inexact". */
@ -108,3 +93,5 @@ __lround (double x)
} }
libm_alias_double (__lround, lround) libm_alias_double (__lround, lround)
#endif

80
sysdeps/ieee754/dbl-64/s_modf.c
View File

@ -1,3 +1,4 @@
/* Rewritten for 64-bit machines by Ulrich Drepper <drepper@gmail.com>. */
/* /*
* ==================================================== * ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
@ -22,63 +23,42 @@
#include <math.h> #include <math.h>
#include <math_private.h> #include <math_private.h>
#include <libm-alias-double.h> #include <libm-alias-double.h>
#include <stdint.h>
static const double one = 1.0; static const double one = 1.0;
double double
__modf (double x, double *iptr) __modf(double x, double *iptr)
{ {
int32_t i0, i1, j0; int64_t i0;
uint32_t i; int32_t j0;
EXTRACT_WORDS (i0, i1, x); EXTRACT_WORDS64(i0,x);
j0 = ((i0 >> 20) & 0x7ff) - 0x3ff; /* exponent of x */ j0 = ((i0>>52)&0x7ff)-0x3ff; /* exponent of x */
if (j0 < 20) /* integer part in high x */ if(j0<52) { /* integer part in x */
{ if(j0<0) { /* |x|<1 */
if (j0 < 0) /* |x|<1 */ /* *iptr = +-0 */
{ INSERT_WORDS64(*iptr,i0&UINT64_C(0x8000000000000000));
INSERT_WORDS (*iptr, i0 & 0x80000000, 0); /* *iptr = +-0 */ return x;
return x; } else {
} uint64_t i = UINT64_C(0x000fffffffffffff)>>j0;
else if((i0&i)==0) { /* x is integral */
{ *iptr = x;
i = (0x000fffff) >> j0; /* return +-0 */
if (((i0 & i) | i1) == 0) /* x is integral */ INSERT_WORDS64(x,i0&UINT64_C(0x8000000000000000));
{ return x;
*iptr = x; } else {
INSERT_WORDS (x, i0 & 0x80000000, 0); /* return +-0 */ INSERT_WORDS64(*iptr,i0&(~i));
return x; return x - *iptr;
} }
else
{
INSERT_WORDS (*iptr, i0 & (~i), 0);
return x - *iptr;
} }
} else { /* no fraction part */
*iptr = x*one;
/* We must handle NaNs separately. */
if (j0 == 0x400 && (i0 & UINT64_C(0xfffffffffffff)))
return x*one;
INSERT_WORDS64(x,i0&UINT64_C(0x8000000000000000)); /* return +-0 */
return x;
} }
}
else if (__glibc_unlikely (j0 > 51)) /* no fraction part */
{
*iptr = x * one;
/* We must handle NaNs separately. */
if (j0 == 0x400 && ((i0 & 0xfffff) | i1))
return x * one;
INSERT_WORDS (x, i0 & 0x80000000, 0); /* return +-0 */
return x;
}
else /* fraction part in low x */
{
i = ((uint32_t) (0xffffffff)) >> (j0 - 20);
if ((i1 & i) == 0) /* x is integral */
{
*iptr = x;
INSERT_WORDS (x, i0 & 0x80000000, 0); /* return +-0 */
return x;
}
else
{
INSERT_WORDS (*iptr, i0, i1 & (~i));
return x - *iptr;
}
}
} }
#ifndef __modf #ifndef __modf
libm_alias_double (__modf, modf) libm_alias_double (__modf, modf)

43
sysdeps/ieee754/dbl-64/s_remquo.c
View File

@ -21,7 +21,7 @@
#include <math_private.h> #include <math_private.h>
#include <libm-alias-double.h> #include <libm-alias-double.h>
#include <stdint.h>
static const double zero = 0.0; static const double zero = 0.0;
@ -29,50 +29,49 @@ static const double zero = 0.0;
double double
__remquo (double x, double y, int *quo) __remquo (double x, double y, int *quo)
{ {
int32_t hx, hy; int64_t hx, hy;
uint32_t sx, lx, ly; uint64_t sx, qs;
int cquo, qs; int cquo;
EXTRACT_WORDS (hx, lx, x); EXTRACT_WORDS64 (hx, x);
EXTRACT_WORDS (hy, ly, y); EXTRACT_WORDS64 (hy, y);
sx = hx & 0x80000000; sx = hx & UINT64_C(0x8000000000000000);
qs = sx ^ (hy & 0x80000000); qs = sx ^ (hy & UINT64_C(0x8000000000000000));
hy &= 0x7fffffff; hy &= UINT64_C(0x7fffffffffffffff);
hx &= 0x7fffffff; hx &= UINT64_C(0x7fffffffffffffff);
/* Purge off exception values. */ /* Purge off exception values. */
if ((hy | ly) == 0) if (__glibc_unlikely (hy == 0))
return (x * y) / (x * y); /* y = 0 */ return (x * y) / (x * y); /* y = 0 */
if ((hx >= 0x7ff00000) /* x not finite */ if (__builtin_expect (hx >= UINT64_C(0x7ff0000000000000) /* x not finite */
|| ((hy >= 0x7ff00000) /* p is NaN */ || hy > UINT64_C(0x7ff0000000000000), 0))/* y is NaN */
&& (((hy - 0x7ff00000) | ly) != 0)))
return (x * y) / (x * y); return (x * y) / (x * y);
if (hy <= 0x7fbfffff) if (hy <= UINT64_C(0x7fbfffffffffffff))
x = __ieee754_fmod (x, 8 * y); /* now x < 8y */ x = __ieee754_fmod (x, 8 * y); /* now x < 8y */
if (((hx - hy) | (lx - ly)) == 0) if (__glibc_unlikely (hx == hy))
{ {
*quo = qs ? -1 : 1; *quo = qs ? -1 : 1;
return zero * x; return zero * x;
} }
x = fabs (x); x = fabs (x);
y = fabs (y); INSERT_WORDS64 (y, hy);
cquo = 0; cquo = 0;
if (hy <= 0x7fcfffff && x >= 4 * y) if (hy <= UINT64_C(0x7fcfffffffffffff) && x >= 4 * y)
{ {
x -= 4 * y; x -= 4 * y;
cquo += 4; cquo += 4;
} }
if (hy <= 0x7fdfffff && x >= 2 * y) if (hy <= UINT64_C(0x7fdfffffffffffff) && x >= 2 * y)
{ {
x -= 2 * y; x -= 2 * y;
cquo += 2; cquo += 2;
} }
if (hy < 0x00200000) if (hy < UINT64_C(0x0020000000000000))
{ {
if (x + x > y) if (x + x > y)
{ {

79
sysdeps/ieee754/dbl-64/s_roundeven.c
View File

@ -1,5 +1,4 @@
/* Round to nearest integer value, rounding halfway cases to even. /* Round to nearest integer value, rounding halfway cases to even.
dbl-64 version.
Copyright (C) 2016-2021 Free Software Foundation, Inc. Copyright (C) 2016-2021 Free Software Foundation, Inc.
This file is part of the GNU C Library. This file is part of the GNU C Library.
@ -29,10 +28,10 @@
double double
__roundeven (double x) __roundeven (double x)
{ {
uint32_t hx, lx, uhx; uint64_t ix, ux;
EXTRACT_WORDS (hx, lx, x); EXTRACT_WORDS64 (ix, x);
uhx = hx & 0x7fffffff; ux = ix & 0x7fffffffffffffffULL;
int exponent = uhx >> (MANT_DIG - 1 - 32); int exponent = ux >> (MANT_DIG - 1);
if (exponent >= BIAS + MANT_DIG - 1) if (exponent >= BIAS + MANT_DIG - 1)
{ {
/* Integer, infinity or NaN. */ /* Integer, infinity or NaN. */
@ -42,63 +41,29 @@ __roundeven (double x)
else else
return x; return x;
} }
else if (exponent >= BIAS + MANT_DIG - 32)
{
/* Not necessarily an integer; integer bit is in low word.
Locate the bits with exponents 0 and -1. */
int int_pos = (BIAS + MANT_DIG - 1) - exponent;
int half_pos = int_pos - 1;
uint32_t half_bit = 1U << half_pos;
uint32_t int_bit = 1U << int_pos;
if ((lx & (int_bit | (half_bit - 1))) != 0)
{
/* Carry into the exponent works correctly. No need to test
whether HALF_BIT is set. */
lx += half_bit;
hx += lx < half_bit;
}
lx &= ~(int_bit - 1);
}
else if (exponent == BIAS + MANT_DIG - 33)
{
/* Not necessarily an integer; integer bit is bottom of high
word, half bit is top of low word. */
if (((hx & 1) | (lx & 0x7fffffff)) != 0)
{
lx += 0x80000000;
hx += lx < 0x80000000;
}
lx = 0;
}
else if (exponent >= BIAS) else if (exponent >= BIAS)
{ {
/* At least 1; not necessarily an integer, integer bit and half /* At least 1; not necessarily an integer. Locate the bits with
bit are in the high word. Locate the bits with exponents 0 exponents 0 and -1 (when the unbiased exponent is 0, the bit
and -1 (when the unbiased exponent is 0, the bit with with exponent 0 is implicit, but as the bias is odd it is OK
exponent 0 is implicit, but as the bias is odd it is OK to to take it from the low bit of the exponent). */
take it from the low bit of the exponent). */ int int_pos = (BIAS + MANT_DIG - 1) - exponent;
int int_pos = (BIAS + MANT_DIG - 33) - exponent;
int half_pos = int_pos - 1; int half_pos = int_pos - 1;
uint32_t half_bit = 1U << half_pos; uint64_t half_bit = 1ULL << half_pos;
uint32_t int_bit = 1U << int_pos; uint64_t int_bit = 1ULL << int_pos;
if (((hx & (int_bit | (half_bit - 1))) | lx) != 0) if ((ix & (int_bit | (half_bit - 1))) != 0)
hx += half_bit; /* Carry into the exponent works correctly. No need to test
hx &= ~(int_bit - 1); whether HALF_BIT is set. */
lx = 0; ix += half_bit;
} ix &= ~(int_bit - 1);
else if (exponent == BIAS - 1 && (uhx > 0x3fe00000 || lx != 0))
{
/* Interval (0.5, 1). */
hx = (hx & 0x80000000) | 0x3ff00000;
lx = 0;
} }
else if (exponent == BIAS - 1 && ux > 0x3fe0000000000000ULL)
/* Interval (0.5, 1). */
ix = (ix & 0x8000000000000000ULL) | 0x3ff0000000000000ULL;
else else
{ /* Rounds to 0. */
/* Rounds to 0. */ ix &= 0x8000000000000000ULL;
hx &= 0x80000000; INSERT_WORDS64 (x, ix);
lx = 0;
}
INSERT_WORDS (x, hx, lx);
return x; return x;
} }
hidden_def (__roundeven) hidden_def (__roundeven)

63
sysdeps/ieee754/dbl-64/s_scalbln.c
View File

@ -20,43 +20,40 @@
#include <math_private.h> #include <math_private.h>
static const double static const double
two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */ two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */ twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */
huge = 1.0e+300, huge = 1.0e+300,
tiny = 1.0e-300; tiny = 1.0e-300;
double double
__scalbln (double x, long int n) __scalbln (double x, long int n)
{ {
int32_t k, hx, lx; int64_t ix;
EXTRACT_WORDS (hx, lx, x); int64_t k;
k = (hx & 0x7ff00000) >> 20; /* extract exponent */ EXTRACT_WORDS64(ix,x);
if (__glibc_unlikely (k == 0)) /* 0 or subnormal x */ k = (ix >> 52) & 0x7ff; /* extract exponent */
{ if (__builtin_expect(k==0, 0)) { /* 0 or subnormal x */
if ((lx | (hx & 0x7fffffff)) == 0) if ((ix & UINT64_C(0xfffffffffffff))==0) return x; /* +-0 */
return x; /* +-0 */ x *= two54;
x *= two54; EXTRACT_WORDS64(ix,x);
GET_HIGH_WORD (hx, x); k = ((ix >> 52) & 0x7ff) - 54;
k = ((hx & 0x7ff00000) >> 20) - 54; }
} if (__builtin_expect(k==0x7ff, 0)) return x+x; /* NaN or Inf */
if (__glibc_unlikely (k == 0x7ff)) if (__builtin_expect(n< -50000, 0))
return x + x; /* NaN or Inf */ return tiny*copysign(tiny,x); /*underflow*/
if (__glibc_unlikely (n < -50000)) if (__builtin_expect(n> 50000 || k+n > 0x7fe, 0))
return tiny * copysign (tiny, x); /*underflow*/ return huge*copysign(huge,x); /* overflow */
if (__glibc_unlikely (n > 50000 || k + n > 0x7fe)) /* Now k and n are bounded we know that k = k+n does not
return huge * copysign (huge, x); /* overflow */ overflow. */
/* Now k and n are bounded we know that k = k+n does not k = k+n;
overflow. */ if (__builtin_expect(k > 0, 1)) /* normal result */
k = k + n; {INSERT_WORDS64(x,(ix&UINT64_C(0x800fffffffffffff))|(k<<52));
if (__glibc_likely (k > 0)) /* normal result */ return x;}
{ if (k <= -54)
SET_HIGH_WORD (x, (hx & 0x800fffff) | (k << 20)); return x; return tiny*copysign(tiny,x); /*underflow*/
} k += 54; /* subnormal result */
if (k <= -54) INSERT_WORDS64(x,(ix&INT64_C(0x800fffffffffffff))|(k<<52));
return tiny * copysign (tiny, x); /*underflow*/ return x*twom54;
k += 54; /* subnormal result */
SET_HIGH_WORD (x, (hx & 0x800fffff) | (k << 20));
return x * twom54;
} }
#ifdef NO_LONG_DOUBLE #ifdef NO_LONG_DOUBLE
strong_alias (__scalbln, __scalblnl) strong_alias (__scalbln, __scalblnl)

63
sysdeps/ieee754/dbl-64/s_scalbn.c
View File

@ -20,43 +20,40 @@
#include <math_private.h> #include <math_private.h>
static const double static const double
two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */ two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */ twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */
huge = 1.0e+300, huge = 1.0e+300,
tiny = 1.0e-300; tiny = 1.0e-300;
double double
__scalbn (double x, int n) __scalbn (double x, int n)
{ {
int32_t k, hx, lx; int64_t ix;
EXTRACT_WORDS (hx, lx, x); int64_t k;
k = (hx & 0x7ff00000) >> 20; /* extract exponent */ EXTRACT_WORDS64(ix,x);
if (__glibc_unlikely (k == 0)) /* 0 or subnormal x */ k = (ix >> 52) & 0x7ff; /* extract exponent */
{ if (__builtin_expect(k==0, 0)) { /* 0 or subnormal x */
if ((lx | (hx & 0x7fffffff)) == 0) if ((ix & UINT64_C(0xfffffffffffff))==0) return x; /* +-0 */
return x; /* +-0 */ x *= two54;
x *= two54; EXTRACT_WORDS64(ix,x);
GET_HIGH_WORD (hx, x); k = ((ix >> 52) & 0x7ff) - 54;
k = ((hx & 0x7ff00000) >> 20) - 54; }
} if (__builtin_expect(k==0x7ff, 0)) return x+x; /* NaN or Inf */
if (__glibc_unlikely (k == 0x7ff)) if (__builtin_expect(n< -50000, 0))
return x + x; /* NaN or Inf */ return tiny*copysign(tiny,x); /*underflow*/
if (__glibc_unlikely (n < -50000)) if (__builtin_expect(n> 50000 || k+n > 0x7fe, 0))
return tiny * copysign (tiny, x); /*underflow*/ return huge*copysign(huge,x); /* overflow */
if (__glibc_unlikely (n > 50000 || k + n > 0x7fe)) /* Now k and n are bounded we know that k = k+n does not
return huge * copysign (huge, x); /* overflow */ overflow. */
/* Now k and n are bounded we know that k = k+n does not k = k+n;
overflow. */ if (__builtin_expect(k > 0, 1)) /* normal result */
k = k + n; {INSERT_WORDS64(x,(ix&UINT64_C(0x800fffffffffffff))|(k<<52));
if (__glibc_likely (k > 0)) /* normal result */ return x;}
{ if (k <= -54)
SET_HIGH_WORD (x, (hx & 0x800fffff) | (k << 20)); return x; return tiny*copysign(tiny,x); /*underflow*/
} k += 54; /* subnormal result */
if (k <= -54) INSERT_WORDS64(x,(ix&INT64_C(0x800fffffffffffff))|(k<<52));
return tiny * copysign (tiny, x); /*underflow*/ return x*twom54;
k += 54; /* subnormal result */
SET_HIGH_WORD (x, (hx & 0x800fffff) | (k << 20));
return x * twom54;
} }
#ifdef NO_LONG_DOUBLE #ifdef NO_LONG_DOUBLE
strong_alias (__scalbn, __scalbnl) strong_alias (__scalbn, __scalbnl)

42
sysdeps/ieee754/dbl-64/s_setpayload_main.c
View File

@ -1,4 +1,4 @@
/* Set NaN payload. dbl-64 version. /* Set NaN payload.
Copyright (C) 2016-2021 Free Software Foundation, Inc. Copyright (C) 2016-2021 Free Software Foundation, Inc.
This file is part of the GNU C Library. This file is part of the GNU C Library.
@ -30,41 +30,25 @@
int int
FUNC (double *x, double payload) FUNC (double *x, double payload)
{ {
uint32_t hx, lx; uint64_t ix;
EXTRACT_WORDS (hx, lx, payload); EXTRACT_WORDS64 (ix, payload);
int exponent = hx >> (EXPLICIT_MANT_DIG - 32); int exponent = ix >> EXPLICIT_MANT_DIG;
/* Test if argument is (a) negative or too large; (b) too small, /* Test if argument is (a) negative or too large; (b) too small,
except for 0 when allowed; (c) not an integer. */ except for 0 when allowed; (c) not an integer. */
if (exponent >= BIAS + PAYLOAD_DIG if (exponent >= BIAS + PAYLOAD_DIG
|| (exponent < BIAS && !(SET_HIGH_BIT && hx == 0 && lx == 0))) || (exponent < BIAS && !(SET_HIGH_BIT && ix == 0))
|| (ix & ((1ULL << (BIAS + EXPLICIT_MANT_DIG - exponent)) - 1)) != 0)
{ {
INSERT_WORDS (*x, 0, 0); INSERT_WORDS64 (*x, 0);
return 1; return 1;
} }
int shift = BIAS + EXPLICIT_MANT_DIG - exponent; if (ix != 0)
if (shift < 32
? (lx & ((1U << shift) - 1)) != 0
: (lx != 0 || (hx & ((1U << (shift - 32)) - 1)) != 0))
{ {
INSERT_WORDS (*x, 0, 0); ix &= (1ULL << EXPLICIT_MANT_DIG) - 1;
return 1; ix |= 1ULL << EXPLICIT_MANT_DIG;
ix >>= BIAS + EXPLICIT_MANT_DIG - exponent;
} }
if (exponent != 0) ix |= 0x7ff0000000000000ULL | (SET_HIGH_BIT ? 0x8000000000000ULL : 0);
{ INSERT_WORDS64 (*x, ix);
hx &= (1U << (EXPLICIT_MANT_DIG - 32)) - 1;
hx |= 1U << (EXPLICIT_MANT_DIG - 32);
if (shift >= 32)
{
lx = hx >> (shift - 32);
hx = 0;
}
else if (shift != 0)
{
lx = (lx >> shift) | (hx << (32 - shift));
hx >>= shift;
}
}
hx |= 0x7ff00000 | (SET_HIGH_BIT ? 0x80000 : 0);
INSERT_WORDS (*x, hx, lx);
return 0; return 0;
} }

32
sysdeps/ieee754/dbl-64/s_totalorder.c
View File

@ -1,4 +1,4 @@
/* Total order operation. dbl-64 version. /* Total order operation.
Copyright (C) 2016-2021 Free Software Foundation, Inc. Copyright (C) 2016-2021 Free Software Foundation, Inc.
This file is part of the GNU C Library. This file is part of the GNU C Library.
@ -18,8 +18,8 @@
#include <math.h> #include <math.h>
#include <math_private.h> #include <math_private.h>
#include <libm-alias-double.h>
#include <nan-high-order-bit.h> #include <nan-high-order-bit.h>
#include <libm-alias-double.h>
#include <stdint.h> #include <stdint.h>
#include <shlib-compat.h> #include <shlib-compat.h>
#include <first-versions.h> #include <first-versions.h>
@ -27,30 +27,26 @@
int int
__totalorder (const double *x, const double *y) __totalorder (const double *x, const double *y)
{ {
int32_t hx, hy; int64_t ix, iy;
uint32_t lx, ly; EXTRACT_WORDS64 (ix, *x);
EXTRACT_WORDS (hx, lx, *x); EXTRACT_WORDS64 (iy, *y);
EXTRACT_WORDS (hy, ly, *y);
#if HIGH_ORDER_BIT_IS_SET_FOR_SNAN #if HIGH_ORDER_BIT_IS_SET_FOR_SNAN
uint32_t uhx = hx & 0x7fffffff, uhy = hy & 0x7fffffff;
/* For the preferred quiet NaN convention, this operation is a /* For the preferred quiet NaN convention, this operation is a
comparison of the representations of the arguments interpreted as comparison of the representations of the arguments interpreted as
sign-magnitude integers. If both arguments are NaNs, invert the sign-magnitude integers. If both arguments are NaNs, invert the
quiet/signaling bit so comparing that way works. */ quiet/signaling bit so comparing that way works. */
if ((uhx > 0x7ff00000 || (uhx == 0x7ff00000 && lx != 0)) if ((ix & 0x7fffffffffffffffULL) > 0x7ff0000000000000ULL
&& (uhy > 0x7ff00000 || (uhy == 0x7ff00000 && ly != 0))) && (iy & 0x7fffffffffffffffULL) > 0x7ff0000000000000ULL)
{ {
hx ^= 0x00080000; ix ^= 0x0008000000000000ULL;
hy ^= 0x00080000; iy ^= 0x0008000000000000ULL;
} }
#endif #endif
uint32_t hx_sign = hx >> 31; uint64_t ix_sign = ix >> 63;
uint32_t hy_sign = hy >> 31; uint64_t iy_sign = iy >> 63;
hx ^= hx_sign >> 1; ix ^= ix_sign >> 1;
lx ^= hx_sign; iy ^= iy_sign >> 1;
hy ^= hy_sign >> 1; return ix <= iy;
ly ^= hy_sign;
return hx < hy || (hx == hy && lx <= ly);
} }
#ifdef SHARED #ifdef SHARED
# define CONCATX(x, y) x ## y # define CONCATX(x, y) x ## y

24
sysdeps/ieee754/dbl-64/s_totalordermag.c
View File

@ -1,4 +1,4 @@
/* Total order operation on absolute values. dbl-64 version. /* Total order operation on absolute values.
Copyright (C) 2016-2021 Free Software Foundation, Inc. Copyright (C) 2016-2021 Free Software Foundation, Inc.
This file is part of the GNU C Library. This file is part of the GNU C Library.
@ -18,8 +18,8 @@
#include <math.h> #include <math.h>
#include <math_private.h> #include <math_private.h>
#include <libm-alias-double.h>
#include <nan-high-order-bit.h> #include <nan-high-order-bit.h>
#include <libm-alias-double.h>
#include <stdint.h> #include <stdint.h>
#include <shlib-compat.h> #include <shlib-compat.h>
#include <first-versions.h> #include <first-versions.h>
@ -27,25 +27,23 @@
int int
__totalordermag (const double *x, const double *y) __totalordermag (const double *x, const double *y)
{ {
uint32_t hx, hy; uint64_t ix, iy;
uint32_t lx, ly; EXTRACT_WORDS64 (ix, *x);
EXTRACT_WORDS (hx, lx, *x); EXTRACT_WORDS64 (iy, *y);
EXTRACT_WORDS (hy, ly, *y); ix &= 0x7fffffffffffffffULL;
hx &= 0x7fffffff; iy &= 0x7fffffffffffffffULL;
hy &= 0x7fffffff;
#if HIGH_ORDER_BIT_IS_SET_FOR_SNAN #if HIGH_ORDER_BIT_IS_SET_FOR_SNAN
/* For the preferred quiet NaN convention, this operation is a /* For the preferred quiet NaN convention, this operation is a
comparison of the representations of the absolute values of the comparison of the representations of the absolute values of the
arguments. If both arguments are NaNs, invert the arguments. If both arguments are NaNs, invert the
quiet/signaling bit so comparing that way works. */ quiet/signaling bit so comparing that way works. */
if ((hx > 0x7ff00000 || (hx == 0x7ff00000 && lx != 0)) if (ix > 0x7ff0000000000000ULL && iy > 0x7ff0000000000000ULL)
&& (hy > 0x7ff00000 || (hy == 0x7ff00000 && ly != 0)))
{ {
hx ^= 0x00080000; ix ^= 0x0008000000000000ULL;
hy ^= 0x00080000; iy ^= 0x0008000000000000ULL;
} }
#endif #endif
return hx < hy || (hx == hy && lx <= ly); return ix <= iy;
} }
#ifdef SHARED #ifdef SHARED
# define CONCATX(x, y) x ## y # define CONCATX(x, y) x ## y

68
sysdeps/ieee754/dbl-64/wordsize-64/e_acosh.c
View File

@ -1,68 +0,0 @@
/* Optimized for 64-bit by Ulrich Drepper <drepper@gmail.com>, 2012 */
/*
* ====================================================
* 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.
* ====================================================
*/
/* __ieee754_acosh(x)
* Method :
* Based on
* acosh(x) = log [ x + sqrt(x*x-1) ]
* we have
* acosh(x) := log(x)+ln2, if x is large; else
* acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x>2; else
* acosh(x) := log1p(t+sqrt(2.0*t+t*t)); where t=x-1.
*
* Special cases:
* acosh(x) is NaN with signal if x<1.
* acosh(NaN) is NaN without signal.
*/
#include <math.h>
#include <math_private.h>
#include <libm-alias-finite.h>
static const double
one = 1.0,
ln2 = 6.93147180559945286227e-01; /* 0x3FE62E42, 0xFEFA39EF */
double
__ieee754_acosh (double x)
{
int64_t hx;
EXTRACT_WORDS64 (hx, x);
if (hx > INT64_C (0x4000000000000000))
{
if (__glibc_unlikely (hx >= INT64_C (0x41b0000000000000)))
{
/* x > 2**28 */
if (hx >= INT64_C (0x7ff0000000000000))
/* x is inf of NaN */
return x + x;
else
return __ieee754_log (x) + ln2;/* acosh(huge)=log(2x) */
}
/* 2**28 > x > 2 */
double t = x * x;
return __ieee754_log (2.0 * x - one / (x + sqrt (t - one)));
}
else if (__glibc_likely (hx > INT64_C (0x3ff0000000000000)))
{
/* 1<x<2 */
double t = x - one;
return __log1p (t + sqrt (2.0 * t + t * t));
}
else if (__glibc_likely (hx == INT64_C (0x3ff0000000000000)))
return 0.0; /* acosh(1) = 0 */
else /* x < 1 */
return (x - x) / (x - x);
}
libm_alias_finite (__ieee754_acosh, __acosh)

85
sysdeps/ieee754/dbl-64/wordsize-64/e_cosh.c
View File

@ -1,85 +0,0 @@
/* Optimized by Ulrich Drepper <drepper@gmail.com>, 2011 */
/*
* ====================================================
* 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.
* ====================================================
*/
/* __ieee754_cosh(x)
* Method :
* mathematically cosh(x) if defined to be (exp(x)+exp(-x))/2
* 1. Replace x by |x| (cosh(x) = cosh(-x)).
* 2.
* [ exp(x) - 1 ]^2
* 0 <= x <= ln2/2 : cosh(x) := 1 + -------------------
* 2*exp(x)
*
* exp(x) + 1/exp(x)
* ln2/2 <= x <= 22 : cosh(x) := -------------------
* 2
* 22 <= x <= lnovft : cosh(x) := exp(x)/2
* lnovft <= x <= ln2ovft: cosh(x) := exp(x/2)/2 * exp(x/2)
* ln2ovft < x : cosh(x) := huge*huge (overflow)
*
* Special cases:
* cosh(x) is |x| if x is +INF, -INF, or NaN.
* only cosh(0)=1 is exact for finite x.
*/
#include <math.h>
#include <math_private.h>
#include <libm-alias-finite.h>
static const double one = 1.0, half=0.5, huge = 1.0e300;
double
__ieee754_cosh (double x)
{
double t,w;
int32_t ix;
/* High word of |x|. */
GET_HIGH_WORD(ix,x);
ix &= 0x7fffffff;
/* |x| in [0,22] */
if (ix < 0x40360000) {
/* |x| in [0,0.5*ln2], return 1+expm1(|x|)^2/(2*exp(|x|)) */
if(ix<0x3fd62e43) {
if (ix<0x3c800000) /* cosh(tiny) = 1 */
return one;
t = __expm1(fabs(x));
w = one+t;
return one+(t*t)/(w+w);
}
/* |x| in [0.5*ln2,22], return (exp(|x|)+1/exp(|x|)/2; */
t = __ieee754_exp(fabs(x));
return half*t+half/t;
}
/* |x| in [22, log(maxdouble)] return half*exp(|x|) */
if (ix < 0x40862e42) return half*__ieee754_exp(fabs(x));
/* |x| in [log(maxdouble), overflowthresold] */
int64_t fix;
EXTRACT_WORDS64(fix, x);
fix &= UINT64_C(0x7fffffffffffffff);
if (fix <= UINT64_C(0x408633ce8fb9f87d)) {
w = __ieee754_exp(half*fabs(x));
t = half*w;
return t*w;
}
/* x is INF or NaN */
if(ix>=0x7ff00000) return x*x;
/* |x| > overflowthresold, cosh(x) overflow */
return huge*huge;
}
libm_alias_finite (__ieee754_cosh, __cosh)

106
sysdeps/ieee754/dbl-64/wordsize-64/e_fmod.c
View File

@ -1,106 +0,0 @@
/* Rewritten for 64-bit machines by Ulrich Drepper <drepper@gmail.com>. */
/*
* ====================================================
* 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.
* ====================================================
*/
/*
* __ieee754_fmod(x,y)
* Return x mod y in exact arithmetic
* Method: shift and subtract
*/
#include <math.h>
#include <math_private.h>
#include <stdint.h>
#include <libm-alias-finite.h>
static const double one = 1.0, Zero[] = {0.0, -0.0,};
double
__ieee754_fmod (double x, double y)
{
int32_t n,ix,iy;
int64_t hx,hy,hz,sx,i;
EXTRACT_WORDS64(hx,x);
EXTRACT_WORDS64(hy,y);
sx = hx&UINT64_C(0x8000000000000000); /* sign of x */
hx ^=sx; /* |x| */
hy &= UINT64_C(0x7fffffffffffffff); /* |y| */
/* purge off exception values */
if(__builtin_expect(hy==0
|| hx >= UINT64_C(0x7ff0000000000000)
|| hy > UINT64_C(0x7ff0000000000000), 0))
/* y=0,or x not finite or y is NaN */
return (x*y)/(x*y);
if(__builtin_expect(hx<=hy, 0)) {
if(hx<hy) return x; /* |x|<|y| return x */
return Zero[(uint64_t)sx>>63]; /* |x|=|y| return x*0*/
}
/* determine ix = ilogb(x) */
if(__builtin_expect(hx<UINT64_C(0x0010000000000000), 0)) {
/* subnormal x */
for (ix = -1022,i=(hx<<11); i>0; i<<=1) ix -=1;
} else ix = (hx>>52)-1023;
/* determine iy = ilogb(y) */
if(__builtin_expect(hy<UINT64_C(0x0010000000000000), 0)) { /* subnormal y */
for (iy = -1022,i=(hy<<11); i>0; i<<=1) iy -=1;
} else iy = (hy>>52)-1023;
/* set up hx, hy and align y to x */
if(__builtin_expect(ix >= -1022, 1))
hx = UINT64_C(0x0010000000000000)|(UINT64_C(0x000fffffffffffff)&hx);
else { /* subnormal x, shift x to normal */
n = -1022-ix;
hx<<=n;
}
if(__builtin_expect(iy >= -1022, 1))
hy = UINT64_C(0x0010000000000000)|(UINT64_C(0x000fffffffffffff)&hy);
else { /* subnormal y, shift y to normal */
n = -1022-iy;
hy<<=n;
}
/* fix point fmod */
n = ix - iy;
while(n--) {
hz=hx-hy;
if(hz<0){hx = hx+hx;}
else {
if(hz==0) /* return sign(x)*0 */
return Zero[(uint64_t)sx>>63];
hx = hz+hz;
}
}
hz=hx-hy;
if(hz>=0) {hx=hz;}
/* convert back to floating value and restore the sign */
if(hx==0) /* return sign(x)*0 */
return Zero[(uint64_t)sx>>63];
while(hx<UINT64_C(0x0010000000000000)) { /* normalize x */
hx = hx+hx;
iy -= 1;
}
if(__builtin_expect(iy>= -1022, 1)) { /* normalize output */
hx = ((hx-UINT64_C(0x0010000000000000))|((uint64_t)(iy+1023)<<52));
INSERT_WORDS64(x,hx|sx);
} else { /* subnormal output */
n = -1022 - iy;
hx>>=n;
INSERT_WORDS64(x,hx|sx);
x *= one; /* create necessary signal */
}
return x; /* exact output */
}
libm_alias_finite (__ieee754_fmod, __fmod)

90
sysdeps/ieee754/dbl-64/wordsize-64/e_log10.c
View File

@ -1,90 +0,0 @@
/* @(#)e_log10.c 5.1 93/09/24 */
/*
* ====================================================
* 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.
* ====================================================
*/
/* __ieee754_log10(x)
* Return the base 10 logarithm of x
*
* Method :
* Let log10_2hi = leading 40 bits of log10(2) and
* log10_2lo = log10(2) - log10_2hi,
* ivln10 = 1/log(10) rounded.
* Then
* n = ilogb(x),
* if(n<0) n = n+1;
* x = scalbn(x,-n);
* log10(x) := n*log10_2hi + (n*log10_2lo + ivln10*log(x))
*
* Note 1:
* To guarantee log10(10**n)=n, where 10**n is normal, the rounding
* mode must set to Round-to-Nearest.
* Note 2:
* [1/log(10)] rounded to 53 bits has error .198 ulps;
* log10 is monotonic at all binary break points.
*
* Special cases:
* log10(x) is NaN with signal if x < 0;
* log10(+INF) is +INF with no signal; log10(0) is -INF with signal;
* log10(NaN) is that NaN with no signal;
* log10(10**N) = N for N=0,1,...,22.
*
* Constants:
* The hexadecimal values are the intended ones for the following constants.
* The decimal values may be used, provided that the compiler will convert
* from decimal to binary accurately enough to produce the hexadecimal values
* shown.
*/
#include <math.h>
#include <fix-int-fp-convert-zero.h>
#include <math_private.h>
#include <stdint.h>
#include <libm-alias-finite.h>
static const double two54 = 1.80143985094819840000e+16; /* 0x4350000000000000 */
static const double ivln10 = 4.34294481903251816668e-01; /* 0x3FDBCB7B1526E50E */
static const double log10_2hi = 3.01029995663611771306e-01; /* 0x3FD34413509F6000 */
static const double log10_2lo = 3.69423907715893078616e-13; /* 0x3D59FEF311F12B36 */
double
__ieee754_log10 (double x)
{
double y, z;
int64_t i, hx;
int32_t k;
EXTRACT_WORDS64 (hx, x);
k = 0;
if (hx < INT64_C(0x0010000000000000))
{ /* x < 2**-1022 */
if (__glibc_unlikely ((hx & UINT64_C(0x7fffffffffffffff)) == 0))
return -two54 / fabs (x); /* log(+-0)=-inf */
if (__glibc_unlikely (hx < 0))
return (x - x) / (x - x); /* log(-#) = NaN */
k -= 54;
x *= two54; /* subnormal number, scale up x */
EXTRACT_WORDS64 (hx, x);
}
/* scale up resulted in a NaN number */
if (__glibc_unlikely (hx >= UINT64_C(0x7ff0000000000000)))
return x + x;
k += (hx >> 52) - 1023;
i = ((uint64_t) k & UINT64_C(0x8000000000000000)) >> 63;
hx = (hx & UINT64_C(0x000fffffffffffff)) | ((0x3ff - i) << 52);
y = (double) (k + i);
if (FIX_INT_FP_CONVERT_ZERO && y == 0.0)
y = 0.0;
INSERT_WORDS64 (x, hx);
z = y * log10_2lo + ivln10 * __ieee754_log (x);
return z + y * log10_2hi;
}
libm_alias_finite (__ieee754_log10, __log10)

66
sysdeps/ieee754/dbl-64/wordsize-64/s_frexp.c
View File

@ -1,66 +0,0 @@
/* Copyright (C) 2011-2021 Free Software Foundation, Inc.
This file is part of the GNU C Library.
Contributed by Ulrich Drepper <drepper@gmail.com>, 2011.
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
<https://www.gnu.org/licenses/>. */
#include <inttypes.h>
#include <math.h>
#include <math_private.h>
#include <libm-alias-double.h>
/*
* for non-zero, finite x
* x = frexp(arg,&exp);
* return a double fp quantity x such that 0.5 <= |x| <1.0
* and the corresponding binary exponent "exp". That is
* arg = x*2^exp.
* If arg is inf, 0.0, or NaN, then frexp(arg,&exp) returns arg
* with *exp=0.
*/
double
__frexp (double x, int *eptr)
{
int64_t ix;
EXTRACT_WORDS64 (ix, x);
int32_t ex = 0x7ff & (ix >> 52);
int e = 0;
if (__glibc_likely (ex != 0x7ff && x != 0.0))
{
/* Not zero and finite. */
e = ex - 1022;
if (__glibc_unlikely (ex == 0))
{
/* Subnormal. */
x *= 0x1p54;
EXTRACT_WORDS64 (ix, x);
ex = 0x7ff & (ix >> 52);
e = ex - 1022 - 54;
}
ix = (ix & INT64_C (0x800fffffffffffff)) | INT64_C (0x3fe0000000000000);
INSERT_WORDS64 (x, ix);
}
else
/* Quiet signaling NaNs. */
x += x;
*eptr = e;
return x;
}
libm_alias_double (__frexp, frexp)

38
sysdeps/ieee754/dbl-64/wordsize-64/s_getpayload.c
View File

@ -1,38 +0,0 @@
/* Get NaN payload. dbl-64/wordsize-64 version.
Copyright (C) 2016-2021 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
<https://www.gnu.org/licenses/>. */
#include <math.h>
#include <math_private.h>
#include <libm-alias-double.h>
#include <stdint.h>
#include <fix-int-fp-convert-zero.h>
double
__getpayload (const double *x)
{
uint64_t ix;
EXTRACT_WORDS64 (ix, *x);
if ((ix & 0x7ff0000000000000ULL) != 0x7ff0000000000000ULL
|| (ix & 0xfffffffffffffULL) == 0)
return -1;
ix &= 0x7ffffffffffffULL;
if (FIX_INT_FP_CONVERT_ZERO && ix == 0)
return 0.0f;
return (double) ix;
}
libm_alias_double (__getpayload, getpayload)

43
sysdeps/ieee754/dbl-64/wordsize-64/s_issignaling.c
View File

@ -1,43 +0,0 @@
/* Test for signaling NaN.
Copyright (C) 2013-2021 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
<https://www.gnu.org/licenses/>. */
#include <math.h>
#include <math_private.h>
#include <nan-high-order-bit.h>
int
__issignaling (double x)
{
uint64_t xi;
EXTRACT_WORDS64 (xi, x);
#if HIGH_ORDER_BIT_IS_SET_FOR_SNAN
/* We only have to care about the high-order bit of x's significand, because
having it set (sNaN) already makes the significand different from that
used to designate infinity. */
return (xi & UINT64_C (0x7ff8000000000000)) == UINT64_C (0x7ff8000000000000);
#else
/* To keep the following comparison simple, toggle the quiet/signaling bit,
so that it is set for sNaNs. This is inverse to IEEE 754-2008 (as well as
common practice for IEEE 754-1985). */
xi ^= UINT64_C (0x0008000000000000);
/* We have to compare for greater (instead of greater or equal), because x's
significand being all-zero designates infinity not NaN. */
return (xi & UINT64_C (0x7fffffffffffffff)) > UINT64_C (0x7ff8000000000000);
#endif
}
libm_hidden_def (__issignaling)

85
sysdeps/ieee754/dbl-64/wordsize-64/s_llround.c
View File

@ -1,85 +0,0 @@
/* Round double value to long long int.
Copyright (C) 1997-2021 Free Software Foundation, Inc.
This file is part of the GNU C Library.
Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
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
<https://www.gnu.org/licenses/>. */
#define lround __hidden_lround
#define __lround __hidden___lround
#include <fenv.h>
#include <limits.h>
#include <math.h>
#include <sysdep.h>
#include <math_private.h>
#include <libm-alias-double.h>
#include <fix-fp-int-convert-overflow.h>
long long int
__llround (double x)
{
int32_t j0;
int64_t i0;
long long int result;
int sign;
EXTRACT_WORDS64 (i0, x);
j0 = ((i0 >> 52) & 0x7ff) - 0x3ff;
sign = i0 < 0 ? -1 : 1;
i0 &= UINT64_C(0xfffffffffffff);
i0 |= UINT64_C(0x10000000000000);
if (j0 < (int32_t) (8 * sizeof (long long int)) - 1)
{
if (j0 < 0)
return j0 < -1 ? 0 : sign;
else if (j0 >= 52)
result = i0 << (j0 - 52);
else
{
i0 += UINT64_C(0x8000000000000) >> j0;
result = i0 >> (52 - j0);
}
}
else
{
#ifdef FE_INVALID
/* The number is too large. Unless it rounds to LLONG_MIN,
FE_INVALID must be raised and the return value is
unspecified. */
if (FIX_DBL_LLONG_CONVERT_OVERFLOW && x != (double) LLONG_MIN)
{
feraiseexcept (FE_INVALID);
return sign == 1 ? LLONG_MAX : LLONG_MIN;
}
#endif
return (long long int) x;
}
return sign * result;
}
libm_alias_double (__llround, llround)
/* long has the same width as long long on LP64 machines, so use an alias. */
#undef lround
#undef __lround
#ifdef _LP64
strong_alias (__llround, __lround)
libm_alias_double (__lround, lround)
#endif

97
sysdeps/ieee754/dbl-64/wordsize-64/s_lround.c
View File

@ -1,97 +0,0 @@
/* Round double value to long int.
Copyright (C) 1997-2021 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
<https://www.gnu.org/licenses/>. */
#include <fenv.h>
#include <limits.h>
#include <math.h>
#include <math_private.h>
#include <libm-alias-double.h>
#include <fix-fp-int-convert-overflow.h>
/* For LP64, lround is an alias for llround. */
#ifndef _LP64
long int
__lround (double x)
{
int32_t j0;
int64_t i0;
long int result;
int sign;
EXTRACT_WORDS64 (i0, x);
j0 = ((i0 >> 52) & 0x7ff) - 0x3ff;
sign = i0 < 0 ? -1 : 1;
i0 &= UINT64_C(0xfffffffffffff);
i0 |= UINT64_C(0x10000000000000);
if (j0 < (int32_t) (8 * sizeof (long int)) - 1)
{
if (j0 < 0)
return j0 < -1 ? 0 : sign;
else if (j0 >= 52)
result = i0 << (j0 - 52);
else
{
i0 += UINT64_C(0x8000000000000) >> j0;
result = i0 >> (52 - j0);
#ifdef FE_INVALID
if (sizeof (long int) == 4
&& sign == 1
&& result == LONG_MIN)
/* Rounding brought the value out of range. */
feraiseexcept (FE_INVALID);
#endif
}
}
else
{
/* The number is too large. Unless it rounds to LONG_MIN,
FE_INVALID must be raised and the return value is
unspecified. */
#ifdef FE_INVALID
if (FIX_DBL_LONG_CONVERT_OVERFLOW
&& !(sign == -1
&& (sizeof (long int) == 4
? x > (double) LONG_MIN - 0.5
: x >= (double) LONG_MIN)))
{
feraiseexcept (FE_INVALID);
return sign == 1 ? LONG_MAX : LONG_MIN;
}
else if (!FIX_DBL_LONG_CONVERT_OVERFLOW
&& sizeof (long int) == 4
&& x <= (double) LONG_MIN - 0.5)
{
/* If truncation produces LONG_MIN, the cast will not raise
the exception, but may raise "inexact". */
feraiseexcept (FE_INVALID);
return LONG_MIN;
}
#endif
return (long int) x;
}
return sign * result;
}
libm_alias_double (__lround, lround)
#endif

65
sysdeps/ieee754/dbl-64/wordsize-64/s_modf.c
View File

@ -1,65 +0,0 @@
/* Rewritten for 64-bit machines by Ulrich Drepper <drepper@gmail.com>. */
/*
* ====================================================
* 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.
* ====================================================
*/
/*
* modf(double x, double *iptr)
* return fraction part of x, and return x's integral part in *iptr.
* Method:
* Bit twiddling.
*
* Exception:
* No exception.
*/
#include <math.h>
#include <math_private.h>
#include <libm-alias-double.h>
#include <stdint.h>
static const double one = 1.0;
double
__modf(double x, double *iptr)
{
int64_t i0;
int32_t j0;
EXTRACT_WORDS64(i0,x);
j0 = ((i0>>52)&0x7ff)-0x3ff; /* exponent of x */
if(j0<52) { /* integer part in x */
if(j0<0) { /* |x|<1 */
/* *iptr = +-0 */
INSERT_WORDS64(*iptr,i0&UINT64_C(0x8000000000000000));
return x;
} else {
uint64_t i = UINT64_C(0x000fffffffffffff)>>j0;
if((i0&i)==0) { /* x is integral */
*iptr = x;
/* return +-0 */
INSERT_WORDS64(x,i0&UINT64_C(0x8000000000000000));
return x;
} else {
INSERT_WORDS64(*iptr,i0&(~i));
return x - *iptr;
}
}
} else { /* no fraction part */
*iptr = x*one;
/* We must handle NaNs separately. */
if (j0 == 0x400 && (i0 & UINT64_C(0xfffffffffffff)))
return x*one;
INSERT_WORDS64(x,i0&UINT64_C(0x8000000000000000)); /* return +-0 */
return x;
}
}
#ifndef __modf
libm_alias_double (__modf, modf)
#endif

111
sysdeps/ieee754/dbl-64/wordsize-64/s_remquo.c
View File

@ -1,111 +0,0 @@
/* Compute remainder and a congruent to the quotient.
Copyright (C) 1997-2021 Free Software Foundation, Inc.
This file is part of the GNU C Library.
Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
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
<https://www.gnu.org/licenses/>. */
#include <math.h>
#include <math_private.h>
#include <libm-alias-double.h>
#include <stdint.h>
static const double zero = 0.0;
double
__remquo (double x, double y, int *quo)
{
int64_t hx, hy;
uint64_t sx, qs;
int cquo;
EXTRACT_WORDS64 (hx, x);
EXTRACT_WORDS64 (hy, y);
sx = hx & UINT64_C(0x8000000000000000);
qs = sx ^ (hy & UINT64_C(0x8000000000000000));
hy &= UINT64_C(0x7fffffffffffffff);
hx &= UINT64_C(0x7fffffffffffffff);
/* Purge off exception values. */
if (__glibc_unlikely (hy == 0))
return (x * y) / (x * y); /* y = 0 */
if (__builtin_expect (hx >= UINT64_C(0x7ff0000000000000) /* x not finite */
|| hy > UINT64_C(0x7ff0000000000000), 0))/* y is NaN */
return (x * y) / (x * y);
if (hy <= UINT64_C(0x7fbfffffffffffff))
x = __ieee754_fmod (x, 8 * y); /* now x < 8y */
if (__glibc_unlikely (hx == hy))
{
*quo = qs ? -1 : 1;
return zero * x;
}
x = fabs (x);
INSERT_WORDS64 (y, hy);
cquo = 0;
if (hy <= UINT64_C(0x7fcfffffffffffff) && x >= 4 * y)
{
x -= 4 * y;
cquo += 4;
}
if (hy <= UINT64_C(0x7fdfffffffffffff) && x >= 2 * y)
{
x -= 2 * y;
cquo += 2;
}
if (hy < UINT64_C(0x0020000000000000))
{
if (x + x > y)
{
x -= y;
++cquo;
if (x + x >= y)
{
x -= y;
++cquo;
}
}
}
else
{
double y_half = 0.5 * y;
if (x > y_half)
{
x -= y;
++cquo;
if (x >= y_half)
{
x -= y;
++cquo;
}
}
}
*quo = qs ? -cquo : cquo;
/* Ensure correct sign of zero result in round-downward mode. */
if (x == 0.0)
x = 0.0;
if (sx)
x = -x;
return x;
}
libm_alias_double (__remquo, remquo)

71
sysdeps/ieee754/dbl-64/wordsize-64/s_roundeven.c
View File

@ -1,71 +0,0 @@
/* Round to nearest integer value, rounding halfway cases to even.
dbl-64/wordsize-64 version.
Copyright (C) 2016-2021 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
<https://www.gnu.org/licenses/>. */
#include <math.h>
#include <math_private.h>
#include <libm-alias-double.h>
#include <stdint.h>
#define BIAS 0x3ff
#define MANT_DIG 53
#define MAX_EXP (2 * BIAS + 1)
double
__roundeven (double x)
{
uint64_t ix, ux;
EXTRACT_WORDS64 (ix, x);
ux = ix & 0x7fffffffffffffffULL;
int exponent = ux >> (MANT_DIG - 1);
if (exponent >= BIAS + MANT_DIG - 1)
{
/* Integer, infinity or NaN. */
if (exponent == MAX_EXP)
/* Infinity or NaN; quiet signaling NaNs. */
return x + x;
else
return x;
}
else if (exponent >= BIAS)
{
/* At least 1; not necessarily an integer. Locate the bits with
exponents 0 and -1 (when the unbiased exponent is 0, the bit
with exponent 0 is implicit, but as the bias is odd it is OK
to take it from the low bit of the exponent). */
int int_pos = (BIAS + MANT_DIG - 1) - exponent;
int half_pos = int_pos - 1;
uint64_t half_bit = 1ULL << half_pos;
uint64_t int_bit = 1ULL << int_pos;
if ((ix & (int_bit | (half_bit - 1))) != 0)
/* Carry into the exponent works correctly. No need to test
whether HALF_BIT is set. */
ix += half_bit;
ix &= ~(int_bit - 1);
}
else if (exponent == BIAS - 1 && ux > 0x3fe0000000000000ULL)
/* Interval (0.5, 1). */
ix = (ix & 0x8000000000000000ULL) | 0x3ff0000000000000ULL;
else
/* Rounds to 0. */
ix &= 0x8000000000000000ULL;
INSERT_WORDS64 (x, ix);
return x;
}
hidden_def (__roundeven)
libm_alias_double (__roundeven, roundeven)

60
sysdeps/ieee754/dbl-64/wordsize-64/s_scalbln.c
View File

@ -1,60 +0,0 @@
/*
* ====================================================
* 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.
* ====================================================
*/
/*
* scalbn (double x, int n)
* scalbn(x,n) returns x* 2**n computed by exponent
* manipulation rather than by actually performing an
* exponentiation or a multiplication.
*/
#include <math.h>
#include <math_private.h>
static const double
two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */
huge = 1.0e+300,
tiny = 1.0e-300;
double
__scalbln (double x, long int n)
{
int64_t ix;
int64_t k;
EXTRACT_WORDS64(ix,x);
k = (ix >> 52) & 0x7ff; /* extract exponent */
if (__builtin_expect(k==0, 0)) { /* 0 or subnormal x */
if ((ix & UINT64_C(0xfffffffffffff))==0) return x; /* +-0 */
x *= two54;
EXTRACT_WORDS64(ix,x);
k = ((ix >> 52) & 0x7ff) - 54;
}
if (__builtin_expect(k==0x7ff, 0)) return x+x; /* NaN or Inf */
if (__builtin_expect(n< -50000, 0))
return tiny*copysign(tiny,x); /*underflow*/
if (__builtin_expect(n> 50000 || k+n > 0x7fe, 0))
return huge*copysign(huge,x); /* overflow */
/* Now k and n are bounded we know that k = k+n does not
overflow. */
k = k+n;
if (__builtin_expect(k > 0, 1)) /* normal result */
{INSERT_WORDS64(x,(ix&UINT64_C(0x800fffffffffffff))|(k<<52));
return x;}
if (k <= -54)
return tiny*copysign(tiny,x); /*underflow*/
k += 54; /* subnormal result */
INSERT_WORDS64(x,(ix&INT64_C(0x800fffffffffffff))|(k<<52));
return x*twom54;
}
#ifdef NO_LONG_DOUBLE
strong_alias (__scalbln, __scalblnl)
#endif

60
sysdeps/ieee754/dbl-64/wordsize-64/s_scalbn.c
View File

@ -1,60 +0,0 @@
/*
* ====================================================
* 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.
* ====================================================
*/
/*
* scalbn (double x, int n)
* scalbn(x,n) returns x* 2**n computed by exponent
* manipulation rather than by actually performing an
* exponentiation or a multiplication.
*/
#include <math.h>
#include <math_private.h>
static const double
two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */
huge = 1.0e+300,
tiny = 1.0e-300;
double
__scalbn (double x, int n)
{
int64_t ix;
int64_t k;
EXTRACT_WORDS64(ix,x);
k = (ix >> 52) & 0x7ff; /* extract exponent */
if (__builtin_expect(k==0, 0)) { /* 0 or subnormal x */
if ((ix & UINT64_C(0xfffffffffffff))==0) return x; /* +-0 */
x *= two54;
EXTRACT_WORDS64(ix,x);
k = ((ix >> 52) & 0x7ff) - 54;
}
if (__builtin_expect(k==0x7ff, 0)) return x+x; /* NaN or Inf */
if (__builtin_expect(n< -50000, 0))
return tiny*copysign(tiny,x); /*underflow*/
if (__builtin_expect(n> 50000 || k+n > 0x7fe, 0))
return huge*copysign(huge,x); /* overflow */
/* Now k and n are bounded we know that k = k+n does not
overflow. */
k = k+n;
if (__builtin_expect(k > 0, 1)) /* normal result */
{INSERT_WORDS64(x,(ix&UINT64_C(0x800fffffffffffff))|(k<<52));
return x;}
if (k <= -54)
return tiny*copysign(tiny,x); /*underflow*/
k += 54; /* subnormal result */
INSERT_WORDS64(x,(ix&INT64_C(0x800fffffffffffff))|(k<<52));
return x*twom54;
}
#ifdef NO_LONG_DOUBLE
strong_alias (__scalbn, __scalbnl)
#endif

54
sysdeps/ieee754/dbl-64/wordsize-64/s_setpayload_main.c
View File

@ -1,54 +0,0 @@
/* Set NaN payload. dbl-64/wordsize-64 version.
Copyright (C) 2016-2021 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
<https://www.gnu.org/licenses/>. */
#include <math.h>
#include <math_private.h>
#include <libm-alias-double.h>
#include <nan-high-order-bit.h>
#include <stdint.h>
#define SET_HIGH_BIT (HIGH_ORDER_BIT_IS_SET_FOR_SNAN ? SIG : !SIG)
#define BIAS 0x3ff
#define PAYLOAD_DIG 51
#define EXPLICIT_MANT_DIG 52
int
FUNC (double *x, double payload)
{
uint64_t ix;
EXTRACT_WORDS64 (ix, payload);
int exponent = ix >> EXPLICIT_MANT_DIG;
/* Test if argument is (a) negative or too large; (b) too small,
except for 0 when allowed; (c) not an integer. */
if (exponent >= BIAS + PAYLOAD_DIG
|| (exponent < BIAS && !(SET_HIGH_BIT && ix == 0))
|| (ix & ((1ULL << (BIAS + EXPLICIT_MANT_DIG - exponent)) - 1)) != 0)
{
INSERT_WORDS64 (*x, 0);
return 1;
}
if (ix != 0)
{
ix &= (1ULL << EXPLICIT_MANT_DIG) - 1;
ix |= 1ULL << EXPLICIT_MANT_DIG;
ix >>= BIAS + EXPLICIT_MANT_DIG - exponent;
}
ix |= 0x7ff0000000000000ULL | (SET_HIGH_BIT ? 0x8000000000000ULL : 0);
INSERT_WORDS64 (*x, ix);
return 0;
}

76
sysdeps/ieee754/dbl-64/wordsize-64/s_totalorder.c
View File

@ -1,76 +0,0 @@
/* Total order operation. dbl-64/wordsize-64 version.
Copyright (C) 2016-2021 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
<https://www.gnu.org/licenses/>. */
#include <math.h>
#include <math_private.h>
#include <nan-high-order-bit.h>
#include <libm-alias-double.h>
#include <stdint.h>
#include <shlib-compat.h>
#include <first-versions.h>
int
__totalorder (const double *x, const double *y)
{
int64_t ix, iy;
EXTRACT_WORDS64 (ix, *x);
EXTRACT_WORDS64 (iy, *y);
#if HIGH_ORDER_BIT_IS_SET_FOR_SNAN
/* For the preferred quiet NaN convention, this operation is a
comparison of the representations of the arguments interpreted as
sign-magnitude integers. If both arguments are NaNs, invert the
quiet/signaling bit so comparing that way works. */
if ((ix & 0x7fffffffffffffffULL) > 0x7ff0000000000000ULL
&& (iy & 0x7fffffffffffffffULL) > 0x7ff0000000000000ULL)
{
ix ^= 0x0008000000000000ULL;
iy ^= 0x0008000000000000ULL;
}
#endif
uint64_t ix_sign = ix >> 63;
uint64_t iy_sign = iy >> 63;
ix ^= ix_sign >> 1;
iy ^= iy_sign >> 1;
return ix <= iy;
}
#ifdef SHARED
# define CONCATX(x, y) x ## y
# define CONCAT(x, y) CONCATX (x, y)
# define UNIQUE_ALIAS(name) CONCAT (name, __COUNTER__)
# define do_symbol(orig_name, name, aliasname) \
strong_alias (orig_name, name) \
versioned_symbol (libm, name, aliasname, GLIBC_2_31)
# undef weak_alias
# define weak_alias(name, aliasname) \
do_symbol (name, UNIQUE_ALIAS (name), aliasname);
#endif
libm_alias_double (__totalorder, totalorder)
#if SHLIB_COMPAT (libm, GLIBC_2_25, GLIBC_2_31)
int
attribute_compat_text_section
__totalorder_compat (double x, double y)
{
return __totalorder (&x, &y);
}
#undef do_symbol
#define do_symbol(orig_name, name, aliasname) \
strong_alias (orig_name, name) \
compat_symbol (libm, name, aliasname, \
CONCAT (FIRST_VERSION_libm_, aliasname))
libm_alias_double (__totalorder_compat, totalorder)
#endif

73
sysdeps/ieee754/dbl-64/wordsize-64/s_totalordermag.c
View File

@ -1,73 +0,0 @@
/* Total order operation on absolute values. dbl-64/wordsize-64 version.
Copyright (C) 2016-2021 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
<https://www.gnu.org/licenses/>. */
#include <math.h>
#include <math_private.h>
#include <nan-high-order-bit.h>
#include <libm-alias-double.h>
#include <stdint.h>
#include <shlib-compat.h>
#include <first-versions.h>
int
__totalordermag (const double *x, const double *y)
{
uint64_t ix, iy;
EXTRACT_WORDS64 (ix, *x);
EXTRACT_WORDS64 (iy, *y);
ix &= 0x7fffffffffffffffULL;
iy &= 0x7fffffffffffffffULL;
#if HIGH_ORDER_BIT_IS_SET_FOR_SNAN
/* For the preferred quiet NaN convention, this operation is a
comparison of the representations of the absolute values of the
arguments. If both arguments are NaNs, invert the
quiet/signaling bit so comparing that way works. */
if (ix > 0x7ff0000000000000ULL && iy > 0x7ff0000000000000ULL)
{
ix ^= 0x0008000000000000ULL;
iy ^= 0x0008000000000000ULL;
}
#endif
return ix <= iy;
}
#ifdef SHARED
# define CONCATX(x, y) x ## y
# define CONCAT(x, y) CONCATX (x, y)
# define UNIQUE_ALIAS(name) CONCAT (name, __COUNTER__)
# define do_symbol(orig_name, name, aliasname) \
strong_alias (orig_name, name) \
versioned_symbol (libm, name, aliasname, GLIBC_2_31)
# undef weak_alias
# define weak_alias(name, aliasname) \
do_symbol (name, UNIQUE_ALIAS (name), aliasname);
#endif
libm_alias_double (__totalordermag, totalordermag)
#if SHLIB_COMPAT (libm, GLIBC_2_25, GLIBC_2_31)
int
attribute_compat_text_section
__totalordermag_compat (double x, double y)
{
return __totalordermag (&x, &y);
}
#undef do_symbol
#define do_symbol(orig_name, name, aliasname) \
strong_alias (orig_name, name) \
compat_symbol (libm, name, aliasname, \
CONCAT (FIRST_VERSION_libm_, aliasname))
libm_alias_double (__totalordermag_compat, totalordermag)
#endif

1
sysdeps/mips/mips64/Implies
View File

@ -1,5 +1,4 @@
# MIPS uses IEEE 754 floating point. # MIPS uses IEEE 754 floating point.
mips/ieee754 mips/ieee754
ieee754/flt-32 ieee754/flt-32
ieee754/dbl-64/wordsize-64
ieee754/dbl-64 ieee754/dbl-64

1
sysdeps/s390/s390-64/Implies
View File

@ -1,2 +1 @@
wordsize-64 wordsize-64
ieee754/dbl-64/wordsize-64

1
sysdeps/sparc/sparc64/Implies
View File

@ -1,6 +1,5 @@
wordsize-64 wordsize-64
# SPARC uses IEEE 754 floating point. # SPARC uses IEEE 754 floating point.
ieee754/ldbl-128 ieee754/ldbl-128
ieee754/dbl-64/wordsize-64
ieee754/dbl-64 ieee754/dbl-64
ieee754/flt-32 ieee754/flt-32

1
sysdeps/x86_64/Implies
View File

@ -1,6 +1,5 @@
x86 x86
ieee754/float128 ieee754/float128
ieee754/ldbl-96 ieee754/ldbl-96
ieee754/dbl-64/wordsize-64
ieee754/dbl-64 ieee754/dbl-64
ieee754/flt-32 ieee754/flt-32