MINLP instead of indicator constraints?

In this post, I want to discuss indicator constraints and how to replace them with simple multiplications. As we shall see, this is a somewhat harebrained but still interesting idea. 

Indicator constraints are implications of the form: \[\delta=0 \implies \text{linear constraint}\] or \[\delta=1 \implies \text{linear constraint}\] where \(\delta \in \{0,1\}\) is a binary decision variable.

There are two aspects of indicator constraints: 

• Indicator constraints help with MIP models where otherwise we would use big-M constraints. This will help address the numerical issues that result from using big-M constraints. These include small (and sometimes not-so-small) solution values where you expect zero, and poor solver performance. With big-M constraints, we need to pay much attention to the size of the big-M constants. This is the algorithmic aspect.
• Indicator constraints provide a convenient modeling construct. They form a useful abstraction that makes MIP modeling easier and more straightforward. This is the modeling aspect.

Another example where a solver feature has such a dual benefit (algorithmic and modeling) are SOS2 sets for piecewise linear functions. 

Sorting: minimize number of swaps

In this post, I want to delve further into sorting. In a question on or.stackexchange.com [1], the subject was minimizing the number of swaps in sorting algorithms. A swap, i.e., an interchange of two items, is a basic operation in sorting. We usually don't pay much attention to this. First, we assume swaps are cheap. If they are not, we can sort not the real data (which can be complex, and thus more complicated and expensive to swap), but rather pointers to the data or keys. But what if we really have to reorder bulky physical things? Then the number of swaps is a more meaningful concept. Can we minimize the work to reorder say \(n\) boxes (or containers to make it more dramatic) by minimizing the number of swaps (assuming interchanging the position of two boxes is the way reordering works).

Clock problem

From [1]:

The hour, minute and second hands of this clock are all the same length and move smoothly in a circle. The dial contains hour and minute markers, but the numbers are missing. Therefore, it’s impossible to tell which one of the 12 hour markers belongs to the 12. The two hands on the left are positioned exactly on hour markers, and the hand on the right is positioned between a minute and an hour marker. What time does the clock show?

It is possible to solve this without really any math, but, of course, here I try to model this as a mathematical programming model.

Graph connectivity as constraints

I was generating, for an example model, a random, sparse, directed graph. Unfortunately, when sparse enough, this is likely to yield a graph that is not connected. Here is my first attempt.

Nodes

I generate simply \(n=25\) nodes with \((x,y)\) coordinates from a uniform distribution \(U(0,100)\). This looks like:

\(n=25\) nodes randomly placed

Revisiting a continuous facility location problem

I am revisiting here a problem from [1]:

We have \(n\) demand points and their locations. How many facilities do we need to service these customers? And where do we place them? We have a restriction: there is a maximum distance between customer and facility.   

Data

We randomly generated 75 demand points inside the \([0,1]\times[0,1]\) square. The maximum allowed distance between a demand point and a facility is 0.25. 

Random Demand Points

Towers of Hanoi: inventory and network formulation

A standard problem with 4 disks requires 15 moves

The towers of Hanoi problem [1] is a famous puzzle demonstrating recursion. The task is to move a stack of disks from one peg to another. We can only move one disk at a time, and we need to obey the rule that larger disks can never be on top of a smaller disk. The moves for the 3 disk problem are:

Click on the picture to enlarge.

Wolf, goat and cabbage problem: MIP and network model

The Wolf, goat, and cabbage problem can be stated as [1]:

A farmer with a wolf, a goat, and a cabbage must cross a river by boat. The boat can carry only the farmer and a single item. If left unattended together, the wolf would eat the goat, or the goat would eat the cabbage. How can they cross the river without anything being eaten?

A long transcontinental flight was a good opportunity to try to attack this problem. Here are some approaches to model this. Not at all a very useful or practical model, but still interesting, I think (although I may be in a small minority here). 

Inventory model

In this model, we keep track of the inventory of items (wolf, goat, cabbage). We assume the starting inventory is: all items are on the left bank of the river. The final inventory should be: all items are on the right bank. In this model, I assume the following numbering scheme:

----     31 SET trip  trips

trip1 ,    trip2 ,    trip3 ,    trip4 ,    trip5 ,    trip6 ,    trip7 ,    trip8 ,    trip9 ,    trip10


----     31 SET dir  direction of trip

L->R,    R->L


----     31 SET tripDir  trip direction combos

              L->R        R->L

trip1          YES
trip2                      YES
trip3          YES
trip4                      YES
trip5          YES
trip6                      YES
trip7          YES
trip8                      YES
trip9          YES
trip10                     YES

Small MIP, proving optimality is difficult

This is a simple, smallish MIP model that is difficult to solve to proven optimality. The problem is stated as [1]:

