2020-10-31 11:28:29 +00:00
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# TP2 - Modeling using docplex"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 1. The `docplex` python package\n",
"\n",
"`DOcplex` is a python package developped by IBM — It provides easy-to-use API for IBM solvers Cplex and Cpoptimizer.\n",
"\n",
"DOcplex documentation for mathematical programming can be found here: http://ibmdecisionoptimization.github.io/docplex-doc/#mathematical-programming-modeling-for-python-using-docplex-mp-docplex-mp"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 2. Solving TSP using `docplex`\n",
"\n",
"### 2.1. TSP model using `docplex`\n",
"\n",
"**Exercice:** Using `docplex`, create a model for the travelling salesman problem using the MTZ or Flow formulation and compare them."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from docplex.mp.model import Model\n",
"import tsp.data as data\n",
"\n",
"distances = data.grid42\n",
"N = len(distances)\n",
"\n",
"tsp = Model(\"TSP\")\n",
"tsp.log_output = True\n",
"\n",
"... # TODO\n",
"\n",
"solution = tsp.solve()\n",
2023-11-16 19:22:47 +00:00
"print(\"z* =\", solution.objective_value)\n"
2020-10-31 11:28:29 +00:00
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The largest set of distances contains 42 nodes, and should be easily solved by `docplex`.\n",
"\n",
"### 2.2. Generating random TSP instances\n",
"\n",
"**Question:** What method could you implement to generate a realistic set of distances for $n$ customers?\n",
"\n",
"**Exercice:** Implement the method proposed above and test it."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"import scipy.spatial.distance\n",
"\n",
"\n",
"def generate_distances(n: int):\n",
" ... # TODO\n",
"\n",
"\n",
"from docplex.mp.model import Model\n",
"\n",
"distances = generate_distances(50)\n",
"print(distances)\n",
"\n",
"N = len(distances)\n",
"\n",
"tsp = Model(\"TSP\")\n",
"tsp.log_output = True\n",
"\n",
"... # TODO: Copy your model from the first question here.\n",
"\n",
"solution = tsp.solve()\n",
2023-11-16 19:22:47 +00:00
"print(\"z* =\", solution.objective_value)\n"
2020-10-31 11:28:29 +00:00
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 3. Solving Warehouse Allocation using Benders decomposition with `docplex`\n",
"\n",
"### 3.1. The warehouse problem\n",
"\n",
"A company needs to supply a set of $n$ clients and needs to open new warehouses (from a\n",
"set of $m$ possible warehouses).\n",
"Opening a warehouse $j$ costs $f_j$ and supplying a client $i$ from a warehouse $j$ costs $c_{ij}$ per supply unit.\n",
"Which warehouses should be opened in order to satisfy all clients while minimizing the total cost?\n",
"\n",
"### 3.2. Solving the warehouse problem with a single ILP\n",
"\n",
"- $y_{j} = 1$ if and only if warehouse $j$ is opened.\n",
"- $x_{ij}$ is the fraction supplied from warehouse $j$ to customer $i$.\n",
"\n",
"$\n",
"\\begin{align}\n",
" \\text{min.} \\quad & \\sum_{i=1}^{n} \\sum_{j=1}^{m} c_{ij} x_{ij} + \\sum_{j=1}^{m} f_{j} y_{j} & \\\\\n",
" \\text{s.t.} \\quad & \\sum_{j=1}^{m} x_{ij} = 1, & \\forall i \\in\\{1,\\ldots,n\\}\\\\\n",
" & x_{ij} \\leq y_{j}, & \\forall i\\in\\{1,\\ldots,n\\},\\forall j\\in\\{1,\\ldots,m\\}\\\\\n",
" & y_{j} \\in \\left\\{0,~1\\right\\}, & \\forall j \\in \\left\\{1,~\\ldots,~m\\right\\}\\\\\n",
" & x_{ij} \\geq 0, & \\forall i \\in \\left\\{1,~\\ldots,~n\\right\\}, \\forall j \\in \\left\\{1,~\\ldots,~m\\right\\}\n",
"\\end{align}\n",
"$\n",
"\n",
"\n",
"**Exercice:** Implement the ILP model for the warehouse allocation problem and test it on the given instance."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from docplex.mp.model import Model\n",
"\n",
"# We will start with a small instances with 3 warehouses and 3 clients:\n",
"N = 3\n",
"M = 3\n",
"\n",
"# Opening and distribution costs:\n",
"f = [20, 20, 20]\n",
"c = [[15, 1, 2], [1, 16, 3], [4, 1, 17]]\n",
"\n",
"wa = Model(\"Warehouse Allocation\")\n",
"wa.log_output = True\n",
"\n",
"... # TODO: Model for the warehouse allocation.\n",
"\n",
"solution = wa.solve()\n",
2023-11-16 19:22:47 +00:00
"print(\"z* =\", solution.objective_value)\n"
2020-10-31 11:28:29 +00:00
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### 3.3. Benders' decomposition for the Warehouse Allocation problem\n",
"\n",
"We are going to use Benders' decomposition to solve the Warehouse Allocation problem. \n",
"\n",
"#### Dual subproblem\n",
"\n",
"$\n",
"\\begin{align*}\n",
"\\text{max.} \\quad & \\sum_{i=1}^{n} v_{i} - \\sum_{i=1}^{n}\\sum_{j=1}^{m} \\bar{y}_{j} u_{ij} & \\\\\n",
"\\text{s.t.} \\quad & v_{i} - u_{ij} \\leq c_{ij}, & \\forall i\\in\\{1,\\ldots,n\\},\\forall j\\in\\{1,\\ldots,m\\}\\\\\n",
" & v_{i} \\in\\mathbb{R},\\ u_{ij} \\geq 0 & \\forall i \\in\\{1,\\ldots,n\\}, \\forall j\\in\\{1,\\ldots,m\\}\n",
"\\end{align*}\n",
"$\n",
"\n",
"#### Master problem\n",
"\n",
"$\n",
"\\begin{align*}\n",
" \\text{min.} \\quad & \\sum_{j=1}^{m} f_j y_j + z & \\\\\n",
" \\text{s.t.} \\quad & z \\geq \\sum_{i=1}^{n}v_i^p - \\sum_{i=1}^{n} \\sum_{j=1}^{m} u_{ij}^p y_j, & \\forall p\\in l_1\\\\\n",
" & 0 \\geq \\sum_{i=1}^{n}v_i^r - \\sum_{i=1}^{n} \\sum_{j=1}^{n} u_{ij}^r y_j, & \\forall r\\in l_2\\\\\n",
" & y_{j} \\in\\{0,1\\}, & \\forall j\\in\\{1,\\ldots,m\\}\n",
"\\end{align*}\n",
"$"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Exercice:** Implement the method `create_master_problem` that creates the initial master problem (without feasibility or optimality constraints) for the warehouse allocation problem.\n",
"\n",
"<div class=\"alert alert-info alert-block\">\n",
"\n",
"You can use `print(m.export_as_lp_string())` to display a textual representation of a `docplex` model `m`.\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from docplex.mp.model import Model\n",
"from docplex.mp.linear import Var\n",
"from docplex.mp.constr import AbstractConstraint\n",
"from typing import List, Sequence, Tuple\n",
"\n",
"\n",
"def create_master_problem(\n",
" N: int, M: int, f: Sequence[float], c: Sequence[Sequence[float]]\n",
") -> Tuple[Model, Var, Sequence[Var]]:\n",
" \"\"\"\n",
" Creates the initial Benders master problem for the Warehouse Allocation problem.\n",
"\n",
" Args:\n",
" N: Number of clients.\n",
" M: Number of warehouses.\n",
" f: Array-like containing costs of opening warehouses.\n",
" c: 2D-array like containing transport costs from client to warehouses.\n",
"\n",
" Returns:\n",
" A 3-tuple containing the docplex problem, the z variable and the y variables.\n",
" \"\"\"\n",
"\n",
" wa = Model(\"Warehouse Allocation - Benders master problem\")\n",
"\n",
" ... # TODO\n",
"\n",
" return wa, z, y\n",
"\n",
"\n",
"# Check your method:\n",
"wa, z, y = create_master_problem(N, M, f, c)\n",
2023-11-16 19:22:47 +00:00
"print(wa.export_as_lp_string())\n"
2020-10-31 11:28:29 +00:00
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Exercice:** Implement the method `add_optimality_constraints` that add optimality constraints to the Benders master problem. "
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def add_optimality_constraint(\n",
" N: int,\n",
" M: int,\n",
" wa: Model,\n",
" z: Var,\n",
" y: Sequence[Var],\n",
" v: Sequence[float],\n",
" u: Sequence[Sequence[float]],\n",
") -> List[AbstractConstraint]:\n",
" \"\"\"\n",
" Adds an optimality constraints to the given Warehouse Allocation model\n",
" using the given optimal values from the Benders dual subproblem.\n",
"\n",
" Args:\n",
" N: Number of clients.\n",
" M: Number of warehouses.\n",
" wa: The Benders master problem (docplex.mp.model.Model).\n",
" z: The z variable of the master problem.\n",
" y: The y variables of the master problem.\n",
" v: The optimal values for the v variables of the Benders dual subproblem.\n",
" u: The optimal values for the u variables of the Benders dual subproblem.\n",
"\n",
" Return: The optimality constraint added.\n",
" \"\"\"\n",
2023-11-16 19:22:47 +00:00
" return ... # TODO\n"
2020-10-31 11:28:29 +00:00
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Exercice:** Implement the method `add_feasibility_constraints` that add feasibility constraints to the Benders master problem. "
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def add_feasibility_constraints(\n",
" N: int,\n",
" M: int,\n",
" wa: Model,\n",
" z: Var,\n",
" y: Sequence[Var],\n",
" v: Sequence[float],\n",
" u: Sequence[Sequence[float]],\n",
") -> List[AbstractConstraint]:\n",
" \"\"\"\n",
" Adds an optimality constraints to the given Warehouse Allocation model\n",
" using the given optimal values from the Benders dual subproblem.\n",
"\n",
" Args:\n",
" - N: Number of clients.\n",
" - M: Number of warehouses.\n",
" - wa: The Benders master problem (docplex.mp.model.Model).\n",
" - z: The z variable of the master problem.\n",
" - y: The y variables of the master problem.\n",
" - v: The extreme rays for the v variables of the Benders dual subproblem.\n",
" - u: The extreme rays for the u variables of the Benders dual subproblem.\n",
"\n",
" Returns:\n",
" The feasibility constraint added.\n",
" \"\"\"\n",
" # TODO:\n",
2023-11-16 19:22:47 +00:00
" return\n"
2020-10-31 11:28:29 +00:00
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Exercice:** Implement the method `create_dual_subproblem` that, given a solution `y` of the master problem, create the corresponding Benders dual subproblem.\n",
"\n",
"$\n",
"\\begin{align*}\n",
"\\text{max.} \\quad & \\sum_{i=1}^{n} v_{i} - \\sum_{i=1}^{n}\\sum_{j=1}^{m} \\bar{y}_{j} u_{ij} & \\\\\n",
"\\text{s.t.} \\quad & v_{i} - u_{ij} \\leq c_{ij}, & \\forall i\\in\\{1,\\ldots,n\\},\\forall j\\in\\{1,\\ldots,m\\}\\\\\n",
" & v_{i} \\in\\mathbb{R},\\ u_{ij} \\geq 0 & \\forall i \\in\\{1,\\ldots,n\\}, \\forall j\\in\\{1,\\ldots,m\\}\n",
"\\end{align*}\n",
"$"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from docplex.mp.model import Model\n",
"\n",
"\n",
"def create_dual_subproblem(\n",
" N: int, M: int, f: Sequence[float], c: Sequence[Sequence[float]], y: Sequence[int]\n",
") -> Tuple[Model, Sequence[Var], Sequence[Sequence[Var]]]:\n",
" \"\"\"\n",
" Creates a Benders dual subproblem for the Warehouse Allocation problem corresponding\n",
" to the given master solution.\n",
"\n",
" Args:\n",
" N: Number of clients.\n",
" M: Number of warehouses.\n",
" f: Array-like containing costs of opening warehouses.\n",
" c: 2D-array like containing transport costs from client to warehouses.\n",
" y: Values of the y variables from the Benders master problem.\n",
"\n",
" Returns:\n",
" A 3-tuple containing the docplex problem, the v variable and the u variables.\n",
" \"\"\"\n",
"\n",
" dsp = Model(\"Warehouse Allocation - Benders dual subproblem\")\n",
"\n",
" # We disable pre-solve to be able to retrieve a meaningful status in the main\n",
" # algorithm:\n",
" dsp.parameters.preprocessing.presolve.set(0)\n",
"\n",
" ... # TODO\n",
"\n",
" return dsp, v, u\n",
"\n",
"\n",
"# Check your method (assuming y = [1 1 1 ... 1]):\n",
"dsp, v, u = create_dual_subproblem(N, M, f, c, [1] * M)\n",
2023-11-16 19:22:47 +00:00
"print(dsp.export_as_lp_string())\n"
2020-10-31 11:28:29 +00:00
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Exercice:** Using the methods you implemented, write the Benders decomposition algorithm for the warehouse allocation problem.\n",
"\n",
"<div class=\"alert alert-block alert-info\">\n",
"\n",
"The `get_extreme_rays` function can be used to retrieve the extreme rays associated with an unbounded solution of the dual subproblem.\n",
" \n",
"</div>\n",
"\n",
"<div class=\"alert alert-block alert-info\">\n",
" \n",
"You can use `model.get_solve_status()` to obtain the status of the resolution and compare it to members of `JobSolveStatus`:\n",
" \n",
"```python\n",
"if model.get_solve_status() == JobSolveStatus.OPTIMAL_SOLUTION:\n",
" pass\n",
"```\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from docplex.mp.model import Model\n",
"from docplex.util.status import JobSolveStatus\n",
"\n",
"\n",
"def get_extreme_rays(\n",
" N: int, M: int, model: Model, v: Sequence[Var], u: Sequence[Sequence[Var]]\n",
") -> Tuple[Sequence[float], Sequence[Sequence[float]]]:\n",
" \"\"\"\n",
" Retrieves the extreme rays associated to the dual subproblem.\n",
"\n",
" Args:\n",
" N: Number of clients.\n",
" M: Number of warehouses.\n",
" model: The Benders dual subproblem model (docplex.mp.model.Model).\n",
" v: 1D array containing the v variables of the subproblem.\n",
" u: Either a 2D array of a tuple-index dictionary containing the u variables\n",
" of the subproblem.\n",
"\n",
" Returns:\n",
" A 2-tuple containing the list of extreme rays correspondig to v,\n",
" and the 2D-list of extreme rays corresponding to u.\n",
" \"\"\"\n",
" ray = model.get_engine().get_cplex().solution.advanced.get_ray()\n",
"\n",
" if isinstance(u, dict):\n",
"\n",
" def get_uij(i, j):\n",
" return u[i, j]\n",
"\n",
" else:\n",
"\n",
" def get_uij(i, j):\n",
" return u[i][j]\n",
"\n",
" return (\n",
" [ray[v[i].index] for i in range(N)],\n",
" [[ray[get_uij(i, j).index] for j in range(M)] for i in range(N)],\n",
" )\n",
"\n",
"\n",
"# We will start with a small instances with 3 warehouses and 3 clients:\n",
"N = 3\n",
"M = 3\n",
"\n",
"# Opening and distribution costs:\n",
"f = [20, 20, 20]\n",
"c = [[15, 1, 2], [1, 16, 3], [4, 1, 17]]\n",
"\n",
"# We stop iterating if the new solution is less than epsilon\n",
"# better than the previous one:\n",
"epsilon = 1e-6\n",
"\n",
"wa, z, y = create_master_problem(N, M, f, c)\n",
"\n",
"n = 0\n",
"while True:\n",
"\n",
" # Print iteration:\n",
" n = n + 1\n",
" print(\"Iteration {}\".format(n))\n",
"\n",
" ... # TODO\n",
"\n",
2023-11-16 19:22:47 +00:00
"print(\"Done.\")\n"
2020-10-31 11:28:29 +00:00
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### 3.4. Generating instances for the Warehouse Allocation problem\n",
"\n",
"**Exercice:** Using the TSP instances contained in `tsp.data` or the `generate_distances` method, create instances for the warehouse allocation problem with randomized opening costs."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"<div class=\"alert alert-block alert-danger\"></div>"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
2023-11-16 19:22:47 +00:00
"version": "3.9.13"
2020-10-31 11:28:29 +00:00
}
},
"nbformat": 4,
"nbformat_minor": 4
}
