Use __ldexp_exp() to simplify things and improve accuracy for x near

the overflow threshold.
This commit is contained in:
David Schultz 2011-10-21 06:28:47 +00:00
parent 12188b77a2
commit d4657ac799
4 changed files with 10 additions and 27 deletions

View File

@ -45,7 +45,6 @@ __ieee754_cosh(double x)
{
double t,w;
int32_t ix;
u_int32_t lx;
/* High word of |x|. */
GET_HIGH_WORD(ix,x);
@ -72,13 +71,8 @@ __ieee754_cosh(double x)
if (ix < 0x40862E42) return half*__ieee754_exp(fabs(x));
/* |x| in [log(maxdouble), overflowthresold] */
GET_LOW_WORD(lx,x);
if (ix<0x408633CE ||
((ix==0x408633ce)&&(lx<=(u_int32_t)0x8fb9f87d))) {
w = __ieee754_exp(half*fabs(x));
t = half*w;
return t*w;
}
if (ix<=0x408633CE)
return __ldexp_exp(fabs(x), -1);
/* |x| > overflowthresold, cosh(x) overflow */
return huge*huge;

View File

@ -51,11 +51,8 @@ __ieee754_coshf(float x)
if (ix < 0x42b17217) return half*__ieee754_expf(fabsf(x));
/* |x| in [log(maxfloat), overflowthresold] */
if (ix<=0x42b2d4fc) {
w = __ieee754_expf(half*fabsf(x));
t = half*w;
return t*w;
}
if (ix<=0x42b2d4fc)
return __ldexp_expf(fabsf(x), -1);
/* |x| > overflowthresold, cosh(x) overflow */
return huge*huge;

View File

@ -40,9 +40,8 @@ static const double one = 1.0, shuge = 1.0e307;
double
__ieee754_sinh(double x)
{
double t,w,h;
double t,h;
int32_t ix,jx;
u_int32_t lx;
/* High word of |x|. */
GET_HIGH_WORD(jx,x);
@ -66,12 +65,8 @@ __ieee754_sinh(double x)
if (ix < 0x40862E42) return h*__ieee754_exp(fabs(x));
/* |x| in [log(maxdouble), overflowthresold] */
GET_LOW_WORD(lx,x);
if (ix<0x408633CE || ((ix==0x408633ce)&&(lx<=(u_int32_t)0x8fb9f87d))) {
w = __ieee754_exp(0.5*fabs(x));
t = h*w;
return t*w;
}
if (ix<=0x408633CE)
return h*2.0*__ldexp_exp(fabs(x), -1);
/* |x| > overflowthresold, sinh(x) overflow */
return x*shuge;

View File

@ -24,7 +24,7 @@ static const float one = 1.0, shuge = 1.0e37;
float
__ieee754_sinhf(float x)
{
float t,w,h;
float t,h;
int32_t ix,jx;
GET_FLOAT_WORD(jx,x);
@ -48,11 +48,8 @@ __ieee754_sinhf(float x)
if (ix < 0x42b17217) return h*__ieee754_expf(fabsf(x));
/* |x| in [logf(maxfloat), overflowthresold] */
if (ix<=0x42b2d4fc) {
w = __ieee754_expf((float)0.5*fabsf(x));
t = h*w;
return t*w;
}
if (ix<=0x42b2d4fc)
return h*2.0F*__ldexp_expf(fabsf(x), -1);
/* |x| > overflowthresold, sinh(x) overflow */
return x*shuge;