[Numpy-svn] r8324 - in branches/1.4.x/numpy/core: src/npymath tests

numpy-svn@scip... numpy-svn@scip...
Sun Apr 4 23:51:33 CDT 2010


Author: charris
Date: 2010-04-04 23:51:33 -0500 (Sun, 04 Apr 2010)
New Revision: 8324

Modified:
   branches/1.4.x/numpy/core/src/npymath/npy_math.c.src
   branches/1.4.x/numpy/core/tests/test_regression.py
Log:
BUG: Backport r8320-r8223 from trunk.


Modified: branches/1.4.x/numpy/core/src/npymath/npy_math.c.src
===================================================================
--- branches/1.4.x/numpy/core/src/npymath/npy_math.c.src	2010-04-05 04:38:49 UTC (rev 8323)
+++ branches/1.4.x/numpy/core/src/npymath/npy_math.c.src	2010-04-05 04:51:33 UTC (rev 8324)
@@ -65,13 +65,14 @@
 #ifndef HAVE_EXPM1
 double npy_expm1(double x)
 {
-    double u = npy_exp(x);
+    const double u = npy_exp(x);
+
     if (u == 1.0) {
         return x;
-    } else if (u-1.0 == -1.0) {
+    } else if (u - 1.0 == -1.0) {
         return -1;
     } else {
-        return (u-1.0) * x/npy_log(u);
+        return (u - 1.0) * x/npy_log(u);
     }
 }
 #endif
@@ -79,11 +80,13 @@
 #ifndef HAVE_LOG1P
 double npy_log1p(double x)
 {
-    double u = 1. + x;
-    if (u == 1.0) {
+    const double u = 1. + x;
+    const double d = u - 1.;
+
+    if (d == 0) {
         return x;
     } else {
-        return npy_log(u) * x / (u - 1);
+        return npy_log(u) * x / d;
     }
 }
 #endif
@@ -154,11 +157,11 @@
     }
 
     /* compute y/x */
