实战指南)
简介MDVRP多配送中心车辆路径规划是VRP的关键扩展涉及多中心、多车辆协同调度在物流配送、供应链优化和智能交通中应用广泛。代码包基于MATLAB实现遗传算法求解MDVRP适合运筹优化方向的科研人员、算法工程师及相关专业学生快速上手。压缩包共22个文件由9个M脚本和13个MAT数据文件组成整体约20KBM脚本覆盖MDVRP主程序、选择/交叉/变异算子、距离计算与初始路径生成等核心环节MAT文件提供多组不同规模的测试实例。已有379人学习使用代码结构模块化易读易改运行主程序即可观察路径寻优过程也能灵活调整适应度函数、遗传参数或约束条件为扩展到带时间窗、容量受限等复杂MDVRP变体提供良好基础。1. 多中心车辆路径问题为什么单中心 VRP 在真实配送中不够用当手里的压缩包叫MDVRP.zip时里面大概率装着多车场车辆路径问题MDVRP的求解代码或数据集。MDVRP 是 VRP 家族里更接近真实运营的一类车队从多个配送中心出发每个中心有各自的车辆每个客户只被一辆车服务一次车辆完成路线后允许回到任意一个中心或必须回到原中心目标是让总行驶距离、总成本或用车数最小。很多人拿到多中心数据的第一反应是把客户按位置硬拆成多个单中心 VRP 单独求解但这样做往往会让车从 A 中心跑很远去送一个本该由 B 中心送的订单。反直觉的结论是多中心问题的难点不在“多几辆车”而在“客户与中心的归属决策”和“路径构造”必须同时被优化拆开做通常不会收敛到全局优解。这篇文章适合计划与调度系统开发者、算法工程师以及想基于 OR-Tools 自建多车场路径规划能力的 IT 从业者。2. 把 MDVRP 拆成可计算模型从距离矩阵到车辆中心归属2.1 多中心 VRP 和单中心 VRP 在建模上的三点差异单中心 VRP 的模型里所有路线共享同一个起点和终点求解器只需要回答“哪些客户排在同一条路线以及路线内的顺序”。MDVRP 至少要多回答一个问题“这个客户由哪个中心的哪辆车服务”。这三点差异会影响建模方式。第一仓库节点不再只有一个。假设有 m 个中心每个中心各有若干车辆车辆列表的起点坐标和终点坐标可以从不同中心选择。第二距离矩阵需要扩展。因为车辆可以从任意中心出发客户之间的距离以及客户到每个中心的距离都要参与计算。如果把原始客户节点数记为 N中心节点数记为 M求解时需要构造一个包含 M N 个节点或者更多看实现方式的代价矩阵。第三约束变量增加。每个中心有自己的车辆数、容量、最大行驶时长不同中心甚至可以有不同的车型。这意味着约束不仅要写在车辆维度上还要能区分中心。2.2 用 OR-Tools 建立多中心 Routing 模型的最小骨架常见做法是用 Google OR-Tools 的pywrapcp模块。OR-Tools 本身只有一个虚拟仓库depot但支持为每辆车单独指定起点和终点索引这给多中心留下了口子。做法是为每个中心在节点列表里预留一个“中心节点索引”让每辆车的起点和终点都指向它所属中心对应的索引。下面的代码演示了数据模型的最小结构只包含坐标、车辆起点终点映射和距离回调。from ortools.constraint_solver import routing_enums_pb2, pywrapcp # 中心坐标和客户坐标 center_coords [(40, 40), (80, 80)] # 两个中心 customer_coords [(10, 30), (70, 10), (30, 70), (90, 50), (50, 60)] all_nodes center_coords customer_coords # 索引0,1是中心, 2..N1是客户 has_center 2 # 每辆车所属的起点中心索引和终点中心索引 starts [0, 1] # 两辆车分别从中心0、中心1出发 ends [0, 1] # 返回原中心 # 距离矩阵要包含中心到中心通常置0或大数以及中心到客户、客户到客户 def distance_matrix(): import math n len(all_nodes) mat [[0] * n for _ in range(n)] for i in range(n): for j in range(n): if i j: mat[i][j] 0 else: mat[i][j] int(math.hypot(all_nodes[i][0] - all_nodes[j][0], all_nodes[i][1] - all_nodes[j][1]) * 10) return mat matrix distance_matrix() manager pywrapcp.RoutingIndexManager(len(matrix), len(starts), starts, ends) routing pywrapcp.RoutingModel(manager) def dist_callback(from_index, to_index): from_node manager.IndexToNode(from_index) to_node manager.IndexToNode(to_index) return matrix[from_node][to_node] transit_callback routing.RegisterTransitCallback(dist_callback) routing.SetArcCostEvaluatorOfAllVehicles(transit_callback)这段代码的关键是RoutingIndexManager的构造参数节点总数、车辆数、起点数组、终点数组。两个中心分别对应索引 0 和 1因此车辆 0 从索引 0 出发并回到 0车辆 1 从索引 1 出发并回到 1。若允许车辆回到任意中心可以让ends指向一个“虚拟汇点”或者让不同车辆的终点不同更简单的办法还是保持原中心回场因为大多业务要求车回原车场。距离矩阵里中心到自身距离取 0客户与客户之间正常计算。2.3 车辆容量、时间窗和最大行驶时长的落地方式在 OR-Tools 里容量和时间窗都通过AddDimensionWithVehicleCapacity或AddDimension实现。以容量为例先定义一个节点需求列表中心节点的需求为 0然后添加维度约束。车辆容量可以按中心设置不同值。demands [0, 0, 5, 8, 3, 6, 2] # 索引0/1为中心, 后面为各客户需求 vehicle_capacities [15, 15] # 两辆车容量都是15 def demand_callback(from_index): node manager.IndexToNode(from_index) return demands[node] demand_callback_index routing.RegisterUnaryTransitCallback(demand_callback) routing.AddDimensionWithVehicleCapacity( demand_callback_index, slack_max0, vehicle_capacityvehicle_capacities, fix_start_cumul_to_zeroTrue, namecapacity) search_parameters pywrapcp.DefaultRoutingSearchParameters() search_parameters.first_solution_strategy ( routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC) solution routing.SolveWithParameters(search_parameters) if solution: print(f总距离: {solution.ObjectiveValue() / 10}) # 因为前面乘了10容量维度的slack_max表示每个节点允许的松弛量0 表示不允许违背容量fix_start_cumul_to_zero让每个车辆在出发中心的累计负载从 0 开始。如果业务里有时间窗需要再加一个时间维度同样可以用AddDimension但必须为每辆车定义速度或服务时间。时间维度的最大跨度不要物理距离混为一谈距离和时长是两个独立的 transit callback。新手最容易错的地方是把demand_callback写成from_index和to_index两个参数导致回调返回值对不上。参数含义可以参照下表参数作用多中心场景建议slack_max每个节点可以等待或延迟的量容量设为 0时间窗设为允许等待上限vehicle_capacity每辆车的容量数组长度必须等于车辆数每个中心可不同fix_start_cumul_to_zero出发时累计值是否归零容量和时间都设为 Truefirst_solution_strategy初始解的构造方法小规模用PATH_CHEAPEST_ARC大规模用SAVINGS3. 用 Python 和 OR-Tools 写出第一个可运行的 MDVRP 脚本3.1 输入数据格式与距离矩阵预处理在实际项目里客户和中心通常以 CSV 传入包含id, x, y, demand四列。第一步是把数据读成列表然后把中心节点拼到客户节点前面。这一步必须严格保证顺序因为后面所有索引都依赖这里。import csv, math centers [] # 每个中心 dict(id,x,y) customers [] # 每个客户 dict(id,x,y,demand) with open(mdvrp_input.csv) as f: reader csv.DictReader(f) for row in reader: if row.get(type) center: centers.append({id: row[id], x: float(row[x]), y: float(row[y])}) else: customers.append({id: row[id], x: float(row[x]), y: float(row[y]), demand: float(row[demand])}) node_coords centers customers node_demands [0] * len(centers) [c[demand] for c in customers]距离矩阵我一般用欧氏距离按米为单位取整因为 OR-Tools 的 arc cost 需要整数。如果数据量在几百个节点以内直接算完整矩阵没问题如果上千建议在回调里做缓存避免重复计算。import functools functools.lru_cache(maxsizeNone) def dist(i, j): if i j: return 0 dx node_coords[i][x] - node_coords[j][x] dy node_coords[i][y] - node_coords[j][y] return int(math.hypot(dx, dy))3.2 完整求解脚本车辆分配、约束与解输出下面是可直接运行的脚本骨架覆盖多中心、多车辆、容量约束。示例中两个中心各配一辆车目的是验证索引设置是否正确。from ortools.constraint_solver import routing_enums_pb2, pywrapcp # 每个中心的车数 vehicles_per_center [1, 1] vehicle_starts, vehicle_ends [], [] for c_idx, count in enumerate(vehicles_per_center): vehicle_starts.extend([c_idx] * count) vehicle_ends.extend([c_idx] * count) num_vehicles len(vehicle_starts) manager pywrapcp.RoutingIndexManager(len(node_coords), num_vehicles, vehicle_starts, vehicle_ends) routing pywrapcp.RoutingModel(manager) def dist_callback(a, b): from_node manager.IndexToNode(a) to_node manager.IndexToNode(b) return dist(from_node, to_node) transit routing.RegisterTransitCallback(dist_callback) routing.SetArcCostEvaluatorOfAllVehicles(transit) # 容量约束 demand_callback lambda idx: node_demands[manager.IndexToNode(idx)] demand_index routing.RegisterUnaryTransitCallback(demand_callback) routing.AddDimensionWithVehicleCapacity( demand_index, 0, vehicle_capacities[15] * num_vehicles, fix_start_cumul_to_zeroTrue, namecapacity) # 求解 search_parameters pywrapcp.DefaultRoutingSearchParameters() search_parameters.first_solution_strategy routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC search_parameters.time_limit.seconds 10 solution routing.SolveWithParameters(search_parameters) if solution: print(f总路程(米): {solution.ObjectiveValue()}) for vehicle_id in range(num_vehicles): index routing.Start(vehicle_id) route [] while not routing.IsEnd(index): node manager.IndexToNode(index) route.append(node) index solution.Value(routing.NextVar(index)) node manager.IndexToNode(index) route.append(node) print(f车辆{vehicle_id} 所属中心{vehicle_starts[vehicle_id]} 路径节点: {route})这个脚本里vehicle_starts每个中心只对应 1 辆车所以车辆 0 起始节点是 0车辆 1 起始节点是 1。当某个中心有多辆车时起始节点仍指向同一个中心节点OR-Tools 允许不同车辆从同一物理节点出发前提是starts数组里可以重复。输出里能看到某条路径是不是从中心 0 出发并在中心 0 结束。如果看到非中心节点作为起终点就说明vehicle_starts和vehicle_ends的索引与node_coords顺序不对应。3.3 运行命令与结果解析在项目根目录执行python mdvrp_solver.py默认输入文件是mdvrp_input.csv输出打印每条路线以及总里程。结果中最先要检查的三件事是每个客户节点是否被访问且仅被访问一次每条路线首尾是不是中心节点总路程数值是否比手工拆分为单中心 VRP 时更低。这里有一个常见误用有人会把客户坐标和中心坐标混在一个数组里节点顺序一变node_demands就对不上。所以在读 CSV 后要立刻打印前几个节点的坐标与需求排序确认无误再建矩阵。4. 多中心路径规划的求解策略搜索参数与全局优化4.1 关键搜索参数对照与选择原则OR-Tools 的默认参数对小规模实例可用但在 MDVRP 中容易陷入局部最优因为客户归属决策带来的组合爆炸比单中心更严重。first_solution_strategy决定初始解构造方式下面是常用策略及其适用场景策略名称构造逻辑适合场景PATH_CHEAPEST_ARC每次扩展当前路径代价最小的弧节点少、单机求解速度快GLOBAL_CHEAPEST_ARC全局挑选代价最小的弧插入多中心、需要均衡路线时SAVINGSClark-Wright 节约法车辆无容量约束或容量宽松PARALLEL_CHEAPEST_INSERTION并行插入客户客户密集、需要并行初始解SWEEP按极角扫描划分区域中心在几何中心附近时MDVRP 中客户到哪个中心最便宜直接影响初始解质量。我常用GLOBAL_CHEAPEST_ARC配合元启发式GUIDED_LOCAL_SEARCH前者保证初始路由不差后者让搜索跳坑。search_parameters.local_search_metaheuristic ( routing_enums_pb2.LocalSearchMetaheuristic.GUIDED_LOCAL_SEARCH) search_parameters.time_limit.seconds 30 search_parameters.solution_limit 100004.2 大规模实例拆分配置、超时和车辆数调优当节点数超过 200单次求解很容易跑满几分钟。常见做法是先按中心对客户做 Voronoi 分割把大问题切成若干个单中心子问题分别求解后再拼接。但千万要注意这种拆分应该以“距离最近中心”为依据但最终归属可能因为容量限制而改变。所以更稳妥的是让 OR-Tools 自己决定归属只把算法时间限制加大。另一个技巧是利用solution_limit控制提前终止条件设定搜索 10000 个解后停止避免死循环。车辆数也是重要调优项。多中心场景下车辆数不是越多越好因为每加一辆车就会多一份固定成本。如果只是求路线可以通过设置一个非常大的vehicle_capacity让求解器自动减少用车。但 MDVRP 通常要求所有携带车辆资源可以被部分闲置此时需要在目标函数中加入车辆启用成本。OR-Tools 里可以用SetFixedCostOfVehicle为每辆车设置固定成本实现“少用车”的优化目标。4.3 常见报错与排查索引越界、维度不匹配、无可行解第一个高频报错是IndexError或Out of range原因多半是manager.IndexToNode返回的节点编号超过了node_coords的长度。这发生在距离回调中起始索引被 OR-Tools 内部的虚拟节点干扰。解决方法是在回调入口强制断言assert from_node len(node_coords)。第二个坑是AddDimensionWithVehicleCapacity的车辆容量数组长度不等于车辆数比如运行时发现车辆数变成 0这通常是因为vehicle_starts列表为空或循环次数写错。第三个坑是无可行解给出的容量小于任何单个客户需求或车辆数少于中心数。排查时先检查node_demands的最大值是否大于vehicle_capacities的最小值再用routing.GetStatus()打印求解状态。OR-Tools 的GetStatus返回 1 是最优解2 是可行解4 是无可行解。当状态为 4 时优先检查每辆车是否必须访问的节点集合以及是否有客户需求超过所有车容量。5. 把 MDVRP 结果落地到业务验证、可视化与动态避障衔接5.1 用 matplotlib 绘制多车场路径图求解后直接看数字不如看到线路清楚。把每条路线按不同颜色画出来中心用方块标记客户用圆点标记。这样能一眼判断有没有线路交叉、有没有车从较远中心跑去送近点。import matplotlib.pyplot as plt colors [blue, green, orange, purple] for v in range(num_vehicles): index routing.Start(v) xs, ys [], [] while not routing.IsEnd(index): n manager.IndexToNode(index) xs.append(node_coords[n][x]) ys.append(node_coords[n][y]) index solution.Value(routing.NextVar(index)) n manager.IndexToNode(index) xs.append(node_coords[n][x]) ys.append(node_coords[n][y]) plt.plot(xs, ys, colorcolors[v % len(colors)], markero, labelfVehicle {v}) for c in centers: plt.scatter(c[x], c[y], markers, s120, cred) plt.legend() plt.savefig(mdvrp_route.png, dpi150)5.2 约束校验脚本自动检查容量和多中心归属可视化之外还要有硬校验。我一般写一个小函数遍历解里的每条路线累加客户需求确认不超过容量同时确认manager.IndexToNode(routing.Start(v))等于vehicle_starts[v]终点也一致。这一步能防止求解器在某些边界条件下给出非法路径。5.3 从静态 MDVRP 到动态避障小车路径规划的衔接MDVRP 得到的是全局静态路线但实际园区、仓储环境里会有临时障碍或动态车辆这时需要在局部做避让。业内常用做法是把 MDVRP 输出的路径按时间窗口拆成一个个中间路点再交给局部路径规划器例如动态避障小车路径规划常用的 DWA 或 ROS2 本地规划器。MDVRP 负责全局多中心、多车辆协同局部规划器负责两个路点之间的避障平滑。两者之间用统一坐标和时间戳对接关键技巧是把车辆速度和服务时间作为条件写入路径文件而不是只在输出里打印顺序这样业务方才能直接把路径喂给执行系统。多重中心路径规划项目里最后再验证一下你拼接的路径文件是否满足两个条件线路上的起点中心索引与车辆 id 一致每个客户的需求总和不超过该车容量。然后剩下的就是重复运行脚本调整搜索策略直到总路线的耗时满足调度窗口。本文还有配套的精品资源点击获取