Assignment 5: Cut Generation and Convex Optimisation COMP4691-8691 2020 School of Computer Science, ANU Guidelines • Answers to go in answers.pdf (replace the file with your pdf). • In order to get full marks for the implementation parts, your code should be clear, well commented, sensibly structured and correct. 1 Cut Generation for TSP [50 marks] The objective of the travelling salesman problem (TSP) is to find the shortest tour of a set of cities C := {1, 2, . . . , n}. The distance between each pair of cities is given by di,j ≥ 0 ∀i, j ∈ C where di,j = dj,i. A valid tour is a sequence of cities, representing the order they are visited, where each city appears exactly once. The total length of a tour is the sum of distances between adjacent cities in the sequence, where we also consider the first and last city in the sequence as being adjacent (the tour must finish in the same city as it started in). A key difficulty in TSP-like problems is simultaneously modelling: • the requirement for a single large tour covering all cities (instead of multiple subtours as shown in figure 1); and • the tour length. 1 2000 1500 1000 500 0 500 1000 1500 2000 2000 1500 1000 500 0 500 1000 Figure 1: Example of an invalid TSP solution, which consists of 10 subtours rather than a single large one. City Sequence Variables One encoding of the problem is to consider decision variables that indicate the position of each city in the sequence xu,i ∈ {0, 1} ∀u, i ∈ C. Here the index u represents the sequence position, and i the city. To ensure we only visit each city exactly once, and that we can only visit one city at a time in the sequence we can impose the following constraints:∑ i∈C xu,i = 1 ∀u ∈ C (1)∑ u∈C xu,i = 1 ∀i ∈ C (2) The downside of this representation is that we cannot easily represent the 2 tour length. Link Variables A second approach is to instead have decision variables for each possible link in the tour (each pair of cities) yi,j ∈ {0, 1} ∀i 6= j ∈ C where yi,j = yj,i. A tour then consists of a set of links, where each city must appear in exactly two links (one for the incoming trip, the other for the outgoing):∑ j∈C\i yi,j = 2 ∀i ∈ C (3) In this representation it is easy to write out the objective:∑ i∈C ∑ j>i∈C di,jyi,j (4) The sticking point with this representation is that our constraints to represent a valid tour are insufficient, as we haven’t prevented subtours from forming. 1.1 Combining Representations [5 marks] One option is to utilise both representations and implement constraints to ensure the two sets of variables are consistent. Mathematically formulate a set of linear constraints (should result in of the order of n3 constraints) to connect the city sequencing variables to the link variables. Unfortunately, because this requires n3 constraints the problem formulations quickly get too large to solve efficiently (we’ll solve a problem with n = 216 shortly), and these constraints tend to have a relatively weak linear relaxation. 1.2 Subtour Elimination and Strength [10 marks] An alternative to combining the representations, is to stick with just the link variables, and directly implement constraints on these variables to eliminate subtours. Let S represent the set of all city subsets with size greater than 2 and less than n: S ⊂ C, 2 < |S| < n ∀S ∈ S. Each of these city subsets 3 S ∈ S represents a potential subtour, that we want to prevent from forming by enforcing the constraint:∑ iyi,j ≤ |S| − 1 ∀S ∈ S (5) Come up with and explain a minimal example that demonstrates that a formulation that utilises link variables (3)–(4) and subtour elimination con- straints (5) results in a stronger linear relaxation than combining the city sequence variables (1)–(2) and link variables (3)–(4) with the constraints you derived in section 1.1. 1.3 Subtour Cut Generation [35 marks] The major downside of the this formulation with link variables and subtour elimination constraints is that the set S grows exponentially with the number of cities. In order to take advantage of the strength of these constraints, without blow- ing out the formulation size, we are going to generate the subtour elimination constraints on demand. The expectation is that we will hopefully only need to generate a small subset of the subtour elimination constraints in order to find the optimal tour. The steps of our overall algorithm are: 1. Initialise the subtour set S ′ to ∅. 2. Solve MILP consisting of (3)–(4), and (5) using S ′. 3. If the solution to the MILP is a single valid tour, return optimal. 4. Otherwise add subtours present in MILP solution to S ′. 5. Repeat from 2 with updated S ′. Explain why we can be confident we have the optimal to the overall TSP problem if we get a single valid tour at step 3 [5 marks]. Implement the above algorithm in tsp.py. Report the results including the resulting plot, the objective, the overall solve time, the number of iterations and the number of subtour cuts generated [30 marks]. 4 Hint: If you are implementing constraints or objectives in PuLP that have expressions with many terms, you will likely want to use the pulp functions pulp.lpSum and pulp.lpDot. These are many times faster that iterating through a loop or calling python’s inbuilt sum function. Hint: You will likely want to set a time limit for the subproblem MILP solves (i.e. using the CBC maxSeconds option through PuLP). The idea is that you don’t need to solve every subproblem optimally (only certain ones), and it might in fact be better to stop early, generate more cuts, and then start working on the next MILP. You sholud make this time limit adaptive, so that the algorithm remains complete (e.g., you don’t want to get stuck returning infeasible results or solutions that don’t lead to useful subtour cuts). Hint: The PuLP status output will be zero when the solver times out, regardless of whether or not an integer feasible solution has been found yet. Depending on how tight the time limits are you choose, you might need to check for yourself whether the returned value is integer feasible and satisfies the important tour requirements, before attempting to generate cuts from it. Hint: You should be able to solve the full 216 city instance in around 2 minutes. If you want small instances to work with, the script allows you to filter out cities below a certain size. E.g., pass --min-population 1e6 to limit the cities to those over 1 million in population. Note: If you are keen to push your implementation even further, you can take a look at some of the problem instances available here: http://comopt. ifi.uni-heidelberg.de/software/TSPLIB95/tsp/. 2 Dual Gradient Ascent [50 marks] Dual gradient ascent is an algorithm for solving strictly convex constrained optimisation problems. In this question you will implement the approach, and use it to solve several problems. The problems.py file contains two constrained optimisation problems hard- coded into classes Problem1 and Problem2. These classes offer functions to evaluate the objective, the gradient of the objective, the constraints and the constraint jacobian. This file also contains a more general class to represent quadratically constrained quadratic programs (QCQP). It can load a specific 5 problem instance from file. Your code should go in the file dual ascent.py. This can be called as a script. If you provide it with a JSON file it will solve the corresponding QCQP, otherwise it will solve the two hardcoded Problem1 and Problem2. 2.1 Gradient Descent [15 marks] As a subcomponent, we need to be able to first solve gradient descent. Im- plement the gradient descent algorithm in the gradient descent function, of the dual ascent.py file. Solve unconstrained versions of Problem1 and Problem2 (just ignore the constraints) using your gradient descent algorithm. Report the number of iterations required and the solution. Hint: You’ll likely want to make use of one of the rule-based step sizes, e.g., one of those by Barzilai and Borwein from the lectures. Hint: Implement a stopping criteria that checks whether the gradient is within a provided tolerance of being zero. 2.2 Dual Gradient Ascent [25 marks] Next implement the dual ascent algorithm in dual ascent, making use of gradient descent. Solve the constrained versions of Problem1 and Problem2 using your algo- rithm. Report the number of iterations required and the solution. Hint: Implement a stopping criteria based on satisfaction of the KKT con- ditions to within a tolerance. 2.3 More General Quadratic Problems [10 marks] Report the results you get by running your dual ascent algorithm on instances qcqp-{1,2,3,4}.json. Based on the dual variables you obtain, explain which constraints are active for each instance [5 marks]. You will likely struggle to solve the instance qcqp-3.json. If so explain why [5 marks]. 6 欢迎咨询51作业君
