Moved the optimization for tiny x from __kernel_{cos,sin}[f](x) to

{cos_sin}[f](x) so that x doesn't need to be reclassified in the
"kernel" functions to determine if it is tiny (it still needs to be
reclassified in the cosine case for other reasons that will go away).

This optimization is quite large for exponentially distributed x, since
x is tiny for almost half of the domain, but it is a pessimization for
uniformally distributed x since it takes a little time for all cases
but rarely applies.  Arg reduction on exponentially distributed x
rarely gives a tiny x unless the reduction is null, so it is best to
only do the optimization if the initial x is tiny, which is what this
commit arranges.  The imediate result is an average optimization of
1.4% relative to the previous version in a case that doesn't favour
the optimization (double cos(x) on all float x) and a large
pessimization for the relatively unimportant cases of lgamma[f][_r](x)
on tiny, negative, exponentially distributed x.  The optimization should
be recovered for lgamma*() as part of fixing lgamma*()'s low-quality
arg reduction.

Fixed various wrong constants for the cutoff for "tiny".  For cosine,
the cutoff is when x**2/2! == {FLT or DBL}_EPSILON/2.  We round down
to an integral power of 2 (and for cos() reduce the power by another
1) because the exact cutoff doesn't matter and would take more work
to determine.  For sine, the exact cutoff is larger due to the ration
of terms being x**2/3! instead of x**2/2!, but we use the same cutoff
as for cosine.  We now use a cutoff of 2**-27 for double precision and
2**-12 for single precision.  2**-27 was used in all cases but was
misspelled 2**27 in comments.  Wrong and sloppy cutoffs just cause
missed optimizations (provided the rounding mode is to nearest --
other modes just aren't supported).
This commit is contained in:
Bruce Evans 2005-10-24 14:08:36 +00:00
parent 4b401243a4
commit 4339c67c48
8 changed files with 22 additions and 22 deletions

View File

@ -69,9 +69,6 @@ __kernel_cos(double x, double y)
int32_t ix;
GET_HIGH_WORD(ix,x);
ix &= 0x7fffffff; /* ix = |x|'s high word*/
if(ix<0x3e400000) { /* if x < 2**27 */
if(((int)x)==0) return one; /* generate inexact */
}
z = x*x;
r = z*(C1+z*(C2+z*(C3+z*(C4+z*(C5+z*C6)))));
if(ix < 0x3FD33333) /* if |x| < 0.3 */

View File

@ -36,9 +36,6 @@ __kernel_cosf(float x, float y)
int32_t ix;
GET_FLOAT_WORD(ix,x);
ix &= 0x7fffffff; /* ix = |x|'s high word*/
if(ix<0x32000000) { /* if x < 2**27 */
if(((int)x)==0) return one; /* generate inexact */
}
z = x*x;
r = z*(C1+z*(C2+z*(C3+z*(C4+z*(C5+z*C6)))));
if(ix < 0x3e99999a) /* if |x| < 0.3 */

View File

@ -59,11 +59,7 @@ double
__kernel_sin(double x, double y, int iy)
{
double z,r,v;
int32_t ix;
GET_HIGH_WORD(ix,x);
ix &= 0x7fffffff; /* high word of x */
if(ix<0x3e400000) /* |x| < 2**-27 */
{if((int)x==0) return x;} /* generate inexact */
z = x*x;
v = z*x;
r = S2+z*(S3+z*(S4+z*(S5+z*S6)));

View File

@ -33,11 +33,7 @@ float
__kernel_sinf(float x, float y, int iy)
{
float z,r,v;
int32_t ix;
GET_FLOAT_WORD(ix,x);
ix &= 0x7fffffff; /* high word of x */
if(ix<0x32000000) /* |x| < 2**-27 */
{if((int)x==0) return x;} /* generate inexact */
z = x*x;
v = z*x;
r = S2+z*(S3+z*(S4+z*(S5+z*S6)));

View File

@ -59,7 +59,11 @@ cos(double x)
/* |x| ~< pi/4 */
ix &= 0x7fffffff;
if(ix <= 0x3fe921fb) return __kernel_cos(x,z);
if(ix <= 0x3fe921fb) {
if(ix<0x3e400000) /* if x < 2**-27 */
if(((int)x)==0) return 1.0; /* generate inexact */
return __kernel_cos(x,z);
}
/* cos(Inf or NaN) is NaN */
else if (ix>=0x7ff00000) return x-x;

View File

@ -20,8 +20,6 @@ static char rcsid[] = "$FreeBSD$";
#include "math.h"
#include "math_private.h"
static const float one=1.0;
float
cosf(float x)
{
@ -32,7 +30,11 @@ cosf(float x)
/* |x| ~< pi/4 */
ix &= 0x7fffffff;
if(ix <= 0x3f490fd8) return __kernel_cosf(x,z);
if(ix <= 0x3f490fd8) {
if(ix<0x39800000) /* if x < 2**-12 */
if(((int)x)==0) return 1.0; /* generate inexact */
return __kernel_cosf(x,z);
}
/* cos(Inf or NaN) is NaN */
else if (ix>=0x7f800000) return x-x;

View File

@ -59,7 +59,11 @@ sin(double x)
/* |x| ~< pi/4 */
ix &= 0x7fffffff;
if(ix <= 0x3fe921fb) return __kernel_sin(x,z,0);
if(ix <= 0x3fe921fb) {
if(ix<0x3e400000) /* |x| < 2**-27 */
{if((int)x==0) return x;} /* generate inexact */
return __kernel_sin(x,z,0);
}
/* sin(Inf or NaN) is NaN */
else if (ix>=0x7ff00000) return x-x;

View File

@ -30,7 +30,11 @@ sinf(float x)
/* |x| ~< pi/4 */
ix &= 0x7fffffff;
if(ix <= 0x3f490fd8) return __kernel_sinf(x,z,0);
if(ix <= 0x3f490fd8) {
if(ix<0x39800000) /* if x < 2**-12 */
if(((int)x)==0) return x; /* generate inexact */
return __kernel_sinf(x,z,0);
}
/* sin(Inf or NaN) is NaN */
else if (ix>=0x7f800000) return x-x;