# Packing squares into a circle

I am currently trying to solve and understand a problem about packing the maximum of squares into a circle. What I am trying to do, is to calculate the maximum number of variable sized squares (with fixed side) into a circle of variable diameter. For example the result you get at WolframAlpha typing 'pack squares of side 50mm into a circle of diameter 300mm':

Densest known packing:

Furthermore, I would like to be able to add a variable for horizontal and vertical spacing between each square. I have found this to be exactly what I am trying to understand and build.

• Hi! and welcome to Mma.SE! It's not clear to me what you're asking and a more specific question might better attract answers. Is there some part of Mathematica you don't understand? Or would you like someone build such a program for you? Or do you have some other specific question? (In case you haven't seen the, here are some guidelines for asking questions.) Commented Aug 17, 2013 at 13:30
• this relates to the problem in general math.stackexchange.com/questions/115735/…. According to www2.stetson.edu/~efriedma/squincir you can fit 29 unit squares inside a circle with 6.825 diameter. wolframalpha.com/input/… wolframalpha says otherwise Commented Aug 17, 2013 at 14:28

(* i, j, -n ...0 ...n*)
(* calculate a square center with vertical and horizontal spacings*)
center[i_, j_, s_, spx_, spy_] := {i (s + spx), j (s + spy)}
(* calculate the vartices of a square given the center*)
verts[{cx_, cy_}, s_] := {cx, cy} + # s/2 & /@ Tuples[{-1, 1}, 2]
(* calculate all squares to consider, given a circle center*)
(* note: it's overestimating, but it doesn't matter *)
allSquares[cirCenterX_, cirCenterY_, circR_, s_, spx_, spy_] :=
N@center[i, j, s, spx, spy] /. Solve[(#.#) &@({cirCenterX, cirCenterY} -
center[i, j, s, spx, spy]) < (circR + s + spx + spy)^2, {i, j},
Integers]
(* calculate all verttices to consider, given a circle center*)
allVerts[cirCenterX_, cirCenterY_, circR_, s_, spx_, spy_] :=
N@verts[#, s] & /@ allSquares[cirCenterX, cirCenterY, circR, s, spx, spy]
(* test if a square is inside a given circle *)
testVertInsideCirc[cirCenterX_, cirCenterY_, circR_, vert_] :=
Norm[{cirCenterX, cirCenterY} - vert] <= circR

(* set a problem*)
s = 10; (* square side *)
r = 50; (* circle radius *)
spx = 2/10; (*horiz spacing*)
spy = 3/10; (*vert spacing*)
(* calculate all possible involved squares for a given problem*)
allvs = allVerts[0, 0, r, s, spx, spy];
(* Solve*)
nm = NMaximize[{Tr[ Boole /@  And @@@ Map[testVertInsideCirc[ccx, ccy, r, #] &, N@allvs, {2}]],
{-s/2 <= ccx <= s/2, -s/2 <= ccy <= s/2}}, {ccx, ccy},
Method -> "DifferentialEvolution"]

{59., {ccx -> -3.33338, ccy -> -0.0896481}}

Edit

Drawing our solution:

With[{cx = ccx /. nm[[2]], cy = ccy /. nm[[2]]},
Graphics[{EdgeForm[Directive[Thick, Black]],
Flatten@Transpose@{And @@@
Map[testVertInsideCirc[cx, cy, r, #] &, N@allvs, {2}] /.
{True -> Blue, False -> Yellow},
Rectangle @@@ (allvs[[All, {1, 4}]])},
Red, Thick, Circle[{cx, cy}, r], PointSize[Large], Point[{cx, cy}]}]]

• This is awesome! But I have tried playing with the spacing settings a little. As I go smaller, calculations fail.
– Alex
Commented Aug 27, 2013 at 18:49
• @Alex Could you please try it again? I corrected a minor bug. If you still have problems please post the example so I can look for the problem Commented Sep 7, 2013 at 19:14
• I seem to get calculation errors, when I use s=1,r=150,spx=spy=1/1000 for example. I am not quite sure, but it seems, that the s variable is causing this?
– Alex
Commented Sep 19, 2013 at 17:07