OR-Tools VRP 问题从入门到升天(一) TSP问题
Ortools的VRP求解器简介
谷歌的Ortools整合了许多对运筹优化问题的求解器，其中最好用的部分就是VRP求解器。
在ortools中，VRP求解器是建立在constraint programming求解器之上的，因此除了一些经典的VRP问题约束，例如最大负载，时间窗以外，还可以通过约束规划求解器，灵活地添加一系列高度定制化的约束，例如要求某辆车经过某个节点，要求一个节点被多台车访问等。
在ortools中求解VRP问题时，最常用的有两类变量，我们可以从抽象层面了解一下其含义：
-
路径变量(Path variable)：
next(i)- 代表节点 i 在路径中的后继节点的下标。用IndexToNode()可以从下表获取节点变量。vehicle(i)- 表示节点 i 所属于的车辆路径编号，即第 i 个节点由哪辆车来服务。active(i)- 布尔变量，如果节点 i 被访问，那么值为True，否则为False。下面的几组关系是等价的：
active(i) == 0 <==> next(i) == i <==> vehicle(i) == -1，next(i) == j ==> vehicle(j) == vehicle(i)
-
维度变量(Dimension variable) - 用于标记经过路径时累加的值，例如车辆的负载、消耗时间、行驶距离等。常用的有两类：
cumul(i, d)- 代表车辆到达节点 i 时维度 d 的值transit(i, d)- 代表车辆经过节点 i 时，需要在维度 d 上累加的值这两者的关系为：
next(i)==j ==> cumul(j,d)=cumul(i,d)+transit(i,d)
因此，从抽象层面上来讲，每台车从depot出发，不断进行next(i)操作，最终回到depot，即可求得该车辆的路径。每次next(i)操作的同时，我们也通过transit(i,d)不断累加我们关注的值，例如车辆的总行驶距离、消耗时间等。这就是ortools对VRP问题的抽象概括。
数据
这里我们使用的测试数据来自TSPlib中的gr120，其中包含120个节点。
数据的读入使用python中的tsplib95包来完成，节省自己编写数据读入部分的时间。
所有节点的可视化如下：
整体代码
我们先看一下整体的代码，了解调用ortools求解VRP的基本框架，然后再在后面去对各个重要的类进行详细分析：
## 包的导入
from ortools.constraint_solver import routing_enums_pb2
from ortools.constraint_solver import pywrapcp
import tsplib95
import matplotlib.pyplot as plt
%matplotlib inline
## 设定数据目录并读入数据
fname: str = "gr120.tsp"
input_fpath: str = r"data/" + fname
data = tsplib95.load(input_fpath).as_name_dict()
## 建立模型
# 创建routing index manager
n_vehicle = 1 # 车辆数，在TSP问题中为1
depot_index = 0 # 车库，设置为0号节点
manager = pywrapcp.RoutingIndexManager(data['dimension'], n_vehicle, depot_index)
# 创建一个Routing model
routing = pywrapcp.RoutingModel(manager)
# 创建一个distance call back，这里使用欧几里得距离
def distance_callback(from_index: int, to_index: int):
from_node = manager.IndexToNode(from_index)
to_node = manager.IndexToNode(to_index)
# +1是因为读入数据时节点下标从1开始
from_coor = data['display_data'][from_node+1]
to_coor = data['display_data'][to_node+1]
return (to_coor[0] - from_coor[0])**2 + (to_coor[1] - from_coor[1])**2
# 因为我们需要在车辆经过节点i是为其在距离维度d上累加一定的值，所以这里我们使用transitCallBack; 之后我们会看到在设置时间窗口时，会用到cumul类的维度变量
transit_callback_index = routing.RegisterTransitCallback(distance_callback)
# 给每条边设置cost
routing.SetArcCostEvaluatorOfAllVehicles(transit_callback_index)
# 设置求解的启发式算法
search_parameters = pywrapcp.DefaultRoutingSearchParameters()
search_parameters.first_solution_strategy = (routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC)
# 求解
solution = routing.SolveWithParameters(search_parameters)
## 输出结果：
if solution:
print("路径总长度: {} miles".format(solution.ObjectiveValue()))
index = routing.Start(0) # routing.Start(0)代表获取0号车辆的起点
path: str = "旅行商的路径：\n"
visiting_sequence: list = []
while not routing.IsEnd(index):
visiting_sequence.append(manager.IndexToNode(index)+1)
path += " {} ->".format(manager.IndexToNode(index)+1)
previous_index = index
index = solution.Value(routing.NextVar(index))
visiting_sequence.append(manager.IndexToNode(index)+1)
path += " {}\n".format(manager.IndexToNode(index)+1)
print(path)
else:
print("没有得到可行解")
得到结果：
路径总长度: 28422 miles
旅行商的路径：
1 -> 76 -> 59 -> 15 -> 30 -> 29 -> 120 -> 32 -> 92 -> 28 -> 45 -> 78 -> 86 -> 94 -> 81 -> 22 -> 66 -> 31 -> 117 -> 85 -> 18 -> 19 -> 25 -> 108 -> 43 -> 79 -> 52 -> 33 -> 100 -> 58 -> 91 -> 68 -> 65 -> 69 -> 113 -> 107 -> 20 -> 46 -> 50 -> 98 -> 17 -> 118 -> 13 -> 49 -> 42 -> 41 -> 56 -> 44 -> 75 -> 14 -> 87 -> 74 -> 105 -> 40 -> 72 -> 7 -> 38 -> 4 -> 119 -> 103 -> 9 -> 23 -> 51 -> 11 -> 115 -> 21 -> 93 -> 2 -> 82 -> 3 -> 80 -> 27 -> 53 -> 64 -> 109 -> 88 -> 97 -> 12 -> 95 -> 77 -> 5 -> 63 -> 73 -> 57 -> 39 -> 83 -> 67 -> 37 -> 62 -> 99 -> 10 -> 35 -> 104 -> 106 -> 114 -> 101 -> 102 -> 48 -> 110 -> 112 -> 36 -> 84 -> 6 -> 89 -> 55 -> 47 -> 71 -> 26 -> 34 -> 116 -> 70 -> 8 -> 54 -> 90 -> 96 -> 111 -> 24 -> 60 -> 16 -> 61 -> 1
我们还可以对路径进行可视化：
plt.figure(figsize=(8, 6))
# 画出depot
depot_coor = data['display_data'][depot_index + 1]
plt.plot(depot_coor[0], depot_coor[1], 'r*', markersize=11)
# 路径可视化
for i, j in zip(visiting_sequence, visiting_sequence[1:]):
start_coor = data['display_data'][i]
end_coor = data['display_data'][j]
plt.arrow(start_coor[0], start_coor[1], end_coor[0] - start_coor[0], end_coor[1] - start_coor[1])
plt.xlabel("X coordinate", fontsize = 14)
plt.ylabel("Y coordinate", fontsize = 14)
plt.title("TSP path for {}".format(fname), fontsize = 16)
结果如下图（红色星星代表depot）：
相关的类
RoutingIndexManager
RoutingIndexManager类是用来管理 NodeIndex <-> variable index之间的映射的。它之所以存在是因为在VRP问题中，节点的编号和变量编号并非是一一对应的，例如我们需要让一个节点被多个车辆访问时，这个节点会被分拆成多个变量，从而产生了单个节点编号和多个变量编号的对应。
它有两类初始化方法：
当所有车辆具有相同起点和终点时，用
RoutingIndexManager(num_nodes: int, num_vehicles: int, depot: NodeIndex)- 其中num_nodes表示问题中的节点数量，num_vehicles表示车辆数，depot表示所有车辆的起点和重点当每个车辆的起终点都需要单独设定时，用
RoutingIndexManager(num_nodes: int, num_vehicles: int, starts: List[NodeIndex], ends: List[NodeIndex])或者RoutingIndexManager(num_nodes: int, num_vehicles: int, starts_ends: List[Tuple[NodeIndex, NodeIndex]]- 其中starts数组给出每辆车的起点下标，ends数组给出每辆车的终点下标，starts_ends以成对的方式给出每辆车的起点和终点。
RoutingModel
在ortools中，这个类的作用相对混杂：
在其中可以对模型的约束进行修改，例如使用
SetPickUpAndDeliveryPolicyOfAllVehicles可以设置pick up & delivery过程中，车辆是FIFO还是FILO。可以对求解算法进行设置，如利用
SolveWithParameters(search_parameters: RoutingSearchParameters)，可以设置求解时的参数。甚至也可以通过一些方法设置求解时用的邻域搜索方法等。在求解后，其中也保存了解的信息，例如我们可以通过
End(vehicle_id:int)获得某辆车的终点`。
它有两种初始化方法：
采用默认的模型参数时：
RoutingModel(manager: RoutingIndexManager)，直接传入一个IndexManager类即可完成初始化；采用特定的模型参数时，可以显式的传入参数：
RoutingModel(manager: RoutingIndexManager, paramters: RoutingModelParameters)
由于它包含的内容实在太多，具体的一些方法遇到的时候我们再行详解。
Assignment
Assignment类中包含了变量->域的映射，用来保存呈现给用户的结果。在上面的例子中，我们求解之后获得的solution就是一个Assignment类实例。
这个类中我们比较常用的方法有：
Value(var)- 返回变量var的值；ObjectiveValue()- 返回当前解的目标函数值
