CODEWITHSUNDEEP
cplexoplmixed-integer-programming
Josh Cho
Josh Cho

Reputation: 159

Vehicle Routing with Pickup and Dropoffs with MIP

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];

Upvotes: 3

Views: 371

Answers (1)

Josh Cho
Josh Cho

Reputation: 159

OK, I have been thinking for a long time to solve this problem. And at last, I have solved it with reasonable amount of running time. I will share it for someone who perhaps wants to dispute what i have done or improve running time. Meanwhile, it will be awesome if someone could verify what i have done is indeed correct logically.

Oh, I have solved it an extension that i have planned with extra several constraints. I hope my comments are enough for understanding what i have done.

// Constants
int TotalTimeLimit      = ...;
int LoadTime            = ...;
int UnloadTime          = ...;
int TransitCost = ...;
int MaxTransit          = ...;

// Cities
int     n      = ...;
int     pickupN= ...;       //4
range   Cities  = 1..n;
range   Cities1 = 2..n-1;   //just for start/end restriction
range   Pickups = 2..pickupN+1; //2-5
range   Dropoffs= pickupN+2..n-1;   //7-17

// Edges -- sparse set
tuple       edge        {int i; int j;}
setof(edge) Edges       = {<i,j> | i,j in Cities};
int         dist[Edges] = ...;
int         time[Edges] = ...;

// Products
int     P        = ...;
range   Products = 1..P;

int Pickup[Cities][Products] = ...;
int Dropoff[Cities][Products] = ...;
int weight[Products] = ...;
int volume[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];

//Trucks
{string} Type = ...;
int     Max[Type] = ...;
int     c   = sum(t in Type) Max[t];
range   Cars= 1..c;
int     VolumeCapacity[Cars] = ...;
int     WeightCapacity[Cars] = ...;
int     NumberCapacity[Cars] = ...;
int     FixedCost[Cars]      = ...;
int     forbid[Cities][Cars] = ...;

// Decision variables
dvar boolean x[Cars][Edges];
dvar boolean y[Cars][Edges][Products];
dvar boolean z[Cars][Cities][Products];
dvar boolean isCarUsed[Cars];

// Cost function
dexpr float cost = sum (c in Cars) (isCarUsed[c]*FixedCost[c] + 
    (sum(e in Edges) (x[c][e]*(dist[e] + TransitCost) + (sum (p in Products) y[c][e][p] * weight[p]))))
    - 30 - sum(p in Products) weight[p];

