Tuples (applied in routing models)

After completing this module the student should:

- Implement and solve a simple TSP model

Video

A simple TSP model

When first implementing a TSP model, many will start with a distance matrix. These use an implicit representation of the nodes to visit, in this case represented by the range N. The decision variables is one for each pair i,j in N, indicating if we visit j after we visit i, ie. if we use the edge going from i to j.

The model will look like this:

int n = 5;
range N = 1..n;
float dists[N][N] = [
  [0.0, 1.0, 1000.0, 4.0, 1000.0],
  [1.0, 0.0, 1.5, 1.0, 1000.0],
  [1000.0, 7.0, 0.0, 3.0, 1.0],
  [2.0, 1000.0, 2.0, 0.0, 1000.0],
  [1000.0, 1.0, 2.0, 5.0, 0.0]
  ];

// Do we go from i to j?
dvar boolean x[N][N];
dvar int u[2..n];

minimize
sum(i,j in N) dists[i][j]*x[i][j];
subject to {
  // Ignore edges to self
  forall(i in N)
    x[i][i] == 0;

  // Ignore 1000 edges
  forall(i,j in N : dists[i][j]>=1000)
    x[i][j] == 0;

  // Flow conservation (For each ingoing edge, one must go out)
  forall(i in N)
    sum(j in N) x[i][j] == sum(j in N) x[j][i];

  // Coverage constraint, ie. one edge going in must be chosen
  forall(i in N)
    sum(j in N) x[j][i] == 1;

  // Subtour elimination using Miller-Tuker-Zemlin formulation
  forall(i,j in N : i>1 && j>1)
    u[i] + 1 <= u[j] + 1000*(1 - x[i][j]);

}

The distance matrix contains entries for a nodes connection to itself. These should of course be ignored. It also contains distances of 1000 representing connections that should not be used. After that we implement flow conservation constraints: All nodes visited must also be left again. And the coverage constraint: One should go to a node exactly once, ie. select one edge going into each node.

Finally we need a subtour elimination constraint. If we try to solve the problem without these the following solution will be found: 1 -> 2 -> 4 (-> 1) and 3 -> 5 (-> 3). Thus forms two subtours, each of which fulfill the other constraints.

The subtour constraint fixes this by assigning each visit except the first a value u where the u value must be increasing when following edges. If two subtours were formed, the one not containing the node 1, would be breaking this constraint. In the example from before we would have u3 < u5 < u3 which is a contradiction invalidating the solution containing the subtour.

Adding the constraint we get the tour: 1 -> 4 -> 3 -> 5 -> 2 (-> 1) which is indeed valid.

Tuples

The distance matrix method is sub optimal in three different ways:

The last point becomes even more important as the size of the graph goes up, and the ration of usable edges to the total number of edges goes down.

What we would like to do is specify the edges that are valid manually. The requires us to couple three different values together. An origin node, a destination node and a distance. This is done using tuples, in this case name Arc.

tuple Arc{
  int fromNode;
  int toNode;
  float dist;
}

We can now specify the valid arcs manually:

{Arc} Arcs = {<1,4,4>,<2,1,1>,<2,3,1.5>,<2,4,1>,<3,2,7>,<3,4,3>,<3,5,1>,<4,1,2>,<4,3,2>,
  <5,2,1>,<5,3,2>,<5,4,5>};

The fist tuple <1,4,4> represents an edge going from node 1 to node 4 having a distance of 4.

The decisions we need to take is which of the edges to select (And we still need the us for subtour elimination):

dvar boolean x[Arcs];
dvar int u[2..n];

The objective function also need to be expressed in terms of the arcs:

minimize
  sum(a in Arcs) a.dist*x[a];

So too is true for the remaining constraints. The final model will look like this:

dvar boolean x[Arcs];
dvar int u[2..n];

minimize
  sum(a in Arcs) a.dist*x[a];
subject to {
  // Flow conservation (For each ingoing edge, one must go out)
  forall(i in N)
    sum(a in Arcs : a.toNode==i) x[a] == sum(a in Arcs : a.fromNode==i) x[a];

  // All nodes must be visited, ie. one edge going in must be choosen
  forall(i in N)
    sum(a in Arcs : a.toNode==i) x[a] == 1;

  // Subtour elemination using Miller-Tuker-Zemlin formulation
  forall(a in Arcs : a.toNode>1 && a.fromNode>1)
    u[a.toNode] + 1 <= u[a.fromNode] + 1000*(1 - x[a]);

}

execute showSolution {
  writeln("Solution");
  for(var a in Arcs){
    if(x[a]>0){
      writeln(a.fromNode+" -> "+a.toNode);
    }
  }
}

Note that some of the first constraints eliminating edges can be omitted, but that an execute environment has been added to facilitate reading the solution.

In case you were wondering if the edges have been typed in manually, they weren't. This code was used to generate them originally:

{Arc} Arcs = {<i,j,dists[i][j]>|i,j in N: j!=i && dists[i][j]<1000};