Mail Archives: djgpp-workers/1998/10/04/07:47:26
Below please find patches for libm.a sources from the last alpha
release. These problems were revealed with a test suite (which I will be
soon uploading for inclusion in djtst202.zip). That test suite
originally was part of the Cygnus library, but now has many improvements,
mainly by K.B. Williams, who also gets most of the credit for the patches
below.
*** src/libm/math/e_atan2.c~1 Fri Sep 25 15:39:38 1998
--- src/libm/math/e_atan2.c Fri Oct 2 18:23:48 1998
*************** pi_lo = 1.2246467991473531772E-16; /*
*** 100,106 ****
} else {
switch(m) {
case 0: return zero ; /* atan(+...,+INF) */
! case 1: return -zero ; /* atan(-...,+INF) */
case 2: return pi+tiny ; /* atan(+...,-INF) */
case 3: return -pi-tiny ; /* atan(-...,-INF) */
}
--- 100,112 ----
} else {
switch(m) {
case 0: return zero ; /* atan(+...,+INF) */
! case 1: /* atan(-...,+INF) */
! {
! /* Make sure we return -0, not +0. If we
! say "return -zero;", GCC might optimize
! that with FLDZ, which we don't want. */
! return copysign(zero, -1);
! }
case 2: return pi+tiny ; /* atan(+...,-INF) */
case 3: return -pi-tiny ; /* atan(-...,-INF) */
}
*** src/libm/math/e_hypot.c~0 Sat Feb 7 13:47:40 1998
--- src/libm/math/e_hypot.c Fri Sep 18 14:37:12 1998
***************
*** 57,62 ****
--- 57,64 ----
double a=x,b=y,t1,t2,y_1,y_2,w;
__int32_t j,k,ha,hb;
+ if (isnan(x) || isnan(y)) return (x-x)/(y-y);
+ if (isinf(x) || isinf(y)) return __kernel_standard(x, y, 104);
GET_HIGH_WORD(ha,x);
ha &= 0x7fffffff;
GET_HIGH_WORD(hb,y);
*** src/libm/math/e_sinh.c~0 Sat Jun 27 15:46:30 1998
--- src/libm/math/e_sinh.c Fri Sep 18 14:37:00 1998
***************
*** 29,36 ****
* only sinh(0)=0 is exact for finite x.
*/
#include "fdlibm.h"
-
#ifndef _DOUBLE_IS_32BITS
#ifdef __STDC__
--- 29,36 ----
* only sinh(0)=0 is exact for finite x.
*/
+ #include <float.h>
#include "fdlibm.h"
#ifndef _DOUBLE_IS_32BITS
#ifdef __STDC__
*************** static double one = 1.0, shuge = 1.0e307
*** 73,83 ****
/* |x| in [log(maxdouble), overflowthresold] */
GET_LOW_WORD(lx,x);
! if ((ix<0x408633CE) || (ix==0x408633ce && lx<=(__uint32_t)0x8fb9f87dU)) {
w = __ieee754_exp(0.5*fabs(x));
t = h*w;
return t*w;
}
/* |x| > overflowthresold, sinh(x) overflow */
return x*shuge;
--- 73,84 ----
/* |x| in [log(maxdouble), overflowthresold] */
GET_LOW_WORD(lx,x);
! if ((ix<0x408633CE) || (ix==0x408633ce && lx<0x8fb9f87eU)) {
w = __ieee754_exp(0.5*fabs(x));
t = h*w;
return t*w;
}
+ if ((ix==0x408633ce) && (lx==0x8fb9f87eU)) return DBL_MAX;
/* |x| > overflowthresold, sinh(x) overflow */
return x*shuge;
*** src/libm/math/ef_atan2.c~1 Sat Feb 7 14:13:26 1998
--- src/libm/math/ef_atan2.c Fri Oct 2 18:23:20 1998
*************** pi_lo = 1.5099578832e-07; /* 0x3422216
*** 72,78 ****
} else {
switch(m) {
case 0: return zero ; /* atan(+...,+INF) */
! case 1: return -zero ; /* atan(-...,+INF) */
case 2: return pi+tiny ; /* atan(+...,-INF) */
case 3: return -pi-tiny ; /* atan(-...,-INF) */
}
--- 72,85 ----
} else {
switch(m) {
case 0: return zero ; /* atan(+...,+INF) */
! case 1: /* atan(-...,+INF) */
! {
! /* Make sure we return -0, not +0. If we
! say "return -zero;", GCC might optimize
! that with FLDZ, which we don't want. */
! return copysignf(zero, -1);
! }
!
case 2: return pi+tiny ; /* atan(+...,-INF) */
case 3: return -pi-tiny ; /* atan(-...,-INF) */
}
*** src/libm/math/ef_exp.c~0 Sat Feb 7 14:13:26 1998
--- src/libm/math/ef_exp.c Fri Sep 18 14:54:12 1998
***************
*** 12,18 ****
* is preserved.
* ====================================================
*/
!
#include "fdlibm.h"
#ifdef __v810__
--- 12,18 ----
* is preserved.
* ====================================================
*/
! #include <float.h>
#include "fdlibm.h"
#ifdef __v810__
*************** one = 1.0,
*** 28,34 ****
halF[2] = {0.5,-0.5,},
huge = 1.0e+30,
twom100 = 7.8886090522e-31, /* 2**-100=0x0d800000 */
! o_threshold= 8.8721679688e+01, /* 0x42b17180 */
u_threshold= -1.0397208405e+02, /* 0xc2cff1b5 */
ln2HI[2] ={ 6.9313812256e-01, /* 0x3f317180 */
-6.9313812256e-01,}, /* 0xbf317180 */
--- 28,34 ----
halF[2] = {0.5,-0.5,},
huge = 1.0e+30,
twom100 = 7.8886090522e-31, /* 2**-100=0x0d800000 */
! o_threshold= 8.872283911e+01, /* 0x42b17218 */
u_threshold= -1.0397208405e+02, /* 0xc2cff1b5 */
ln2HI[2] ={ 6.9313812256e-01, /* 0x3f317180 */
-6.9313812256e-01,}, /* 0xbf317180 */
*************** P5 = 4.1381369442e-08; /* 0x3331bb4c
*** 58,63 ****
--- 58,67 ----
/* filter out non-finite argument */
if(hx >= 0x42b17218) { /* if |x|>=88.721... */
+ if(hx == 0x42b17218) { /* if |x| == 88.7228... */
+ if (xsb == 0)
+ return FLT_MAX;
+ }
if(hx>0x7f800000)
return x+x; /* NaN */
if(hx==0x7f800000)
*** src/libm/math/ef_hypot.c~0 Sat Feb 7 14:13:26 1998
--- src/libm/math/ef_hypot.c Fri Sep 18 14:55:34 1998
***************
*** 13,18 ****
--- 13,21 ----
* ====================================================
*/
+ #define SQRT_FLT_MAX 1.84467429742e+19 /* 0x5f7fffff */
+
+ #include <float.h>
#include "fdlibm.h"
#ifdef __STDC__
***************
*** 22,82 ****
float x, y;
#endif
{
! float a=x,b=y,t1,t2,y_1,y_2,w;
! __int32_t j,k,ha,hb;
GET_FLOAT_WORD(ha,x);
ha &= 0x7fffffff;
GET_FLOAT_WORD(hb,y);
hb &= 0x7fffffff;
! if(hb > ha) {a=y;b=x;j=ha; ha=hb;hb=j;} else {a=x;b=y;}
SET_FLOAT_WORD(a,ha); /* a <- |a| */
SET_FLOAT_WORD(b,hb); /* b <- |b| */
! if((ha-hb)>0xf000000) {return a+b;} /* x/y > 2**30 */
! k=0;
! if(ha > 0x58800000) { /* a>2**50 */
! if(ha >= 0x7f800000) { /* Inf or NaN */
! w = a+b; /* for sNaN */
! if(ha == 0x7f800000) w = a;
! if(hb == 0x7f800000) w = b;
return w;
- }
- /* scale a and b by 2**-60 */
- ha -= 0x5d800000; hb -= 0x5d800000; k += 60;
- SET_FLOAT_WORD(a,ha);
- SET_FLOAT_WORD(b,hb);
- }
- if(hb < 0x26800000) { /* b < 2**-50 */
- if(hb <= 0x007fffff) { /* subnormal b or 0 */
- if(hb==0) return a;
- SET_FLOAT_WORD(t1,0x3f000000); /* t1=2^126 */
- b *= t1;
- a *= t1;
- k -= 126;
- } else { /* scale a and b by 2^60 */
- ha += 0x5d800000; /* a *= 2^60 */
- hb += 0x5d800000; /* b *= 2^60 */
- k -= 60;
- SET_FLOAT_WORD(a,ha);
- SET_FLOAT_WORD(b,hb);
- }
- }
- /* medium size a and b */
- w = a-b;
- if (w>b) {
- SET_FLOAT_WORD(t1,ha&0xfffff000U);
- t2 = a-t1;
- w = __ieee754_sqrtf(t1*t1-(b*(-b)-t2*(a+t1)));
- } else {
- a = a+a;
- SET_FLOAT_WORD(y_1,hb&0xfffff000U);
- y_2 = b - y_1;
- SET_FLOAT_WORD(t1,ha+0x00800000);
- t2 = a - t1;
- w = __ieee754_sqrtf(t1*y_1-(w*(-w)-(t1*y_2+t2*b)));
- }
- if(k!=0) {
- SET_FLOAT_WORD(t1,0x3f800000+(k<<23));
- return t1*w;
- } else return w;
}
--- 25,52 ----
float x, y;
#endif
{
! float a=x,b=y,t1,t2,w;
! __int32_t j,ha,hb;
+ if (isnanf(x) || isnanf(y))
+ return (x-x)/(y-y);
+ if (isinff(x) || isinff(y))
+ return __kernel_standard(x, y, 104);
GET_FLOAT_WORD(ha,x);
ha &= 0x7fffffff;
GET_FLOAT_WORD(hb,y);
hb &= 0x7fffffff;
! if(hb < ha) {a=y;b=x;j=ha; ha=hb;hb=j;} else {a=x;b=y;}
SET_FLOAT_WORD(a,ha); /* a <- |a| */
SET_FLOAT_WORD(b,hb); /* b <- |b| */
! t1 = (a > 0) ? (a/b) : 0;
! if ((t1 >= SQRT_FLT_MAX) && finitef(x) && finitef(y))
! return __kernel_standard(x, y, 104);
! t1 *= t1;
! t2 = sqrtf(++t1);
! if ((t2 > 1) && (b >= FLT_MAX))
! return __kernel_standard(x, y, 104);
! w = t2 * b;
!
return w;
}
*** src/libm/math/fdlibm.h~0 Sat Jun 27 15:27:44 1998
--- src/libm/math/fdlibm.h Fri Oct 2 17:55:28 1998
***************
*** 11,16 ****
--- 11,19 ----
* ====================================================
*/
+ #ifndef __fdlibm_h
+ #define __fdlibm_h
+
/* CYGNUS LOCAL: Include files. */
#include <math.h>
#ifdef __DJGPP__
***************
*** 157,168 ****
typedef union
{
- double value;
struct
{
__uint32_t msw;
__uint32_t lsw;
} parts;
} ieee_double_shape_type;
#endif
--- 160,171 ----
typedef union
{
struct
{
__uint32_t msw;
__uint32_t lsw;
} parts;
+ double value;
} ieee_double_shape_type;
#endif
***************
*** 171,182 ****
typedef union
{
- double value;
struct
{
__uint32_t lsw;
__uint32_t msw;
} parts;
} ieee_double_shape_type;
#endif
--- 174,185 ----
typedef union
{
struct
{
__uint32_t lsw;
__uint32_t msw;
} parts;
+ double value;
} ieee_double_shape_type;
#endif
***************
*** 244,251 ****
typedef union
{
- float value;
__uint32_t word;
} ieee_float_shape_type;
/* Get a 32 bit int from a float. */
--- 247,254 ----
typedef union
{
__uint32_t word;
+ float value;
} ieee_float_shape_type;
/* Get a 32 bit int from a float. */
***************
*** 265,267 ****
--- 268,272 ----
sf_u.word = (i); \
(d) = sf_u.value; \
} while (0)
+
+ #endif /* __fdlibm_h */
*** src/libm/math/k_standard.c~0 Tue Apr 15 08:39:48 1997
--- src/libm/math/k_standard.c Fri Oct 2 18:59:44 1998
*************** static double zero = 0.0; /* used as con
*** 105,111 ****
/* acos(|x|>1) */
exc.type = DOMAIN;
exc.name = type < 100 ? "acos" : "acosf";
! exc.retval = zero;
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
--- 105,111 ----
/* acos(|x|>1) */
exc.type = DOMAIN;
exc.name = type < 100 ? "acos" : "acosf";
! exc.retval = nan();
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
*************** static double zero = 0.0; /* used as con
*** 120,126 ****
/* asin(|x|>1) */
exc.type = DOMAIN;
exc.name = type < 100 ? "asin" : "asinf";
! exc.retval = zero;
if(_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
--- 120,126 ----
/* asin(|x|>1) */
exc.type = DOMAIN;
exc.name = type < 100 ? "asin" : "asinf";
! exc.retval = nan();
if(_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
*************** static double zero = 0.0; /* used as con
*** 440,448 ****
y *= 0.5;
if(x<zero&&rint(y)!=y) exc.retval = -HUGE;
} else {
! exc.retval = HUGE_VAL;
y *= 0.5;
! if(x<zero&&rint(y)!=y) exc.retval = -HUGE_VAL;
}
if (_LIB_VERSION == _POSIX_)
errno = ERANGE;
--- 440,448 ----
y *= 0.5;
if(x<zero&&rint(y)!=y) exc.retval = -HUGE;
} else {
! exc.retval = infinity();
y *= 0.5;
! if(x<zero&&rint(y)!=y) exc.retval = -infinity();
}
if (_LIB_VERSION == _POSIX_)
errno = ERANGE;
*************** static double zero = 0.0; /* used as con
*** 482,488 ****
break;
case 24:
case 124:
! /* neg**non-integral */
exc.type = DOMAIN;
exc.name = type < 100 ? "pow" : "powf";
if (_LIB_VERSION == _SVID_)
--- 482,488 ----
break;
case 24:
case 124:
! /* pow(neg**non-integral) */
exc.type = DOMAIN;
exc.name = type < 100 ? "pow" : "powf";
if (_LIB_VERSION == _SVID_)
*************** static double zero = 0.0; /* used as con
*** 539,545 ****
if (_LIB_VERSION == _SVID_)
exc.retval = x;
else
! exc.retval = zero/zero;
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
--- 539,545 ----
if (_LIB_VERSION == _SVID_)
exc.retval = x;
else
! exc.retval = nan();
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
*************** static double zero = 0.0; /* used as con
*** 554,560 ****
/* remainder(x,0) */
exc.type = DOMAIN;
exc.name = type < 100 ? "remainder" : "remainderf";
! exc.retval = zero/zero;
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
--- 554,560 ----
/* remainder(x,0) */
exc.type = DOMAIN;
exc.name = type < 100 ? "remainder" : "remainderf";
! exc.retval = nan();
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
*************** static double zero = 0.0; /* used as con
*** 569,575 ****
/* acosh(x<1) */
exc.type = DOMAIN;
exc.name = type < 100 ? "acosh" : "acoshf";
! exc.retval = zero/zero;
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
--- 569,575 ----
/* acosh(x<1) */
exc.type = DOMAIN;
exc.name = type < 100 ? "acosh" : "acoshf";
! exc.retval = nan();
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
*************** static double zero = 0.0; /* used as con
*** 584,590 ****
/* atanh(|x|>1) */
exc.type = DOMAIN;
exc.name = type < 100 ? "atanh" : "atanhf";
! exc.retval = zero/zero;
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
--- 584,590 ----
/* atanh(|x|>1) */
exc.type = DOMAIN;
exc.name = type < 100 ? "atanh" : "atanhf";
! exc.retval = nan();
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
*************** static double zero = 0.0; /* used as con
*** 599,605 ****
/* atanh(|x|=1) */
exc.type = SING;
exc.name = type < 100 ? "atanh" : "atanhf";
! exc.retval = x/zero; /* sign(x)*inf */
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
--- 599,606 ----
/* atanh(|x|=1) */
exc.type = SING;
exc.name = type < 100 ? "atanh" : "atanhf";
! /* return sign(x)*inf */
! exc.retval = x < zero ? -infinity() : infinity();
if (_LIB_VERSION == _POSIX_)
errno = EDOM;
else if (!matherr(&exc)) {
*** src/libm/math/makefile.~0 Mon Jul 6 16:37:10 1998
--- src/libm/math/makefile Sat Sep 19 11:29:36 1998
***************
*** 25,31 ****
SRC += e_pow.c
SRC += e_rem_pio2.c
SRC += e_remainder.c
- SRC +=
SRC += e_scalb.c
SRC += e_sinh.c
SRC += e_sqrt.c
--- 25,30 ----
***************
*** 199,205 ****
$(GCC) -D_HAVE_STDC -O2 -o $@ $<
clean ::
! -$(MISC) rm *.def chew.exe
# Additional dependencies
$(OBJS): fdlibm.h $(TOP)/../../include/libm/math.h
--- 198,204 ----
$(GCC) -D_HAVE_STDC -O2 -o $@ $<
clean ::
! -$(MISC) rm *.def chew.exe targetdep.texi
# Additional dependencies
$(OBJS): fdlibm.h $(TOP)/../../include/libm/math.h
*** src/libm/math/s_log1p.c~0 Mon Jul 6 16:01:58 1998
--- src/libm/math/s_log1p.c Fri Sep 18 16:44:32 1998
***************
*** 161,170 ****
k = 1;
if (hx < 0x3FDA827A) { /* x < 0.41422 */
! if(ax>=0x3ff00000) { /* x <= -1.0 */
! if(x==-1.0) return -two54/zero; /* log1p(-1)=+inf */
! else return (x-x)/(x-x); /* log1p(x<-1)=NaN */
! }
if(ax<0x3e200000) { /* |x| < 2**-29 */
if(two54+x>zero /* raise inexact */
&&ax<0x3c900000) /* |x| < 2**-54 */
--- 161,168 ----
k = 1;
if (hx < 0x3FDA827A) { /* x < 0.41422 */
! if(ax>=0x3ff00000) /* x <= -1.0 */
! return log(1+x); /* log1p(-1) */
if(ax<0x3e200000) { /* |x| < 2**-29 */
if(two54+x>zero /* raise inexact */
&&ax<0x3c900000) /* |x| < 2**-54 */
*** src/libm/math/s_scalbn.c~0 Mon Jul 6 16:02:50 1998
--- src/libm/math/s_scalbn.c Sat Sep 19 22:12:52 1998
***************
*** 26,37 ****
TRAD_SYNOPSIS
#include <math.h>
! double scalbn(<[x]>,<[y]>)
double <[x]>;
! int <[y]>;
! float scalbnf(<[x]>,<[y]>)
float <[x]>;
! int <[y]>;
DESCRIPTION
<<scalbn>> and <<scalbnf>> scale <[x]> by <[n]>, returning <[x]> times
--- 26,37 ----
TRAD_SYNOPSIS
#include <math.h>
! double scalbn(<[x]>,<[n]>)
double <[x]>;
! int <[n]>;
! float scalbnf(<[x]>,<[n]>)
float <[x]>;
! int <[n]>;
DESCRIPTION
<<scalbn>> and <<scalbnf>> scale <[x]> by <[n]>, returning <[x]> times
***************
*** 88,99 ****
}
if (k==0x7ff) return x+x; /* NaN or Inf */
k = k+n;
! if (k > 0x7fe) return huge*copysign(huge,x); /* overflow */
if (k > 0) /* normal result */
{SET_HIGH_WORD(x,(hx&0x800fffffU)|(k<<20)); return x;}
if (k <= -54) {
if (n > 50000) /* in case integer overflow in n+k */
! return huge*copysign(huge,x); /*overflow*/
else return tiny*copysign(tiny,x); /*underflow*/
}
k += 54; /* subnormal result */
--- 88,99 ----
}
if (k==0x7ff) return x+x; /* NaN or Inf */
k = k+n;
! if (k > 0x7fe) return copysign(HUGE_VAL,x); /* overflow */
if (k > 0) /* normal result */
{SET_HIGH_WORD(x,(hx&0x800fffffU)|(k<<20)); return x;}
if (k <= -54) {
if (n > 50000) /* in case integer overflow in n+k */
! return copysign(HUGE_VAL,x); /*overflow*/
else return tiny*copysign(tiny,x); /*underflow*/
}
k += 54; /* subnormal result */
*** src/libm/math/sf_ldexp.c~0 Sat Feb 7 14:13:20 1998
--- src/libm/math/sf_ldexp.c Sat Sep 19 23:02:34 1998
***************
*** 15,20 ****
--- 15,21 ----
#include "fdlibm.h"
#include <errno.h>
+ #include <float.h>
#ifdef __STDC__
float ldexpf(float value, int expon)
***************
*** 25,31 ****
{
if(!finitef(value)||value==(float)0.0) return value;
value = scalbnf(value,expon);
! if(!finitef(value)||value==(float)0.0) errno = ERANGE;
return value;
}
--- 26,33 ----
{
if(!finitef(value)||value==(float)0.0) return value;
value = scalbnf(value,expon);
! if(!finitef(value)||
! (value < FLT_MIN && value > -FLT_MIN)) errno = ERANGE;
return value;
}
*** src/libm/math/sf_log1p.c~0 Mon Jul 6 16:04:38 1998
--- src/libm/math/sf_log1p.c Fri Sep 18 16:52:44 1998
***************
*** 52,61 ****
k = 1;
if (hx < 0x3ed413d7) { /* x < 0.41422 */
! if(ax>=0x3f800000) { /* x <= -1.0 */
! if(x==(float)-1.0) return -two25/zero; /* log1p(-1)=+inf */
! else return (x-x)/(x-x); /* log1p(x<-1)=NaN */
! }
if(ax<0x31000000) { /* |x| < 2**-29 */
if(two25+x>zero /* raise inexact */
&&ax<0x24800000) /* |x| < 2**-54 */
--- 52,59 ----
k = 1;
if (hx < 0x3ed413d7) { /* x < 0.41422 */
! if(ax>=0x3f800000) /* x <= -1.0 */
! return logf(1+x); /* log1p(x<=-1) */
if(ax<0x31000000) { /* |x| < 2**-29 */
if(two25+x>zero /* raise inexact */
&&ax<0x24800000) /* |x| < 2**-54 */
*** src/libm/math/sf_scalbn.c~0 Mon Jul 6 16:05:00 1998
--- src/libm/math/sf_scalbn.c Fri Sep 18 16:57:52 1998
***************
*** 51,62 ****
}
if (k==0xff) return x+x; /* NaN or Inf */
k = k+n;
! if (k > 0xfe) return huge*copysignf(huge,x); /* overflow */
if (k > 0) /* normal result */
{SET_FLOAT_WORD(x,(ix&0x807fffffU)|(k<<23)); return x;}
if (k <= -25) {
if (n > OVERFLOW_INT) /* in case integer overflow in n+k */
! return huge*copysignf(huge,x); /*overflow*/
else return tiny*copysignf(tiny,x); /*underflow*/
}
k += 25; /* subnormal result */
--- 51,62 ----
}
if (k==0xff) return x+x; /* NaN or Inf */
k = k+n;
! if (k > 0xfe) return copysignf(infinityf(),x); /* overflow */
if (k > 0) /* normal result */
{SET_FLOAT_WORD(x,(ix&0x807fffffU)|(k<<23)); return x;}
if (k <= -25) {
if (n > OVERFLOW_INT) /* in case integer overflow in n+k */
! return copysignf(infinityf(),x);/*overflow*/
else return tiny*copysignf(tiny,x); /*underflow*/
}
k += 25; /* subnormal result */
*** src/libm/math/wf_exp.c~0 Tue Apr 15 08:39:56 1997
--- src/libm/math/wf_exp.c Fri Sep 18 16:58:16 1998
***************
*** 24,30 ****
#else
static float
#endif
! o_threshold= 8.8721679688e+01, /* 0x42b17180 */
u_threshold= -1.0397208405e+02; /* 0xc2cff1b5 */
#ifdef __STDC__
--- 24,30 ----
#else
static float
#endif
! o_threshold= 8.872283911e+01, /* 0x42b17218 */
u_threshold= -1.0397208405e+02; /* 0xc2cff1b5 */
#ifdef __STDC__
*** src/libm/math/wf_pow.c~0 Tue Apr 15 08:39:56 1997
--- src/libm/math/wf_pow.c Fri Sep 18 16:59:06 1998
***************
*** 16,22 ****
/*
* wrapper powf(x,y) return x**y
*/
!
#include "fdlibm.h"
#ifdef __STDC__
--- 16,22 ----
/*
* wrapper powf(x,y) return x**y
*/
! #include <float.h>
#include "fdlibm.h"
#ifdef __STDC__
***************
*** 58,64 ****
return (float)__kernel_standard((double)x,(double)y,121);
}
}
! if(z==(float)0.0&&finitef(x)&&finitef(y))
/* powf underflow */
return (float)__kernel_standard((double)x,(double)y,122);
return z;
--- 58,64 ----
return (float)__kernel_standard((double)x,(double)y,121);
}
}
! if(z < FLT_MIN &&finitef(x)&&finitef(y))
/* powf underflow */
return (float)__kernel_standard((double)x,(double)y,122);
return z;
*** src/libm/math/s_nan.c~0 Sat Feb 7 14:13:28 1998
--- src/libm/math/s_nan.c Fri Oct 2 18:11:16 1998
*************** QUICKREF
*** 36,51 ****
#ifndef _DOUBLE_IS_32BITS
#ifdef __STDC__
double nan(void)
#else
double nan()
#endif
{
! double x;
!
! INSERT_WORDS(x,0x7ff80000,0);
! return x;
}
#endif /* _DOUBLE_IS_32BITS */
--- 36,58 ----
#ifndef _DOUBLE_IS_32BITS
+ /* WARNING: Some versions of GCC optimize expressions which
+ involve NaNs by emitting constant bit patterns hard-wired
+ into the compiler, and ignore the bit pattern in the source
+ code. For example, 7FF80000h below could be ignored and
+ FFF80000h used instead. The trick with initializing a union
+ used below was suggested by K.B. Williams and seems to work
+ for now, but you better watch future versions of compiler
+ to not mess this up. */
+ static const ieee_double_shape_type a_nan = { {0, 0x7ff80000} };
+
#ifdef __STDC__
double nan(void)
#else
double nan()
#endif
{
! return a_nan.value;
}
#endif /* _DOUBLE_IS_32BITS */
*** src/libm/math/sf_nan.c~0 Tue Apr 15 08:39:52 1997
--- src/libm/math/sf_nan.c Fri Oct 2 18:11:42 1998
***************
*** 5,16 ****
#include "fdlibm.h"
float nanf()
{
! float x;
!
! SET_FLOAT_WORD(x,0x7fc00000);
! return x;
}
#ifdef _DOUBLE_IS_32BITS
--- 5,23 ----
#include "fdlibm.h"
+ /* WARNING: Some versions of GCC optimize expressions which
+ involve NaNs by emitting constant bit patterns hard-wired
+ into the compiler, and ignore the bit pattern in the source
+ code. For example, 7FC00000h below could be ignored and
+ FFC00000h used instead. The trick with initializing a union
+ used below was suggested by K.B. Williams and seems to work
+ for now, but you better watch future versions of compiler
+ to not mess this up. */
+ static const ieee_float_shape_type a_nan = { 0x7fc00000 };
+
float nanf()
{
! return a_nan.value;
}
#ifdef _DOUBLE_IS_32BITS
*** src/libm/math/s_ldexp.c~0 Sat Feb 7 13:47:34 1998
--- src/libm/math/s_ldexp.c Sat Sep 19 23:00:32 1998
*************** PORTABILITY
*** 62,67 ****
--- 62,68 ----
#include "fdlibm.h"
#include <errno.h>
+ #include <float.h>
#ifndef _DOUBLE_IS_32BITS
*************** PORTABILITY
*** 74,80 ****
{
if(!finite(value)||value==0.0) return value;
value = scalbn(value,expon);
! if(!finite(value)||value==0.0) errno = ERANGE;
return value;
}
--- 75,82 ----
{
if(!finite(value)||value==0.0) return value;
value = scalbn(value,expon);
! if(!finite(value)||
! (value < DBL_MIN && value > -DBL_MIN)) errno = ERANGE;
return value;
}
- Raw text -