Another better lower bound for n=17 square packing

The idea is to improve a recent result and prove that 4.5058(?)≤s(17) using these weights:
But let’s first define s(17). Quoting the old article about the topic
Let s(n) be the side of the smallest square into which we can pack n unit squares.
For n=16, the best is obviously a 4x4 array, so s(16)=4.
For n=15, the 15 unit squares can also obviously be enclosed in a 4x4 square so s(15)≤4. Proving that it’s the smaller square is not obvious at all. Anyway, Erich Friedman proved that in 1999, so s(15)=4
For n=17, the obvious enclosing square is the 5x5, but in 1998 John Bidwell found an example that shows that a square of 4.6756… is enough, so s(17)≤4.6756… It’s a very interesting arrangement of the squares, so it’s worth visiting the collection to see it and the versions for other numbers.
On the other hand, Trevor Green proved in 2000 that 4.4452…≤s(17), (more details later). So there was a huge gap 4.4452…≤s(17)≤4.6756…
A few weeks ago, Sam Burns with ChapGPT 5.6 Sol improved(?) the lower bound. The new bound is still not reviewed by the community. I took a look and it makes a lot of sense and I think it’s correct, but I may be missing a small corner case in the proof or the accompanying program, or I may be missing a huge hole. I’ll add a small (?) to the number just in case, but I’m quite optimistic and confident it’s correct so I’ll use only a half font size. So the current bound is 4.4452…≤4.4811(?)≤s(17)≤4.6756…
My main objection to Sam Burns is that it really deserved a nice graphic! So my first step will be to add a nice graphic here. Also, making a few improvements to the program, I found a new lower bound that is 4.5058 So now we have 4.4452…≤4.4811(?)≤4.5058(?)≤s(17)≤4.6756…
My new example and the modification of the code are here, but the more technical details about finding the new bound are part of a second post.
Trevor Green’s bound
The idea of the old proof (19+40*sqrt(2))/17≅4.4452…≤s(17) of Trevor Green is to pick 16 very interesting "unavoidable" points in a square of side 4.4452… and then he uses a lot of geometry to prove that any unit square must include at least one of them. So if we try to fit 17 unit squares there, at least two unit squares must share one of the 16 interesting points. The construction chooses 16 points out of a 4x6 grid.
I only found an image of the points in the old article, but I couldn't find the analytical definition. Looking at the formula for the side of the square, and using a rule, and some guessing, I think that the empty left/right margin is 0.5 and the empty top/bottom margin is sqrt(2)-1/2≅0.9142… With these choices, the diagonal segment in the original graphic has length 1, which is a very useful number to make triangles that have vertices that are unavoidable points. (I’d be glad to hear a confirmation.)
It uses a 6x4 grid with an empty margin of 0.9142… and 0.5000, and the total size of the grid is 2.6168… and 2.4452…
To compare the construction to the newer constructions, it’s better to symmetrize it. In this symmetrized version each unit square includes at least 4 points, but some points are thicker, and they count as double points (more details later).
Sam Burns’ bound
The idea to prove 4.4811(?)≤s(17) posted by Sam Burns using ChatGPT picks 268 somewhat interesting points in a square of side 4.4811 The points have different weights, and the total weight is only 16.9476. After some reductions, it’s only necessary to test a finite number of directions and they use a program in Python to test “all” the possible “almost-unit” (actually .9973) squares and verify that the sum of weight inside each one of them is at least 1 (actually 1.0003). So if we try to fit 17 unit squares there, at least two unit squares must share at least one of the 268 somewhat interesting points. (More details in the second post.)
This method has false negatives. If it verifies a solution then it’s surely correct, but if the program fails there is a tiny chance that it’s a mistake. This is fine to ensure the weight proves a lower bound.
It’s not clear how the weights were selected. Comparing this solution to all the examples in the old article, the 0.5 margin is too narrow because most examples use ~1.0 or ~9.1 or something like that. The selections of weight agree with me, and all the weights in the first/last row/column of the grid are zero. In my handwaving opinion, the second/penultimate row/columns should be empty too, but there is a non-zero weight in (1, 11) of the grid and the symmetric images, I hope it is not necessary in a better example. The third/penpenultimate row/column is quite full. It’s closer to the border than in the old examples, so it looks like adding more points near the border may be a good idea to improve the bound.
I draw the images using Racket with the Metapict package. The radius of each circle is calculated from the weight as
r = sqrt(weight^(1/gamma)) * scale
With gamma = 1.0 the area is proportional to the weight, but the small weights are too small in the image. After some tweaking, gamma=2.0 looks nice because the smaller weights are easier to see. The scale is not so mysterious, and I should have used pi somewhere in it, but scale=0.07 looks nice in my machine. The circles are semi-transparent, so it’s possible to see when they overlap if you ever increase the scale. The code is at the bottom, and it divides the weight by 1.0003 that is the actual minimal sum.
New bound
My idea was to try different combinations of the margin and internal grid size. As I said, it’s not clear how the weights were selected in the example of Sam Burns. So for each fixed size, I decided to use linear programming to find them.
Then I used a combination of brute force search and luck to get the best grid I could find. After that, I rounded the weight so they look nice and are nice fractions. (More details in the second post.)
After a lot of time, the best I got is 4.5058(?)≤s(17). The new solution uses 168 somewhat interesting points in a square of side 4.5058 in a 29x29 grid. They sum only 16.9166… Each unit square includes at least a total weight of 1.There is an empty margin of 0.77565 and the internal grid has a total side of 3.9545.
As I said, the weights are closer to the border than what I expected looking at the old examples, close to the second/penpenultimate row/column of the previous one. It also uses fewer weights, so I hope it’s easier to prove that it’s correct without a computer. I’d like to make a non symmetric version, that may be even better.
The program published by Sam Burns assumes that the empty margin is 0.5, so I had to modify it slightly to allow arbitrary borders with a variable M that is the double of the margin. The version with that modification, the new sizes and the new table of weight is at the bottom. Running that program and making the obvious changes to the explanation posted by Sam Burns proves(?) the new bound.
Conclusion and Future Work
- The distributions look quite discrete in the corner, but it has some strange bars near the center. It would be nice to increase the grid size and take a look. Also, the narrow empty margins appear to be useful.
- My search program in the second article is too slow (like 1 hour), so I avoided changing the size of the grid. It may be useful to explore other grid sizes in case there are some interesting coincidences.
- Adding more digits takes only a few minutes, I didn't bother because it looks like refining the grid or using more directions for the rotations would make bigger changes.
- This result also automatically improves the lower bound of s(18), s(19) and s(20). But a more deep search for those values should provide even better bounds. I’ve seen too many cases where the total sum of the weight is 18. There is something interesting about 18.
- I’d like to find the non-symmetrical version. I have some ideas to try, so check again in a few days. A non-symmetrical version hopefully has like 1/8 of the weight and hopefully shows the almost equilateral triangles and is easier to understand without a computer.
You may like to read the second post with details about how I got the new weights.
Program to verify the bound
from __future__ import annotations
from bisect import bisect_left, bisect_right
from fractions import Fraction as F
import numpy as np
# Original version posted by Sam Burns 2026
# Modified by Gustavo Massaccesi 2026
# Proposed exact lower-bound certificate for packing 17 unit squares in a square.
#
# All geometric quantities and predicates are rational. NumPy is used only for
# integer range-addition and cumulative sums; no floating-point geometry is used.
L = F(45058, 10000) # side of the square
M = F(15513, 10000) # both empty borders
B = F(9973, 10000)
T = F(207107, 500000)
KMAX = 180
D = T / KMAX
WEIGHT_SCALE = 576 # min weight
NGRID = 29
LAST = NGRID - 1
# (i, j, w): every distinct D4 image of grid point (i,j) receives weight w/WEIGHT_SCALE.
CERT = [
(0, 2, 165),
(0, 11, 129),
(1, 8, 36),
(1, 10, 21),
(1, 11, 15),
(2, 2, 246),
(2, 8, 129),
(2, 9, 105),
(2, 10, 36),
(2, 11, 105),
(5, 10, 36),
(6, 10, 63),
(6, 11, 12),
(7, 10, 21),
(8, 9, 33),
(8, 11, 15),
(9, 11, 75),
(9, 14, 39),
(10, 11, 25),
(10, 12, 21),
(10, 13, 24),
(10, 14, 3),
(11, 11, 16)
]
def orbit(i: int, j: int) -> set[tuple[int, int]]:
n = LAST
return {
(i, j), (n - i, j), (i, n - j), (n - i, n - j),
(j, i), (n - j, i), (j, n - i), (n - j, n - i),
}
def build_atoms() -> list[tuple[F, F, int]]:
step = (L - M) / LAST
coord = [M / 2 + step * i for i in range(NGRID)]
by_index: dict[tuple[int, int], int] = {}
for i, j, w in CERT:
for ij in orbit(i, j):
if ij in by_index:
raise ValueError(f"duplicate orbit assignment at {ij}")
by_index[ij] = w
return [
(coord[i], coord[j], w)
for (i, j), w in sorted(by_index.items())
]
# Clip a convex rational polygon against U >= bound or U <= bound.
def clip_u(
poly: list[tuple[F, F]],
bound: F,
keep_ge: bool,
) -> list[tuple[F, F]]:
if not poly:
return []
out: list[tuple[F, F]] = []
def inside(p: tuple[F, F]) -> bool:
return p[0] >= bound if keep_ge else p[0] <= bound
prev = poly[-1]
prev_in = inside(prev)
for cur in poly:
cur_in = inside(cur)
if cur_in != prev_in:
u1, v1 = prev
u2, v2 = cur
if u2 == u1:
v = v1
else:
lam = (bound - u1) / (u2 - u1)
v = v1 + lam * (v2 - v1)
out.append((bound, v))
if cur_in:
out.append(cur)
prev, prev_in = cur, cur_in
return out
def center_domain(c: F, s: F) -> list[tuple[F, F]]:
# A B-square at orientation (c,s) lies in [0,L]^2 exactly when its
# center lies in [h,L-h]^2, with h=B(c+s)/2.
# Transform that square to the B-square's (U,V) frame.
h = B * (c + s) / 2
lo, hi = h, L - h
corners_xy = [(lo, lo), (hi, lo), (hi, hi), (lo, hi)]
return [(c * x + s * y, -s * x + c * y) for x, y in corners_xy]
def verify_orientation(
c: F,
s: F,
atoms: list[tuple[F, F, int]],
) -> int:
"""Return the exact minimum integer score for one rational orientation."""
half = B / 2
dom = center_domain(c, s)
u_dom_min = min(u for u, _ in dom)
u_dom_max = max(u for u, _ in dom)
v_dom_min = min(v for _, v in dom)
v_dom_max = max(v for _, v in dom)
rects: list[tuple[F, F, F, F, int]] = []
u_events = {u_dom_min, u_dom_max}
v_events = {v_dom_min, v_dom_max}
# In center coordinates, atom membership is an axis-aligned rectangle.
for x, y, w in atoms:
pu = c * x + s * y
pv = -s * x + c * y
u1, u2 = pu - half, pu + half
v1, v2 = pv - half, pv + half
rects.append((u1, u2, v1, v2, w))
u_events.add(u1)
u_events.add(u2)
v_events.add(v1)
v_events.add(v2)
ue = sorted(u_events)
ve = sorted(v_events)
ui = {x: i for i, x in enumerate(ue)}
vi = {x: i for i, x in enumerate(ve)}
# Exact integer 2D difference array. Scores are constant in every open
# event cell. NumPy performs only integer arithmetic here.
diff = np.zeros((len(ue), len(ve)), dtype=np.int64)
for u1, u2, v1, v2, w in rects:
a, b = ui[u1], ui[u2]
p, q = vi[v1], vi[v2]
diff[a, p] += w
diff[b, p] -= w
diff[a, q] -= w
diff[b, q] += w
scores = diff.cumsum(axis=0).cumsum(axis=1)
nu, nv = len(ue) - 1, len(ve) - 1
best = 10**18
for i in range(nu):
u0, u1 = ue[i], ue[i + 1]
if u1 <= u_dom_min or u0 >= u_dom_max:
continue
slab = clip_u(dom, u0, True)
slab = clip_u(slab, u1, False)
if not slab:
continue
vlo = min(v for _, v in slab)
vhi = max(v for _, v in slab)
if vhi <= vlo:
continue
# This may examine a superset of feasible event cells, which is
# conservative for a lower-bound verification.
j0 = max(0, bisect_right(ve, vlo) - 1)
j1 = min(nv - 1, bisect_left(ve, vhi) - 1)
if j0 <= j1:
row_min = int(scores[i, j0:j1 + 1].min())
best = min(best, row_min)
if best == 10**18:
raise RuntimeError("center domain was not enumerated")
return best
def angle_net() -> list[tuple[F, F]]:
out: list[tuple[F, F]] = []
for k in range(KMAX + 1):
t = T * k / KMAX
den = 1 + t * t
c = (1 - t * t) / den
s = 2 * t / den
assert c * c + s * s == 1
out.append((c, s))
# The final adjacent pair brackets pi/4.
assert out[-2][1] < out[-2][0]
assert out[-1][1] >= out[-1][0]
# If psi_k=2 arctan(t_k), half an adjacent angular gap is
# arctan(t_{k+1})-arctan(t_k), whose tangent is
# D/(1+t_k*t_{k+1}) <= D. Therefore every angle in [0,pi/4]
# is within an error epsilon < D of a net direction.
for k in range(KMAX):
t0 = T * k / KMAX
t1 = T * (k + 1) / KMAX
tan_half_gap = (t1 - t0) / (1 + t0 * t1)
assert tan_half_gap <= D
return out
def main() -> None:
atoms = build_atoms()
total = sum(w for _, _, w in atoms)
print(f"atoms = {len(atoms)}")
print(
f"total_weight = {total}/{WEIGHT_SCALE}"
f" = {total / WEIGHT_SCALE:.4f}"
)
#assert len(atoms) == 268
#assert total == 169476
assert total < 17 * WEIGHT_SCALE
net = angle_net()
# For an orientation error epsilon <= D,
# cos(epsilon)+sin(epsilon) <= 1+epsilon <= 1+D.
contain = B * (1 + D)
print(f"angle_net_size = {len(net)}")
print(f"b*(1+d) = {contain} = {float(contain):.12f} < 1")
assert contain < 1
global_min = 10**18
argmin = -1
for k, (c, s) in enumerate(net):
m = verify_orientation(c, s, atoms)
if m < global_min:
global_min, argmin = m, k
if k % 30 == 0 or k == KMAX:
print(
f"orientation {k:3d}/{KMAX}: "
f"min={m}/{WEIGHT_SCALE}, "
f"global={global_min}/{WEIGHT_SCALE}"
)
print(
f"minimum_score = {global_min}/{WEIGHT_SCALE}"
f" = {global_min / WEIGHT_SCALE:.4f} at k={argmin}"
)
assert global_min >= WEIGHT_SCALE
print("CERTIFICATE CONDITIONS VERIFIED.")
print(f"By the scaling argument: s(17) >= {L} = {L:.4f}.")
if __name__ == "__main__":
main() Program to draw the images
#lang racket
(require racket/list)
(require metapict)
{define-syntax-rule (for/append clauses body ...)
; Todo: Add support for #:breack and #:final
(append* (for/list clauses (begin body ...)))}
{define (mirror-x N atoms)
(for/append ([a (in-list atoms)])
(match a
[(list x y w)
(list (list x y w) (list (- N 1 x) y w))]))}
{define (mirror-y N atoms)
(for/append ([a (in-list atoms)])
(match a
[(list x y w)
(list (list x y w) (list x (- N 1 y) w))]))}
{define (mirror-d atoms)
(for/append ([a (in-list atoms)])
(match a
[(list x y w)
(list (list x y w) (list y x w))]))}
{define-values (L-Green atoms-Green)
(let ()
; Todo: Confirm these are the correct lenghts.
(define grid-nx 6)
(define grid-ny 4)
(define border-size-y (- (sqrt 2) 1/2))
(define grid-size-y (/ (+ 12 (sqrt 8)) 17))
(define border-size-x 1)
(define L
(+ (* border-size-y 2) (* grid-size-y 3)))
(define grid-size-x (/ (- L (* border-size-x 2)) 5))
{define atoms/int '(#;()
(0 3 1) (1 3 1) (3 3 1) (5 3 1)
(0 2 1) (2 2 1) (4 2 1) (5 2 1)
(0 1 1) (1 1 1) (3 1 1) (5 1 1)
(0 0 1) (2 0 1) (4 0 1) (5 0 1))}
(define min-weight 1)
{define atoms (for/list ([a (in-list atoms/int)])
(match a
[(list x y w)
(list (+ border-size-x (* grid-size-x x))
(+ border-size-y (* grid-size-y y))
(/ w min-weight))]))}
(values L atoms))}
{define-values (L-Green/S atoms-Green/S)
(let ()
; Todo: Confirm these are the correct lenghts.
(define grid-nx 6)
(define grid-ny 4)
(define border-size-y (- (sqrt 2) 1/2))
(define grid-size-y (/ (+ 12 (sqrt 8)) 17))
(define border-size-x 1)
(define L
(+ (* border-size-y 2) (* grid-size-y (- grid-ny 1))))
(define grid-size-x (/ (- L (* border-size-x 2)) (- grid-nx 1)))
; It's easier to calculate the overlaps by hand
(define atoms/gen '(#;()
(0 1 2) (1 1 1) (2 1 1)
(0 0 2) (1 0 1) (2 0 1)))
(define min-weight 4)
(define atoms/int (remove-duplicates
(mirror-x grid-nx
(mirror-y grid-ny
atoms/gen))))
(define atoms/one-dir (for/list ([a (in-list atoms/int)])
(match a
[(list x y w)
(list (+ border-size-x (* grid-size-x x))
(+ border-size-y (* grid-size-y y))
(/ w min-weight))])))
(define atoms (mirror-d atoms/one-dir))
(values L atoms))}
{define-values (L-Burns atoms-Burns)
(let ()
(define grid-n 29)
(define L 44811/10000)
(define M 1)
(define grid-size (/ (- L M) grid-n))
(define border-size (/ M 2))
(define min-weight 10003)
{define atoms/gen '(#;()
(1 11 107) (2 4 137) (2 9 214) (2 11 107) (2 12 137)
(3 4 3884) (3 7 214) (3 8 913) (3 9 214)
(3 10 214) (3 11 1234) (3 12 2189) (3 14 384)
(4 4 1961) (4 7 520) (4 8 214) (4 9 1413) (4 10 1234)
(4 11 1083) (4 13 137) (4 14 292)
(7 11 529) (7 12 33) (8 10 906) (8 11 384) (8 12 351)
(9 9 340) (9 10 180) (9 11 204) (9 12 549)
(10 12 879) (10 13 201) (10 14 378)
(11 11 396) (11 12 622) (11 13 204) (11 14 204))}
(define atoms/int (remove-duplicates
(mirror-x grid-n
(mirror-y grid-n
(mirror-d
atoms/gen)))))
(define atoms (for/list ([a (in-list atoms/int)])
(match a
[(list x y w)
(list (+ border-size (* grid-size x))
(+ border-size (* grid-size y))
(/ w min-weight))])))
(values L atoms))}
{define-values (L-Massaccesi atoms-Massaccesi)
(let ()
(define grid-n 29)
(define L 45058/10000)
(define M 15513/10000)
(define grid-size (/ (- L M) grid-n))
(define border-size (/ M 2))
(define min-weight 576)
{define atoms/gen '(#;()
(0 2 165) (0 11 129) (1 8 36) (1 10 21) (1 11 15)
(2 2 246) (2 8 129) (2 9 105) (2 10 36) (2 11 105)
(5 10 36) (6 10 63) (6 11 12) (7 10 21)
(8 9 33) (8 11 15) (9 11 75) (9 14 39)
(10 11 25) (10 12 21) (10 13 24) (10 14 3) (11 11 16))}
(define atoms/int (remove-duplicates
(mirror-x grid-n
(mirror-y grid-n
(mirror-d
atoms/gen)))))
(define atoms (for/list ([a (in-list atoms/int)])
(match a
[(list x y w)
(list (+ border-size (* grid-size x))
(+ border-size (* grid-size y))
(/ w min-weight))])))
(values L atoms))}
{define (draw-example L atoms #:gamma [gamma 2.0] #:scale [scale 0.07])
[with-window (window -.1 (+ L .1) -.1 (+ L .1))
(define big-fill-color "whitesmoke")
(define big-border-color "black")
(define dots-color (change-alpha "darkred" 0.75))
(define big-square (curve (pt 0 0) -- (pt 0 L) -- (pt L L) -- (pt L 0) -- cycle))
(draw (color big-fill-color (fill big-square))
(penscale .1 (color big-border-color (draw big-square)))
(draw* (for/list ([a (in-list atoms)])
(match a
[(list x y w)
(define s (* (sqrt (expt w (/ 1. gamma))) scale))
(penstyle 'transparent (color dots-color (filldraw (circle (pt x y) s))))]))))
]}
(scale 4 (draw-example L-Green atoms-Green #:gamma 2.0 #:scale .07))
(scale 4 (draw-example L-Green/S atoms-Green/S #:gamma 2.0 #:scale .07))
(scale 4 (draw-example L-Burns atoms-Burns #:gamma 2.0 #:scale .07))
(scale 4 (draw-example L-Massaccesi atoms-Massaccesi #:gamma 2.0 #:scale .07))