Knapsack + packing: a difficult MIQCP

Packing Circles in a Rectangle

The problem of placing the maximum number of circles (or spheres) with radius \(r\) in a container (e.g. rectangle or box) such that there is no overlap is a well-known packing problem. Another version is: minimize the size of the container, such that it can contain \(n\) spheres. In [1] some nice pictures are produced for different type of containers:

Smallest containers with 100 balls. Picture from [1]

Smallest container problem

The problem of finding the smallest enclosing box to fit \(n\) circles (spheres) with radius \(r_i\), can be stated as a quadratic problem: \[\begin{align} \min\>&S\\ &(x_i - x_j)^2 +(y_i - y_j)^2 \ge (r_i+r_j)^2 &\forall i \lt j\\ & r_i \le x_i \le S- r_i \\ & r_i \le y_i \le S - r_i\end{align}\] For the 3D case, just extend this model with \(z_i\).

This looks like a simple model, but it is not. The constraints are non-convex, so solvers like Cplex and Gurobi are out the window.  Global solvers like Baron, Couenne and Antigone have a very tough time on this model. Even for \(n=10\) points they have troubles finding a proven optimal solution.

A simple heuristic that works well, is to solve this with a local solver (I used SNOPT here) using a multi-start approach. To make things easy, I just ran \(N=25\) copies of the problem simultaneously, each with a different random starting point.

Notes:

• My solutions are marginally worse than reported in [3]. One reason is that SNOPT is using tighter tolerances than [3]: their solutions are declared infeasible by SNOPT. The circles in [3] just slightly overlap.
• I simply used random starting point (to be precise: random locations in a larger box than I would expect the final solution to be). It may be interesting to see if we can do better, such as trying to generate a feasible starting point. Alternatively, we could start with circle organized in a grid and then jitter them a bit by a random perturbation.
• I  ran all 25 solves in parallel. My laptop has 4 real cores, so it may be slightly better to run 4 problems in parallel. I doubt this will make much difference. 
• Circles and spheres are somewhat easy. If the items are boxes, we need to worry about possible rotations. 
• In [2] a slightly different formulation is used. Instead of having \((x_i - x_j)^2 +(y_i - y_j)^2 \ge (r_i+r_j)^2\) for each \(i\lt j\), we can do \[\min_{i \lt j} \left[ (x_i - x_j)^2 +(y_i - y_j)^2 -(r_i+r_j)^2\right] \ge  0\] Instead [2] proposes:  \[\sum_{i\lt j} \left[ \max\{0,(r_i+r_j)^2- (x_i - x_j)^2 -(y_i - y_j)^2 \}   \right]^2 = 0 \] The outer squares are used to alleviate that \(\max\{0,x\}\) is non-differentiable: the function \(\left( \max\{0,x\}\right)^2 \) is differentiable.

A problem with \(n=50\) circles with radius \(r_i=1\) can be packed in a square of size \(S=14.021\).

A second problem with \(n=50\) and \(r_i \sim U(0,2)\) (i.e. uniform between 0 and 2) needs a slightly larger square with \(S=14.598\). Now we get the more psychedelic plot:

Knapsack problem

We can make this problem a lot more difficult. Suppose we have a collection of circles with a value \(v_i\) and a radius \(r_i\). Given a rectangular container of size \(W \times H\) try to design a payload with maximum value. Our random data looks like:

----     19 PARAMETER v  value

i1  4.237,    i2  4.163,    i3  2.183,    i4  2.351,    i5  6.302,    i6  8.478,    i7  3.077,    i8  6.992
i9  7.983,    i10 3.733,    i11 1.994,    i12 5.521,    i13 2.442,    i14 8.852,    i15 3.386,    i16 3.572
i17 6.346,    i18 7.504,    i19 6.654,    i20 5.174


----     19 PARAMETER r  radius

i1  1.273,    i2  4.295,    i3  2.977,    i4  1.855,    i5  1.815,    i6  1.508,    i7  2.074,    i8  4.353
i9  0.802,    i10 2.751,    i11 4.992,    i12 3.104,    i13 4.960,    i14 3.930,    i15 1.088,    i16 3.379
i17 1.218,    i18 1.625,    i19 3.510,    i20 2.459


----     19 PARAMETER d  diameter

i1  2.546,    i2  8.589,    i3  5.953,    i4  3.710,    i5  3.630,    i6  3.016,    i7  4.148,    i8  8.706
i9  1.604,    i10 5.502,    i11 9.983,    i12 6.209,    i13 9.920,    i14 7.860,    i15 2.176,    i16 6.757
i17 2.436,    i18 3.251,    i19 7.020,    i20 4.918

