Return to Answer

np.random.seed(42)
c = np.random.rand(10000, 2) * 800
dists = squareform(pdist(c, metric = 'chebyshev'))     # distance matrix, manhattan here since you seem to want squares
indices = np.arange(c.shape[0])  # indices that haven't been dropped (all to start)
out = [0]                        # always want the first index
while True:
    try:
        indices = indices[dists[indices[0], indices] > 400] #drop indices that are inside threshhold
        out.append(indices[0])   # add the next index that hasn't been dropped to the output
    except:
        break   # once you run out of indices, you'll get an IndexError and you're done
print(out, pdist(c[out], metric = 'chebyshev'))
[0, 2, 6, 17] [635.77582886 590.70015659 472.87353138 541.13920029 647.69071411
 476.84658995]

np.random.seed(42)
c = np.random.rand(10000, 2) * 800
dists = squareform(pdist(c, metric = 'chebyshev'))     # distance matrix, checbyshev here since you seem to want squares
indices = np.arange(c.shape[0])  # indices that haven't been dropped (all to start)
out = [0]                        # always want the first index
while True:
    try:
        indices = indices[dists[indices[0], indices] > 400] #drop indices that are inside threshhold
        out.append(indices[0])   # add the next index that hasn't been dropped to the output
    except:
        break   # once you run out of indices, you'll get an IndexError and you're done
print(out, pdist(c[out], metric = 'chebyshev'))
[0, 2, 6, 17] [635.77582886 590.70015659 472.87353138 541.13920029 647.69071411
 476.84658995]

c = np.random.rand(10000, 2) * 800
dists = squareform(pdist(c, metric = 'chebyshev'))     # distance matrix, manhattan here since you seem to want squares
indices = np.arange(c.shape[0])  # indices that haven't been dropped (all to start)
out = [0]                        # always want the first index
while True:
    try:
        indices = indices[dists[indices[0], indices] > 400] #drop indices that are inside threshhold
        out.append(indices[0])   # add the next index that hasn't been dropped to the output
    except:
        break   # once you run out of indices, you'll get an IndexError and you're done
print(out, pdist(c[out], metric = 'chebyshev'))
[0, 1, 3, 50] [423.75088712 557.53552808 475.17691811 471.15111887 403.16282495
 781.35464998]

So, 4 points (makes sense since 4 400x400 blocks tile a 800x800 space with 4 tiles), mostly low values (50 << 10000) and distance between kept points is always > 400

np.random.seed(42)
c = np.random.rand(10000, 2) * 800
dists = squareform(pdist(c, metric = 'chebyshev'))     # distance matrix, manhattan here since you seem to want squares
indices = np.arange(c.shape[0])  # indices that haven't been dropped (all to start)
out = [0]                        # always want the first index
while True:
    try:
        indices = indices[dists[indices[0], indices] > 400] #drop indices that are inside threshhold
        out.append(indices[0])   # add the next index that hasn't been dropped to the output
    except:
        break   # once you run out of indices, you'll get an IndexError and you're done
print(out, pdist(c[out], metric = 'chebyshev'))
[0, 2, 6, 17] [635.77582886 590.70015659 472.87353138 541.13920029 647.69071411
 476.84658995]

So, 4 points (makes sense since 4 400x400 blocks tile a 800x800 space with 4 tiles), mostly low values (17 << 10000) and distance between kept points is always > 400

let's try with a whole bunch of points:

c = np.random.rand(10000, 2) * 800
dists = squareform(pdist(c, metric = 'chebyshev'))     # distance matrix, manhattan here since you seem to want squares
indices = np.arange(c.shape[0])  # indices that haven't been dropped (all to start)
out = [0]                        # always want the first index
while True:
    try:
        indices = indices[dists[indices[0], indices] > 400] #drop indices that are inside threshhold
        out.append(indices[0])   # add the next index that hasn't been dropped to the output
    except:
        break   # once you run out of indices, you'll get an IndexError and you're done
print(out, pdist(c[out], metric = 'chebyshev'))
[0, 1, 3, 50] [423.75088712 557.53552808 475.17691811 471.15111887 403.16282495
 781.35464998]

So, 4 points (makes sense since 4 400x400 blocks tile a 800x800 space with 4 tiles), mostly low values (50 << 10000) and distance between kept points is always > 400