518 lines
17 KiB
Plaintext
518 lines
17 KiB
Plaintext
{
|
|
"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",
|
|
"print(\"z* =\", solution.objective_value)\n"
|
|
]
|
|
},
|
|
{
|
|
"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",
|
|
"print(\"z* =\", solution.objective_value)\n"
|
|
]
|
|
},
|
|
{
|
|
"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",
|
|
"print(\"z* =\", solution.objective_value)\n"
|
|
]
|
|
},
|
|
{
|
|
"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",
|
|
"print(wa.export_as_lp_string())\n"
|
|
]
|
|
},
|
|
{
|
|
"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",
|
|
" return ... # TODO\n"
|
|
]
|
|
},
|
|
{
|
|
"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",
|
|
" return\n"
|
|
]
|
|
},
|
|
{
|
|
"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",
|
|
"print(dsp.export_as_lp_string())\n"
|
|
]
|
|
},
|
|
{
|
|
"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",
|
|
"print(\"Done.\")\n"
|
|
]
|
|
},
|
|
{
|
|
"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",
|
|
"version": "3.9.13"
|
|
}
|
|
},
|
|
"nbformat": 4,
|
|
"nbformat_minor": 4
|
|
}
|