Can't realize CVRPTW (Capacitated Vehicle Routing Problem with Time Windows) in IBM ILOG CPLEX with adequate solution.
My .mod code:
int n = ...; // Set of vertices.
int nVehicles = ...; // Number of vehicles.
range N = 0..n-1; // Range of all points (including depot).
range C = 1..n-1; // Range of all clients (excluding depot).
range V = 0..nVehicles-1; // Range of all vehicles.
int c[N][N] = ...; // Transportation cost.
dvar boolean x[N][N][V]; // 1 if the k-th vehicle goes from i to j, otherwise 0.
dvar int+ Z[V];
int Q = ...; // Vehicle capacity.
int d[C] = ...; // Client demand.
dvar int+ s[N][V]; // Arrival time.
int t[N][N] = ...; // Transportation time.
int a[N] = ...; // Start of time windows.
int b[N] = ...; // End of time windows.
// 1. Objective function: minimize total transportation cost.
minimize sum(k in V) Z[k];
subject to {
// All vehicles depart from the depot at time 0.
//forall(k in V)
//s[0][k] == 0;
// 2. Each client is visited exactly once.
forall(i in C)
sum(k in V, j in N) x[i][j][k] == 1;
// 3. Capacity constraint.
forall(k in V)
sum(i in C, j in N) (d[i] * x[i][j][k]) <= Q;
// 4. Each vehicle starts its route from the depot.
forall(k in V)
sum(j in N) x[0][j][k] == 1;
// 5. Movement between clients.
forall(k in V, h in C)
sum(i in N) x[i][h][k] == sum(j in N) x[h][j][k];
// 6. Each vehicle returns to the depot.
forall(k in V)
sum(i in N) x[i][0][k] == 1; // Possibly 0 instead of n-1.
// 7. Departure time from one point and arrival time at another point.
forall(i in N, j in N, k in V)
s[i][k] + t[i][j] <= s[j][k] + (1 - x[i][j][k]) * 10000; // Large number to disable constraint when x[i][j][k] = 0.
// 8. Time window constraints.
forall(i in N, k in V)
a[i] <= s[i][k] <= b[i];
// 10. Constraint on using a specified number of vehicles.
sum(k in V, j in N) x[0][j][k] <= nVehicles;
// Constraint linking variable Z[k] with total transportation cost
forall(k in V)
Z[k] == sum(i in N, j in N) (c[i][j] * x[i][j][k]);
}
.dat code:
n = 6; // Set of vertices.
nVehicles = 4; // Number of vehicles.
// Transportation cost.
c = [
[0, 10, 15, 20, 30, 25],
[10, 0, 9, 14, 23, 18],
[15, 9, 0, 10, 12, 25],
[20, 14, 10, 0, 11, 15],
[30, 23, 12, 11, 0, 12],
[25, 18, 25, 15, 12, 0]
];
// Transportation time.
t = [
[0, 2, 3, 4, 5, 6],
[2, 0, 1, 3, 5, 7],
[3, 1, 0, 2, 4, 6],
[4, 3, 2, 0, 3, 5],
[5, 5, 4, 3, 0, 2],
[6, 7, 6, 5, 2, 0]
];
Q = 100; // Vehicle capacity.
d = [2, 3, 4, 2, 1]; // Demand.
// Time windows
a = [0, 3, 5, 7, 9, 12];
b = [100, 10, 12, 14, 16, 18];
In solution I get zeros in Z. When I change the c matrix diagonally on big numbers (e. g. 100000), I get 600000 100000 100000 100000.
The results I would like to get should look like this:
Sample Output
Suppose we have 4 vehicles and 5 clients. The results might look like this:
Route for Vehicle 1: Depot → Client 3 → Client 1 → Depot
Route for Vehicle 2: Depot → Client 2 → Client 4 → Client 5 → Depot
Route for Vehicle 3: Depot (unused)
Route for Vehicle 4: Depot (unused)
Total cost: 150 (e.g., total distance or route time)
Arrival times:
Vehicle 1: Client 3 — 5, Client 1 — 9
Vehicle 2: Client 2 — 4, Client 4 — 8, Client 5 — 12
starting point
.mod
using CP;
tuple Node {
key string nodeID;
int x;
int y;
}
{Node} nodes = ...;
tuple Visit {
key int visitID;
string nodeID;
int quantity;
int minTime;
int maxTime;
int dropTime;
};
{Visit} clientVisits = ...;
tuple Vehicle {
key string vehicleID;
string firstVisitID;
string lastVisitID;
int capacity;
int start;
int end;
}
{Vehicle} vehicles = ...;
int firstDepotVisitID = min(v in clientVisits) v.visitID - 1;
int lastDepotVisitID = max(v in clientVisits) v.visitID + 1;
Visit firstDepotVisit = <firstDepotVisitID,"depot",0,0,230,0>;
Visit lastDepotVisit = <lastDepotVisitID,"depot",0,0,230,0>;
{Visit} allVisits = clientVisits union
// Needs to be generalized for multiple depots
{firstDepotVisit, lastDepotVisit };
int xPerVisit[i in allVisits]=item(nodes,<i.nodeID>).x;
int yPerVisit[i in allVisits]=item(nodes,<i.nodeID>).y;
int horizon = max (v in allVisits) v.maxTime;
// Create transition distance
tuple triplet { int c1; int c2; int d; };
{triplet} Dist = {
<ord(allVisits,v1), ord(allVisits,v2),
ftoi(round(sqrt(pow(n2.x-n1.x,2)+pow(n2.y-n1.y,2))))>
| v1,v2 in allVisits, n1, n2 in nodes : v1.nodeID == n1.nodeID && v2.nodeID == n2.nodeID };
dvar interval visitInterval[v in clientVisits] in v.minTime..(v.maxTime+v.dropTime) size v.dropTime;
dvar interval wtvisitInterval [v in clientVisits] size v.dropTime..horizon;
dvar interval tvisitInterval [v in allVisits][veh in vehicles]
optional(v.visitID!=firstDepotVisitID && v.visitID!=lastDepotVisitID);
dvar sequence route[veh in vehicles] in all(v in allVisits) tvisitInterval[v][veh]
types all(v in allVisits) ord(allVisits,v);
dvar interval truck [veh in vehicles] optional;
execute {
cp.param.TimeLimit = 10
}
dexpr int nonTravelTime = sum(v in clientVisits) sizeOf(wtvisitInterval[v]);
dexpr float travelTime = sum(veh in vehicles) endOf(tvisitInterval[lastDepotVisit][veh]) - nonTravelTime;
dexpr int nbUsed = sum(veh in vehicles) presenceOf(truck[veh]);
dexpr int load[veh in vehicles] = sum(v in clientVisits) presenceOf(tvisitInterval[v][veh])*v.quantity;
minimize staticLex(nbUsed,travelTime);
subject to {
forall(v in clientVisits ) {
endAtEnd(wtvisitInterval[v], visitInterval[v]);
startBeforeStart(wtvisitInterval[v], visitInterval[v]);
}
forall(veh in vehicles) {
span (truck[veh], all(v in clientVisits) tvisitInterval[v][veh]);
noOverlap(route[veh], Dist); // Travel time
startOf(tvisitInterval[firstDepotVisit][veh])==0; // Truck t starts at time 0 from depot
last (route[veh],tvisitInterval[lastDepotVisit] [veh]); // Truck t returns at depot
load[veh] <= veh.capacity; // Truck capacity
}
forall(v in clientVisits)
alternative(wtvisitInterval[v], all(t in vehicles) tvisitInterval[v][t]); // Truck selection
}
execute {
writeln(nbUsed + " vehicles are used");
writeln("Total travelled distance is " + travelTime);
}
.dat
nodes = {<"depot" 35 35> <"node1" 41 49> <"node2" 35 17>
<"node3" 55 45> <"node4" 55 20> <"node5" 15 30>
<"node6" 25 39> <"node7" 20 50> <"node8" 10 43>
<"node9" 55 60> <"node10" 30 60> <"node11" 20 65>
<"node12" 50 35> <"node13" 30 25> <"node14" 15 10>
<"node15" 30 5> <"node16" 10 20> <"node17" 5 30>
<"node18" 20 40> <"node19" 15 60> <"node20" 45 65>};
clientVisits = {<1 "node1" 10 0 204 10>
<2 "node2" 7 0 202 10>
<3 "node3" 13 0 197 10>
<4 "node4" 19 149 159 10>
<5 "node5" 26 0 199 10>
<6 "node6" 3 0 208 10>
<7 "node7" 5 0 198 10>
<8 "node8" 9 95 105 10>
<9 "node9" 16 97 107 10>
<10 "node10" 16 0 194 10>
<11 "node11" 12 67 77 10>
<12 "node12" 19 0 205 10>
<13 "node13" 23 159 169 10>
<14 "node14" 20 0 187 10>
<15 "node15" 8 61 71 10>
<16 "node16" 19 0 190 10>
<17 "node17" 2 0 189 10>
<18 "node18" 12 0 204 10>
<19 "node19" 17 0 187 10>
<20 "node20" 9 0 188 10>};
vehicles = {<"vehicle1" "depot" "depot" 200 0 230>
<"vehicle2" "depot" "depot" 200 0 230>
<"vehicle3" "depot" "depot" 200 0 230>
<"vehicle4" "depot" "depot" 200 0 230>
<"vehicle5" "depot" "depot" 200 0 230>};