mirror of
https://github.com/fluencelabs/musl
synced 2025-06-30 15:11:55 +00:00
math: remove STRICT_ASSIGN macro
gcc did not always drop excess precision according to c99 at assignments before version 4.5 even if -std=c99 was requested which caused badly broken mathematical functions on i386 when FLT_EVAL_METHOD!=0 but STRICT_ASSIGN was not used consistently and it is worked around for old compilers with -ffloat-store so it is no longer needed the new convention is to get the compiler respect c99 semantics and when excess precision is not harmful use float_t or double_t or to specialize code using FLT_EVAL_METHOD
This commit is contained in:
@ -155,15 +155,4 @@ long double __tanl(long double, long double, int);
|
|||||||
long double __polevll(long double, const long double *, int);
|
long double __polevll(long double, const long double *, int);
|
||||||
long double __p1evll(long double, const long double *, int);
|
long double __p1evll(long double, const long double *, int);
|
||||||
|
|
||||||
#if 0
|
|
||||||
/* Attempt to get strict C99 semantics for assignment with non-C99 compilers. */
|
|
||||||
#define STRICT_ASSIGN(type, lval, rval) do { \
|
|
||||||
volatile type __v = (rval); \
|
|
||||||
(lval) = __v; \
|
|
||||||
} while (0)
|
|
||||||
#else
|
|
||||||
/* Should work with -fexcess-precision=standard (>=gcc-4.5) or -ffloat-store */
|
|
||||||
#define STRICT_ASSIGN(type, lval, rval) ((lval) = (type)(rval))
|
|
||||||
#endif
|
|
||||||
|
|
||||||
#endif
|
#endif
|
||||||
|
@ -112,7 +112,7 @@ int __rem_pio2(double x, double *y)
|
|||||||
uint32_t high;
|
uint32_t high;
|
||||||
medium:
|
medium:
|
||||||
/* Use a specialized rint() to get fn. Assume round-to-nearest. */
|
/* Use a specialized rint() to get fn. Assume round-to-nearest. */
|
||||||
STRICT_ASSIGN(double, fn, x*invpio2 + 0x1.8p52);
|
fn = x*invpio2 + 0x1.8p52;
|
||||||
fn = fn - 0x1.8p52;
|
fn = fn - 0x1.8p52;
|
||||||
// FIXME
|
// FIXME
|
||||||
#ifdef HAVE_EFFICIENT_IRINT
|
#ifdef HAVE_EFFICIENT_IRINT
|
||||||
|
@ -415,7 +415,8 @@ recompute:
|
|||||||
fw = 0.0;
|
fw = 0.0;
|
||||||
for (i=jz; i>=0; i--)
|
for (i=jz; i>=0; i--)
|
||||||
fw += fq[i];
|
fw += fq[i];
|
||||||
STRICT_ASSIGN(double,fw,fw);
|
// TODO: drop excess precision here once double_t is used
|
||||||
|
fw = (double)fw;
|
||||||
y[0] = ih==0 ? fw : -fw;
|
y[0] = ih==0 ? fw : -fw;
|
||||||
fw = fq[0]-fw;
|
fw = fq[0]-fw;
|
||||||
for (i=1; i<=jz; i++)
|
for (i=1; i<=jz; i++)
|
||||||
|
@ -44,7 +44,7 @@ int __rem_pio2f(float x, double *y)
|
|||||||
/* 25+53 bit pi is good enough for medium size */
|
/* 25+53 bit pi is good enough for medium size */
|
||||||
if (ix < 0x4dc90fdb) { /* |x| ~< 2^28*(pi/2), medium size */
|
if (ix < 0x4dc90fdb) { /* |x| ~< 2^28*(pi/2), medium size */
|
||||||
/* Use a specialized rint() to get fn. Assume round-to-nearest. */
|
/* Use a specialized rint() to get fn. Assume round-to-nearest. */
|
||||||
STRICT_ASSIGN(double, fn, x*invpio2 + 0x1.8p52);
|
fn = x*invpio2 + 0x1.8p52;
|
||||||
fn = fn - 0x1.8p52;
|
fn = fn - 0x1.8p52;
|
||||||
// FIXME
|
// FIXME
|
||||||
#ifdef HAVE_EFFICIENT_IRINT
|
#ifdef HAVE_EFFICIENT_IRINT
|
||||||
|
@ -94,7 +94,7 @@ double exp(double x)
|
|||||||
return x;
|
return x;
|
||||||
if (x > 709.782712893383973096) {
|
if (x > 709.782712893383973096) {
|
||||||
/* overflow if x!=inf */
|
/* overflow if x!=inf */
|
||||||
STRICT_ASSIGN(double, x, 0x1p1023 * x);
|
x *= 0x1p1023;
|
||||||
return x;
|
return x;
|
||||||
}
|
}
|
||||||
if (x < -708.39641853226410622) {
|
if (x < -708.39641853226410622) {
|
||||||
@ -113,7 +113,7 @@ double exp(double x)
|
|||||||
k = 1 - sign - sign;
|
k = 1 - sign - sign;
|
||||||
hi = x - k*ln2hi; /* k*ln2hi is exact here */
|
hi = x - k*ln2hi; /* k*ln2hi is exact here */
|
||||||
lo = k*ln2lo;
|
lo = k*ln2lo;
|
||||||
STRICT_ASSIGN(double, x, hi - lo);
|
x = hi - lo;
|
||||||
} else if (hx > 0x3e300000) { /* if |x| > 2**-28 */
|
} else if (hx > 0x3e300000) { /* if |x| > 2**-28 */
|
||||||
k = 0;
|
k = 0;
|
||||||
hi = x;
|
hi = x;
|
||||||
|
@ -340,7 +340,7 @@ double exp2(double x)
|
|||||||
if (ix >= 0x408ff000) { /* |x| >= 1022 or nan */
|
if (ix >= 0x408ff000) { /* |x| >= 1022 or nan */
|
||||||
if (ix >= 0x40900000 && u.i>>63 == 0) { /* x >= 1024 or nan */
|
if (ix >= 0x40900000 && u.i>>63 == 0) { /* x >= 1024 or nan */
|
||||||
/* overflow */
|
/* overflow */
|
||||||
STRICT_ASSIGN(double, x, x * 0x1p1023);
|
x *= 0x1p1023;
|
||||||
return x;
|
return x;
|
||||||
}
|
}
|
||||||
if (ix >= 0x7ff00000) /* -inf or -nan */
|
if (ix >= 0x7ff00000) /* -inf or -nan */
|
||||||
|
@ -41,7 +41,7 @@ float expf(float x)
|
|||||||
if (hx >= 0x42aeac50) { /* if |x| >= -87.33655f or NaN */
|
if (hx >= 0x42aeac50) { /* if |x| >= -87.33655f or NaN */
|
||||||
if (hx >= 0x42b17218 && !sign) { /* x >= 88.722839f */
|
if (hx >= 0x42b17218 && !sign) { /* x >= 88.722839f */
|
||||||
/* overflow */
|
/* overflow */
|
||||||
STRICT_ASSIGN(float, x, x * 0x1p127f);
|
x *= 0x1p127f;
|
||||||
return x;
|
return x;
|
||||||
}
|
}
|
||||||
if (sign) {
|
if (sign) {
|
||||||
@ -60,7 +60,7 @@ float expf(float x)
|
|||||||
k = 1 - sign - sign;
|
k = 1 - sign - sign;
|
||||||
hi = x - k*ln2hi; /* k*ln2hi is exact here */
|
hi = x - k*ln2hi; /* k*ln2hi is exact here */
|
||||||
lo = k*ln2lo;
|
lo = k*ln2lo;
|
||||||
STRICT_ASSIGN(float, x, hi - lo);
|
x = hi - lo;
|
||||||
} else if (hx > 0x39000000) { /* |x| > 2**-14 */
|
} else if (hx > 0x39000000) { /* |x| > 2**-14 */
|
||||||
k = 0;
|
k = 0;
|
||||||
hi = x;
|
hi = x;
|
||||||
|
@ -155,7 +155,7 @@ double expm1(double x)
|
|||||||
hi = x - t*ln2_hi; /* t*ln2_hi is exact here */
|
hi = x - t*ln2_hi; /* t*ln2_hi is exact here */
|
||||||
lo = t*ln2_lo;
|
lo = t*ln2_lo;
|
||||||
}
|
}
|
||||||
STRICT_ASSIGN(double, x, hi - lo);
|
x = hi-lo;
|
||||||
c = (hi-x)-lo;
|
c = (hi-x)-lo;
|
||||||
} else if (hx < 0x3c900000) { /* |x| < 2**-54, return x */
|
} else if (hx < 0x3c900000) { /* |x| < 2**-54, return x */
|
||||||
if (hx < 0x00100000)
|
if (hx < 0x00100000)
|
||||||
|
@ -65,7 +65,7 @@ float expm1f(float x)
|
|||||||
hi = x - t*ln2_hi; /* t*ln2_hi is exact here */
|
hi = x - t*ln2_hi; /* t*ln2_hi is exact here */
|
||||||
lo = t*ln2_lo;
|
lo = t*ln2_lo;
|
||||||
}
|
}
|
||||||
STRICT_ASSIGN(float, x, hi - lo);
|
x = hi-lo;
|
||||||
c = (hi-x)-lo;
|
c = (hi-x)-lo;
|
||||||
} else if (hx < 0x33000000) { /* when |x|<2**-25, return x */
|
} else if (hx < 0x33000000) { /* when |x|<2**-25, return x */
|
||||||
if (hx < 0x00800000)
|
if (hx < 0x00800000)
|
||||||
|
@ -122,7 +122,7 @@ double log1p(double x)
|
|||||||
return x+x;
|
return x+x;
|
||||||
if (k != 0) {
|
if (k != 0) {
|
||||||
if (hx < 0x43400000) {
|
if (hx < 0x43400000) {
|
||||||
STRICT_ASSIGN(double, u, 1.0 + x);
|
u = 1 + x;
|
||||||
GET_HIGH_WORD(hu, u);
|
GET_HIGH_WORD(hu, u);
|
||||||
k = (hu>>20) - 1023;
|
k = (hu>>20) - 1023;
|
||||||
c = k > 0 ? 1.0-(u-x) : x-(u-1.0); /* correction term */
|
c = k > 0 ? 1.0-(u-x) : x-(u-1.0); /* correction term */
|
||||||
|
@ -61,7 +61,7 @@ float log1pf(float x)
|
|||||||
return x+x;
|
return x+x;
|
||||||
if (k != 0) {
|
if (k != 0) {
|
||||||
if (hx < 0x5a000000) {
|
if (hx < 0x5a000000) {
|
||||||
STRICT_ASSIGN(float, u, 1.0f + x);
|
u = 1 + x;
|
||||||
GET_FLOAT_WORD(hu, u);
|
GET_FLOAT_WORD(hu, u);
|
||||||
k = (hu>>23) - 127;
|
k = (hu>>23) - 127;
|
||||||
/* correction term */
|
/* correction term */
|
||||||
|
Reference in New Issue
Block a user