trjhtr

This is a follow up to a previously asked Finding the closest pair of points divide-and-conquer question . The original code changed significantly based on answers from @AJNeufeld therefore I created a new question per codereview guidelines for additional help.

After incorporating suggestions from answers to the original question linked above, the run-time improved by 30% which is significant. Nevertheless, the code is still slow and I think there is room for further improvements. My hunch is that find_min_distance_in_rec is slowing the code down. Below is the current code. Again, the code has been stress tested so I am confident that it's correct, however, it's slow.

#Uses python3
import math
import statistics as stats

# helper functions:
def two_point_distance(p0,p1):
 # returns distance between two (x,y) pairs
 return math.sqrt( ((p0[0]-p1[0])*(p0[0]-p1[0])) + 
 ((p0[1] - p1[1])*(p0[1] - p1[1])) )

def combine_xy(x_arr,y_arr):
 # combine x_arr and y_arr to combined list of (x,y) tuples 
 return list(zip(x_arr,y_arr))

def find_closest_distance_brute(xy_arr):
 # brute force approach to find closest distance 
 dmin = math.inf
 for i, pnt_i in enumerate(xy_arr[:-1]): 
 dis_storage_min = min( two_point_distance(pnt_i, pnt_j) for pnt_j in xy_arr[i+1:]) 
 if dis_storage_min < dmin:
 dmin = dis_storage_min 
 return dmin

def calc_median_x(xy_arr):
 # return median of x values in list of (x,y) points
 return stats.median( val[0] for val in xy_arr )

def filter_set(xy_arr_y_sorted, median, distance):
# filter initial set such than |x-median|<= distance
 return [ val for val in xy_arr_y_sorted if abs(val[0] - median) <= distance ]

def x_sort(xy_arr):
 # sort array according to x value
 return sorted(xy_arr, key=lambda val: val[0])

def y_sort(xy_arr):
 # sort array according to y value
 return sorted(xy_arr, key=lambda val: val[1])


def split_array(arr_x_sorted, arr_y_sorted,median):
 # split array of size n to two arrays of n/2
 # input is the same array twice, one sorted wrt x, the other wrt y
 leq_arr_x_sorted = [ val for val in arr_x_sorted if val[0] < median ]
 geq_arr_x_sorted = [ val for val in arr_x_sorted if val[0] > median ]
 eq_arr_x = [ val for val in arr_x_sorted if val[0] == median ]

 n = len(eq_arr_x)//2
 leq_arr_x_sorted = leq_arr_x_sorted + eq_arr_x[:n]
 geq_arr_x_sorted = eq_arr_x[n:] + geq_arr_x_sorted

 leq_arr_y_sorted = [ val for val in arr_y_sorted if val[0] < median ]
 geq_arr_y_sorted = [ val for val in arr_y_sorted if val[0] > median ]
 eq_arr_y = [ val for val in arr_y_sorted if val[0] == median ]

 n = len(eq_arr_y)//2
 leq_arr_y_sorted = leq_arr_y_sorted + eq_arr_y[:n]
 geq_arr_y_sorted = eq_arr_y[n:] + geq_arr_y_sorted

 return leq_arr_x_sorted, leq_arr_y_sorted, geq_arr_x_sorted, geq_arr_y_sorted

def find_min_distance_in_rec(xy_arr_y_sorted,dmin):
 # takes in array sorted in y, and minimum distance of n/2 halves
 # for each point it computes distance to 7 subsequent points
 # output min distance encountered

 dmin_rec = dmin

 if len(xy_arr_y_sorted) == 1:
 return math.inf

 if len(xy_arr_y_sorted) > 7: 
 for i, pnt_i in enumerate(xy_arr_y_sorted[:-7]):
 dis_storage_min = min(two_point_distance(pnt_i, pnt_j) 
 for pnt_j in xy_arr_y_sorted[i+1:i+1+7])
 if dis_storage_min < dmin_rec:
 dmin_rec = dis_storage_min

 dis_storage_min = find_closest_distance_brute(xy_arr_y_sorted[-7:])
 if dis_storage_min < dmin_rec:
 dmin_rec = dis_storage_min
 else:
 for k, pnt_k in enumerate(xy_arr_y_sorted[:-1]): 
 dis_storage_min = min( two_point_distance(pnt_k, pnt_l) 
 for pnt_l in xy_arr_y_sorted[k+1:]) 
 if dis_storage_min < dmin_rec:
 dmin_rec = dis_storage_min 

 return dmin_rec 

def find_closest_distance_recur(xy_arr_x_sorted, xy_arr_y_sorted):
 # recursive function to find closest distance between points
 if len(xy_arr_x_sorted) <=3 :
 return find_closest_distance_brute(xy_arr_x_sorted)

 median = calc_median_x(xy_arr_x_sorted)
 leq_arr_x_sorted, leq_arr_y_sorted , grt_arr_x_sorted, grt_arr_y_sorted = split_array(xy_arr_x_sorted, xy_arr_y_sorted, median)

 distance_left = find_closest_distance_recur(leq_arr_x_sorted, leq_arr_y_sorted)
 distance_right = find_closest_distance_recur(grt_arr_x_sorted, grt_arr_y_sorted)
 distance_min = min(distance_left, distance_right)

 filt_out = filter_set(xy_arr_y_sorted, median, distance_min)
 distance_filt = find_min_distance_in_rec(filt_out, distance_min)

 return min(distance_min, distance_filt)