// Objective
minimize cost;
subject to {
// Total pickups for each p equal the total sum of each colum of z in pickups
forall (p in Products)
    sum(i in Pickups, c in Cars) z[c][i][p] == sum(i in Pickups) Pickup[i][p];

// For each node, each car can pick at most the node's stock
forall(i in Pickups, p in Products, c in Cars)
    z[c][i][p] <= Pickup[i][p];

// For each car, picked up products should be smaller than car's capacities
forall (c in Cars, e in Edges) {
    sum(p in Products) y[c][e][p] <= NumberCapacity[c];
    sum(p in Products) y[c][e][p]*weight[p] <= WeightCapacity[c];
    sum(p in Products) y[c][e][p]*volume[p] <= VolumeCapacity[c];
}

// For each car and product, its picked up amount should equal dropped off amount
forall (c in Cars, p in Products)
    sum(i in Pickups) z[c][i][p] == sum(i in Dropoffs) z[c][i][p];

// Total dropoffs for each p equal the total sum of each column of z in dropoffs
forall (p in Products)
    sum(i in Dropoffs, c in Cars) z[c][i][p] == sum(i in Dropoffs) Dropoff[i][p];

// For each node, each car can drop at most its needs
forall(i in Dropoffs, p in Products, c in Cars)
    z[c][i][p] <= Dropoff[i][p];
    
// For each car, it cannot stay at one place
forall (c in Cars)
    sum(i in Cities) x[c][<i,i>] == 0;

// Prevents looping, must start at node 1 and end at node n
forall (c in Cars) {
    sum (e in Edges : e.i==n) x[c][e] == 0;
    sum (e in Edges : e.j==1) x[c][e] == 0;
}

// For each car, it cannot go back and forth
forall (c in Cars, i,j in Cities)
    x[c][<i,j>]+x[c][<j,i>] <= 1;

// Each city is linked with two other cities
forall (j in Cities1, c in Cars)
    sum (<i,j> in Edges) x[c][<i,j>] - sum (<j,k> in Edges) x[c][<j,k>] == 0;

// If car is used, must start at node 1 and end at node n
forall(c in Cars) {
    100*sum(e in Edges : e.i==1) x[c][e] >= sum (p in Products, i in Cities1) z[c][i][p];
    100*sum(e in Edges : e.j==n) x[c][e] >= sum (p in Products, i in Cities1) z[c][i][p];
    100*isCarUsed[c] >= sum(p in Products, i in Cities1) z[c][i][p];
}
forall (c in Cars) {
    sum(e in Edges : e.i==1) x[c][e] <= 1;
    sum(e in Edges : e.j==n) x[c][e] <= 1;
}
    
// For each car, it needs to cross nodes that picks or drops something
forall (c in Cars, i in Cities)
    100*sum (j in Cities) x[c][<i,j>] >= sum (p in Products) z[c][i][p];

// For each car, its transit count <= maxTransit
forall (c in Cars)
    sum(e in Edges) x[c][e] <= MaxTransit;

// For each car, and for each node, we need to check whether time took is <= totalTimeLimit
forall(c in Cars)
    sum(e in Edges : e.i in Pickups || e.i in Dropoffs) x[c][e]*time[e] 
        + LoadTime*sum(i in Pickups, p in Products) z[c][i][p] + UnloadTime*sum(i in Dropoffs, p in Products) z[c][i][p]
            <= TotalTimeLimit;

// Certain type of cars cannot go to certain city
forall (c in Cars, i in Cities)
    if (forbid[i][c] == 1)
        sum(j,k in Cities) (x[c][<i,j>] + x[c][<k,i>]) == 0;

// Keeps culmulated pacakges through each arcs
forall (c in Cars, i in Pickups, p in Products)
    sum(k in Cities) y[c][<i,k>][p] == sum(j in Cities) y[c][<j,i>][p] + z[c][i][p];
forall (c in Cars, i in Pickups, j in Cities)
    100*x[c][<i,j>] >= sum(p in Products) y[c][<i,j>][p];

// MTZ subtour elimination constraints
forall (c in Cars, j in Dropoffs, i,k in Cities, p in Products : j != i && j != k)
    y[c][<i,j>][p] - y[c][<j,k>][p] >= z[c][j][p] - 100*(1-x[c][<i,j>]);
};

Here is an example data. I made the setting such that we first pickup, and dropoff later. I also made 3 dummy variables for start and end, and transit position between pickups and dropoffs.

n = 12;
pickupN = 3;
dist = [
9999999
0
0
0
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
390
661
0
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
390
9999999
228
0
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
661
228
9999999
0
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
63
34
264
360
208
329
9999999
9999999
9999999
9999999
9999999
9999999
9999999
232
444
292
297
47
0
9999999
9999999
9999999
9999999
9999999
232
9999999
 29
232
444
292
0
9999999
9999999
9999999
9999999
9999999
444
29
9999999
249
402
250
0
9999999
9999999
9999999
9999999
9999999
292
232
249
9999999
495
352
0
9999999
9999999
9999999
9999999
9999999
297
444
402
495
9999999
154
0
9999999
9999999
9999999
9999999
9999999
47
292
250
352
154
9999999
0
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
9999999
];
time = [
9999999 
0 
0 
0 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
20 
52 
0 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
14 
9999999 
26 
0 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
43 
31 
9999999 
0 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
12 
31 
10 
17 
26 
26 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
21 
19 
20 
18 
7 
0 
9999999 
9999999 
9999999 
9999999 
9999999 
35 
9999999 
6 
13 
30 
28 
0 
9999999 
9999999 
9999999 
9999999 
9999999 
39 
23 
9999999 
38 
23 
38 
0 
9999999 
9999999 
9999999 
9999999 
9999999 
23 
32 
20 
9999999 
45 
17 
0 
9999999 
9999999 
9999999 
9999999 
9999999 
28 
35 
37 
24 
9999999 
24 
0 
9999999 
9999999 
9999999 
9999999 
9999999 
2 
39 
8 
37 
21 
9999999 
0 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999 
9999999
9999999 
9999999 
];

