A nonlinear programming model with implicit variables for packing ellipsoids

Abstract

The problem of packing ellipsoids is considered in the present work. Usually, the computational effort associated with numerical optimization methods devoted to packing ellipsoids grows quadratically with respect to the number of ellipsoids being packed. The reason is that the number of variables and constraints of ellipsoids’ packing models is associated with the requirement that every pair of ellipsoids must not overlap. As a consequence, it is hard to solve the problem when the number of ellipsoids is large. In this paper, we present a nonlinear programming model for packing ellipsoids that contains a linear number of variables and constraints. The proposed model finds its basis in a transformation-based non-overlapping model recently introduced by Birgin et al. (J Glob Optim 65(4):709–743, 2016). For solving large-sized instances of ellipsoids’ packing problems with up to 1000 ellipsoids, a multi-start strategy that combines clever initial random guesses with a state-of-the-art (local) nonlinear optimization solver is presented. Numerical experiments show the efficiency and effectiveness of the proposed model and methodology.

Appendix: Derivatives

The computation of the derivatives of the function defined in (33) is nontrivial. That is because this function depends on the functions $\mathcal {X}$ and $\mathcal {U}$ whose values are given by the solution of an optimization problem. Firstly, we will show the derivatives of the terms that compose the function defined in (33) in terms of the derivatives of the functions $\mathcal {X}$ and $\mathcal {U}$. Next, we will show how to compute the derivatives of the functions $\mathcal {X}$ and $\mathcal {U}$. To simplify the notation, we will denote by $\mathcal {X}_{ij}$ the value $\mathcal {X}(c_i,c_j,\Omega _i,\Omega _j;P_i,P_j)$ and by $\mathcal {U}_{ij}$ the value $\mathcal {U}(c_i,c_j,\Omega _i,\Omega _j;P_i,P_j)$.

1.1 First order derivatives

Let $i,j \in \{1,\dots ,m\}$ such that $i < j$. We have that $\mathcal {X}(c_i,c_j,\Omega _i,\Omega _j;P_i,P_j)$ is a solution to the problem

$$\begin{aligned} \begin{array}{ll} \text {minimize} &{} \frac{1}{2}\Big \Vert x - c_{i}^{ij}\Big \Vert ^2\\ \text {subject to} &{} x^{\top } S_{ij} x = 1, \end{array} \end{aligned}$$

(38)

where $c_{i}^{ij} = P_i^{-\frac{1}{2}}Q_i^{\top }(c_i-c_j)$, and $\mathcal {U}(c_i,c_j,\Omega _i,\Omega _j;P_i,P_j)$ is the corresponding Lagrange multiplier. According to the Karush–Kuhn–Tucker first-order necessary conditions for problem (38), we have

$$\begin{aligned}&\mathcal {X}_{ij} + \mathcal {U}_{ij} S_{ij}\mathcal {X}_{ij} - c_{i}^{ij} = 0\\&{\mathcal {X}_{ij}}^{\top }S_{ij}\mathcal {X}_{ij} - 1 = 0. \end{aligned}$$

Thus, by defining the function $F: \mathbb {R}^{n} \times \mathbb {R}^{n} \times \mathbb {R}^{q} \times \mathbb {R}^{q} \rightarrow \mathbb {R}^{n+1}$ as

$$\begin{aligned} F(c_i,c_j,\Omega _i,\Omega _j) = \left[ \begin{array}{c} \mathcal {X}_{ij} + \mathcal {U}_{ij} S_{ij}\mathcal {X}_{ij} - c_{i}^{ij}\\ \frac{1}{2}\left( {\mathcal {X}_{ij}}^{\top }S_{ij}\mathcal {X}_{ij} - 1\right) \end{array} \right] , \end{aligned}$$

(39)

we have that $F(c_i,c_j,\Omega _i,\Omega _j) = 0$ for all $c_i,c_j \in \mathbb {R}^n$ and for all $\Omega _i,\Omega _j \in \mathbb {R}^q$. That is, F is an identically zero function. Therefore, we have that the derivative of function F is also an identically zero function. Hence, for each variable v of the function F and for each component $\ell \in \{1,\dots ,n+1\}$ of F, we have

