# OptCorner: Christmas edition

$ \newcommand{\R}{\mathbb{R}} \newcommand{\costR}{\mathcal{R}} $

Merry Christmas! Welcome to the Christmas special edition of The NAG Optimization Corner series. For the occasion, we will show you how to solve a crucial problem in these festive times: Finding out where best to hide Santa's candy!

Let's say that after an extensive survey of family members, friends, coworkers and Christmas elves, you miraculously manage to create a nonlinear model of the candy attractivity that satisfies everybody. You still need to be able to optimize it over the intricate region of the tree and the present: \begin{equation*} \min_{(x,y)\in Tree \ \cup\ Present} F(x,y) \end{equation*} The following demonstrates one way to model such a feasible region.

## Let's build a Christmas tree!

### Starting with the top

If you limit yourself to the top of the tree, the problem is fairly simple. It is just a nonlinear minimization problem subject to some linear constraints: \begin{equation*} \min_{(x,y)\in \mathbb{R}^n} F(x,y)\\ \text{s.t. } \left\{ \begin{aligned} x + y &\leq 6\\ -x + y &\leq 6\\ y &\geq 5\\ \end{aligned}\right. \end{equation*} which define the region as follows:

It can be solved easily with a nonlinear optimization solver. We used the NAG e04wd solver to find the optimum, using the top of the tree as the starting point.

Similarly, you could limit yourself only to the middle of the tree and solve a very similar problem \begin{equation*} \min_{(x,y)\in \mathbb{R}^n} F(x,y)\\ \text{s.t. } \left\{ \begin{aligned} x + y &\leq 5.5\\ -x + y &\leq 5.5\\ y &\geq 4\\ \end{aligned}\right. \end{equation*} However, in each of these models, you focus only on one part of the region and you have to ignore the rest. This is clearly not a good solution to our problem!

### Joining two domains

It is possible to join both the middle and the top of the tree to the same model by introducing additional decision variables to the problem. The goal is to introduce a binary variable $\delta_K$ for each subdomain $K$ such that $\delta_K=1 \Rightarrow (x,y) \in K$. In general, for a constraint of the form $ax+by \leq c$, it is possible to create such a decision variable $\delta$ by replacing the original inequality by: $$ax+by+M\delta \leq M + c$$ where $M$ is an upper bound on the quantity $ax+by-c$. You then have a constraint where if $\delta=1$, the original inequality applies, and if $\delta=0$, the inequality becomes $ax+by-c \leq M$ which is always satisfied so has no influence.

As each of the subdomain is defined by three linear inequalities, we can apply this technique to all of our constraints and introduce two binary variables $\delta_{top}$ and $\delta_{mid} \in \{0,1\}$. The model can then be written as: \begin{equation*} \min_{(x,y)\in \mathbb{R}^n} F(x,y)\\ \text{s.t. } \left\{ \begin{aligned} x + y + \delta_{top}M &\leq M + 6\\ -x + y + \delta_{top}M &\leq M + 6\\ y - \delta_{top}M &\geq -M + 5\\ x + y + \delta_{mid}M &\leq M + 5.5\\ -x + y + \delta_{mid}M &\leq M + 5.5\\ y - \delta_{mid}M &\geq -M + 4\\ \end{aligned}\right. \end{equation*} For simplicity we chose a unique $M=9$ which works as a bound for all of these inequalities.

As it stands, the model is still incomplete. The solver could indeed choose $\delta_{mid}=\delta_{top}=0$ and ignore all constraints. To force the model to choose a point that is in at least one of the domains, we need to add a binding constraint: $$ \delta_{mid}+\delta_{top} \geq 1 $$ specifying that at least one of the sets of constraints has to be satisfied.

As the model now involves binary variables, a classical NLP solver is no longer suitable. We solved this model with NAG's MINLP solver h02da.

As you can see, the solver explores both subdomains and converges towards a minimum. Note that this is a local solver and that it could declare convergence at any of the local minima in the domain depending on the starting point.

### Finishing touches

It is possible to generalize the technique described above to create an arbitrary number of domains. For example, the full tree model can be solved with the same solver and two additional binary variables as shown in the picture below

It is even possible to create a completely disjointed domain with this method. Let's for example add a present under the tree! The full model becomes: \begin{equation*} \min_{(x,y)\in \mathbb{R}^n} F(x,y)\\ \text{s.t. } \left\{ \begin{aligned} \left. \delta_{top} + \delta_{mid} + \delta_{bot} + \delta_{tr} + \delta_{pr} \geq 1 \right\}\text{binding constraint}\\ \left.\begin{aligned} x + y + \delta_{top}M &\leq M + 6\\ -x + y + \delta_{top}M &\leq M + 6\\ y - \delta_{top}M &\geq -M + 5\\ \end{aligned}\right\}\text{top of the tree}\\ \left.\begin{aligned} x + y + \delta_{mid}M &\leq M + 5.5\\ -x + y + \delta_{mid}M &\leq M + 5.5\\ y - \delta_{mid}M &\geq -M + 4\\ \end{aligned}\right\}\text{middle of the tree}\\ \left.\begin{aligned} x + y + \delta_{bot}M &\leq M + 5\\ -x + y + \delta_{bot}M &\leq M + 5\\ y - \delta_{bot}M &\geq -M + 2.5\\ \end{aligned}\right\}\text{bottom of the tree}\\ \left.\begin{aligned} x - \delta_{tr}M &\geq -M - 0.75\\ x + \delta_{tr}M &\leq M + 0.75\\ y - \delta_{tr}M &\geq -M + 2\\ y + \delta_{tr}M &\leq M + 2.5\\ \end{aligned}\right\}\text{trunk}\\ \left.\begin{aligned} x - \delta_{pr}M &\geq -M - 1\\ x + \delta_{pr}M &\leq M - 0.25\\ y - \delta_{pr}M &\geq -M + 1.8\\ y + \delta_{pr}M &\leq M + 1.3\\ \end{aligned}\right\}\text{present}\\ \end{aligned}\right. \end{equation*}

If we solve this model with h02da, we finally obtain the results we were all waiting for, the candy should be in the present!

## Things to remember

While we described a fairly simple example, it is quite representative of the kind of model that is possible to create when modelling with binary variables. In particular, disjoint sets of constraints can appear quite frequently in practical applications. Just remember not to abuse integer variables as the complexity of your model will increase greatly as you increase their number.

Have a great holiday and stay tuned for more posts from The NAG Optimization Corner in January!

The experiment presented above was performed using the NAG Toolbox *for MATLAB®* with the following files (3,124 bytes).