[Scipy-svn] r6813 - trunk/scipy/spatial

[scipy-svn at scipy.org]

Author: ptvirtan
Date: 2010-09-22 04:00:38 -0500 (Wed, 22 Sep 2010)
New Revision: 6813

Modified:
   trunk/scipy/spatial/qhull.pyx
Log:
BUG: spatial/qhull: ensure that _find_simplex always terminates, even when compiled with optimizations on

Modified: trunk/scipy/spatial/qhull.pyx
===================================================================
--- trunk/scipy/spatial/qhull.pyx	2010-09-18 04:58:40 UTC (rev 6812)
+++ trunk/scipy/spatial/qhull.pyx	2010-09-22 09:00:38 UTC (rev 6813)
@@ -27,6 +27,9 @@
     extern void *stderr
     extern void *stdout
 
+cdef extern from "math.h":
+    double fabs(double x) nogil
+
 cdef extern from "qhull/src/qset.h":
     ctypedef union setelemT:
         void *p
@@ -832,7 +835,14 @@
             if ineigh == -1:
                 continue
             dist = _distplane(d, ineigh, z)
-            if dist > best_dist:
+
+            # Note addition of eps -- otherwise, this code does not
+            # necessarily terminate! The compiler may use extended
+            # accuracy of the FPU so that (dist > best_dist), but
+            # after storing to double size, dist == best_dist,
+            # resulting to non-terminating loop
+
+            if dist > best_dist + eps*(1 + fabs(best_dist)):
                 # Note: this is intentional: we jump in the middle of the cycle,
                 #       and continue the cycle from the next k.
                 #