【VRPCB】Python+Gurobi求解运输问题建模实践三

采用Python+Gurobi求解带有集群回程需求的VRPCB问题

1. 模型

1.1 VRPB问题介绍

带有回程需求的VRP问题（VRP with Backhauls,VRPB)最初由 Deif, . Bodin, LD. 于1984年出。相关学者将其进一步细分为四类子问题：
第一类VRPCB：

第一类VRPB(VRPCB)：
– 客户只能是linehaul客户和backhaul客户中的一类
– 车辆在访问backhaul客户集群前，必须先向linehaul客户集群交付货物
– 客户只能被访问一次
第二类VRPB（VRPBM)
– 不考虑客户集群限制，车辆可以交叉服务linehaul客户和backhaul客户
– 客户只能被访问一次
第三类VRPB（VRPPD）
– 客户可以同时是linehaul客户和backhaul客户
– 客户可以2次被访问，即：先访问几个客户交付货物，以便部分清空车辆装载。然后，再返回访问该客户提货
第四类VRPB（VRPPD)
– 客户可以同时是linehaul客户和backhaul客户
– 客户只能被访问一次

1.2 数学模型

这里暂时只研究VRPCB问题，借用 Toth, P., Vigo, D. 学者于1997给出的整数线性规模。

1.2.1 模型参数

在这里插入图片描述

1.2.2 数学模型

在这里插入图片描述

1.2.3 模型分解

上述模型在复现时有一定困难，尤其是约束（6）和约束（7）中的 $\mathcal{S, F}$
这里根据帖子（其实有一点小bug）提供的分解模型进行复现，模型如下：

在这里插入图片描述

2. 数据结构

原网址暂时想不起来了，回头再补。这里直接贴上代码读取的文件结构，各列依次是：节点id，节点横坐标，节点纵坐标，linehaul客户需求和backhaul客户需求

在这里插入图片描述

3. Gurobi源码