$$\begin{aligned} \frac{\text {d}F_{\ell }}{\text {d}v} = \frac{\partial F_{\ell }}{\partial v}\frac{\text {d}v}{\text {d}v} + \frac{\partial F_{\ell }}{\partial \mathcal {U}_{ij}}\frac{\text {d}\mathcal {U}_{ij}}{\text {d}v} + \sum _{k=1}^{n} \frac{\partial F_{\ell }}{\partial [\mathcal {X}_{ij}]_k}\frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}v} = 0. \end{aligned}$$

(40)

Once the values of $\mathcal {X}_{ij}$ and $\mathcal {U}_{ij}$ are known, we have, for each $\ell \in \{1,\dots ,n+1\}$, analytical expressions for $\frac{\partial F_{\ell }}{\partial v}$, $\frac{\partial F_{\ell }}{\partial \mathcal {U}_{ij}}$, and $\frac{\partial F_{\ell }}{\partial [\mathcal {X}_{ij}]_k}$ for each $k \in \{1,\dots ,n\}$. On the other hand, the values of $\frac{\text {d}\mathcal {U}_{ij}}{\text {d}v}$ and $\frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}v}$ for each $k \in \{1,\dots ,n\}$ are unknown, but can be computed by solving the linear system provided by (40):

$$\begin{aligned} \begin{bmatrix} \frac{\partial F_{1}}{\partial [\mathcal {X}_{ij}]_1}&\quad \cdots&\quad \frac{\partial F_{1}}{\partial [\mathcal {X}_{ij}]_n}&\quad \frac{\partial F_{1}}{\partial \mathcal {U}_{ij}}\\ \frac{\partial F_{2}}{\partial [\mathcal {X}_{ij}]_1}&\quad \cdots&\quad \frac{\partial F_{2}}{\partial [\mathcal {X}_{ij}]_n}&\quad \frac{\partial F_{2}}{\partial \mathcal {U}_{ij}}\\ \vdots&\quad \vdots&\quad \vdots&\quad \vdots \\ \frac{\partial F_{n+1}}{\partial [\mathcal {X}_{ij}]_1}&\quad \cdots&\quad \frac{\partial F_{n+1}}{\partial [\mathcal {X}_{ij}]_n}&\quad \frac{\partial F_{n+1}}{\partial \mathcal {U}_{ij}}\\ \end{bmatrix} \begin{bmatrix} \frac{\text {d}[\mathcal {X}_{ij}]_1}{\text {d}v}\\ \vdots \\ \frac{\text {d}[\mathcal {X}_{ij}]_n}{\text {d}v}\\ \frac{\text {d}\mathcal {U}_{ij}}{\text {d}v} \end{bmatrix} = - \begin{bmatrix} \frac{\partial F_{1}}{\partial v}\\ \frac{\partial F_{2}}{\partial v}\\ \vdots \\ \frac{\partial F_{n+1}}{\partial v} \end{bmatrix}. \end{aligned}$$

Then, for each $i, j \in \{1,\dots ,m\}$ such that $i < j$, we need to solve $2(n+q)$ linear systems with $n+1$ equations and $n+1$ variables (one linear system for each variable among $c_i$, $c_j$, $\Omega _i$ and $\Omega _j$).

Once i and j are fixed, observe that the $2(n+q)$ linear systems have the same coefficient matrix. The only difference between these systems are their right-hand sides. Thus, in order to solve these linear systems, we can factorize the coefficient matrix only once and then, for each right-hand side, solve the linear system with the coefficient matrix already factorized.

1.2 Second order derivatives

For each variable v of the function F defined in (39) and for each component $\ell$ of F, we define the function $G^v_\ell$ as the total derivative of the function $F_\ell$ with respect to v:

$$\begin{aligned} G^v_\ell (c_i,c_j,\Omega _i,\Omega _j) = \frac{\text {d}F_\ell }{\text {d}v}(c_i,c_j,\Omega _i,\Omega _j). \end{aligned}$$

Since the function F is identically zero, its derivative is also identically zero. Then, for each variable u of the function $G^v_\ell$, we have

