# Select Points

*
*

Given multiple sets of data points. Select one point per set so that the distances between the points is minimal.

^{1}

## Model

### Non-convex MIQP Model

$$ \min \sum_{i,j| ok_{i,j}} dist_{i,j} \cdot x_i \cdot x_j \\ \sum_{i|group_{i,g}} x_i = 1 \forall g \\ x_i \in {0,1} $$

### Linear MIP Model

To linearize the model a standard procedure can be applied:

Let’s introduce a new binary variable, $z_{i,j}$, that represents the product of the binary variables $\color{darkred}x_i$ and $\color{darkred}x_j$. This new variable will replace the nonlinear term in the objective function, making the model linear. We will also add constraints to link the new variable $z_{i,j}$ with the original variables $\color{darkred}x_i$ and $\color{darkred}x_j$.

\begin{align} \min&\sum_{i,j|ok_{i,j}} dist_{i,j}\cdot z_{i,j} \\ & \sum_{i|group_{i,g}} x_i = 1 && \forall g\\ & x_i \in {0,1} \\ & z_{i,j} \leq x_i && \forall i,j|ok_{i,j} \\ & z_{i,j} \leq x_j && \forall i,j|ok_{i,j} \\ & z_{i,j} \geq x_i + x_j - 1 && \forall i,j|ok_{i,j} \\ & z_{i,j} \in {0,1} && \forall i,j|ok_{i,j} \end{align}

## Implementation

Calculating the sparse $ok_{i,j}$ set is key for a simple Pyomo implementation of the model. With this set and using Pyomo’s decorator syntax the model looks compact:

```
m.I = pyo.RangeSet(config['N'])
m.G = pyo.RangeSet(config['G'])
m.Ok = pyo.Set(initialize=Ok(config).ok)
@m.Param(m.Ok)
def distance(model, i, j):
p = Point(*self.points[i])
q = Point(*self.points[j])
return sqrt((p.x - q.x) ** 2 + (p.y - q.y) ** 2)
m.groups = pyo.Param(m.G, initialize=group(config))
m.x = pyo.Var(m.I, domain=pyo.Binary, initialize=0)
m.pair = pyo.Var(m.Ok, domain=pyo.Binary, initialize=0)
@m.Constraint(m.G)
def one_point_per_group(m, g):
return sum(m.x[i] for i in m.groups[g]) == 1
@m.Constraint(m.Ok)
def both_selected(m, i, j):
return m.pair[i, j] >= m.x[i] + m.x[j] - 1
@m.Objective(sense=pyo.minimize)
def total_distance(m):
return sum(
m.pair[i, j] * m.distance[i, j] for (i, j) in m.Ok)
```

## Result

The CBC solver has no problems with this model:

- feasible and optimal solution found
- Solution 794.15391509
- Number of constraints : 4476
- Number of variables : 4566
- Duration: 00:00:22

The edges of the 10 solution points are visualized in the graph. The complete code can be found here.

Interesting to see is the solution time between CPLEX ^{1} 0.672s and CBC 22s.
There is a huge gap between commercial solvers and OSS solvers.

### Maximize

Maximizing the same model seems to be much harder for the solver. The same model size as above did not solve within 6000s.

Halve the size solved:

- feasible and optimal solution found
- Solution 777.58370723
- Number of constraints : 2933
- Number of variables : 1026
- Duration: 00:00:33