def find_closest_point(x_arr, y_arr):
 # input is x,y points in two arrays, all x's in x_arr, all y's in y_arr
 xy_arr = combine_xy(x_arr,y_arr)
 xy_arr_x_sorted = x_sort(xy_arr)
 xy_arr_y_sored = y_sort(xy_arr)

 min_distance = find_closest_distance_recur(xy_arr_x_sorted, xy_arr_y_sored)

 return min_distance

The original code changed significantly based on answers

I think that the plural is inappropriate there: the changes seem to be solely based on AJNeufeld's answer. Not a single one of the issues I pointed out has been addressed, not even the bug:

1. My understanding is that the end result should be four arrays such that leq_arr_x_sorted and leq_arr_y_sorted contain the same points
   sorted on different projections. But the dmy_x / dmy_y separation
   of points whose x coordinate is median doesn't seem to guarantee
   that those points will be sent to the same half in the x-projection
   and the y-projection. In order to guarantee that an additional
   pre-condition is needed: most obviously that in arr_x_sorted ties
   are broken by y; but that pre-condition is not produced by x_sort.

This bug can be demonstrated easily by adding

 print(leq_arr_x_sorted, leq_arr_y_sorted, geq_arr_x_sorted, geq_arr_y_sorted)

before the return of split_array and calling with the simple test case find_closest_point([1, 2, 2, 3], [4, 3, 2, 1]).

I guess I didn't make this clear in my answer on the original question.

The minimum of the square-root of the square distances is the square-root of the minimum of the square distances. You can save time by calculating the square-root only when needed.

Replace this (and all calls to it):

def two_point_distance(p0,p1):
 # returns distance between two (x,y) pairs
 return math.sqrt( ((p0[0]-p1[0])*(p0[0]-p1[0])) + 
 ((p0[1] - p1[1])*(p0[1] - p1[1])) )

with this:

def square_distance(p0, p1):
 dx = p0[0] - p1[0]
 dy = p0[1] - p1[1]
 return dx * dx + dy * dy

For example, find_closest_distance_brute() becomes:

def find_closest_distance_brute(xy_arr):
 # brute force approach to find closest distance
 dist_sqr_min = math.inf
 for i, pnt_i in enumerate(xy_arr[:-1]): 
 d_sqr_min = min( square_distance(pnt_i, pnt_j) for pnt_j in xy_arr[i+1:]) 
 if d_sqr_min < dist_sqr_min:
 dist_sqr_min = d_sqr_min

 return math.sqrt(dist_sqr_min) # Only calculate square-root of final value

Or more pythonically:

def find_closest_distance_brute(xy_arr):
 # brute force approach to find closest distance
 dist_sqr_min = min((square_distance(pnt_i, pnt_j) for i, pnt_i in enumerate(xy_arr[:-1])
 for pnt_j in xy_arr[i+1:]),
 default=math.inf)
 return math.sqrt(dist_sqr_min) # Only calculate square-root of final value

Regarding find_min_distance_in_rect():

def find_min_distance_in_rec(xy_arr_y_sorted,dmin):

 dmin_rec = dmin

 if len(xy_arr_y_sorted) == 1:
 return math.inf

What happens if len(xy_arr_y_sorted) == 0? A better test would be <= 1.

if len(xy_arr_y_sorted) > 7: 
 # ... complicated code involving lots of 7's.
else:
 for k, pnt_k in enumerate(xy_arr_y_sorted[:-1]): 
 dis_storage_min = min( two_point_distance(pnt_k, pnt_l) 
 for pnt_l in xy_arr_y_sorted[k+1:]) 
 if dis_storage_min < dmin_rec:
 dmin_rec = dis_storage_min 

This else: code looks exactly like find_closest_distance_brute() code. In fact, you can replace it with a call to that function.

if len(xy_arr_y_sorted) > 7: 
 # ... complicated code involving lots of 7's.
else:
 dmin_rec = find_closest_distance_brute(xy_arr_y_sorted) 

I don't like # ... complicated code involving lots of 7's.. I'll rework it in a later edit.

See Peter Taylor's answer ... especially the part about median = calc_median_x(xy_arr_x_sorted), for an O(n) to O(1) improvement.

As pointed out a second time by Peter Taylor, your code has a bug; xy_arr_x_sorted and xy_arr_y_sorted are supposed to have the same content, but due to the bug, may not. He has suggested a change to fix the bug. I have a different approach: eliminate xy_arr_x_sorted entirely.

In find_closest_distance_recur(), you have:

 if len(xy_arr_x_sorted) <=3 :
 return find_closest_distance_brute(xy_arr_x_sorted)

Since the arrays are supposed to contain the same content, you could replace this with:

 if len(xy_arr_y_sorted) <= 3:
 return find_closest_distance_brute(xy_arr_y_sorted)

Once you do that, you no longer have any useful references to the [1] members of the contents of xy_arr_x_sorted. Instead, you could replace it with a sorted list of just the x-coordinates.

Once that is done, the partitioning code should be revisited. There is no difference between duplicate median values, so a simple slice operation can divide the x-coordinate list in two with zero copying required.