I have grid of dimensions H and W of boolean variables. The only constraint is that if a variable is false then at least one of the adjacent variables (top, right, left, bottom, diagonals don't count) must be true. The goal is to minimize the number of true values in the grid.

A high-level mathematical representation of this can be:

High-level model
\[\begin{align}\min&\sum_{i,j}\color{darkred}x_{i,j}\\ & \color{darkred}x_{i,j}=0 \implies \color{darkred}x_{i-1,j} + \color{darkred}x_{i+1,j} + \color{darkred}x_{i,j-1} + \color{darkred}x_{i,j+1} \ge 1 \\ & \color{darkred}x_{i,j}\in \{0,1\} \end{align}\]

Sorting using a MIP model

This is not, per se, very useful, but sorting a parameter inside a MIP model is not very easy for MIP solvers. Obviously, if you need a sorted parameter in the model, it is better to use a sorting algorithm. But useless models can still be interesting.

I use: 

    Input: a 1-dimensional parameter with values.
    Output: a 1-dimensional variable with the values sorted in ascending order.

We can implement this with a permutation matrix \(\color{darkred}X\), which is a permuted identity matrix. In a MIP context, this becomes a binary variable with some assignment constraints.

MIP Model for sorting \(p_i\)
\[\begin{align}\min\>&\color{darkred}z=0\\& \sum_i \color{darkred}x_{i,j} = 1&&\forall j\\ & \sum_j \color{darkred}x_{i,j} = 1&&\forall i\\ & \color{darkred}y_i = \sum_j \color{darkred}x_{i,j}\cdot\color{darkblue}p_j\\& \color{darkred}y_i \ge \color{darkred}y_{i-1}\\ & \color{darkred}x_{i,j} \in \{0,1\}\end{align}\]

N-queens and solution pool

In [1], I described some chess-related problems. Here, I want to reproduce the \(n\)-queens problem. The single solution problem, placing as many queens on the chess board as possible so they don't attack each other, is pretty standard. I want to focus on the more complex question: How many different ways can we place those queens? In other words: what are all the optimal solutions? We can do this by adding a no-good constraint that forbids the previously found solution. However, as this problem has more than a handful of different solutions, I want to use the Cplex solution pool.

Single Solution Model

We define the decision variables as: \[\color{darkred}x_{i,j} = \begin{cases} 1 & \text{if we place a queen on the square $(i,j)$} \\ 0 & \text{otherwise}\end{cases}\] 

Chess Board

Rounding inside an optimization model

In [1], the question was asked: how can I round to two decimal places inside an optimization model? I.e., \[\color{darkred}y_{i,j} = \mathbf{round}(\color{darkred}x_{i,j},2)\] To get this off my chest first: I have never encountered a situation like this. Rounding to two decimal places is more for reporting than something we want inside model equations. Given that, let me look into this modeling problem a bit more as an exercise. 

One nonzero in set of free variables

In [1] the following question is posed:

I have free variables \(\color{darkred}x_i\). How can I impose the constraint that at least one of the variables is nonzero: \(\color{darkred}x_i\ne 0\).

Grouping items: a difficult combinatorial problem

In [1], a simple problem is described:

• We have \(n\) items (or orders) with a certain width. 
• We need to combine these items in groups (called patterns) with rather tight limits on the total width. The total length of a pattern (the sum of the lengths of the items assigned to this pattern) must be between 335 and 340.
• As a result, we may not be able to assign all items. The remaining items cannot be formed into valid patterns.
• The objective is to try to place as many items as possible into patterns.
• An indication of the size of the problem: \(n \approx 500\).  

Data

Instead of immediately working on a full-known \(n=500\) problem, I generated a random data set with a very manageable \(n=50\) items. The widths were drawn from a discrete uniform distribution between 30 and 300. The data looks like:

----     15 PARAMETER w  item widths

order1   76.000,    order2  258.000,    order3  179.000,    order4  111.000,    order5  109.000,    order6   90.000
order7  124.000,    order8  262.000,    order9   48.000,    order10 165.000,    order11 300.000,    order12 186.000
order13 298.000,    order14 236.000,    order15  65.000,    order16 203.000,    order17  73.000,    order18  97.000
order19 211.000,    order20 147.000,    order21 127.000,    order22 125.000,    order23  65.000,    order24  70.000
order25 189.000,    order26 255.000,    order27  92.000,    order28 210.000,    order29 240.000,    order30 112.000
order31  59.000,    order32 166.000,    order33  73.000,    order34 266.000,    order35 101.000,    order36 107.000
order37 190.000,    order38 225.000,    order39 200.000,    order40 155.000,    order41 142.000,    order42  61.000
order43 115.000,    order44  42.000,    order45 121.000,    order46  79.000,    order47 204.000,    order48 181.000
order49 238.000,    order50 110.000

I stick to the pattern limits \(\color{darkblue}L=335\) and \(\color{darkblue}U=340\).

We need some estimate of the number of patterns to use. We could just guess. But a better approach is the following. An upper bound for the number patterns can be established quite easily: \[{\mathit{maxj}} = \left\lfloor \frac{\sum_i \color{darkblue}w_i}{\color{darkblue}L}\right\rfloor\] For our data set this number is:

----     29 PARAMETER maxj                 =       22.000  max number of patterns we can fill

This means we can safely use this number as the number of bins (patterns). 

Scheduling Team Meetings

This is a simple problem from [1]:

I'm rusty on constraint optimization and am looking for help in this particular

use case. There are individuals who are each member to several teams. This is a

fixed many-to-many relationship and is determined a-priori. There are 3 time

slots where the teams can be scheduled to conduct a business meeting, but if an

individual is a member of more than one team which are both meeting at a given

time slot, they'll only be able to attend one. The objective is to schedule the

teams into the time slots, minimizing the number of overlaps of individuals.

For beginners, it is often a good idea to split the task in two:

1. Formulate a mathematical model (on a piece of paper)
2. Implement the model in code

Generate all solutions that sum up to one

In a post, the following question was posed:

We can select unique values  \(\displaystyle\frac{1}{i}\) for \(i=1,\dots,n\). Find all combinations that add up to 1.

A complete enumeration scheme was slow even for \(n=10\). Can we use a MIP model for this or something related?

A single solution is easily found using the model:

Mathematical Model
\[ \begin{align} & \sum_{i=1}^n \frac{1}{i} \cdot \color{darkred}x_i = 1 \\ & \color{darkred}x_i \in \{0,1\} \end{align}\]

Finding common patterns

In [1], the following problem is stated:

Given a boolean matrix, with \(m\) rows and \(n\) columns, find the largest pattern of ones that is found in at least \(\color{darkblue}K\) rows. We can ignore cells where the pattern has a zero value: they don't count.

A small example [1] is given: 

row 1:[0, 1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 1, 1, 0, 1, 0, 0, 0, 1, 0, 1]
row 2:[0, 1, 0, 1, 1, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 1, 1, 1, 0, 1, 1]
row 3:[0, 0, 0, 1, 1, 0, 0, 1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 1, 0, 0, 1]
row 4:[1, 1, 0, 1, 1, 1, 0, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1]
row 5:[1, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 0, 1, 1]
row 6:[1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1]

 With \(K=3\), we can form a pattern with 10 nonzero elements:

 

row 1:  [0, 1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 1, 1, 0, 1, 0, 0, 0, 1, 0, 1]
row 2:  [0, 1, 0, 1, 1, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 1, 1, 1, 0, 1, 1]
row 3:  [0, 0, 0, 1, 1, 0, 0, 1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 1, 0, 0, 1]
row 4:  [1, 1, 0, 1, 1, 1, 0, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1]
row 5:  [1, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 0, 1, 1]
row 6:  [1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1] 
pattern:[1, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1]

The pattern is shared by rows 4, 5, and 6.

Can we formulate a MIP model for this? My first attempt is as follows.

A large MIP model that should be solved as LP: the root node

In this post, I want to discuss an observation about the root node when solving a MIP model.

Problem Description

The model I want to tackle here is from [1]:

We have a large matrix \(\color{darkblue}A=\color{darkblue}a_{i,j}\) with values 0 and 1. In addition, we have a minimum on the row and column totals. These are called \(\color{darkblue}r_i\) and \(\color{darkblue}c_j\). The goal is to remove as many 1's in the matrix \(\color{darkblue}A\) subject to these minimum row and column totals.

Another fast MIP model: covering

In [1], the following problem is stated:

• There is a collection of \(n=1,000\) test questions.
• Each question covers a number of skills.
• Given is a requirement for a number of questions for each required skill (e.g., 4 questions about skill 1, 3 questions about skill 2, etc.).
• Create a test with the minimum number of questions that fulfills the requirements.

Of course, we see some remarks about NP-Hardness. Complexity theory is about bounds, about worst-case and asymptotic behavior. Even if we can't prove that a particular problem with a particular dataset given to a particular solver must be fast, it does not mean the opposite: that it must be slow. It is not at all unusual that we can solve things very quickly, even though the complexity bounds are terrible. Of course, this depends on the model, on the data set, and on the solver. Here, we see an open-source MIP solver can solve this particular problem in less than 0.01 seconds. Many well-educated people misunderstand complexity theory – to their detriment. They miss out on great ways to solve actual problems! Obviously, this is a problem with how complexity theory is taught. 

Assigning jobs to machines without overlap

 Here we consider the following problem from [1]:

• We have jobs with a given start time and completion time
• Jobs can be repeated on given days (e.g. job 1 needs to run on Monday, Wednesday, and Friday)
• We want to assign jobs to machines in such a way that there is no overlap
• The objective is to minimize the number of machines needed to execute all jobs

The planning horizon is a week, and the problems has 50 jobs.

Supplier selection: an easy MIP

This is a fairly standard supplier selection model. It was posted in [1]. It is interesting as it is a good fit for a MIP model: it is simple to model, it solves fast, and it is not obvious how to solve without an optimization model.

The problem is:

• We want to order items in different quantities from suppliers.
• Suppliers have an available inventory for these items. This can be zero.
• We can split the ordering over different suppliers.
• The cost structure is as follows:
• Shipping cost is a fixed cost per supplier.
• Item cost is a variable per-unit cost.