something about the changes I made (abs -> np.abs) made my stuff succeed more often

lib/baker/tools.py

@@ -71,12 +71,12 @@ def exact_func_3D(X):
 def improved_answer(answer, exact, verbose=False):
   if not answer['error']:
     return True
-  smblog.debug('exact: %s' % exact)
-  smblog.debug('qlin:  %s' % answer['qlin'])
   smblog.debug('error: %s' % answer['error'])
+  smblog.debug('qlin:  %s' % answer['qlin'])
   smblog.debug('final: %s' % answer['final'])
+  smblog.debug('exact: %s' % exact)
 
-  if abs(answer['final'] - exact) <= abs(answer['qlin'] - exact):
+  if np.abs(answer['final'] - exact) <= np.abs(answer['qlin'] - exact):
     smblog.debug(":) improved result")
     return True
   else:

lib/grid/DD.py

@@ -66,14 +66,14 @@ class random_grid(rect_grid):
     smblog.debug("appx_side_res: %d" % appx_side_res)
     delta = 1.0 / float(appx_side_res)
 
-    for x in xrange(appx_side_res):
+    for x in xrange(appx_side_res + 1):
       cur_x = x * delta
       for cur_y in (0, 1):
         new_point = [cur_x, cur_y]
         points.append(new_point)
         q.append(exact_func(new_point))
 
-    for y in xrange(appx_side_res):
+    for y in xrange(appx_side_res + 1):
       cur_y = y * delta
       for cur_x in (0, 1):
         new_point = [cur_x, cur_y]

lib/grid/__init__.py

@@ -1,6 +1,7 @@
 import sys
 import re
 from   collections import defaultdict
+import inspect
 
 import numpy as np
 import scipy.spatial
@@ -44,13 +45,18 @@ class grid(object):
     self.facets_for_point = defaultdict(list)
 
   def create_mesh(self, indicies):
+    """
+      this function takes a list of indicies, and then creates
+      and returns a grid object (collection of points and q).
+
+      note: the input is indicies, the grid contains points
+    """
     p = [self.points[i] for i in indicies]
     q = [self.q[i]      for i in indicies]
     return grid(p, q)
 
   def get_containing_simplex(self, X):
     if not self.faces:
-      smblog.debug('setting up connectivity')
       self.construct_connectivity()
 
     # get closest point
@@ -67,16 +73,21 @@ class grid(object):
     facets_to_check = self.facets_for_point[closest_point]
 
     attempts = 0
-    while not simplex:
+    while not simplex and facets_to_check:
       attempts += 1
       # if attempts > 20:
       #   raise smberror("probably recursing to many times")
-      cur_facet = facets_to_check.pop()
+      cur_facet = facets_to_check.pop(0)
       checked_facets.append(cur_facet)
-      facets_to_check.extend([i for i in cur_facet.neighbors if i not in checked_facets])
 
       if cur_facet.contains(X, self):
         simplex = cur_facet
+        continue
+
+      new_facest = []
+      for neighbor in cur_facet.neighbors:
+        if (neighbor not in checked_facets) and (neighbor not in facets_to_check):
+          facets_to_check.append(neighbor)
 
     if not simplex:
       raise AssertionError('no containing simplex found')
@@ -96,12 +107,31 @@ class grid(object):
 
       S is S_j from baker's paper : some points from all point that are not the simplex
     """
+    smblog.debug(inspect.stack()[1][3])
+    smblog.debug("extra points: %d" % extra_points)
+    smblog.debug("simplex size: %d" % simplex_size)
+
+    r_mesh = self.get_containing_simplex(X)
+    # smblog.debug("R:\n%s" % r_mesh)
+
+    # and some UNIQUE extra points
     (dist, indicies) = self.tree.query(X, simplex_size + extra_points)
 
-    # r_mesh = self.create_mesh(indicies[:simplex_size])
+    unique_indicies = []
-    r_mesh = self.get_containing_simplex(X)
+    for index in indicies:
-    # and some extra points
+      if self.points[index] not in r_mesh.points:
-    s_mesh = self.create_mesh(indicies[simplex_size:])
+        unique_indicies.append(index)
+
+    smblog.debug("indicies: %s" % ",".join([str(i) for i in indicies]))
+    smblog.debug("indicies: %s" % ",".join([str(i) for i in unique_indicies]))
+    s_mesh = self.create_mesh(unique_indicies)# indicies[simplex_size:])
+
+    # TODO: eventually remove this test:
+    for point in s_mesh.points:
+      if point in r_mesh.points:
+        smblog.error("ERROR")
+        smblog.error("\n%s\nin\n%s" % (point, r_mesh))
+        raise smberror("repeating point S and R")
 
     return (r_mesh, s_mesh)
 
@@ -136,8 +166,8 @@ class grid(object):
     if not simplex:
       raise AssertionError('no containing simplex found')
 
-    R = self.get_containing_simplex(X)# self.create_mesh(simplex.verts)
+    # self.create_mesh(simplex.verts)
-
+    R = self.get_containing_simplex(X)
 
     s = []
     for c,i in enumerate(simplex.neighbors):
@@ -169,7 +199,7 @@ class grid(object):
 
       this is part of the __init__ for a rect_grid, but can be called from any grid object
     """
-    smblog.debug()
+    smblog.debug('start')
     qdelaunay_string = get_qdelaunay_dump_str(self)
     facet_to_facets = []
     for matcher in grid.facet_re.finditer(qdelaunay_string):
@@ -194,13 +224,8 @@ class grid(object):
       if rel[1] in self.faces:
         self.faces[rel[0]].add_neighbor(self.faces[rel[1]])
 
-    # for matcher in grid.point_re.finditer(qdelaunay_string):
+    smblog.debug('end')
-    #   d = matcher.groupdict()
 
-    #   point               = d['point']
-    #   neighboring_facets  = d['neigh']
-
-    #   self.facets_for_point[int(point[1:])] = [i for i in neighboring_facets.split() if i in self.faces]
 
   def __str__(self):
     r = ''

setup.py

@@ -12,4 +12,6 @@ setup(
                             'numpy',
                             'scipy',
                            ],
+
+        entry_points = {'console_scripts': ['booze = baker.booze:main',]},
      )

tools/mpl_plot.py

@@ -0,0 +1,14 @@
+#!/usr/bin/env python
+
+good = [82, 164, 247, 328, 411, 491, 565, 649, 729, 808]
+bad = [18, 36, 53, 72, 89, 109, 135, 151, 171, 192]
+import pylab
+
+pylab.xlabel('total runs')
+pylab.ylabel('count')
+pylab.grid(True)
+pylab.plot(bad)
+pylab.plot(good)
+pylab.legend(('bad', 'good'))
+
+pylab.show()