From FreeBSD:
authorPeter Avalos <pavalos@dragonflybsd.org>
Sun, 24 Jun 2007 05:17:51 +0000 (05:17 +0000)
committerPeter Avalos <pavalos@dragonflybsd.org>
Sun, 24 Jun 2007 05:17:51 +0000 (05:17 +0000)
Fixed the threshold for using the simple Taylor approximation.

In e_log.c, there was just a off-by-1 (1 ulp) error in the comment
about the threshold.  The precision of the threshold is unimportant,
but the magic numbers in the code are easier to understand when the
threshold is described precisely.

In e_logf.c, mistranslation of the magic numbers gave an off-by-1
(1 * 16 ulps) error in the intended negative bound for the threshold
and an off-by-7 (7 * 16 ulps) error in the intended positive bound for
the threshold, and the intended bounds were not translated from the
double precision bounds so they were unnecessarily small by a factor
of about 2048.

The optimization of using the simple Taylor approximation for args
near a power of 2 is dubious since it only applies to a relatively
small proportion of args, but if it is done then doing it 2048 times
as often _may_ be more efficient.  (My benchmarks show unexplained
dependencies on the data that increase with further optimizations
in this area.)

lib/libm/src/e_log.c
lib/libm/src/e_logf.c

index 2c304e0..917b2e1 100644 (file)
@@ -10,7 +10,7 @@
  * ====================================================
  *
  * $NetBSD: e_log.c,v 1.12 2002/05/26 22:01:51 wiz Exp $
- * $DragonFly: src/lib/libm/src/e_log.c,v 1.1 2005/07/26 21:15:20 joerg Exp $
+ * $DragonFly: src/lib/libm/src/e_log.c,v 1.2 2007/06/24 05:17:51 pavalos Exp $
  */
 
 /* log(x)
@@ -105,7 +105,7 @@ log(double x)
        SET_HIGH_WORD(x,hx|(i^0x3ff00000));     /* normalize x or x/2 */
        k += (i>>20);
        f = x-1.0;
-       if((0x000fffff&(2+hx))<3) {     /* |f| < 2**-20 */
+       if((0x000fffff&(2+hx))<3) {     /* -2**-20 <= f < 2**-20 */
            if(f==zero) { if(k==0) return zero;  else {dk=(double)k;
                                   return dk*ln2_hi+dk*ln2_lo;}
            }
index 38282d1..1302fea 100644 (file)
@@ -13,7 +13,7 @@
  * ====================================================
  *
  * $NetBSD: e_logf.c,v 1.8 2002/05/26 22:01:51 wiz Exp $
- * $DragonFly: src/lib/libm/src/e_logf.c,v 1.1 2005/07/26 21:15:20 joerg Exp $
+ * $DragonFly: src/lib/libm/src/e_logf.c,v 1.2 2007/06/24 05:17:51 pavalos Exp $
  */
 
 #include <math.h>
@@ -56,7 +56,7 @@ logf(float x)
        SET_FLOAT_WORD(x,ix|(i^0x3f800000));    /* normalize x or x/2 */
        k += (i>>23);
        f = x-(float)1.0;
-       if((0x007fffff&(15+ix))<16) {   /* |f| < 2**-20 */
+       if((0x007fffff&(0x8000+ix))<0xc000) {   /* -2**-9 <= f < 2**-9 */
            if(f==zero) { if(k==0) return zero;  else {dk=(float)k;
                                   return dk*ln2_hi+dk*ln2_lo;}
            }