sysid blog

Dispersion

June 21, 2020 ☕️ 1 min read

Find 10 points out of 50 with maximum distance from each other.1

points

This problem can modelled as a max-min model where we maximize the smallest distance of 10 points.

Model

maxΔs.t.Δdi,j+M(1xixj),i<jixi=kxi{0,1}\max \Delta\\ s.t.\\ \Delta \le d_{i,j} + M(1- x_i x_j), \forall i \lt j \\ \sum_i x_i = k\\ x_i \in \{0,1\}\\

M is an upper bound on the distances:

M=maxijdi,jM = \max_{i\le j} d_{i,j}\\

The quadratic part xixjx_i x_j of the constraint can be linearized:

Δdi,j+M(1xixj)Δdi,j+M(1xi)+M(1xj)\Delta \le d_{i,j} + M(1-x_i x_j)\\ \Leftrightarrow\\ \Delta \le d_{i,j} + M(1-x_i) + M(1- x_j)\\

Implementation

The Pyomo implementation of the constraints follows the mathematical formulation closely:

model.k_c = Constraint(
    model.I, rule=lambda model, i: sum(model.x[i] for i in model.I) == K
)

M = max(model.distance[i, j] for i in model.I for j in model.I if i <= j)
def delta_c(model, i, j):
    if i < j:
        return model.delta <= model.distance[i, j] + M * (1 - model.x[i] + M * (1 - model.x[j]))
    else:
        return Constraint.Skip
model.delta_c = Constraint(model.I, model.I, rule=delta_c)

Result

The CBC solver has no problems with this model:

Solution: 3.2049493
Number of constraints : 1275
Number of variables : 51
Duration: 00:00:03

solution

This model ensures, that the 10 solution points are at least distanced by 3.2 units, which is the maximal possible common distance for 10 points of the given set of 50.

The takeaway here is the linearization of the quadratic constraint. Another trick in our belt.


always be building, by sysid.

© 2022