I am trying to solve a Vehicle Routing Problem with multiple pickups and dropoffs with multiple products carried by just one car. After solving this problem I am going to extend to multiple types of cars as well.
One special setting is that it has a starting point and an ending point which needs not be same. I assumed to be different and setted 1 and n to be the dummy nodes for start and end.
I partially used an example TSP code provided by IBM to solve subtour problem, and got help from internet to print out the optimal tour.
Since I need to find a optimal path that go through all the points. This is NP-hard. But as a first time using ILog, I would like to solve this problem using MIP for practice purpose.
I am having trouble with keeping track of picked up products and droped off products in each arc.
I am trying to minimize the cost of transportation which I set to be
// Decision variables
dvar boolean x[Edges]; //car goes through this arc?
dvar boolean y[Edges][Products]; //at each e, currently loaded products in the car
dvar boolean z[Cities][Products]; //at each cities, what products to load or unload?
// Cost function
// at each arc that car goes through (distance + sum(products in the car*their weights))
dexpr float cost = sum (e in Edges) x[e] *
(dist[e] + (sum (p in Products) y[e][p] * weight[p]));
y is a variable that associates each arc with currently loaded products. z accounts for what to load or unload in each nodes. Since there is only one car I don't think z is actually needed, but for a extension with multiple cars, I think this will be a good thing to have.
If some of these dvars are not necessary, please give me some insight! Below are setups.
// Cities
int n = ...;
range Cities = 1..n;
range Cities1 = 2..n-1; //just for start/end restriction
range Pickups = 2..3;
range Dropoffs= 4..n-1;
// Edges -- sparse set
tuple edge {int i; int j;}
setof(edge) Edges = {<i,j> | ordered i,j in Cities};
int dist[Edges] = ...;
// Products
int p = ...;
range Products = 1..p;
float Pickup[Cities][Products] = ...;
float Dropoff[Cities][Products] = ...;
float weight[Products] = ...;
// Products in pickup and dropoff sums up to equal amount
assert
forall(p in Products)
sum(o in Cities) Pickup[o][p] == sum(d in Cities) Dropoff[d][p];
tuple Subtour { int size; int subtour[Cities]; }
{Subtour} subtours = ...;
Any help on below restrictions will be very helpful. Especially on keeping track of loaded products along the route.
// Objective
minimize cost;
subject to {
// Each city is linked with two other cities
forall (j in Cities1)
sum (<i,j> in Edges) x[<i,j>] + sum (<j,k> in Edges) x[<j,k>] == 2;
// Must start at node 1 and end at node n
sum (e in Edges : e.i==1) x[e] == 1;
sum (e in Edges : e.j==n) x[e] == 1;
// no product remains at 1,n (not necessary?)
sum (p in Products, e in Edges : e.i==1) y[e][p] == 0;
sum (p in Products, e in Edges : e.j==n) y[e][p] == 0;
sum (p in Products) z[1][p] == 0;
sum (p in Products) z[n][p] == 0;
// must pickup all
forall (p in Products) {
sum(i in Pickups) z[i][p] == sum(i in Cities) Pickup[i][p];
sum(i in Dropoffs) z[i][p] == sum(i in Cities) Dropoff[i][p];
}
forall (i in Pickups, p in Products)
z[i][p] <= Pickup[i][p];
//tried to keep track of picked ups, but it is not working
forall (i in Pickups, j,k in Cities, p in Products : k < i < j)
y[<i,j>][p] == y[<k,i>][p] + z[i][p];
// forall (j in Cities, p in Products)
// ctDemand:
// sum(<i,j> in Edges) y[<i,j>][p] + sum(<j,i> in Edges) y[<j,i>][p] == z[j][p];
// tried keeping track of dropoffs. It is partially working but not sure of it
forall (i, k in Cities, j in Dropoffs, p in Products : i < j < k) (y[<i,j>][p] == 1)
=> y[<j,k>][p] == y[<i,j>][p] - Dropoff[j][p];
// Subtour elimination constraints.
forall (s in subtours)
sum (i in Cities : s.subtour[i] != 0)
x[<minl(i, s.subtour[i]), maxl(i, s.subtour[i])>] <= s.size-1;
};
And post processing to find the subtours
// POST-PROCESSING to find the subtours
// Solution information
int thisSubtour[Cities];
int newSubtourSize;
int newSubtour[Cities];
// Auxiliar information
int visited[i in Cities] = 0;
setof(int) adj[j in Cities] = {i | <i,j> in Edges : x[<i,j>] == 1} union
{k | <j,k> in Edges : x[<j,k>] == 1};
execute {
newSubtourSize = n;
for (var i in Cities) { // Find an unexplored node
if (visited[i]==1) continue;
var start = i;
var node = i;
var thisSubtourSize = 0;
for (var j in Cities)
thisSubtour[j] = 0;
while (node!=start || thisSubtourSize==0) {
visited[node] = 1;
var succ = start;
for (i in adj[node])
if (visited[i] == 0) {
succ = i;
break;
}
thisSubtour[node] = succ;
node = succ;
++thisSubtourSize;
}
writeln("Found subtour of size : ", thisSubtourSize);
if (thisSubtourSize < newSubtourSize) {
for (i in Cities)
newSubtour[i] = thisSubtour[i];
newSubtourSize = thisSubtourSize;
}
}
if (newSubtourSize != n)
writeln("Best subtour of size ", newSubtourSize);
}
main {
var opl = thisOplModel
var mod = opl.modelDefinition;
var dat = opl.dataElements;
var status = 0;
var it =0;
while (1) {
var cplex1 = new IloCplex();
opl = new IloOplModel(mod,cplex1);
opl.addDataSource(dat);
opl.generate();
it++;
writeln("Iteration ",it, " with ", opl.subtours.size, " subtours.");
if (!cplex1.solve()) {
writeln("ERROR: could not solve");
status = 1;
opl.end();
break;
}
opl.postProcess();
writeln("Current solution : ", cplex1.getObjValue());
if (opl.newSubtourSize == opl.n) {
// This prints the tour as a cycle
var c = 1; // current city
var lastc = -1; // city visited right before C
write(c);
while (true) {
var nextc = -1; // next city to visit
// Find the next city to visit. To this end we
// find the edge that leaves city C and does not
// end in city LASTC. We know that exactly one such
// edge exists, otherwise the solution would be infeasible.
for (var e in opl.Edges) {
if (opl.x[e] > 0.5) {
if (e.i == c && e.j != lastc) {
nextc = e.j;
break;
}
else if (e.j == c && e.i != lastc) {
nextc = e.i;
break;
}
}
}
// Stop if we are back at the origin.
if (nextc == -1) {
break;
}
// Write next city and update current and last city.
write(" -> ", nextc);
lastc = c;
c = nextc;
}
opl.end();
cplex1.end();
break; // not found
}
dat.subtours.add(opl.newSubtourSize, opl.newSubtour);
opl.end();
cplex1.end();
}
status;
}
Here is a sample dataset that I created. I hope my explanation makes sence to everyone! Thank you very much!!
n = 10;
dist = [
633
257
91
412
150
80
134
259
505
390
661
227
488
572
530
555
289
228
169
112
196
154
372
262
383
120
77
105
175
476
267
351
309
338
196
63
34
264
360
29
232
444
249
402
495
];
// Products
p = 8;
Pickup = [
// 1,2,3,4,5,6,7,8 products
[0 0 0 0 0 0 0 0],//city1
[0 1 0 1 0 1 1 0],//city2
[1 0 1 0 1 0 0 1],//city3
[0 0 0 0 0 0 0 0],//city4
[0 0 0 0 0 0 0 0],//city5
[0 0 0 0 0 0 0 0],//city6
[0 0 0 0 0 0 0 0],//city7
[0 0 0 0 0 0 0 0],//city8
[0 0 0 0 0 0 0 0],//city9
[0 0 0 0 0 0 0 0] //city10
];
Dropoff = [
[0 0 0 0 0 0 0 0],
[0 0 0 0 0 0 0 0],
[0 0 0 0 1 0 0 0],
[0 0 0 0 0 1 0 0],
[0 0 0 0 0 0 1 0],
[1 0 0 0 0 0 0 1],//city6
[0 1 0 0 0 0 0 0],//city7
[0 0 1 0 0 0 0 0],//city8
[0 0 0 1 0 0 0 0],//city9
[0 0 0 0 0 0 0 0] //city10
];
weight = [1, 2, 3, 4, 5, 6, 7, 8];
