The Vehicle Routing Problem with Simultaneous Pickup and Delivery (VRPSPD) is probably the most studied and most general variant of the 1-M-1 problems, and it is also known as the Multiple-Vehicle Hamiltonian 1-M-1-PDP with Combined. Demands

1. 数学模型

见参考文献[1]

2. python调用SCIP求解 VRPSPD

• 完整代码

import os
import re
import math
import pandas as pd
import numpy as np
import itertools
from collections import defaultdict

import pyscipopt as opt
import matplotlib.pyplot as plt
import networkx as nx


class VRPSPD(object):
    def __init__(self):
        self.vehicleNum  = 5         # 车辆数
        self.nodeNum     = None      # 顶点数量
        self.capacity    = 40        # 车容量
        self.cor_X       = []        # x轴坐标
        self.cor_Y       = []        # y轴坐标
        self.disMatrix   = [[]]      # 距离成本
        
    def prepare_data(self):
        data = pd.DataFrame([[  0, 250, 250,   0,   0],
                             [  1, 387, 297,  10,   5],
                             [  2,   5, 297,   0,  20],
                             [  3, 355, 177,  20,  10],
                             [  4,  78, 346,   0,   0],
                             [  5, 286, 159,   5,  10],
                             [  6, 322, 465,  20,   0],
                             [  7, 393, 408,  10,  20],
                             [  8,  89, 216,   0,   5],
                             [  9,  76, 345,   5,   5],
                             [ 10, 410, 285,  20,   0]], columns=['NO', 'XCOORD', 'YCOORD', 'PICKUP', 'DELIVERY'])

        self.nodeNum = len(data)
        self.cor_X = data['XCOORD'].values
        self.cor_Y = data['YCOORD'].values
        self.PICKUP = data['PICKUP'].values
        self.DELIVERY = data['DELIVERY'].values
        
        # 计算欧式距离矩阵
        self.disMatrix = np.full((self.nodeNum, self.nodeNum), 0)
        for i in range(self.nodeNum):
            for j in range(self.nodeNum):
                self.disMatrix[i][j] = math.sqrt((self.cor_X[i] - self.cor_X[j]) ** 2 + (self.cor_Y[i] - self.cor_Y[j]) ** 2)
    
    def print_vehicle_route(self, pdp_route):
        # 使用defaultdict来收集属于同一车辆的所有边
        routes_by_vehicle = defaultdict(list)
        for start, end, vehicle in pdp_route:
            routes_by_vehicle[vehicle].append((start, end))
    
        # 构建车辆及其路径的表示
        vehicles_and_paths = {}

        for vehicle_id, edges in routes_by_vehicle.items():
            routes = []
            for i in range(len(edges)):
                if i == 0:
                    routes.append(edges[0][:2])
                    edges.pop(0)
                else:
                    pre = routes[-1]
                    connected_node = pre[1] 
                    for j in range(len(edges)):
                        aft = edges[j]
                        if connected_node in aft:
                            if aft[0] == connected_node:
                                chosen_node = (aft[0], aft[1])
                            elif aft[1] == connected_node:
                                chosen_node = (aft[1], aft[0])
                            routes.append(chosen_node)
                            edges.pop(j)
                            break
                i = i+1
                        
            path = " -> ".join(str(edge[0]) for edge in routes)
            path += '->' + ' 0'
            vehicles_and_paths[vehicle_id] = path
            
            print(f'车辆 {vehicle_id} 路径 {path}')

    def plt_route_pic(self, vrp_route):
        plt.figure(figsize=(10, 10), dpi=200)
        # 绘制坐标的散点图
        for i in range(0, self.nodeNum):
            if i == 0:
                plt.scatter(self.cor_X[i], self.cor_Y[i], marker='^', c='r', s=50)
                plt.annotate('depot', (self.cor_X[i], self.cor_Y[i]))
            else:
                plt.scatter(self.cor_X[i], self.cor_Y[i], marker='o', c='black')
                mark = str((i, self.PICKUP[i], self.DELIVERY[i]))
                plt.annotate(mark, (self.cor_X[i], self.cor_Y[i] + 5) )
        # 绘制箭头
        for idx in vrp_route:
            i,j,k = idx
            plt.arrow(self.cor_X[i], self.cor_Y[i], self.cor_X[j] - self.cor_X[i], self.cor_Y[j] - self.cor_Y[i],
                      length_includes_head=True, head_width=3.5, head_length=5.8, fc='blue', ec='blue', linewidth=1, shape='full')

        plt.xlabel('x')
        plt.ylabel('y')
        plt.title('VRPSPD')
        plt.savefig('./VRPSPD.png')

    def solve_pdp(self):
        model = opt.Model()
        # ==========定义变量==========
        x = {}  # x_i_j_k: 0-1变量，表示车辆k经过弧（i，j） i != j
        y = {}  # 弧（i, j）上的pickup数量
        z = {}  # 弧（i, j）上的delivery数量
        
        for i in range(self.nodeNum):
            for j in range(self.nodeNum):
                if i != j:
                    y[i, j] = model.addVar(vtype='C', name='p_' + str(i) + '_' + str(j), lb=0)
                    z[i, j] = model.addVar(vtype='C', name='d_' + str(i) + '_' + str(j), lb=0)
                    for k in range(self.vehicleNum):
                        x[i, j, k] = model.addVar(vtype='B', name='x_' + str(i) + '_' + str(j) + '_' + str(k))

        # ==========定义约束==========
        # 约束6.9：对于客户点: 保证每个客户点都被服务一次
        for i in range(1, self.nodeNum):
            model.addCons(opt.quicksum(x[i, j, k] for j in range(self.nodeNum) for k in range(self.vehicleNum) if i!=j) == 1,
                          name='customer_visited_' + str(i))

        # 约束6.10: 对于客户点: 进出是同一辆车，流平衡，服务之后需离开         
        for i in range(self.nodeNum):
            for k in range(self.vehicleNum):
                model.addCons(opt.quicksum(x[i, j, k] for j in range(self.nodeNum) if i!=j) ==
                              opt.quicksum(x[j, i, k] for j in range(self.nodeNum) if i!=j),
                              name='flow_' + str(k) + '_' + str(i))
        
        # 约束6.11: 车辆最多使用一次
        for k in range(self.vehicleNum):
            model.addCons(opt.quicksum(x[0, j, k] for j in range(self.nodeNum) if j!=0) <= 1)
                    
        
        # 约束6.12: 限制了每个弧上的流量最多为Q 
        for i in range(0, self.nodeNum):
            for j in range(0, self.nodeNum):
                if i != j:
                    model.addCons(y[i, j] + z[i, j] <= opt.quicksum(x[i, j, k] for k in range(self.vehicleNum)) * self.capacity)
        
        # 约束6.13: 提货能量守恒
        for i in range(1, self.nodeNum):
            model.addCons(opt.quicksum(y[i, j] for j in range(self.nodeNum) if i!=j) -
                          opt.quicksum(y[j, i] for j in range(self.nodeNum) if i!=j) == self.PICKUP[i])
        
        # 约束6.14: 交货能量守恒  
        for j in range(1, self.nodeNum):
            model.addCons(opt.quicksum(z[i, j] for j in range(self.nodeNum) if i!=j) -
                          opt.quicksum(z[j, i] for j in range(self.nodeNum) if i!=j) == self.DELIVERY[i])

        # ==========定义目标==========
        model.setObjective(opt.quicksum(x[i, j, k] * self.disMatrix[i, j] for (i, j, k) in x))
        model.setMinimize()
        model.setParam('limits/gap', 0.0001)
        model.optimize()

        # ==========输出结果==========
        if model.getStatus() != 'infeasible':
            print('model_status = ', model.getStatus())
            print('model_gap =', model.getGap())
            print('model_obj =', model.getObjVal())
    
            # model.writeProblem('pdp.lp')
            pdp_route = []
            for k in range(self.vehicleNum):
                for i in range(self.nodeNum):
                    for j in range(self.nodeNum):
                        if(i != j and model.getVal(x[i,j,k]) > 0):
                            pdp_route.append((i,j,k))      
            # print('pdp_route: ', pdp_route) 
        
            # 打印车辆路径
            self.print_vehicle_route(pdp_route)
            # 绘制结果路径
            self.plt_route_pic(pdp_route)
        else:
            print('model is infeasible')
        
if __name__ == '__main__':
    VRPSPD = VRPSPD()
    # 预处理数据
    VRPSPD.prepare_data()
    # 求解模型并绘制结果
    VRPSPD.solve_pdp()
        

• 运行日志 & 结果

参考文献

[1] Toth P , Vigo D .Vehicle Routing: Problems, Methods, and Applications, Second Edition[J].Society for Industrial and Applied Mathematics, 2014.DOI:10.1137/1.9781611973594.ch9.