Using the .tree attribute of cKDTree gives you direct access to the kdtree. You can recursively navigate it like this:
from scipy.spatial import cKDTree
tree = cKDTree(points) # let points be an array of shape (n,2)
groups = 
depth = 5 # depth of the tree, ie: depth=5 results in 64 bins
def recurse(node, i=0):
g = node.greater
l = node.lesser
has about a dozen local derivative-free optimizers, including SLSQP in C.
Only COBYLA currently supports arbitrary nonlinear inequality and equality constraints;
the rest of them support bound-constrained or unconstrained problems only.
(However, any of them can be applied to nonlinearly constrained problems by combining them with the
You could compute the Voronoi diagram. Each cell in this diagram contains exactly one point and taken together the cells cover the domain. If you need, say, $n$ points in each bin, it is then a matter of picking a cell, merging it with its $n-1$ neighbors, and doing the same for the remaining cells.