Making a non-overlapping bubble chart in Matplotlib (circle packing)
The following would be a brute-force approach.
You can first place all circles on a grid, with a gridspacing as large as twice the maximum radius of any of the circles.
Then you let the circles do a random walk and check in each step if the "potential energy" of the bunch of cicles has become smaller and if the positions obtained are valid (i.e. no overlaps).
if (e < self.E and self.isvalid(i)):
As a "potential" we can simply use a square radial function.
self.p = lambda x,y: np.sum((x**2+y**2)**2)
The code:
import numpy as np
import matplotlib.pyplot as plt
# create 10 circles with different radii
r = np.random.randint(5,15, size=10)
class C():
def __init__(self,r):
self.N = len(r)
self.x = np.ones((self.N,3))
self.x[:,2] = r
maxstep = 2*self.x[:,2].max()
length = np.ceil(np.sqrt(self.N))
grid = np.arange(0,length*maxstep,maxstep)
gx,gy = np.meshgrid(grid,grid)
self.x[:,0] = gx.flatten()[:self.N]
self.x[:,1] = gy.flatten()[:self.N]
self.x[:,:2] = self.x[:,:2] - np.mean(self.x[:,:2], axis=0)
self.step = self.x[:,2].min()
self.p = lambda x,y: np.sum((x**2+y**2)**2)
self.E = self.energy()
self.iter = 1.
def minimize(self):
while self.iter < 1000*self.N:
for i in range(self.N):
rand = np.random.randn(2)*self.step/self.iter
self.x[i,:2] += rand
e = self.energy()
if (e < self.E and self.isvalid(i)):
self.E = e
self.iter = 1.
else:
self.x[i,:2] -= rand
self.iter += 1.
def energy(self):
return self.p(self.x[:,0], self.x[:,1])
def distance(self,x1,x2):
return np.sqrt((x1[0]-x2[0])**2+(x1[1]-x2[1])**2)-x1[2]-x2[2]
def isvalid(self, i):
for j in range(self.N):
if i!=j:
if self.distance(self.x[i,:], self.x[j,:]) < 0:
return False
return True
def plot(self, ax):
for i in range(self.N):
circ = plt.Circle(self.x[i,:2],self.x[i,2] )
ax.add_patch(circ)
c = C(r)
fig, ax = plt.subplots(subplot_kw=dict(aspect="equal"))
ax.axis("off")
c.minimize()
c.plot(ax)
ax.relim()
ax.autoscale_view()
plt.show()
Because of the random walk nature of this, finding the solution will take a little time (~10 seconds in this case); you may of course play with the parameters (mainly the number of steps 1000*self.N
until a solution is fixed) and see what suits your needs.