$$\begin{aligned} \frac{\text {d}G^v_{\ell }}{\text {d}u} = \frac{\partial G^v_{\ell }}{\partial u}\frac{\text {d}u}{\text {d}u} + \frac{\partial G^v_{\ell }}{\partial \mathcal {U}_{ij}}\frac{\text {d}\mathcal {U}_{ij}}{\text {d}u} + \sum _{k=1}^{n} \frac{\partial G^v_{\ell }}{\partial [\mathcal {X}_{ij}]_k}\frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}u} = 0. \end{aligned}$$

(41)

Next, we present the partial derivatives of the function $G^v_{\ell }$, that appear in the expression (41). The partial derivative of $G^v_{\ell }$ with respect to the variable u is given by

$$\begin{aligned} \frac{\partial G^v_{\ell }}{\partial u}= & {} \frac{\partial ^2 F_\ell }{\partial u \partial v} + \frac{\partial ^2 F_\ell }{\partial u \partial \mathcal {U}_{ij}} \frac{\text {d}\mathcal {U}_{ij}}{\text {d}v} + \frac{\partial F_\ell }{\partial \mathcal {U}_{ij}} \frac{\partial }{\partial u}\left( \frac{\text {d}\mathcal {U}_{ij}}{\text {d}v}\right) + \sum _{k=1}^{n} {\frac{\partial ^2 F_\ell }{\partial u \partial [\mathcal {X}_{ij}]_k}} \frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}v}\\&+ \sum _{k=1}^{n} \frac{\partial F_\ell }{\partial [\mathcal {X}_{ij}]_k} \frac{\partial }{\partial u}\left( \frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}v}\right) . \end{aligned}$$

The partial derivative of $G^v_{\ell }$ with respect to $\mathcal {U}_{ij}$ is given by

$$\begin{aligned} \frac{\partial G^v_{\ell }}{\partial \mathcal {U}_{ij}}&= \frac{\partial ^2 F_\ell }{\partial \mathcal {U}_{ij} \partial v} + \frac{\partial ^2 F_\ell }{\partial \mathcal {U}_{ij} \partial \mathcal {U}_{ij}} \frac{\text {d}\mathcal {U}_{ij}}{\text {d}v} + \frac{\partial F_\ell }{\partial \mathcal {U}_{ij}} \frac{\partial }{\partial \mathcal {U}_{ij}}\left( \frac{\text {d}\mathcal {U}_{ij}}{\text {d}v}\right) \\&\quad +\sum _{k=1}^{n} {\frac{\partial ^2 F_\ell }{\partial \mathcal {U}_{ij} \partial [\mathcal {X}_{ij}]_k}} \frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}v} + \sum _{k=1}^{n} \frac{\partial F_\ell }{\partial [\mathcal {X}_{ij}]_k} \frac{\partial }{\partial \mathcal {U}_{ij}}\left( \frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}v}\right) \\&= \frac{\partial ^2 F_\ell }{\partial \mathcal {U}_{ij} \partial v} + \sum _{k=1}^{n} {\frac{\partial ^2 F_\ell }{\partial \mathcal {U}_{ij} \partial [\mathcal {X}_{ij}]_k}} \frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}v}. \end{aligned}$$

Finally, the partial derivative of $G^v_{\ell }$ with respect to $\mathcal {X}_{ij}$ is given by

$$\begin{aligned} \frac{\partial G^v_{\ell }}{\partial [\mathcal {X}_{ij}]_t}&= \frac{\partial ^2 F_\ell }{\partial [\mathcal {X}_{ij}]_t \partial v} + \frac{\partial ^2 F_\ell }{\partial [\mathcal {X}_{ij}]_t \partial \mathcal {U}_{ij}} \frac{\text {d}\mathcal {U}_{ij}}{\text {d}v} + \frac{\partial F_\ell }{\partial \mathcal {U}_{ij}} \frac{\partial }{\partial [\mathcal {X}_{ij}]_t}\left( \frac{\text {d}\mathcal {U}_{ij}}{\text {d}v}\right) \\&\quad + \sum _{k=1}^{n} {\frac{\partial ^2 F_\ell }{\partial [\mathcal {X}_{ij}]_t \partial [\mathcal {X}_{ij}]_k}} \frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}v} + \sum _{k=1}^{n} \frac{\partial F_\ell }{\partial [\mathcal {X}_{ij}]_k} \frac{\partial }{\partial [\mathcal {X}_{ij}]_t}\left( \frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}v}\right) \\&= \frac{\partial ^2 F_\ell }{\partial [\mathcal {X}_{ij}]_t \partial v} + \frac{\partial ^2 F_\ell }{\partial [\mathcal {X}_{ij}]_t \partial \mathcal {U}_{ij}} \frac{\text {d}\mathcal {U}_{ij}}{\text {d}v} + \sum _{k=1}^{n} {\frac{\partial ^2 F_\ell }{\partial [\mathcal {X}_{ij}]_t \partial [\mathcal {X}_{ij}]_k}} \frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}v}. \end{aligned}$$

