aboutsummaryrefslogtreecommitdiff
path: root/lib/msun
diff options
context:
space:
mode:
Diffstat (limited to 'lib/msun')
-rw-r--r--lib/msun/bsdsrc/b_tgamma.c2
-rw-r--r--lib/msun/ld128/e_rem_pio2l.h2
-rw-r--r--lib/msun/ld128/s_logl.c2
-rw-r--r--lib/msun/ld80/e_rem_pio2l.h2
-rw-r--r--lib/msun/ld80/s_logl.c2
-rw-r--r--lib/msun/src/e_fmod.c30
-rw-r--r--lib/msun/src/e_fmodf.c16
-rw-r--r--lib/msun/src/e_rem_pio2.c2
-rw-r--r--lib/msun/src/e_rem_pio2f.c2
-rw-r--r--lib/msun/src/e_remainder.c4
-rw-r--r--lib/msun/src/math_private.h35
-rw-r--r--lib/msun/src/s_ccosh.c23
-rw-r--r--lib/msun/src/s_ccoshf.c18
-rw-r--r--lib/msun/src/s_csinh.c23
-rw-r--r--lib/msun/src/s_csinhf.c18
-rw-r--r--lib/msun/src/s_ilogb.c13
-rw-r--r--lib/msun/src/s_ilogbf.c5
-rw-r--r--lib/msun/src/s_remquo.c30
-rw-r--r--lib/msun/src/s_remquof.c16
19 files changed, 148 insertions, 97 deletions
diff --git a/lib/msun/bsdsrc/b_tgamma.c b/lib/msun/bsdsrc/b_tgamma.c
index cc9f8f70297e..5ed69a28156f 100644
--- a/lib/msun/bsdsrc/b_tgamma.c
+++ b/lib/msun/bsdsrc/b_tgamma.c
@@ -261,7 +261,7 @@ small_gam(double x)
static double
smaller_gam(double x)
{
- double d, rhi, rlo, t, xhi, xlo;
+ double d, t, xhi, xlo;
struct Double r;
if (x < x0 + left) {
diff --git a/lib/msun/ld128/e_rem_pio2l.h b/lib/msun/ld128/e_rem_pio2l.h
index 35ed0b865a7c..d38a027dc215 100644
--- a/lib/msun/ld128/e_rem_pio2l.h
+++ b/lib/msun/ld128/e_rem_pio2l.h
@@ -56,7 +56,7 @@ pio2_2t = 2.0670321098263988236496903051604844e-43L, /* 0x127044533e63a0105df5
pio2_3 = 2.0670321098263988236499468110329591e-43L, /* 0x127044533e63a0105e00000000000.0p-254 */
pio2_3t = -2.5650587247459238361625433492959285e-65L; /* -0x159c4ec64ddaeb5f78671cbfb2210.0p-327 */
-static inline __always_inline int
+static __always_inline int
__ieee754_rem_pio2l(long double x, long double *y)
{
union IEEEl2bits u,u1;
diff --git a/lib/msun/ld128/s_logl.c b/lib/msun/ld128/s_logl.c
index 8961dd0c96db..734e0a8c2d4d 100644
--- a/lib/msun/ld128/s_logl.c
+++ b/lib/msun/ld128/s_logl.c
@@ -445,7 +445,7 @@ struct ld {
#endif
#ifdef STRUCT_RETURN
-static inline __always_inline void
+static __always_inline void
k_logl(long double x, struct ld *rp)
#else
long double
diff --git a/lib/msun/ld80/e_rem_pio2l.h b/lib/msun/ld80/e_rem_pio2l.h
index d9e4d3a0941f..ee35b28c83fd 100644
--- a/lib/msun/ld80/e_rem_pio2l.h
+++ b/lib/msun/ld80/e_rem_pio2l.h
@@ -68,7 +68,7 @@ pio2_2t = 6.36831716351095013979e-25L, /* 0xc51701b839a25205.0p-144 */
pio2_3t = -2.75299651904407171810e-37L; /* -0xbb5bf6c7ddd660ce.0p-185 */
#endif
-static inline __always_inline int
+static __always_inline int
__ieee754_rem_pio2l(long double x, long double *y)
{
union IEEEl2bits u,u1;
diff --git a/lib/msun/ld80/s_logl.c b/lib/msun/ld80/s_logl.c
index 459374d7d164..2305a906bd4e 100644
--- a/lib/msun/ld80/s_logl.c
+++ b/lib/msun/ld80/s_logl.c
@@ -445,7 +445,7 @@ struct ld {
#endif
#ifdef STRUCT_RETURN
-static inline __always_inline void
+static __always_inline void
k_logl(long double x, struct ld *rp)
#else
long double
diff --git a/lib/msun/src/e_fmod.c b/lib/msun/src/e_fmod.c
index fdfb56c2a475..973ab2ef30c4 100644
--- a/lib/msun/src/e_fmod.c
+++ b/lib/msun/src/e_fmod.c
@@ -28,14 +28,14 @@ static const double one = 1.0, Zero[] = {0.0, -0.0,};
double
fmod(double x, double y)
{
- int32_t n,hx,hy,hz,ix,iy,sx,i;
- u_int32_t lx,ly,lz;
+ int32_t hx, hy, hz, ix, iy, n, sx;
+ u_int32_t lx, ly, lz;
EXTRACT_WORDS(hx,lx,x);
EXTRACT_WORDS(hy,ly,y);
sx = hx&0x80000000; /* sign of x */
- hx ^=sx; /* |x| */
- hy &= 0x7fffffff; /* |y| */
+ hx ^= sx; /* |x| */
+ hy &= 0x7fffffff; /* |y| */
/* purge off exception values */
if((hy|ly)==0||(hx>=0x7ff00000)|| /* y=0,or x not finite */
@@ -48,22 +48,16 @@ fmod(double x, double y)
}
/* determine ix = ilogb(x) */
- if(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;
+ if(hx<0x00100000)
+ ix = subnormal_ilogb(hx, lx);
+ else
+ ix = (hx>>20)-1023;
/* determine iy = ilogb(y) */
- if(hy<0x00100000) { /* subnormal y */
- if(hy==0) {
- for (iy = -1043, i=ly; i>0; i<<=1) iy -=1;
- } else {
- for (iy = -1022,i=(hy<<11); i>0; i<<=1) iy -=1;
- }
- } else iy = (hy>>20)-1023;
+ if(hy<0x00100000)
+ iy = subnormal_ilogb(hy, ly);
+ else
+ iy = (hy>>20)-1023;
/* set up {hx,lx}, {hy,ly} and align y to x */
if(ix >= -1022)
diff --git a/lib/msun/src/e_fmodf.c b/lib/msun/src/e_fmodf.c
index 0e6633fbe739..fc6fee879f52 100644
--- a/lib/msun/src/e_fmodf.c
+++ b/lib/msun/src/e_fmodf.c
@@ -28,7 +28,7 @@ static const float one = 1.0, Zero[] = {0.0, -0.0,};
float
fmodf(float x, float y)
{
- int32_t n,hx,hy,hz,ix,iy,sx,i;
+ int32_t hx, hy, hz, ix, iy, n, sx;
GET_FLOAT_WORD(hx,x);
GET_FLOAT_WORD(hy,y);
@@ -45,14 +45,16 @@ fmodf(float x, float y)
return Zero[(u_int32_t)sx>>31]; /* |x|=|y| return x*0*/
/* determine ix = ilogb(x) */
- if(hx<0x00800000) { /* subnormal x */
- for (ix = -126,i=(hx<<8); i>0; i<<=1) ix -=1;
- } else ix = (hx>>23)-127;
+ if(hx<0x00800000)
+ ix = subnormal_ilogbf(hx);
+ else
+ ix = (hx>>23)-127;
/* determine iy = ilogb(y) */
- if(hy<0x00800000) { /* subnormal y */
- for (iy = -126,i=(hy<<8); i>=0; i<<=1) iy -=1;
- } else iy = (hy>>23)-127;
+ if(hy<0x00800000)
+ iy = subnormal_ilogbf(hy);
+ else
+ iy = (hy>>23)-127;
/* set up {hx,lx}, {hy,ly} and align y to x */
if(ix >= -126)
diff --git a/lib/msun/src/e_rem_pio2.c b/lib/msun/src/e_rem_pio2.c
index ef4107af94cb..d31cce4d6912 100644
--- a/lib/msun/src/e_rem_pio2.c
+++ b/lib/msun/src/e_rem_pio2.c
@@ -47,7 +47,7 @@ pio2_3 = 2.02226624871116645580e-21, /* 0x3BA3198A, 0x2E000000 */
pio2_3t = 8.47842766036889956997e-32; /* 0x397B839A, 0x252049C1 */
#ifdef INLINE_REM_PIO2
-static __inline __always_inline
+static __always_inline
#endif
int
__ieee754_rem_pio2(double x, double *y)
diff --git a/lib/msun/src/e_rem_pio2f.c b/lib/msun/src/e_rem_pio2f.c
index 26f6bc85791b..1740b4f49a58 100644
--- a/lib/msun/src/e_rem_pio2f.c
+++ b/lib/msun/src/e_rem_pio2f.c
@@ -39,7 +39,7 @@ pio2_1 = 1.57079631090164184570e+00, /* 0x3FF921FB, 0x50000000 */
pio2_1t = 1.58932547735281966916e-08; /* 0x3E5110b4, 0x611A6263 */
#ifdef INLINE_REM_PIO2F
-static __inline __always_inline
+static __always_inline
#endif
int
__ieee754_rem_pio2f(float x, double *y)
diff --git a/lib/msun/src/e_remainder.c b/lib/msun/src/e_remainder.c
index bd1ce8950619..5b71088911a2 100644
--- a/lib/msun/src/e_remainder.c
+++ b/lib/msun/src/e_remainder.c
@@ -66,8 +66,8 @@ remainder(double x, double p)
if(x>=p_half) x -= p;
}
}
- GET_HIGH_WORD(hx,x);
- if ((hx&0x7fffffff)==0) hx = 0;
+ EXTRACT_WORDS(hx, lx, x);
+ if (((hx&0x7fffffff)|lx) == 0) hx = 0;
SET_HIGH_WORD(x,hx^sx);
return x;
}
diff --git a/lib/msun/src/math_private.h b/lib/msun/src/math_private.h
index b71a8d03c0a1..6513c2ce14fc 100644
--- a/lib/msun/src/math_private.h
+++ b/lib/msun/src/math_private.h
@@ -740,6 +740,41 @@ irintl(long double x)
(ar) = (x) - (ai); \
} while (0)
+/*
+ * For a subnormal double entity split into high and low parts, compute ilogb.
+ */
+static inline int32_t
+subnormal_ilogb(int32_t hi, int32_t lo)
+{
+ int32_t j;
+ uint32_t i;
+
+ j = -1022;
+ if (hi == 0) {
+ j -= 21;
+ i = (uint32_t)lo;
+ } else
+ i = (uint32_t)hi << 11;
+
+ for (; i < 0x7fffffff; i <<= 1) j -= 1;
+
+ return (j);
+}
+
+/*
+ * For a subnormal float entity represented as an int32_t, compute ilogb.
+ */
+static inline int32_t
+subnormal_ilogbf(int32_t hx)
+{
+ int32_t j;
+ uint32_t i;
+ i = (uint32_t) hx << 8;
+ for (j = -126; i < 0x7fffffff; i <<= 1) j -=1;
+
+ return (j);
+}
+
#ifdef DEBUG
#if defined(__amd64__) || defined(__i386__)
#define breakpoint() asm("int $3")
diff --git a/lib/msun/src/s_ccosh.c b/lib/msun/src/s_ccosh.c
index 95ed3a32ddd7..145d9673f392 100644
--- a/lib/msun/src/s_ccosh.c
+++ b/lib/msun/src/s_ccosh.c
@@ -1,7 +1,7 @@
/*-
* SPDX-License-Identifier: BSD-2-Clause
*
- * Copyright (c) 2005 Bruce D. Evans and Steven G. Kargl
+ * Copyright (c) 2005-2025 Bruce D. Evans and Steven G. Kargl
* All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
@@ -49,7 +49,7 @@ static const double huge = 0x1p1023;
double complex
ccosh(double complex z)
{
- double x, y, h;
+ double c, h, s, x, y;
int32_t hx, hy, ix, iy, lx, ly;
x = creal(z);
@@ -65,14 +65,16 @@ ccosh(double complex z)
if (ix < 0x7ff00000 && iy < 0x7ff00000) {
if ((iy | ly) == 0)
return (CMPLX(cosh(x), x * y));
+
+ sincos(y, &s, &c);
if (ix < 0x40360000) /* |x| < 22: normal case */
- return (CMPLX(cosh(x) * cos(y), sinh(x) * sin(y)));
+ return (CMPLX(cosh(x) * c, sinh(x) * s));
/* |x| >= 22, so cosh(x) ~= exp(|x|) */
if (ix < 0x40862e42) {
/* x < 710: exp(|x|) won't overflow */
- h = exp(fabs(x)) * 0.5;
- return (CMPLX(h * cos(y), copysign(h, x) * sin(y)));
+ h = exp(fabs(x)) / 2;
+ return (CMPLX(h * c, copysign(h, x) * s));
} else if (ix < 0x4096bbaa) {
/* x < 1455: scale to avoid overflow */
z = __ldexp_cexp(CMPLX(fabs(x), y), -1);
@@ -80,7 +82,7 @@ ccosh(double complex z)
} else {
/* x >= 1455: the result always overflows */
h = huge * x;
- return (CMPLX(h * h * cos(y), h * sin(y)));
+ return (CMPLX(h * h * c, h * s));
}
}
@@ -130,7 +132,9 @@ ccosh(double complex z)
if (ix == 0x7ff00000 && lx == 0) {
if (iy >= 0x7ff00000)
return (CMPLX(INFINITY, x * (y - y)));
- return (CMPLX(INFINITY * cos(y), x * sin(y)));
+
+ sincos(y, &s, &c);
+ return (CMPLX(INFINITY * c, x * s));
}
/*
@@ -155,3 +159,8 @@ ccos(double complex z)
/* ccos(z) = ccosh(I * z) */
return (ccosh(CMPLX(-cimag(z), creal(z))));
}
+
+#if (LDBL_MANT_DIG == 53)
+__weak_reference(ccosh, ccoshl);
+__weak_reference(ccos, ccosl);
+#endif
diff --git a/lib/msun/src/s_ccoshf.c b/lib/msun/src/s_ccoshf.c
index ba97a390c832..2851157ead42 100644
--- a/lib/msun/src/s_ccoshf.c
+++ b/lib/msun/src/s_ccoshf.c
@@ -1,7 +1,7 @@
/*-
* SPDX-License-Identifier: BSD-2-Clause
*
- * Copyright (c) 2005 Bruce D. Evans and Steven G. Kargl
+ * Copyright (c) 2005-2025 Bruce D. Evans and Steven G. Kargl
* All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
@@ -41,7 +41,7 @@ static const float huge = 0x1p127;
float complex
ccoshf(float complex z)
{
- float x, y, h;
+ float c, h, s, x, y;
int32_t hx, hy, ix, iy;
x = crealf(z);
@@ -56,14 +56,16 @@ ccoshf(float complex z)
if (ix < 0x7f800000 && iy < 0x7f800000) {
if (iy == 0)
return (CMPLXF(coshf(x), x * y));
+
+ sincosf(y, &s, &c);
if (ix < 0x41100000) /* |x| < 9: normal case */
- return (CMPLXF(coshf(x) * cosf(y), sinhf(x) * sinf(y)));
+ return (CMPLXF(coshf(x) * c, sinhf(x) * s));
/* |x| >= 9, so cosh(x) ~= exp(|x|) */
if (ix < 0x42b17218) {
/* x < 88.7: expf(|x|) won't overflow */
- h = expf(fabsf(x)) * 0.5F;
- return (CMPLXF(h * cosf(y), copysignf(h, x) * sinf(y)));
+ h = expf(fabsf(x)) / 2;
+ return (CMPLXF(h * c, copysignf(h, x) * s));
} else if (ix < 0x4340b1e7) {
/* x < 192.7: scale to avoid overflow */
z = __ldexp_cexpf(CMPLXF(fabsf(x), y), -1);
@@ -71,7 +73,7 @@ ccoshf(float complex z)
} else {
/* x >= 192.7: the result always overflows */
h = huge * x;
- return (CMPLXF(h * h * cosf(y), h * sinf(y)));
+ return (CMPLXF(h * h * c, h * s));
}
}
@@ -87,7 +89,9 @@ ccoshf(float complex z)
if (ix == 0x7f800000) {
if (iy >= 0x7f800000)
return (CMPLXF(INFINITY, x * (y - y)));
- return (CMPLXF(INFINITY * cosf(y), x * sinf(y)));
+
+ sincosf(y, &s, &c);
+ return (CMPLXF(INFINITY * c, x * s));
}
return (CMPLXF(((long double)x * x) * (y - y),
diff --git a/lib/msun/src/s_csinh.c b/lib/msun/src/s_csinh.c
index 1bd78b1e49bf..8ea97970af7f 100644
--- a/lib/msun/src/s_csinh.c
+++ b/lib/msun/src/s_csinh.c
@@ -1,7 +1,7 @@
/*-
* SPDX-License-Identifier: BSD-2-Clause
*
- * Copyright (c) 2005 Bruce D. Evans and Steven G. Kargl
+ * Copyright (c) 2005-2025 Bruce D. Evans and Steven G. Kargl
* All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
@@ -49,7 +49,7 @@ static const double huge = 0x1p1023;
double complex
csinh(double complex z)
{
- double x, y, h;
+ double c, h, s, x, y;
int32_t hx, hy, ix, iy, lx, ly;
x = creal(z);
@@ -65,14 +65,16 @@ csinh(double complex z)
if (ix < 0x7ff00000 && iy < 0x7ff00000) {
if ((iy | ly) == 0)
return (CMPLX(sinh(x), y));
+
+ sincos(y, &s, &c);
if (ix < 0x40360000) /* |x| < 22: normal case */
- return (CMPLX(sinh(x) * cos(y), cosh(x) * sin(y)));
+ return (CMPLX(sinh(x) * c, cosh(x) * s));
/* |x| >= 22, so cosh(x) ~= exp(|x|) */
if (ix < 0x40862e42) {
/* x < 710: exp(|x|) won't overflow */
- h = exp(fabs(x)) * 0.5;
- return (CMPLX(copysign(h, x) * cos(y), h * sin(y)));
+ h = exp(fabs(x)) / 2;
+ return (CMPLX(copysign(h, x) * c, h * s));
} else if (ix < 0x4096bbaa) {
/* x < 1455: scale to avoid overflow */
z = __ldexp_cexp(CMPLX(fabs(x), y), -1);
@@ -80,7 +82,7 @@ csinh(double complex z)
} else {
/* x >= 1455: the result always overflows */
h = huge * x;
- return (CMPLX(h * cos(y), h * h * sin(y)));
+ return (CMPLX(h * c, h * h * s));
}
}
@@ -129,7 +131,9 @@ csinh(double complex z)
if (ix == 0x7ff00000 && lx == 0) {
if (iy >= 0x7ff00000)
return (CMPLX(x, y - y));
- return (CMPLX(x * cos(y), INFINITY * sin(y)));
+
+ sincos(y, &s, &c);
+ return (CMPLX(x * c, INFINITY * s));
}
/*
@@ -155,3 +159,8 @@ csin(double complex z)
z = csinh(CMPLX(cimag(z), creal(z)));
return (CMPLX(cimag(z), creal(z)));
}
+
+#if (LDBL_MANT_DIG == 53)
+__weak_reference(csinh, csinhl);
+__weak_reference(csin, csinl);
+#endif
diff --git a/lib/msun/src/s_csinhf.c b/lib/msun/src/s_csinhf.c
index b1f333955e53..530ac63b6a27 100644
--- a/lib/msun/src/s_csinhf.c
+++ b/lib/msun/src/s_csinhf.c
@@ -1,7 +1,7 @@
/*-
* SPDX-License-Identifier: BSD-2-Clause
*
- * Copyright (c) 2005 Bruce D. Evans and Steven G. Kargl
+ * Copyright (c) 2005-2025 Bruce D. Evans and Steven G. Kargl
* All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
@@ -41,7 +41,7 @@ static const float huge = 0x1p127;
float complex
csinhf(float complex z)
{
- float x, y, h;
+ float c, h, s, x, y;
int32_t hx, hy, ix, iy;
x = crealf(z);
@@ -56,14 +56,16 @@ csinhf(float complex z)
if (ix < 0x7f800000 && iy < 0x7f800000) {
if (iy == 0)
return (CMPLXF(sinhf(x), y));
+
+ sincosf(y, &s, &c);
if (ix < 0x41100000) /* |x| < 9: normal case */
- return (CMPLXF(sinhf(x) * cosf(y), coshf(x) * sinf(y)));
+ return (CMPLXF(sinhf(x) * c, coshf(x) * s));
/* |x| >= 9, so cosh(x) ~= exp(|x|) */
if (ix < 0x42b17218) {
/* x < 88.7: expf(|x|) won't overflow */
- h = expf(fabsf(x)) * 0.5F;
- return (CMPLXF(copysignf(h, x) * cosf(y), h * sinf(y)));
+ h = expf(fabsf(x)) / 2;
+ return (CMPLXF(copysignf(h, x) * c, h * s));
} else if (ix < 0x4340b1e7) {
/* x < 192.7: scale to avoid overflow */
z = __ldexp_cexpf(CMPLXF(fabsf(x), y), -1);
@@ -71,7 +73,7 @@ csinhf(float complex z)
} else {
/* x >= 192.7: the result always overflows */
h = huge * x;
- return (CMPLXF(h * cosf(y), h * h * sinf(y)));
+ return (CMPLXF(h * c, h * h * s));
}
}
@@ -87,7 +89,9 @@ csinhf(float complex z)
if (ix == 0x7f800000) {
if (iy >= 0x7f800000)
return (CMPLXF(x, y - y));
- return (CMPLXF(x * cosf(y), INFINITY * sinf(y)));
+
+ sincosf(y, &s, &c);
+ return (CMPLXF(x * c, INFINITY * s));
}
return (CMPLXF(((long double)x + x) * (y - y),
diff --git a/lib/msun/src/s_ilogb.c b/lib/msun/src/s_ilogb.c
index 0b076edbd9b5..a958182a2a35 100644
--- a/lib/msun/src/s_ilogb.c
+++ b/lib/msun/src/s_ilogb.c
@@ -23,21 +23,18 @@
#include "math.h"
#include "math_private.h"
- int ilogb(double x)
+int
+ilogb(double x)
{
- int32_t hx,lx,ix;
+ int32_t hx, ix, lx;
EXTRACT_WORDS(hx,lx,x);
hx &= 0x7fffffff;
if(hx<0x00100000) {
if((hx|lx)==0)
return FP_ILOGB0;
- else /* subnormal x */
- if(hx==0) {
- for (ix = -1043; lx>0; lx<<=1) ix -=1;
- } else {
- for (ix = -1022,hx<<=11; hx>0; hx<<=1) ix -=1;
- }
+ else
+ ix = subnormal_ilogb(hx, lx);
return ix;
}
else if (hx<0x7ff00000) return (hx>>20)-1023;
diff --git a/lib/msun/src/s_ilogbf.c b/lib/msun/src/s_ilogbf.c
index ff3df1fc5b90..75508ab77120 100644
--- a/lib/msun/src/s_ilogbf.c
+++ b/lib/msun/src/s_ilogbf.c
@@ -19,7 +19,8 @@
#include "math.h"
#include "math_private.h"
- int ilogbf(float x)
+int
+ilogbf(float x)
{
int32_t hx,ix;
@@ -29,7 +30,7 @@
if(hx==0)
return FP_ILOGB0;
else /* subnormal x */
- for (ix = -126,hx<<=8; hx>0; hx<<=1) ix -=1;
+ ix = subnormal_ilogbf(hx);
return ix;
}
else if (hx<0x7f800000) return (hx>>23)-127;
diff --git a/lib/msun/src/s_remquo.c b/lib/msun/src/s_remquo.c
index e3aac25230e0..fba788ed5104 100644
--- a/lib/msun/src/s_remquo.c
+++ b/lib/msun/src/s_remquo.c
@@ -5,7 +5,7 @@
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
- * software is freely granted, provided that this notice
+ * software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
@@ -29,7 +29,7 @@ static const double Zero[] = {0.0, -0.0,};
double
remquo(double x, double y, int *quo)
{
- int32_t n,hx,hy,hz,ix,iy,sx,i;
+ int32_t hx,hy,hz,ix,iy,n,sx;
u_int32_t lx,ly,lz,q,sxy;
EXTRACT_WORDS(hx,lx,x);
@@ -55,25 +55,19 @@ remquo(double x, double y, int *quo)
}
/* determine ix = ilogb(x) */
- if(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;
+ if(hx<0x00100000)
+ ix = subnormal_ilogb(hx, lx);
+ else
+ ix = (hx>>20)-1023;
/* determine iy = ilogb(y) */
- if(hy<0x00100000) { /* subnormal y */
- if(hy==0) {
- for (iy = -1043, i=ly; i>0; i<<=1) iy -=1;
- } else {
- for (iy = -1022,i=(hy<<11); i>0; i<<=1) iy -=1;
- }
- } else iy = (hy>>20)-1023;
+ if(hy<0x00100000)
+ iy = subnormal_ilogb(hy, ly);
+ else
+ iy = (hy>>20)-1023;
/* set up {hx,lx}, {hy,ly} and align y to x */
- if(ix >= -1022)
+ if(ix >= -1022)
hx = 0x00100000|(0x000fffff&hx);
else { /* subnormal x, shift x to normal */
n = -1022-ix;
@@ -85,7 +79,7 @@ remquo(double x, double y, int *quo)
lx = 0;
}
}
- if(iy >= -1022)
+ if(iy >= -1022)
hy = 0x00100000|(0x000fffff&hy);
else { /* subnormal y, shift y to normal */
n = -1022-iy;
diff --git a/lib/msun/src/s_remquof.c b/lib/msun/src/s_remquof.c
index c42bd8c4320d..39b87a91845f 100644
--- a/lib/msun/src/s_remquof.c
+++ b/lib/msun/src/s_remquof.c
@@ -27,7 +27,7 @@ static const float Zero[] = {0.0, -0.0,};
float
remquof(float x, float y, int *quo)
{
- int32_t n,hx,hy,hz,ix,iy,sx,i;
+ int32_t hx, hy, hz, ix, iy, n, sx;
u_int32_t q,sxy;
GET_FLOAT_WORD(hx,x);
@@ -49,14 +49,16 @@ remquof(float x, float y, int *quo)
}
/* determine ix = ilogb(x) */
- if(hx<0x00800000) { /* subnormal x */
- for (ix = -126,i=(hx<<8); i>0; i<<=1) ix -=1;
- } else ix = (hx>>23)-127;
+ if(hx<0x00800000)
+ ix = subnormal_ilogbf(hx);
+ else
+ ix = (hx>>23)-127;
/* determine iy = ilogb(y) */
- if(hy<0x00800000) { /* subnormal y */
- for (iy = -126,i=(hy<<8); i>0; i<<=1) iy -=1;
- } else iy = (hy>>23)-127;
+ if(hy<0x00800000)
+ iy = subnormal_ilogbf(hy);
+ else
+ iy = (hy>>23)-127;
/* set up {hx,lx}, {hy,ly} and align y to x */
if(ix >= -126)