---
jupyter:
  description: 'Circle packing in a square with planex2d: equal circles spread by SciPy, checked as exact discs inside a
    unit square without overlap, and their radii against the proven optimal packings for 2 to 10 circles.'
---

# Circle packing in a square

How large can $n$ equal circles be when they must fit into a unit square without overlapping? The question is old and
hard. The best packings are proven optimal for every $n$ up to 33, and beyond that the best known ones are collected,
as on Eckard Specht's Packomania. For two circles, the answer lies on the diagonal, $r = 1/(2 + \sqrt{2})$,
four sit in a grid of $r = 1/4$, and five put the fifth in the middle, $r = (\sqrt{2} - 1)/2$.

Here, SciPy searches for the packings, and planex2d checks them with exact discs: the discs must lie inside the square,
and their union must have the area of $n$ separate discs, so none of them overlaps another. Their radii are compared
with the proven ones.

![The best packings of 5 to 10 equal circles in a square](../_static/gallery/circle-packing.png){.result}

```python
import math

import numpy as np
import planex as px
from scipy.optimize import minimize
```

```python nbsphinx="hidden"
from docs_figures import Row
```

## The search

Packing circles of radius $r$ in the unit square is the same as spreading $n$ points in a square of side $1 - 2r$ as
far apart as possible: points at least $d$ apart in the unit square give circles of radius $r = d / (2(1 + d))$ once
the square is scaled down by $1 - 2r$. `spread` maximizes the smallest distance $d$ with SciPy's SLSQP from 100 random
starts, since the problem has many local optima, and returns the circles' radius and centers:

```python
def spread(n, starts=100, seed=1):
    rng = np.random.default_rng(seed)
    apart = [{"type": "ineq", "fun": lambda x, i=i, j=j: np.sum((x[2 * i:2 * i + 2] - x[2 * j:2 * j + 2]) ** 2) - x[-1] ** 2}
             for i in range(n) for j in range(i + 1, n)]
    best = None
    for _ in range(starts):
        found = minimize(lambda x: -x[-1], np.r_[rng.random(2 * n), 0.1], method="SLSQP", constraints=apart,
                         bounds=[(0, 1)] * (2 * n) + [(0, 2)], options={"ftol": 1e-15, "maxiter": 500})
        if found.success and (best is None or found.x[-1] > best.x[-1]):
            best = found
    d = best.x[-1]
    r = d / (2 * (1 + d))
    return r, r + best.x[:-1].reshape(n, 2) * (1 - 2 * r)
```

## The checks

For each $n$, the discs are drawn exactly, a hair smaller than the radius found so that circles that touch stay apart.
Nothing of them may lie outside the square, and their union must have the area $n\pi r^2$. The radii agree with the
proven optimal ones, as Packomania lists them, to the twelve digits printed:

```python
proven = {2: 0.292893218813, 3: 0.254333095030, 4: 0.25, 5: 0.207106781187, 6: 0.187680601147,
          7: 0.174457630187, 8: 0.170540688701, 9: 0.166666666667, 10: 0.148204322565}
square = px.QuickSketch().Rectangle(1, 1).Face()
packings = {}
for n, best in proven.items():
    r, centers = spread(n)
    discs = px.Unify([px.Circle(tuple(c), r * (1 - 1e-12)).Face() for c in centers])
    packings[n] = discs
    print("n", n, "| radius", round(r, 12), "| proven", best, "| outside", (discs - square).mass,
          "| union against n π r²", discs.mass / (n * math.pi * (r * (1 - 1e-12)) ** 2) - 1)
```

```python nbsphinx="hide-input"
Row(*[(px.Compound([square.outer, packings[n]]), f"{n} circles") for n in (5, 7, 10)],
    caption="The best packings of 5, 7 and 10 circles")
```

```python nbsphinx="hidden"
import sys
sys.path.insert(0, "..")
import matplotlib.pyplot as plt
from docs_figures import FACE, INK
from gallery import save_figure
fig, axes = plt.subplots(2, 3, figsize=(6, 4))
for ax, n in zip(axes.ravel(), range(5, 11)):
    square.Plot(ax=ax, face_col=FACE, line_col=INK, legend=False)
    packings[n].Plot(ax=ax, face_col="#9fd3e6", line_col=INK, arc_col=INK, spline_col=INK, legend=False)
    ax.set_aspect("equal")
    ax.axis("off")
save_figure(fig, "circle-packing", bare=True)
```