The simplifications in the expressions of the derivatives of $G^v_{\ell }$ with respect to $\mathcal {U}_{ij}$ and $\mathcal {X}_{ij}$ come from the removal of null elements.

Considering that the values of the first order derivatives of the function F are known, the equation (41) provides the following linear system

$$\begin{aligned} \begin{bmatrix} \frac{\partial F_{1}}{\partial [\mathcal {X}_{ij}]_1}&\quad \cdots&\quad \frac{\partial F_{1}}{\partial [\mathcal {X}_{ij}]_n}&\quad \frac{\partial F_{1}}{\partial \mathcal {U}_{ij}}\\ \frac{\partial F_{2}}{\partial [\mathcal {X}_{ij}]_1}&\quad \cdots&\quad \frac{\partial F_{2}}{\partial [\mathcal {X}_{ij}]_n}&\quad \frac{\partial F_{2}}{\partial \mathcal {U}_{ij}}\\ \vdots&\quad \vdots&\quad \vdots&\quad \vdots \\ \frac{\partial F_{n+1}}{\partial [\mathcal {X}_{ij}]_1}&\quad \cdots&\quad \frac{\partial F_{n+1}}{\partial [\mathcal {X}_{ij}]_n}&\quad \frac{\partial F_{n+1}}{\partial \mathcal {U}_{ij}}\\ \end{bmatrix} \begin{bmatrix} \frac{\partial }{\partial u}\left( \frac{\text {d}[\mathcal {X}_{ij}]_1}{\text {d}v}\right) \\ \vdots \\ \frac{\partial }{\partial u}\left( \frac{\text {d}[\mathcal {X}_{ij}]_n}{\text {d}v}\right) \\ \frac{\partial }{\partial u}\left( \frac{\text {d}\mathcal {U}_{ij}}{\text {d}v}\right) \end{bmatrix} = - \begin{bmatrix} b_{1}^{u,v}\\ b_{2}^{u,v}\\ \vdots \\ b_{n+1}^{u,v} \end{bmatrix} \end{aligned}$$

for each variable u and for each variable v of the function F. The components of the right-hand side of this system are given by

$$\begin{aligned} b_{\ell }^{u,v} = \frac{\partial ^2 F_\ell }{\partial u \partial v} + \frac{\partial ^2 F_\ell }{\partial u \partial \mathcal {U}_{ij}} \frac{\text {d}\mathcal {U}_{ij}}{\text {d}v} + \sum _{k=1}^{n} {\frac{\partial ^2 F_\ell }{\partial u \partial [\mathcal {X}_{ij}]_k}} \frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}v} + \frac{\partial G^v_{\ell }}{\partial \mathcal {U}_{ij}}\frac{\text {d}\mathcal {U}_{ij}}{\text {d}u} + \sum _{k=1}^{n} \frac{\partial G^v_{\ell }}{\partial [\mathcal {X}_{ij}]_k}\frac{\text {d}[\mathcal {X}_{ij}]_k}{\text {d}u} \end{aligned}$$

for each $\ell \in \{1,\dots ,n+1\}$. Notice that the coefficient matrix of this linear system does not depend on the variables u and v. Therefore, for each $i,j \in \{1,\dots ,m\}$ such that $i < j$, we have $(n+q)(2(n+q) + 1)$ linear systems (one for each pair of variables u and v of the function F) with $n+1$ variables each one, and all of them have the same coefficient matrix.