We want to place as many high-valued and small circles in our rectangle:

Color coded items we can place in our rectangular knapsack

The rectangular container for this problem has \(W=15, H=10\).

The model to select the most valuable payload can be formulated as:\[\begin{align}\max\>&\sum_i v_i \delta_i\\ & (x_i-x_j)^2 + (y_i-y_j)^2 \ge (r_i+r_j)^2 \delta_i \delta_j & \forall i\lt j\\ & r_i \le x_i \le W - r_i\\ & r_i \le y_i \le H - r_i\\ &\delta_i \in \{0,1\}\end{align}\] Here the binary variable \(\delta_i\) indicates whether we select circle \(i\) for inclusion in the container.

Instead of \((x_i-x_j)^2 + (y_i-y_j)^2 \ge (r_i+r_j)^2 \delta_i \delta_j\), we can use a big-M constraint: \[(x_i-x_j)^2 + (y_i-y_j)^2 \ge (r_i+r_j)^2 - (1-\delta_i) M_{i,j} - (1-\delta_j) M_{i,j} \] with \(M_{i,j}=(r_i+r_j)^2\). Rewriting this a bit leads to: \[(x_i-x_j)^2 + (y_i-y_j)^2 \ge M_{i,j} (\delta_i + \delta_i - 1)\]

This is a MIQCP (Mixed Integer Quadratically Constrained Problem). Again, commercial solvers like Cplex and Gurobi don't help us here. We can try to solve this with a local MINLP solver or be ambitious and use a global solver. Even though we can not prove optimality, we can try to get a good solution with a local solver and then improve the solution with a global solver.

The local solvers (SBB, DICOPT, BONMIN) seem to like formulation \( (x_i-x_j)^2 + (y_i-y_j)^2 \ge (r_i+r_j)^2 \delta_i \delta_j \) better. BONMIN finds a solution of 58.57. LINDOGlobal is able to improve this to 60.36. The solution looks like:

----     92 VARIABLE delta.L  select

i1  1.000,    i4  1.000,    i5  1.000,    i6  1.000,    i7  1.000,    i9  1.000,    i12 1.000,    i15 1.000
i17 1.000,    i18 1.000,    i20 1.000


----     92 VARIABLE x.L  

i1   1.273,    i2   7.726,    i3   3.777,    i4   7.983,    i5   4.313,    i6  13.492,    i7   2.074,    i8  10.647
i9   3.980,    i10  5.454,    i11  5.769,    i12  7.721,    i13  7.856,    i14 11.070,    i15 13.912,    i16 10.971
i17 11.610,    i18  1.625,    i19  4.373,    i20 12.541


----     92 VARIABLE y.L  

i1  1.273,    i2  5.705,    i3  4.005,    i4  1.855,    i5  1.815,    i6  6.350,    i7  7.926,    i8  4.991
i9  5.771,    i10 4.268,    i11 4.992,    i12 6.896,    i13 5.040,    i14 5.713,    i15 8.912,    i16 5.312
i17 8.782,    i18 4.253,    i19 5.705,    i20 2.459

Plotting leads to:

Circles selected for 15 x 10 container. Unselected circles are plotted with transparent color.

Interestingly we don't select the most valuable circle i14. Note that this is not a proven globally optimal solution. To do a sanity check, I fixed \(\delta_{i14}=1\) and retried. Indeed I was not able to beat the best total value of  60.36.

References

1. Ernesto G. Birgin, Applications of nonlinear programming to packing problems, in Applications + Practical Conceptualization + Mathematics = fruitful Innovation, Proceedings of the Forum of Mathematics for Industry 2014, R. S. Anderssen, P. Broadbridge, Y. Fukumoto, K. Kajiwara, T. Takagi, E. Verbitskiy, and M. Wakayama (Eds.), Springer Series Mathematics for Industry 11, pp. 31-39, 2016
2. E. G. Birgin and F. N. C. Sobral, "Minimizing the object dimensions in circle and sphere packing problems", Computers & Operations Research 35, pp. 2357-2375, 2008.
3. Packing Circles and Spheres using nonlinear programming, https://www.ime.usp.br/~egbirgin/packing/packing_by_nlp
4. Mhand Hifi and Rym M'Hallah, “A Literature Review on Circle and Sphere Packing Problems: Models and Methodologies,” Advances in Operations Research, vol. 2009, Article ID 150624, 22 pages, 2009. https://doi.org/10.1155/2009/150624.
5. Giorgio Fasano, Solving Non-standard Packing Problems by Global Optimization and Heuristics, Springer, 2014.