INSA
/
ILP5A
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity
ILP5A/tp1-branch-and-bound.ipynb

1127 lines
38 KiB
Plaintext
Raw Permalink Normal View History

Initial commit.
2020-10-31 11:28:29 +00:00
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# TP1 - Branch & Bound, Cut-Generation for the TSP"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 1. The `minilp` python package\n",
"\n",
"The `minilp` python package is a small python package that allows us to easily model (integer) linear program. The package comes with an interface to common linear programming solvers (`cplex`, `scipy`) but no integer linear programming solver.\n",
"\n",
"<div class=\"alert alert-info\">\n",
"\n",
"The `minilp` module has been implemented for these practical sessions so you will not find any relevant documentation on Google. \n",
"To get help on the module and its components, you can use the built-in `help` function:\n",
" \n",
" help(minilp)\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import minilp\n",
"\n",
"lp = minilp.problem(\"My first LP problem\")\n",
"\n",
"# Create two continuous variables within [0, 4]:\n",
"x1, x2 = lp.continuous_var_list(2, 0, 4)\n",
"\n",
"# Add constraints:\n",
"lp.add_constraint(-3 * x1 + 4 * x2 <= 7)\n",
"lp.add_constraint(2 * x2 <= 5)\n",
"lp.add_constraint(6 * x1 + 4 * x2 <= 25)\n",
"lp.add_constraint(2 * x1 - x2 <= 6)\n",
"\n",
"# Set the objective function:\n",
"lp.set_objective(\"max\", x1 + 2 * x2)\n",
"\n",
"# Solve the problem:\n",
"res = lp.lp_solve()\n",
"print(res)\n",
"print(\"x1 = {:.4f}, x2 = {:.4f}\".format(res.get_value(x1), res.get_value(x2)))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The `minilp` package also allows you to modelise simple integer linear programs.\n",
"\n",
"<div class=\"alert alert-info alert-block\">\n",
" \n",
"There are also `binary_var` and `binary_var_list` method to create binary variable (integer variable constrained to 0 and 1).\n",
"The `_list` methods returns standard python `list` object, so you can combine them to create multi-dimensional lists of\n",
"minilp variables.\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import minilp\n",
"\n",
"lp = minilp.problem(\"My first ILP problem\")\n",
"\n",
"# Create two integer variables within [0, 4]:\n",
"x1, x2 = lp.integer_var_list(2, 0, 4)\n",
"\n",
"# Add constraints:\n",
"lp.add_constraint(-3 * x1 + 4 * x2 <= 7)\n",
"lp.add_constraint(2 * x2 <= 5)\n",
"lp.add_constraint(6 * x1 + 4 * x2 <= 25)\n",
"lp.add_constraint(2 * x1 - x2 <= 6)\n",
"\n",
"# Set the objective function:\n",
"lp.maximize(x1 + 2 * x2)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"While `minilp` allows you to model (mixed-)integer linear program, it does not provide a integer linear program solver &mdash; the `lp_solve` method will always solve the linear relaxation of the problem.\n",
"\n",
"<div class=\"alert alert-info alert-block\">\n",
" \n",
"You can use a different solver by passing a different object to the `lp_solve` method. Available solvers are:\n",
" \n",
"- `minilp.solvers.scipy` &mdash; A solver based on the `scipy` module, wildly available.\n",
"- `minilp.solvers.docplex` &mdash; A solver based on the `docplex` module, which requires a valid CPLEX installation.\n",
" \n",
"The default solver used is `docplex` if available, otherwize `scipy`.\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"res = lp.lp_solve()\n",
"print(res)\n",
"print(\"x1 = {:.4f}, x2 = {:.4f}\".format(res.get_value(x1), res.get_value(x2)))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The `minilp` package allows you to modelise `<=`, `>=` or `==` (in)equalities. You can create linear expression by simply adding, substracting or multiplying values (`int` or `float`) and variables or existing expressions. You can use the standard python `sum` to sum a bunch of expressions or variables, and the `minilp.dot` function to compute the [dot product](https://en.wikipedia.org/wiki/Dot_product) of two vectors.\n",
"\n",
"**Exercice:** Complete the following code to create a simple model for the knapsack problem.\n",
"\n",
"<div class=\"alert alert-block alert-info\">\n",
" \n",
"Make your program as generic as possible, i.e., if `N` or `K` changes, you should not have to\n",
"modify the code of your problem.\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import minilp\n",
"\n",
"N = 5\n",
"p = [1, 4, 5, 3, 5] # profits\n",
"w = [3, 4, 3, 5, 9] # weights\n",
"K = 10 # capacity\n",
"\n",
"assert N == len(w) and N == len(p)\n",
"\n",
"# A simple knapsack\n",
"kp = minilp.problem(\"Simple knapsack\")\n",
"\n",
"# TODO: Create variables, add constraints and set the objective.\n",
"x = ...\n",
"\n",
"# We can solve the linear relaxation:\n",
"res = kp.lp_solve()\n",
"print(res)\n",
"print(res.get_values(x))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 2. Generic Branch & Bound\n",
"\n",
"The purpose of the following section is to implement a generic branch-and-bound procedure based on a `minilp.problem` instance.\n",
"\n",
"### 2.1. Node structure\n",
"\n",
"<div class=\"alert alert-block alert-danger\">\n",
" \n",
"Do not overthink this section! \n",
"The three methods you have to implement in this section (`create_root_node`, `create_root_node` and `is_empty`) are\n",
"all **one-liner** and can be implemented with **a single `return` statement**.\n",
" \n",
"</div>\n",
"\n",
"We are going to use a simple list to represent the current set of leaf nodes in the branch-and-bound tree.\n",
"For each node, we are only going to store the list of separation constraints of the node (and not an actual \n",
"problem).\n",
"\n",
"For instance, assume a problem with variables $x_1$ and $x_2$, we two separation constraints, our node\n",
"could be created with:\n",
"\n",
"```python\n",
"node = [x1 <= 2, x2 >= 4]\n",
"```"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import minilp\n",
"from typing import List, Optional\n",
"\n",
"# minilp.cons is the Python type representing minilp constraints:\n",
"Node = List[minilp.cons]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Question:** What does the root node of the branching tree contain?\n",
"\n",
"**Exercice:** Implement the function `create_root_node()` that creates the root node of the branching tree."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def create_root_node() -> Node:\n",
" \"\"\"\n",
" Creates the root node for a branch-and-bound procedure.\n",
"\n",
" Returns:\n",
" The root node for a branch-and-bound procedure.\n",
" \"\"\"\n",
" return ... # TODO"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"For the sake of simplicity, we are going to process node in a first-in/first-out order.\n",
"\n",
"**Exercice:** Implement the function `get_next_node(nodes)` that extract the next node to process from the tree.\n",
"\n",
"<div class=\"alert alert-info alert-block\">\n",
" \n",
"You can check the `list.pop` method from python.\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def extract_next_node(nodes: List[Node]) -> Node:\n",
" \"\"\"\n",
" Extracts the next node to process from the given list of nodes. The node\n",
" is removed from the given list before being returned.\n",
"\n",
" Args:\n",
" nodes: Current list of nodes to extract the next node from.\n",
"\n",
" Return:\n",
" The next node process.\n",
" \"\"\"\n",
" return ... # TODO"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Exercice:** Implement the function `is_empty(nodes)` that returns `True` if there are no more nodes to process in the list."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def is_empty(nodes: List[Node]):\n",
" \"\"\"\n",
" Checks if there are any nodes remaining to process in the given list.\n",
"\n",
" Args:\n",
" - nodes: The list of nodes to check.\n",
"\n",
" Returns:\n",
" True if there are no more nodes, False otherwise.\n",
" \"\"\"\n",
" return ... # TODO"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### 2.2. Comparing `minilp.result`\n",
"\n",
"The `minilp.problem.lp_solve` method returns an instance of `minilp.result`. The following method compare two `minilp.result` objects for a given \n",
"problem and returns `True` if the left one is better.\n",
"\n",
"<div class=\"alert alert-block alert-warning\">\n",
" \n",
"Don't forget to execute the cell even if you do not have to modify it in order to have access to the `compare_solution` function.\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def compare_solutions(\n",
" problem: minilp.problem, l: minilp.result, r: minilp.result\n",
") -> bool:\n",
" \"\"\"\n",
" Compares the two given solutions, returning True if the left\n",
" one is better than the right one for the current problem.\n",
"\n",
" Args:\n",
" problem: The problem for which the solution are compared.\n",
" l, r: The two solutions (minilp.result) to compare.\n",
"\n",
" Returns:\n",
" True if the left solution is better than the right one, or\n",
" if the right solution has no solution.\n",
" \"\"\"\n",
" if problem.isnan(r.objective):\n",
" return True\n",
" if problem.sense == \"min\":\n",
" return l.objective < r.objective\n",
" return l.objective > r.objective"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### 2.3. Finding bound variable\n",
"\n",
"In order to execute the branch-and-bound procedure, we must be able to find a variable to separate the problem on.\n",
"For this notebook, we are always going to separate on the first non-integer variable.\n",
"\n",
"<div class=\"alert alert-block alert-info\">\n",
" \n",
"You can access the list of variables in the problem with `minilp.problem.variables` and their values in the\n",
"solution with `minilp.result.get_value` or `minilp.result.get_values`.\n",
" \n",
"</div>\n",
"\n",
"<div class=\"alert alert-block alert-warning\">\n",
" \n",
"There are variables in the problem that **do not have** to be integers in the solution (e.g., transport\n",
"variables in the warehouse allocation problem).\n",
"These variables should not be checked against when looking for a non-integer value.\n",
"\n",
"You can retrieve the type of a variable (`int` or `float`) using `minilp.expr.var.category`.\n",
" \n",
"</div>\n",
"\n",
"**Exercire:** Implement the `get_first_non_integral(problem, result)` method that, given a problem and solution, returns\n",
"the first variable of the problem that is not integral in the result (and should be), or `None` if there is no such variable\n",
"\n",
"<div class=\"alert alert-block alert-info\">\n",
" \n",
"In order to solve the linear relaxation of the problem, multiple matrix operations have to be performed. These \n",
"operations are numerically imprecise, thus it is common to find near-integral values (`1.0000001` or `2.9999999`).\n",
"Such values should be considered integral for the purpose of our algorithms, as long as their distance to the\n",
"nearest integer is less than `eps` (default to `1e-6`).\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def find_first_non_integral(\n",
" problem: minilp.problem, result: minilp.result, eps: float = 1e-6\n",
") -> Optional[minilp.var]:\n",
" \"\"\"\n",
" Retrieves the first integer variable in the given problem\n",
" whose value is not integral in the given solution.\n",
"\n",
" Args:\n",
" problem: The problem to find a variable from.\n",
" result: A solution of the problem.\n",
" eps: The maximum allowed distance to consider a variable integral.\n",
"\n",
" Returns:\n",
" The first variable (minilp.expr.var) whose value is not integral,\n",
" or None if no such variable exists.\n",
" \"\"\"\n",
" ... # TODO"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### 2.4. Relaxation, separation, iteration\n",
"\n",
"**Exercice:** Implement the `iterate` method below that performs a single iteration of the branch-and-bound algorithm, i.e., extract a node, \n",
"solve its relaxation, and then update the current best solution (return) or separate the problem. \n",
"\n",
"<div class=\"alert alert-block alert-info\">\n",
" \n",
"You can use the `minilp.problem.add_constraints` and `minilp.problem.del_constraints` methods to add or remove constraints from\n",
"a `minilp.problem`.\n",
" \n",
"Do not forget to return the `minilp.result` you found if it is your new best integer solution!\n",
" \n",
"</div>\n",
"\n",
"<div class=\"alert alert-block alert-danger\">\n",
" \n",
"Again, **do not overthink** this section! While `iterate` is the core function of the branch-and-bound algorith, it is a \n",
"pretty simply and short function that can be implemented in about 10 lines.\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import math\n",
"\n",
"\n",
"def iterate(\n",
" problem: minilp.problem,\n",
" current_best: minilp.result,\n",
" nodes: List[Node],\n",
" solver: minilp.solver = minilp.solvers.get_default_solver(),\n",
" eps: float = 1e-6,\n",
") -> Optional[minilp.result]:\n",
" \"\"\"\n",
" Performs an iteration of the branch-and-bound algorithm.\n",
"\n",
" Args:\n",
" problem: Problem to perform an iteration for.\n",
" current_best: Current best known feasible (integral) solution. Instance of minilp.result.\n",
" nodes: Current list of nodes to update.\n",
" solver: Solver to use to solve the linear relaxation.\n",
" eps: The maximum allowed distance to consider a variable integral.\n",
"\n",
" Returns:\n",
" The new best solution if the current one was improved, or None if no new solution\n",
" was found.\n",
" \"\"\"\n",
"\n",
" ... # TODO"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### 2.5. The branch-and-bound algorithm\n",
"\n",
"The cell below defines the global `branch_and_bound` procedure."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import datetime as dt\n",
"\n",
"\n",
"def log_solution(\n",
" niterations: int, nodes: List[Node], res: minilp.result, new_best: bool = False\n",
"):\n",
" \"\"\"\n",
" Logs the given solution with time information.\n",
"\n",
" Args:\n",
" niterations: Number of iterations.\n",
" nodes: List of nodes.\n",
" res: Solution (minilp.result) to log.\n",
" new_best: Indicates if this solution is new best solution (integer).\n",
" \"\"\"\n",
" print(\n",
" \"{} {:5d} {:5d} {:9g}{}\".format(\n",
" dt.datetime.now().strftime(\"%T\"),\n",
" niterations,\n",
" len(nodes),\n",
" res.objective,\n",
" \"*\" if new_best else \"\",\n",
" )\n",
" )\n",
"\n",
"\n",
"def branch_and_bound(\n",
" problem: minilp.problem,\n",
" lp_solver: minilp.solver = minilp.solvers.get_default_solver(),\n",
" eps: float = 1e-6,\n",
" log_frequency: int = 10,\n",
") -> minilp.result:\n",
" \"\"\"\n",
" Applies a branch-and-bound algorithm to solve the given problem.\n",
"\n",
" Args:\n",
" problem: A minilp.problem instance corresponding to the problem to solve.\n",
" lp_solver: Solver to use to solve the linear relaxation.\n",
" eps: The maximum allowed distance to consider a variable integral.\n",
" log_frequency: Number of iterations between two log (not including exceptional log).\n",
"\n",
" Returns:\n",
" A minilp.result containing the result of running the branch-and-bound\n",
" algorithm.\n",
" \"\"\"\n",
"\n",
" print(\n",
" \"B&B using {} to solve linear relaxation\".format(lp_solver.__class__.__name__)\n",
" )\n",
"\n",
" # Insert the first node in the list:\n",
" nodes = []\n",
" nodes.append(create_root_node())\n",
"\n",
" # Current best result (unknown):\n",
" current_best = minilp.result()\n",
"\n",
" # Counter for iterations:\n",
" nb_iterations = 0\n",
"\n",
" while not is_empty(nodes):\n",
"\n",
" # Increment counter and log.\n",
" if nb_iterations % log_frequency == 0:\n",
" log_solution(nb_iterations, nodes, current_best)\n",
"\n",
" nb_iterations += 1\n",
"\n",
" # Iterate:\n",
" new_best = iterate(problem, current_best, nodes, lp_solver, eps)\n",
"\n",
" # Check if we have a new current best:\n",
" if new_best is not None:\n",
" current_best = new_best\n",
" log_solution(nb_iterations, nodes, current_best, True)\n",
"\n",
" # Return the best solution found (if any).\n",
" return current_best"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Exercice:** Use the `branch_and_bound` method to solve the knapsack instance defined at the beginning of the notebook."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import minilp\n",
"\n",
"N = 5\n",
"p = [1, 4, 5, 3, 5] # profits\n",
"w = [3, 4, 3, 5, 9] # weights\n",
"K = 10 # capacity\n",
"\n",
"assert N == len(w) and N == len(p)\n",
"\n",
"# A simple knapsack\n",
"kp = minilp.problem(\"Simple knapsack\")\n",
"\n",
"# TODO: Create variables, add constraints and set the objective.\n",
"x = kp.binary_var_list(N)\n",
"kp.add_constraint(kp.dot(x, w) <= K)\n",
"\n",
"kp.maximize(kp.dot(x, p))\n",
"\n",
"# We can solve the linear relaxation:\n",
"res = branch_and_bound(kp)\n",
"print(res)\n",
"print(res.get_values(x))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Exercice:** Create other instances of the knapsack problem to reach the \"limits\" of your implementation &mdash; What is the largest instance you can solve in e.g. less than 5 seconds?\n",
"\n",
"<div class=\"alert alert-block alert-info\">\n",
" \n",
"You can use the `numpy.random` module to generate arrays of random integer or floating point values. \n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import minilp\n",
"import numpy as np\n",
"\n",
"N = 20\n",
"p = np.random.randint(50, 100, size=N)\n",
"w = np.random.randint(20, 40, size=N)\n",
"K = np.sum(w) // (1.5 + np.random.rand() * 1)\n",
"\n",
"assert N == len(w) and N == len(p)\n",
"\n",
"print(\"Knapsack problem with {} items and a capacity of {}.\".format(N, K))\n",
"print(\" Profits: {}\".format(p))\n",
"print(\" Weights: {}\".format(w))\n",
"\n",
"# A simple knapsack\n",
"kp = minilp.problem(\"Simple knapsack\")\n",
"\n",
"# TODO: Create variables, add constraints and set the objective.\n",
"x = kp.binary_var_list(N)\n",
"kp.add_constraint(kp.dot(x, w) <= K)\n",
"\n",
"kp.maximize(kp.dot(x, p))\n",
"\n",
"# We can solve the linear relaxation:\n",
"res = branch_and_bound(kp, log_frequency=5)\n",
"print(res)\n",
"print(res.get_values(x))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 3. The Travelling Salesman Problemn (TSP)\n",
"\n",
"\n",
"Given a list of $n$ cities and the distances $c_{ij}$ between each pair of cities, you want to find the shortest circuit that visits each city **exactly once** and comes back to the first visited city.\n",
"\n",
"<div class=\"alert alert-info alert-block\">\n",
"\n",
"The `tsp.data` packages contains grid of distances of various sizes (5, 6, 7, 8, 9, 10, 15, 17, 26, 42).\n",
"\n",
"</div>\n",
"\n",
"The goal of this section is to implement a cut-generation algorithm for the travelling salesman problem, using the `minilp` python package."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### 3.1. Creating a model for the TSP\n",
"\n",
"### 3.1.1. Relaxation of the TSP\n",
"\n",
"We call \"TSP relax\" the relaxation of the TSP problem that do not include constraints to eliminate subtours. \n",
"The model is given below:\n",
"\n",
" - $x_{ij}\\in\\{0, 1\\}$ &mdash; Binary variable indicating if we go directly from city $i$ to city $j$.\n",
"\n",
"$\n",
"\\begin{align}\n",
" \\text{min.} \\quad & \\sum_{i = 1}^{n} \\sum_{j=1}^{n} c_{ij}x_{ij} & \\\\\n",
" \\text{s.t.} \\quad & \\sum_{\\substack{j = 1\\\\ j \\neq i}}^{n} x_{ij} = 1, & \\forall i\\in \\left\\{1,~\\ldots,~n\\right\\}\\\\\n",
" & \\sum_{\\substack{i = 1\\\\i \\neq j}}^{n} x_{ij} = 1, & \\forall j\\in \\left\\{1,~\\ldots,~n\\right\\}\\\\\n",
" & x_{ij} \\in\\{0,1\\}, & \\forall i,j \\in \\{1,\\ldots,n\\}\n",
"\\end{align}\n",
"$\n",
" \n",
"**Exercice:** Create the `tsp_relax` function that creates a `minilp.problem` instance corresponding to\n",
"the relaxed TSP problem."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from typing import Sequence, Tuple\n",
"import minilp\n",
"\n",
"\n",
"def tsp_relax(\n",
" distances: Sequence[Sequence[float]], name: str = \"TSP relax\"\n",
") -> Tuple[minilp.problem, Sequence[Sequence[minilp.var]]]:\n",
" \"\"\"\n",
" Create a 'relaxed' model for the TSP. A relaxed includes all the standard\n",
" constraints of the TPS, but not the subtours constraints.\n",
"\n",
" Args:\n",
" distances: The matrix of distances (diagonal is zeros).\n",
" name: The name of the model.\n",
"\n",
" Returns:\n",
" A tuple (model, vars) where model is the TSP model and vars is the matrix\n",
" of x variables.\n",
" \"\"\"\n",
"\n",
" N = len(distances)\n",
"\n",
" tsp = minilp.problem(name)\n",
"\n",
" ... # TODO\n",
"\n",
" # Returns both the problem and the variables:\n",
" return tsp, x"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"<div class=\"alert alert-info alert-block\">\n",
" \n",
"You can visualize a `minilp` model by using `repr(model)` or simply writting the variable name at the end of a cell.\n",
"\n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import tsp.data as data\n",
"\n",
"tsp, x = tsp_relax(data.grid5)\n",
"tsp"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### 3.1.2. MTZ formulation for the TSP\n",
"\n",
"We call \"TSP relax\" the relaxation of the TSP problem that do not include constraints to eliminate subtours. \n",
"The model is given below:\n",
"\n",
" - $x_{ij}\\in\\{0, 1\\}$ &mdash; Binary variable indicating if we go directly from city $i$ to city $j$.\n",
" - $u_{i} \\in\\{1,\\ldots,n\\}$ &mdash; Subtour elimination variables: $u_i$ is the position of\n",
" city $i$ in the tour.\n",
"\n",
"$\n",
"\\begin{align}\n",
" \\text{min.} \\quad & \\sum_{i = 1}^{n} \\sum_{j=1}^{n} c_{ij}x_{ij} & \\\\\n",
" \\text{s.t.} \\quad & \\sum_{\\substack{j = 1\\\\ j \\neq i}}^{n} x_{ij} = 1, & \\forall i\\in \\left\\{1,~\\ldots,~n\\right\\}\\\\\n",
" & \\sum_{\\substack{i = 1\\\\i \\neq j}}^{n} x_{ij} = 1, & \\forall j\\in \\left\\{1,~\\ldots,~n\\right\\}\\\\\n",
" & u_1 = 1 & \\\\\n",
" & 2 \\leq u_i \\leq n, & \\forall i\\in\\{2,\\ldots,n\\} \\\\\n",
" & u_i - u_j +1 \\leq (n-1)(1 - x_{ij}), &\\forall i, j \\in \\{2,\\ldots,n\\}\\\\\n",
" & x_{ij} \\in\\{0,1\\},\\ u_{i}\\in\\mathbb{N} & \\forall i,j \\in \\{1,\\ldots,n\\}\n",
"\\end{align}\n",
"$\n",
" \n",
"**Exercice:** Create the `tsp_mtz` function that creates a `minilp.problem` instance corresponding to\n",
"the MTZ formulation of the TSP.\n",
"\n",
"<div class=\"alert alert-info alert-block\">\n",
" \n",
"Use the `tsp_relax` method to initialize the problem with the basic TSP constraints and retrieve the $x$ variables.\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import minilp\n",
"\n",
"\n",
"def tsp_mtz(\n",
" distances: Sequence[Sequence[float]], name: str = \"TSP MTZ\"\n",
") -> Tuple[minilp.problem, Tuple[Sequence[Sequence[minilp.var]], Sequence[minilp.var]]]:\n",
" \"\"\"\n",
" Create a MTZ model for the TSP.\n",
"\n",
" Args:\n",
" distances: The matrix of distances (diagonal is zeros).\n",
" name: The name of the model.\n",
"\n",
" Returns:\n",
" A tuple (model, (x, u)) where model is the TSP model and (x, u) are the\n",
" variables of the model (x is a 2D array, u is a 1D array).\n",
" \"\"\"\n",
"\n",
" N = len(distances)\n",
"\n",
" # Use tsp_relax to create the basic model:\n",
" tsp, x = tsp_relax(distances, name)\n",
"\n",
" ... # TODO\n",
"\n",
" # Returns both the problem and the variables:\n",
" return tsp, (x, u)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### 3.1.3. Flow formulation for the TSP\n",
"\n",
"We call \"TSP relax\" the relaxation of the TSP problem that do not include constraints to eliminate subtours. \n",
"The model is given below:\n",
"\n",
" - $x_{ij}\\in\\{0, 1\\}$ &mdash; Binary variable indicating if we go directly from city $i$ to city $j$.\n",
" - $y_{ij} \\in\\mathbb{R}_*^+$ &mdash; Subtour elimination variables: $y_{ij}$ is the *flow* on arc $(i,j)$, each city\n",
" produces one unit of flow.\n",
"\n",
"$\n",
"\\begin{align}\n",
" \\text{min.} \\quad & \\sum_{i = 1}^{n} \\sum_{j=1}^{n} c_{ij}x_{ij} & \\\\\n",
" \\text{s.t.} \\quad & \\sum_{\\substack{j = 1\\\\ j \\neq i}}^{n} x_{ij} = 1, & \\forall i\\in \\left\\{1,~\\ldots,~n\\right\\}\\\\\n",
" & \\sum_{\\substack{i = 1\\\\i \\neq j}}^{n} x_{ij} = 1, & \\forall j\\in \\left\\{1,~\\ldots,~n\\right\\}\\\\\n",
" & \\sum_{j=2}^{n} y_{1j} = 1 & \\\\\n",
" & \\sum_{j=1}^{n} y_{ij} = \\sum_{j=1}^{n} y_{ji} + 1, & \\forall i\\in\\{2,\\ldots,n\\} \\\\\n",
" & y_{ij} \\leq n x_{ij}, &\\forall i, j \\in \\{1,\\ldots,n\\}\\\\\n",
" & x_{ij} \\in\\{0,1\\}, y_{ij}\\in\\mathbb{R}_*^{+} & \\forall i,j \\in \\{1,\\ldots,n\\}\n",
"\\end{align}\n",
"$\n",
" \n",
"**Exercice:** Create the `tsp_flow` function that creates a `minilp.problem` instance corresponding to\n",
"the flow formulation of the TSP.\n",
"\n",
"<div class=\"alert alert-info alert-block\">\n",
" \n",
"Use the `tsp_relax` method to initialize the problem with the basic TSP constraints and retrieve the $x$ variables.\n",
" \n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import minilp\n",
"\n",
"\n",
"def tsp_flow(\n",
" distances: Sequence[Sequence[float]], name: str = \"TSP Flow\"\n",
") -> Tuple[\n",
" minilp.problem,\n",
" Tuple[Sequence[Sequence[minilp.var]], Sequence[Sequence[minilp.var]]],\n",
"]:\n",
" \"\"\"\n",
" Create a Flow model for the TSP.\n",
"\n",
" Args:\n",
" distances: The matrix of distances (diagonal is zeros).\n",
" name: The name of the model.\n",
"\n",
" Returns:\n",
" A tuple (model, (x, y)) where model is the TSP model and (x, y) are the\n",
" variables of the model (x and y are 2D arrays).\n",
" \"\"\"\n",
"\n",
" N = len(distances)\n",
"\n",
" # Use tsp_relax to create the basic model:\n",
" tsp, x = tsp_relax(distances, name)\n",
"\n",
" ... # TODO\n",
"\n",
" # Returns both the problem and the variables:\n",
" return tsp, (x, y)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### 3.2. Solving small instances of the TSP problem\n",
"\n",
"**Exercice:** Using the `tsp_mtz`, `tsp_flow` and `branch_and_bound` functions you implemented, solve the **small** TSP instances found in `tsp.data`.\n",
"\n",
"**Question:** How large are the instances you are able to solve in a reasonable amount of time?"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import tsp.data as data\n",
"\n",
"distances = data.grid7\n",
"\n",
"print(\"=== MTZ ===\")\n",
"tsp, (x, u) = tsp_mtz(distances)\n",
"r = branch_and_bound(tsp)\n",
"print(r)\n",
"\n",
"print(\"=== Flow ===\")\n",
"tsp, (x, y) = tsp_flow(distances)\n",
"r = branch_and_bound(tsp, eps=1e-8, log_freqeuncy=2)\n",
"print(r)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### 3.3. Cut generation algorithm for the TSP\n",
"\n",
"We are going to implement a cut-generation procedure using the following TSP formulation:\n",
"\n",
" - $x_{ij}\\in\\{0, 1\\}$ &mdash; Binary variable indicating if we go directly from city $i$ to city $j$.\n",
"\n",
"$\n",
"\\begin{align}\n",
" \\text{min.} \\quad & \\sum_{i = 1}^{n} \\sum_{j=1}^{n} c_{ij}x_{ij} & \\\\\n",
" \\text{s.t.} \\quad & \\sum_{\\substack{j = 1\\\\ j \\neq i}}^{n} x_{ij} = 1, & \\forall i\\in \\left\\{1,~\\ldots,~n\\right\\}\\label{tsp2:leave}\\\\\n",
" & \\sum_{\\substack{i = 1\\\\i \\neq j}}^{n} x_{ij} = 1, & \\forall j\\in \\left\\{1,~\\ldots,~n\\right\\}\\label{tsp2:enter}\\\\\n",
" & \\sum_{i\\in S} \\sum_{j \\notin S} x_{ij} \\geq 1 & \\mathcal{S} \\subset \\{1, \\ldots, n\\},\\ S \\ne \\emptyset \\label{tsp4:subtour:0} \\\\\n",
" & x_{ij} \\in\\{0,1\\}, & \\forall i,j \\in \\{1,\\ldots,n\\}\n",
"\\end{align}\n",
"$\n",
"\n",
"Where $\\mathcal{S}$ is the set of subtours in the graph. \n",
"The generation procedure is as follow:\n",
"\n",
"1. We start with an empty set $\\mathcal{S} = \\emptyset$.\n",
"2. We solve the problem (completely) using our `branch_and_bound` procedure.\n",
"3. We find all the subtours in the solution:\n",
" - If there is a single subtours, we have found the optimal solution.\n",
" - Otherwize, we add the subtours to the set $\\mathcal{S}$ and we got back to 2."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### 3.3.1. Generating subtour constraints\n",
"\n",
"**Exercice:** Implement the `add_subtour_constraints` method that, given a list of subtours (a list of list of integers, e.g., `[1, 2, 4]` is a subtour going through nodes 1, 2 and 4),\n",
"create subtour constraints and add them to the given problem."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def add_subtour_constraints(\n",
" tsp: minilp.problem,\n",
" x: Sequence[Sequence[minilp.var]],\n",
" subtours: Sequence[Sequence[int]],\n",
") -> Sequence[minilp.cons]:\n",
" \"\"\"\n",
" Generates and adds subtours constraints for the given problem.\n",
"\n",
" Args:\n",
" tsp: The current TSP problem.\n",
" x: The x variables (2D-array) of the TSP formulation.\n",
" subtours: The list of subtours. A subtour is a list of int containing the ID of the node.\n",
"\n",
" Returns:\n",
" The generated constraints.\n",
" \"\"\"\n",
" ... # TODO\n",
"\n",
"\n",
"# We can test the method by generated a relaxed-TSP and asking constraints for random subtours:\n",
"tsp, x = tsp_relax(data.grid5)\n",
"add_subtour_constraints(tsp, x, [[1, 2, 3], [0, 4]])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### 3.3.2. Finding subtours in a solution\n",
"\n",
"**Exercice:** Implement the method `find_subtours` that given a 2D-array `x` of integer (0 or 1), returns a list of subtours in `x`."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def find_subtours(x: Sequence[Sequence[bool]]) -> Sequence[Sequence[int]]:\n",
" \"\"\"\n",
" Extracts subtours from the given 2D-array.\n",
"\n",
" Args:\n",
" x: A two-dimensional array corresponding to the x variable in the TSP formulation, where\n",
" x[i][j] is 1 if arc (i, j) is used.\n",
"\n",
" Returns:\n",
" A list of subtours, where each subtour is a list.\n",
" \"\"\"\n",
" ... # TODO\n",
"\n",
"\n",
"# We can check the method by using a custom x array corresponding\n",
"# to two subtours (0 -> 2 -> 0) and (1 -> 4 -> 3 -> 1):\n",
"find_subtours(\n",
" [\n",
" [0, 0, 1, 0, 0], # 0 -> 2\n",
" [0, 0, 0, 0, 1], # 1 -> 4\n",
" [1, 0, 0, 0, 0], # 2 -> 0\n",
" [0, 1, 0, 0, 0], # 3 -> 1\n",
" [0, 0, 0, 1, 0], # 4 -> 3\n",
" ]\n",
")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### 3.3.3. Branch-and-cut for the TSP\n",
"\n",
"**Exercice:** Complete the `tsp_branch_and_cut` method below."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from typing import Callable\n",
"\n",
"\n",
"def tsp_branch_and_cut(\n",
" distances: Sequence[Sequence[float]],\n",
" ilp_solver: Callable[[minilp.problem], minilp.result] = branch_and_bound,\n",
"):\n",
" \"\"\"\n",
" Solves the given TSP instance using a branch-and-cut with the given solver.\n",
"\n",
" Args:\n",
" distances: Transport costs for the TSP.\n",
" ilp_solver: Function that can be called as ilp_solver(ilp) to solve integer linear program.\n",
"\n",
" Returns:\n",
" A minilp.result solution for the given TSP instance.\n",
" \"\"\"\n",
" # Create the relaxation:\n",
" N = len(distances)\n",
" tsp, x = tsp_relax(distances, \"TSP Branch & Cut\")\n",
"\n",
" while True:\n",
"\n",
" # Solve the problem:\n",
" res = ilp_solver(tsp)\n",
"\n",
" ... # TODO\n",
"\n",
" return res"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Exercice:** Test your `tsp_branch_and_cut` implementation."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import tsp.data as data\n",
"\n",
"res = tsp_branch_and_cut(data.grid5)"
]
},
{
"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.7.2"
}
},
"nbformat": 4,
"nbformat_minor": 4
}