AMPL Capacitated p-Median Problem with GCG

cpmp.ipynb Open In Colab Kaggle Gradient Open In SageMaker Studio Lab Hits

Description: Dantzig-Wolfe decomposition for Capacitated p-Median Problem with GCG

Tags: GCG, cpmp, amplpy, dantzig-wolfe decomposition, branch-price-and-cut, highlights

Notebook author: Jurgen Lentz <jurgenlentz26@gmail.com>

# Install dependencies
%pip install -q amplpy
# Google Colab & Kaggle integration
from amplpy import AMPL, ampl_notebook

ampl = ampl_notebook(
    modules=["gcg"],  # modules to install
    license_uuid="default",  # license to use
)  # instantiate AMPL object and register magics

Capacitated $p$-median problem

The goal of the capacitated $p$-median problem (CPMP) is to fix $p \in \mathbb{N}$ locations as medians from a given number of $n \in \mathbb{N}$ locations. Furthermore, we are given demands $q_{i} \geq 0$ and capacities $Q_{i} \geq 0$ for every location $i \in {1,\ldots,n}$ as well as distances $d_{ij} \geq 0$ for each pair of locations $i,j \in {1,\ldots,n}$. While fixing the medians, we also assign every location $i \in {1,\ldots,n}$ to exactly one median such that the demands of all locations assigned to the median can be satisfied by the capacity of the median. An assignment of locations to a median is called a cluster. The overall objective is to minimize the total distance between the medians and their assigned locations.

CPMP can be modeled as follows (see Ceselli et al.):

$$ \begin{aligned} \text{minimize} \quad &\sum_{i = 1}^{n} \sum_{j = 1}^{n} d_{i j} x_{i j} \ \text{subject to} \quad &\sum_{j = 1}^{n} x_{i j} = 1 \quad \forall i \in {1,…,n} \ &\sum_{i = 1}^{n} q_{i} x_{i j} \leq Q_{j} y_{j} \quad \forall j \in {1,…,n} \ &\sum_{j = 1}^{n} y_{j} = p \ &x_{i j} \in {0,1} \quad \forall i,j \in {1,…,n} \ &y_{j} \in {0,1} \quad \forall j \in {1,…,n} \end{aligned} $$

We will solve an CPMP instance using Dantzig-Wolfe decomposition with the GCG solver in AMPL. Therefore, we first model CPMP in AMPL (shown below).

%%ampl_eval
param n;
param p;

set I = 1..n ordered;
set J = 1..n ordered;

param d {I,J};
param w {I};
param C {J};

var x {I,J} binary;
var y {J} binary;

minimize Cost:  sum {i in I} sum {j in J} d[i,j] * x[i,j];

subject to Allocate {i in I}:
   sum {j in J} x[i,j] = 1;

subject to Capacity {j in J}:
   sum {i in I} w[i] * x[i,j] <= C[j] * y[j];

subject to NFacilities:
  sum{j in J} y[j] <= p;

We generate a small instance with 5 locations and 2 locations as medians.

ampl.param["n"] = 5
ampl.param["p"] = 2

d = [
    [0, 6, 54, 52, 19],
    [6, 0, 28, 75, 61],
    [54, 28, 0, 91, 40],
    [52, 75, 91, 0, 28],
    [19, 61, 40, 28, 0],
]
ampl.param["d"] = {(i, j): d[i - 1][j - 1] for i in range(1, 6) for j in range(1, 6)}

ampl.param["w"] = [14, 13, 9, 15, 6]
ampl.param["C"] = [39, 39, 39, 39, 39]

Automatic Mode in GCG with AMPL

We use AMPL to call the solver GCG to solve our CPMP instance automatically without providing any information about the Dantzig-Wolfe decomposition. Here, GCG detects different decompositions and chooses heuristically the best decomposition. Afterwards, the solver uses a branch-price-and-cut algorithm to solve it to optimality.

ampl.option["solver"] = "gcg"
ampl.option["gcg_options"] = "outlev=1"
ampl.solve()

Using a custom decomposition in GCG with AMPL

AMPL allows the users to create their own decomposition and forwards it to GCG using suffixes. Here, we assign the allocation and pmedian constraints to the master/linking constraints and each capacity constraint to a different pricing problem.

%%ampl_eval

suffix master IN, binary;
suffix block IN, integer;

let {i in I}
   Allocate[i].master := 1;

let NFacilities.master := 1;

let {j in J}
   Capacity[j].block := j;
ampl.solve()

Analogously, the user can fix constraints as pricing block or master constraint and variables as pricing block, master or linking variables. It is not needed to provide a complete decomposition. If the user provides a partial decomposition, GCG completes the decomposition by fixing only the left constraints and variables using its detection loop.

NOTE: The index of pricing block constraints/variables starts at 1.