Two-type VRP
I am trying to formulate a VRP for two types of vehicles and a single depot. V1 represents, say EVs, and V2 represents gas-fueled vehicles. V1-type vehicles consume less fuel/energy in the first M time units. I have the following Gurobipy programming. The problem is that I get solutions with incomplete/disconnected tours, mainly because a vehicle traverses i to j and the reverse to satisfy the connectivity constraints. Although I added the MTZ sub-tour elimination, I still end up receiving such solutions. When I introduce the constraint: x[v,i,j] + x[v,j,i] == 1, I get an infeasible solution. What is wrong in this formulation?
Here are the set, parameter, and variable definitions. Set I represents all customer locations, Depot_loc is the depot location. All locations:= I \cup Depot_loc = J. V1 is the set of EV vehicles, V2 is the set of gas vehicles, and V1 \cup V2 = V. T[i,j] represents the travel time from i to j, and is asymmetric. U is the UB of travel time for each vehicle. Variables are x, y, t, t_xi, t_kappa. x[v,i,j] = 1 if v takes i,j. y[v,j] is the sub-tour elimination variable. Variable t[v] records the travel time for V2 type vehicle. t_xi[v] + t_kappa[v] records the travel time for V1 type vehicle, where t_xi[v] accounts for the cheaper first M time units and t_kappa[v] for the rest.
m = Model("Problem")
x = {}; y = {}
for v in V:
for i in J:
for j in J:
x[v,i,j] = m.addVar(vtype=GRB.BINARY, name="x%s"%str([v,i,j]))for v in V:
for j in I:
y[v,j] = m.addVar(vtype=GRB.CONTINUOUS, name = "y%s"%str([v,j]))
t = {v: m.addVar(vtype=GRB.CONTINUOUS, name="t%s"%str([v])) for v in V2}
t_xi = {v: m.addVar(vtype=GRB.CONTINUOUS, name="t_xi%s"%str([v])) for v in V1}
t_kappa = {v: m.addVar(vtype=GRB.CONTINUOUS, name="t_kappa%s"%str([v])) for v in V1}
m.setObjective(quicksum(First_cost*t_xi[v] +
Second_cost*t_kappa[v] for v in V1) +
quicksum(Third_cost*t[v] for v in V2), GRB.MINIMIZE)
m.addConstrs((quicksum(x[v,Depot_loc,j] for j in I) == 1 for v in V), name = "every_vehicle_starts_at_depot")
m.addConstrs((quicksum(x[v,i,Depot_loc] for i in I) == 1 for v in V), name = "every_vehicle_ends_at_depot")
m.addConstrs((quicksum(x[v,i,j] for v in V for j in J if i != j) >= 1 for i in I),name = "every_node_to_be_visited")
m.addConstrs((quicksum(T[i,j]*x[v,i,j] for i in J for j in J if i != j) <= U for v in V), name = "limited_travel_time")
m.addConstrs((t_xi[v]+t_kappa[v] == quicksum(T[i,j]*x[v,i,j]
for i in J for j in J if i != j) for v in V1), name = "travel_time_for_V1")
m.addConstrs((t[v] == quicksum(T[i,j]*x[v,i,j] for i in J for j in J if i != j) for v in V2), name = "travel_time_for_V2")
m.addConstrs((t_xi[v] <= M for v in V1), name = "limit_electric_travel")
m.addConstrs((quicksum(x[v,i,j] for i in I if i != j) == quicksum(x[v,j,j_prime] for j_prime in I if j_prime != j) for v in V for j in I), name = "connectivity")
m.addConstrs((y[v,j] >= y[v,i]-len(J)*(1-x[v,i,j])
for v in V for i in I for j in I if i != j), name = "subtour_elimination")
### BEGIN SOLVING ###
m.setParam('TimeLimit', 720000); m.setParam('MIPgap', 0.0005)
m.update(); m.optimize(); print("\n\n")
### END SOLVING ###
-
Official comment
This post is more than three years old. Some information may not be up to date. For current information, please check the Gurobi Documentation or Knowledge Base. If you need more help, please create a new post in the community forum. Or why not try our AI Gurobot?.
Post is closed for comments.
Comments
1 comment