-    k = (iy - ix)>>20;
-    if(k > 60) {            /* |y/x| >  2**60 */
+    k = (iy - ix) >> 20;
+    if (k > 60) {            /* |y/x| >  2**60 */
         z = NPY_PI_2 + 0.5 * NPY_DBL_EPSILON;
         m &= 1;
-    } else if(hx < 0 && k < -60) {
+    } else if (hx < 0 && k < -60) {
         z = 0.0;    /* 0 > |y|/x > -2**-60 */
     } else {
         z = npy_atan(npy_fabs(y/x));        /* safe to do y/x */
@@ -209,7 +212,7 @@
 #ifndef HAVE_ACOSH
 double npy_acosh(double x)
 {
-    return 2*npy_log(npy_sqrt((x+1.0)/2)+npy_sqrt((x-1.0)/2));
+    return 2*npy_log(npy_sqrt((x + 1.0)/2) + npy_sqrt((x - 1.0)/2));
 }
 #endif
 
@@ -255,14 +258,15 @@
     y = npy_floor(x);
     r = x - y;
 
-    if (r > 0.5) goto rndup;
+    if (r > 0.5) {
+        y += 1.0;
+    }
 
     /* Round to nearest even */
-    if (r==0.5) {
+    if (r == 0.5) {
         r = y - 2.0*npy_floor(0.5*y);
-        if (r==1.0) {
-        rndup:
-            y+=1.0;
+        if (r == 1.0) {
+            y += 1.0;
         }
     }
     return y;
@@ -277,21 +281,17 @@
 #endif
 
 #ifndef HAVE_EXP2
-#define LOG2 0.69314718055994530943
 double npy_exp2(double x)
 {
-    return npy_exp(LOG2*x);
+    return npy_exp(NPY_LOGE2*x);
 }
-#undef LOG2
 #endif
 
 #ifndef HAVE_LOG2
-#define INVLOG2 1.4426950408889634074
 double npy_log2(double x)
 {
-    return INVLOG2*npy_log(x);
+    return NPY_LOG2E*npy_log(x);
 }
-#undef INVLOG2
 #endif
 
 /*
@@ -370,7 +370,54 @@
 
 /**end repeat**/
 
+
 /*
+ * Decorate all the math functions which are available on the current platform
+ */
+
+/**begin repeat
+ * #type = npy_longdouble,double,float#
+ * #c = l,,f#
+ * #C = L,,F#
+ */
+/**begin repeat1
+ * #kind = sin,cos,tan,sinh,cosh,tanh,fabs,floor,ceil,rint,trunc,sqrt,log10,
+ *         log,exp,expm1,asin,acos,atan,asinh,acosh,atanh,log1p,exp2,log2#
+ * #KIND = SIN,COS,TAN,SINH,COSH,TANH,FABS,FLOOR,CEIL,RINT,TRUNC,SQRT,LOG10,
+ *         LOG,EXP,EXPM1,ASIN,ACOS,ATAN,ASINH,ACOSH,ATANH,LOG1P,EXP2,LOG2#
+ */
+#ifdef HAVE_@KIND@@C@
+@type@ npy_@kind@@c@(@type@ x)
+{
+    return @kind@@c@(x);
+}
+#endif
+
+/**end repeat1**/
+
+/**begin repeat1
+ * #kind = atan2,hypot,pow,fmod,copysign#
+ * #KIND = ATAN2,HYPOT,POW,FMOD,COPYSIGN#
+ */
+#ifdef HAVE_@KIND@@C@
+@type@ npy_@kind@@c@(@type@ x, @type@ y)
+{
+    return @kind@@c@(x, y);
+}
+#endif
+/**end repeat1**/
+
+#ifdef HAVE_MODF@C@
+@type@ npy_modf@c@(@type@ x, @type@ *iptr)
+{
+    return modf@c@(x, iptr);
+}
+#endif
+
+/**end repeat**/
+
+
+/*
  * Non standard functions
  */
 
@@ -397,24 +444,12 @@
 
 @type@ npy_log2_1p@c@(@type@ x)
 {
-    @type@ u = 1 + x;
-    if (u == 1) {
-        return LOG2E*x;
-    } else {
-        return npy_log2@c@(u) * x / (u - 1);
-    }
+    return LOG2E*npy_log1p@c@(x);
 }
 
-@type@ npy_exp2_1m@c@(@type@ x)
+@type@ npy_exp2_m1@c@(@type@ x)
 {
-    @type@ u = npy_exp@c@(x);
-    if (u == 1.0) {
-        return LOGE2*x;
-    } else if (u - 1 == -1) {
-        return -LOGE2;
-    } else {
-        return (u - 1) * x/npy_log2@c@(u);
-    }
+    return npy_expm1@c@(LOGE2*x);
 }
 
 @type@ npy_logaddexp@c@(@type@ x, @type@ y)
@@ -453,48 +488,3 @@
 #undef DEG2RAD
 
 /**end repeat**/
-
-/*
- * Decorate all the math functions which are available on the current platform
- */
-
-/**begin repeat
- * #type = npy_longdouble,double,float#
- * #c = l,,f#
- * #C = L,,F#
- */
-/**begin repeat1
- * #kind = sin,cos,tan,sinh,cosh,tanh,fabs,floor,ceil,rint,trunc,sqrt,log10,
- *         log,exp,expm1,asin,acos,atan,asinh,acosh,atanh,log1p,exp2,log2#
- * #KIND = SIN,COS,TAN,SINH,COSH,TANH,FABS,FLOOR,CEIL,RINT,TRUNC,SQRT,LOG10,
- *         LOG,EXP,EXPM1,ASIN,ACOS,ATAN,ASINH,ACOSH,ATANH,LOG1P,EXP2,LOG2#
- */
-#ifdef HAVE_@KIND@@C@
-@type@ npy_@kind@@c@(@type@ x)
-{
-    return @kind@@c@(x);
-}
-#endif
-
-/**end repeat1**/
-
-/**begin repeat1
- * #kind = atan2,hypot,pow,fmod,copysign#
- * #KIND = ATAN2,HYPOT,POW,FMOD,COPYSIGN#
- */
-#ifdef HAVE_@KIND@@C@
-@type@ npy_@kind@@c@(@type@ x, @type@ y)
-{
-    return @kind@@c@(x, y);
-}
-#endif
-/**end repeat1**/
-
-#ifdef HAVE_MODF@C@
-@type@ npy_modf@c@(@type@ x, @type@ *iptr)
-{
-    return modf@c@(x, iptr);
-}
-#endif
-
-/**end repeat**/

Modified: branches/1.4.x/numpy/core/tests/test_regression.py
===================================================================
--- branches/1.4.x/numpy/core/tests/test_regression.py	2010-04-05 04:38:49 UTC (rev 8323)
+++ branches/1.4.x/numpy/core/tests/test_regression.py	2010-04-05 04:51:33 UTC (rev 8324)
@@ -1220,7 +1220,7 @@
         x[x.nonzero()] = x.ravel()[:1]
         assert sys.getrefcount(strb) == numb
         assert sys.getrefcount(stra) == numa + 2
-        
+
     def test_signed_integer_division_overflow(self):
         """Ticket #1317."""
         def test_type(t):
@@ -1230,7 +1230,11 @@
         for t in (np.int8, np.int16, np.int32, np.int64, np.int, np.long):
             test_type(t)
 
+    def test_log1p_compiler_shenanigans(self):
+        # Check if log1p is behaving on 32 bit intel systems.
+        assert_(np.isfinite(np.log1p(np.exp2(-53))))
 
+
 if __name__ == "__main__":
     run_module_suite()
     def test_duplicate_title_and_name(self):



More information about the Numpy-svn mailing list