TotalTimeLimit  = 150;
LoadTime        = 10;
UnloadTime      = 10;

TransitCost = 10;
//DropoffTransitCost    = 10;
MaxTransit          = 10;

// Products
P = 12;
Pickup = [
//   1,2,3,4,5,6,7,8 9 101112 products
[0 0 0 0 0 0 0 0 0 0 0 0],//city1 pickstart
[0 0 0 1 0 0 1 0 1 1 0 0],//city2
[1 1 0 0 0 0 0 1 0 0 1 0],//city3
[0 0 1 0 1 1 0 0 0 0 0 1],//city4
[0 0 0 0 0 0 0 0 0 0 0 0],//city5 pickend & dropstart
[0 0 0 0 0 0 0 0 0 0 0 0],//city6 
[0 0 0 0 0 0 0 0 0 0 0 0],//city7
[0 0 0 0 0 0 0 0 0 0 0 0],//city8
[0 0 0 0 0 0 0 0 0 0 0 0],//city9
[0 0 0 0 0 0 0 0 0 0 0 0],//city10
[0 0 0 0 0 0 0 0 0 0 0 0],//city11
[0 0 0 0 0 0 0 0 0 0 0 0]//city12
];
Dropoff = [
[0 0 0 0 0 0 0 0 0 0 0 0],
[0 0 0 0 0 0 0 0 0 0 0 0],
[0 0 0 0 0 0 0 0 0 0 0 0],
[0 0 0 0 0 0 0 0 0 0 0 0],
[0 0 0 0 0 0 0 0 0 0 0 0],//city5 pickend & dropstart
[0 0 0 0 1 0 0 0 0 0 0 1],//city6 
[1 0 0 0 0 0 0 0 0 0 0 0],//city7
[0 1 0 0 0 0 1 0 0 0 0 0],//city8
[0 0 0 0 0 1 0 0 1 0 0 0],//city9
[0 0 0 1 0 0 0 0 0 1 0 0],//city10
[0 0 1 0 0 0 0 1 0 0 1 0],//city11
[0 0 0 0 0 0 0 0 0 0 0 0]//city12 dropend
];

weight = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12];
volume = [20, 500, 30, 20, 60, 100, 300, 50, 10, 400, 1, 15];

//Cars
Type ={"SmallTruck", "MiddleTruck", "BigTruck"};
Max = #[
SmallTruck : 2,
MiddleTruck : 2,
BigTruck : 2
]#;
VolumeCapacity = [100, 100, 500, 500, 1000, 1000];
WeightCapacity = [10, 10, 30, 30, 50, 50];
NumberCapacity = [4, 4, 6, 6, 12, 12];
FixedCost = [100, 100, 200, 200, 300, 300];

forbid = [ 
[0 0 0 0 0 0],//city1
[0 0 0 0 0 0],//city2
[0 0 0 0 0 0],//city3
[0 0 0 0 0 0],//city4 pickend & dropstart
[0 0 0 0 0 0],//city5 
[0 0 0 0 0 0],//city6
[0 0 0 0 0 0],//city7
[0 0 0 0 0 0],//city8
[0 0 0 0 0 0],//city9
[0 0 0 0 1 1],//city10
[0 0 0 0 0 0],//city11
[0 0 0 0 0 0]//city12
];

Upvotes: 1

Related Questions