import copy
import csv
import math
import xlsxwriter
import matplotlib.pyplot as plt
from gurobipy import GRB,Model,quicksum# 读取文件
def read_data(filename):Q = {}L = []B = []V = []XY = {}depot = None# 读取网络节点、需求with open(filename, 'r') as f:node_reader = csv.DictReader(f)for row in node_reader:if float(row['linehaul']) == 0 and  float(row['backhaul']) == 0:depot = row['id']Q[depot] = 0elif float(row['linehaul']) > 0 and  float(row['backhaul']) == 0:L.append(row['id'])Q[row['id']] = float(row['linehaul'])elif float(row['linehaul']) == 0 and  float(row['backhaul']) > 0:B.append(row['id'])Q[row['id']] = float(row['backhaul'])V.append(row['id'])XY[row['id']] = (float(row['x_coord']), float(row['y_coord']))# 计算网络弧L0 = L + [depot]B0 = B + [depot]AL = [ (i, j) for i in L0 for j in L if i != j ]AB = [ (i, j) for i in B for j in B0 if i != j]AC = [ (i, j) for i in L for j in B0 if i != j]# 计算网络弧距离Cost = {}for i in V:x1, y1 = XY[i][0], XY[i][0]for j in V:x2, y2 = XY[j][0], XY[j][0]Cost[i,j] = math.sqrt( (x1-x2)**2 + (y1-y2)**2 )return depot,L,L0,B,B0,AL,AB,AC,Q,Cost,XY
# 提取结果，形成车辆路径
def extract_routes(depot,L,B,B0,X,Y,Z):L = copy.deepcopy(L)B = copy.deepcopy(B)B0 = copy.deepcopy(B0)route_list = []V = []while len(L):# 提取 派送阶段路径route = [depot]cur_node = depotfor j in L:if X[depot, j].x > 0:cur_node = jroute.append(j)L.remove(j)breakstop = Truewhile len(L) > 0:for j in L:if X[cur_node, j].x > 0:cur_node = jroute.append(j)if j != depot:L.remove(j)stop = Falsebreakif stop:breakelse:stop = True# 提取 取货阶段路径for j in B:if Z[cur_node,j].x > 0:cur_node = jroute.append(j)B0.remove(j)breakif cur_node in B:while cur_node != depot:for j in B0:if Y[cur_node, j].x > 0:cur_node = jroute.append(j)if j != depot:B0.remove(j)breakelse:route.append(depot)route_list.append(route)V.extend(route[1:-1])print(len(V))return route_list
# 绘制车辆路径
def draw_routes(route_list,XY,L,B):for route in route_list:path_x = []path_y = []for n in route:path_x.append(XY[n][0])path_y.append(XY[n][1])plt.plot(path_x, path_y, linewidth=0.5, ms=5,color='black')linehual_point_x = [XY[n][0] for n in L]linehual_point_y = [XY[n][1] for n in L]backhual_point_x = [XY[n][0] for n in B]backhual_point_y = [XY[n][1] for n in B]plt.scatter(linehual_point_x, linehual_point_y, marker='s', c='b', s=5)plt.scatter(backhual_point_x, backhual_point_y, marker='o', c='r', s=5)plt.show()
# 保存结果
def save_file(route_list,total_cost,Cost):wb = xlsxwriter.Workbook('路径方案.xlsx')ws = wb.add_worksheet()ws.write(0,0,'总费用')ws.write(0,1,total_cost)ws.write(1,0,'车辆')ws.write(1,1,'路径')ws.write(1,2,'距离')for row,route in enumerate(route_list):route_str = [str(i) for i in route]dist = sum(Cost[route[i], route[i + 1]] for i in range(len(route) - 1))ws.write(row + 2, 0, f'{row + 1}')ws.write(row+2,1,'-'.join(route_str))ws.write(row + 2, 2, dist)row += 1wb.close()
# 建模和求解
def solve_model(depot,L,L0,B,B0,AL,AB,AC,Q,Cost,K,CAP,XY):""":param depot:车场id:param L:linehaul节点集合:param B:backhaul节点集合:param AL:linehaul节点衔接弧集合:param AB:backhaul节点衔接弧集合:param AC:linehaul节点和backhaul节点衔接弧集合:param Q:节点需求集合:param Cost:网络弧费用:return:"""model = Model()# 添加变量X = model.addVars(AL,vtype=GRB.BINARY,name='X[i,j]')Y = model.addVars(AB,vtype=GRB.BINARY,name='Y[i,j]')Z = model.addVars(AC,vtype=GRB.BINARY,name='Z[i,j]')U1 = model.addVars(L0,vtype=GRB.CONTINUOUS,name='U[i]')U2 = model.addVars(B0, vtype=GRB.CONTINUOUS, name='U[i]')# 目标函数obj = (quicksum(X[i,j]*Cost[i,j] for i,j in AL) + quicksum(Y[i,j]*Cost[i,j] for i,j in AB) +quicksum(Z[i,j]*Cost[i,j] for i,j in AC))model.setObjective(obj,GRB.MINIMIZE)# linebaul 相关约束model.addConstr( quicksum(X[depot,j] for j in L) == K ) # 车辆数约束model.addConstrs( (quicksum(X[i,j] for i in L0 if i != j) == 1 for j in L) ) # 派送需求约束model.addConstrs( (U1[i] - U1[j] + CAP*X[i,j] <= CAP - Q[j] for i,j in AL) ) # 破圈约束# backbaul 相关约束model.addConstr( quicksum([Y[j,depot] for j in B]) + quicksum(Z[i,depot] for i in L) == K ) # 车辆数约束model.addConstrs( (quicksum(Y[i,j] for j in B0 if i != j) == 1 for i in B) ) # 取货需求约束model.addConstrs( (U2[i] - U2[j] + CAP*Y[i,j] <= CAP - Q[j] for i,j in AB) ) # 破圈约束# connection 相关约束model.addConstr( quicksum(Z[i,j] for i,j in AC) == K ) # 车辆数约束model.addConstrs( quicksum(X[i,j] for j in L if i != j ) + quicksum(Z[i,j] for j in B0) <= 1 for i in L ) # 接续约束model.addConstrs( quicksum(Y[i,j] for i in B if i != j) + quicksum(Z[i,j] for i in L) <= 1 for j in B) # 接续约束# 模型求解model.Params.TimeLimit = 300  # 规模较大时可设置求解时间限制model.optimize()# 判断求解状态if model.status == GRB.Status.OPTIMAL or model.status == GRB.Status.TIME_LIMIT:route_list = extract_routes(depot,L,B,B0,X,Y,Z)draw_routes(route_list, XY, L, B)save_file(route_list, model.objVal,Cost)
if __name__=='__main__':filename=r'demand-X-n120-50-k3.csv'depot,L,L0,B,B0,AL,AB,AC,Q,Cost, XY = read_data(filename)solve_model(depot=depot,L=L,L0=L0,B=B,B0=B0,AL=AL,AB=AB,AC=AC,Q=Q,Cost=Cost,K=3,CAP=21,XY=XY)

4. 求解结果

在这里插入图片描述

参考

https://blog.csdn.net/qq_44149045/article/details/128941919
Koc,Cagri,Laporte,et al.Vehicle routing with backhauls: Review and research perspectives[J].Computers, 2018.
Parragh, S.N., Doerner, K.F. & Hartl, R.F. A survey on pickup and delivery problems . Journal für Betriebswirtschaft 58, 21–51 (2008). https://doi.org/10.1007/s11301-008-0033-7