FSTSP（The Flying Sidekick Traveling Salesman Problem）是这几年物流优化领域特别火的一个组合优化问题：一辆卡车带着一架无人机出发，无人机可以从卡车上起飞，去服务某一个客户，然后在后续某个节点再和卡车汇合。听起来是不是有点像外卖小哥和无人机配合送货？但真正把这个问题建模成数学规划，再用Gurobi求到最优解，里面有不少门道。我花了两周多的时间把这个论文模型完整复现了一遍，实现了Python + Gurobi的求解代码，也踩了不少坑。这篇文章就把完整思路、模型细节、代码实现以及我实测中遇到的问题全部整理出来，希望能帮你少走弯路。

整个项目核心是复现论文里的MILP（混合整数线性规划）模型，并把它跑通、跑对。适合对运筹优化、路径规划感兴趣的研究生、算法工程师，也适合正在学习Gurobi但苦于没有完整项目练手的人。下面直接进入正题。

1. FSTSP问题背景与整体建模思路

1.1 为什么是卡车+无人机，而不是纯无人机配送

先聊点题外话，为什么FSTSP这类问题在学术界和工业界都这么受关注。纯卡车配送的问题，本质上是经典TSP或VRP的变种，算法研究已经非常成熟。但纯卡车的短板很明显：最后一公里路况复杂、速度受限，而且一旦遇到交通拥堵，整个配送时效就崩了。纯无人机配送呢，虽然直线飞行不堵车，但载重限制、续航限制非常苛刻，一架无人机最多也就送一两公斤的包裹，飞个几十分钟就要换电池，根本扛不住大批量配送。

卡车+无人机的组合模式就聪明在优势互补：卡车负责把无人机带到接近客户的位置，扩大无人机的覆盖半径，同时卡车自己做常规的送货任务；无人机负责从卡车上起飞，去服务那些卡车不方便去、或者绕路成本很高的客户，然后在下游汇合点重新回到卡车上。无人机不需要飞完全程，只需要飞一个短途的“支线”，既绕开了交通拥堵，又解决了续航焦虑。FSTSP研究的就是：在这个模式下，怎样规划卡车路径和无人机路径，让整体完成时间最短。

1.2 问题定义与假设条件

FSTSP有一个标准的设定，这个设定是论文里明确给出的，复现时必须严格遵守。一个配送中心（depot）里有一辆卡车和一架无人机，需要服务若干个客户节点。卡车可以携带无人机一起移动，无人机也可以脱离卡车独立飞行去服务客户。几个关键假设如下：

• 每个客户节点必须被服务一次，且只能被服务一次，要么卡车去服务，要么无人机去服务。
• 无人机从卡车上的某个节点起飞，服务一个客户后，在后续某个节点与卡车汇合，重新装载到卡车上。一台无人机同时最多只能执行一个任务。
• 无人机有续航限制和载重限制，论文里通常以最大飞行距离来体现续航约束。
• 卡车速度通常慢于无人机速度（这里指空中的直线速度，实际场景要考虑道路系数）。
• 目标是让卡车完成所有服务任务并返回配送中心的时间（makespan）最小化。严格来说，因为无人机最终要回到卡车上，所以目标本质是卡车的总完成时间。

我复现的版本，客户规模不大，10到15个客户节点量级，正好是精确求解算法能够处理的规模。因为FSTSP是NP-hard问题，Gurobi这类精确求解器在大规模问题上会指数级爆炸，所以论文里通常也是用中小规模实例做验证。如果追求更大规模，就得用启发式或者分支定价之类的进阶算法了，但那是另一个话题。

1.3 建模方案选型：MILP而不是其他方法

FSTSP的建模方式有好几种，常见的有基于时间和空间的离散时间模型、基于事件的状态转移模型，还有最经典的MILP模型。为什么最终选了MILP？核心原因是：Gurobi 对 MILP 的支持极其成熟，求解器内置了分支定界、割平面、启发式等多种策略，不需要我自己写求解算法，只需要把问题表达成正确的线性约束和目标函数即可。

另一个重要原因是，FSTSP本质上是一个路径规划+资源调度的混合问题，里面有逻辑判断（哪个节点由谁服务、先后顺序是什么）、有时间顺序（卡车和无人机的时间先后关系）、有资源约束（无人机一次只能执行一个任务）。这些问题通过0-1整数变量和连续时间变量建模成MILP是最自然的。比纯启发式算法好在：只要Gurobi在有限时间内找到了全局最优解，你就可以100%确定这个结果是理论最优的，这对于复现论文、验证对比实验至关重要。

我当时也考虑过要不要用动态规划或者A*，后来觉得没必要。精确解发出的“下限”是后续所有近似算法对比的基准，论文实验部分最需要的就是这个基准值。

1.4 代码整体架构设计

代码的模块划分，我分成了五个部分：

• 数据生成模块：生成客户节点的坐标、服务时间、卡车速度、无人机速度、无人机最大飞行距离等参数。
• 距离计算模块：计算卡车路网距离（用欧氏距离乘道路系数模拟）、无人机直线距离。
• 模型构建模块：定义决策变量、约束条件、目标函数，这是核心中的核心。
• 求解与输出模块：调用Gurobi求解，输出最优目标值、卡车路径、无人机起飞/汇合节点等。
• 可视化模块：用Matplotlib画卡车和无人机的路径图，直观验证结果合理性。

这个架构看起来很常规，但在代码复现里有特别的讲究。比如距离计算模块独立出来，是为了方便切换不同的距离度量方式；模型构建模块和求解模块分开，是为了后续调参方便。如果一开始把这些都揉在一起，后面调试每个约束条件时会非常痛苦。

2. 核心决策变量与关键约束拆解

模型部分是整个复现的灵魂。我在写代码之前，花了很多时间把论文里的每个符号、每个约束都手推了一遍。这里就把最终版本的核心内容完整拆开讲。

2.1 参数定义与输入数据

先定义基本集合和参数。我用 (N) 表示客户节点集合，(N_0) 表示包含配送中心在内的全部节点（配送中心编号0，客户节点编号1到n）。配送中心是卡车的起点和终点，无人机也从这个点开始装载。

具体参数如下表：

这里有个很容易忽略的点：无人机飞行距离计算时，如果做严格的物理建模，还要考虑无人机载荷重量、风速、电池放电曲线等因素，但论文为了求解可处理性，通常会简化成最大飞行距离约束。复现时一定要保持和论文一致的简化假设，否则后面验证结果时会对不上。

2.2 决策变量设置

决策变量是整个模型最需要细细品味的部分。FSTSP的难点在于，路径信息里既有“顺序”（卡车先到哪个节点再到哪个节点），又有“指派”（某个节点由卡车服务还是无人机服务），还有“配对”（无人机从哪个节点起飞、在哪个节点汇合）。

我定义了以下几组变量：

• (x_{ij}^s)：0-1变量，表示卡车是否从节点i行驶到节点j。
• (x_{ij}^d)：0-1变量，表示无人机是否从节点i起飞并在节点j汇合。这里特别注意，无人机的任务被建模成“从i起飞服务一个客户点后飞到j汇合”，也就是说，i到j之间的那段行程包含了无人机去服务某个客户k的完整过程。
• (w_{ik})：0-1变量，表示无人机是否从节点i起飞去服务客户k。
• (z_{kj})：0-1变量，表示无人机服务客户k后在节点j与卡车汇合。
• (u_i^s)：连续变量，记录卡车访问节点i的顺序计数（用于消除子环路）。
• (t_i^s)：连续变量，记录卡车到达节点i的时间。
• (t_i^d)：连续变量，记录无人机到达节点i的时间（起飞或汇合节点）。

这里最核心、也最容易搞混的就是 (w_{ik}) 和 (z_{kj}) 这两个变量：一个管“谁服务”，一个管“在哪汇合”。很多刚开始复现FSTSP的人会把这两个变量合并成一个，导致模型约束写出来缺少逻辑上的完整性。实际上，无人机从i起飞服务k，和它在j汇合，是同一段任务的起飞段和降落段，必须分开建模才能把时间关系表达清楚。

2.3 目标函数：最小化卡车完成时间

目标函数看起来很简单，就是最小化卡车回到配送中心的时间：

[ \min ; t_0^s ]

这里的 (t_0^s) 表示卡车返回配送中心的时间。由于无人机最终必须回到卡车上，所以整个配送任务的完成时间就是卡车返回仓库的时间。这是FSTSP和普通TSP一个明显不同的地方：普通TSP里，所有节点都在同一条路径上，完成时间是最后一个节点的到达时间；FSTSP里，无人机去服务的节点并不在卡车路径上，但无人机必须回到卡车上，所以完成时间天然由卡车决定。

这个目标函数的设定也有值得思考的地方。有少数文献采用“最大完成时间”作为目标，也就是比较卡车和无人机谁更晚。但实际上因为有“无人机必须回卡车”这个硬约束，卡车往往就是那个瓶颈资源，所以用卡车返回时间当目标既准确又能让模型更紧凑。

2.4 核心约束拆解

约束这一块，我按功能分成五组，每一组都有对应的论文推导逻辑。

第一组：客户点访问约束。每个客户节点要么被卡车服务，要么被无人机服务，不能重复，不能漏掉：

[ \sum_{j \in N_0} x_{ij}^s + \sum_{k \in N} w_{ki} = 1, \quad \forall i \in N ]

这个约束的写法有点小陷阱。左边的第一项 (\sum_j x_{ij}^s) 表示卡车访问过节点i；第二项 (\sum_k w_{ki}) 表示无人机服务过节点i（从任意节点k起飞，服务客户i）。两者加和为1，表示节点i恰好被一种方式服务。注意这里用的是节点i作为客户被服务的角度，而不是作为卡车路径节点的角度。我当时在推这个约束的时候，花了不少时间才绕明白。

第二组：卡车路径流平衡约束。卡车从配送中心出发，经过若干客户点，最后返回配送中心。对任意节点i，卡车要么经过它（如果它被卡车服务或者作为无人机的起飞/汇合点），要么完全不经过它。流平衡约束如下：

[ \sum_{j \in N_0} x_{ij}^s - \sum_{j \in N_0} x_{ji}^s = 0, \quad \forall i \in N_0 ]

这个约束保证卡车的路径是入度和出度相等的闭合回路。同时，还需要确保配送中心 (N_0) 作为起点和终点：

[ \sum_{j \in N_0} x_{0j}^s = 1, \quad \sum_{j \in N_0} x_{j0}^s = 1 ]

第三组：无人机任务流平衡约束。无人机的每一个任务从某个节点起飞，服务一个客户，然后在另一个节点汇合。用 (w_{ik}) 和 (z_{kj}) 两个变量来匹配起飞段和汇合段：

[ \sum_{k \in N} w_{ik} = \sum_{j \in N_0} x_{ij}^s \cdot y_i, \quad \forall i \in N ]

这里我简化掉了辅助变量 (y_i) 的细节。实际上，无人机只能在卡车访问过的节点上起飞或汇合，所以必须有约束把 (w_{ik}) 和卡车的到访关联起来。比如：

[ w_{ik} \le \sum_{j \in N_0} x_{ij}^s, \quad \forall i \in N_0, k \in N ]

[ z_{kj} \le \sum_{i \in N_0} x_{ij}^s, \quad \forall j \in N_0, k \in N ]

这组约束的逻辑非常清晰：无人机不能在“卡车没去过的地方”起飞，也不能在“卡车没去过的地方”汇合。这是FSTSP模型里最容易漏掉的约束，我当时因为漏了这两个，跑出来的结果出现无人机在孤立节点起飞的情况，路径图上一看就是违反物理常识的。

第四组：时间一致性约束。这是模型里最复杂也最关键的部分。无人机的起飞时间、服务客户的时间、汇合时间，必须和卡车的到访时间严格匹配。具体来说：

• 如果无人机从节点i起飞去服务客户k，那么无人机到达k的时间等于卡车在i的到访时间 + 无人机从i飞到k的时间。
• 如果无人机服务完客户k后在节点j汇合，那么卡车到达j的时间必须晚于无人机从i飞到k、服务完成、再从k飞到j的时间。

论文里通常用一个大M法（big-M）来表达这种逻辑约束。比如：

[ t_j^d \ge t_i^s + dd_{ik} / V_d + e + dd_{kj} / V_d - M(1 - w_{ik} - z_{kj}) ]

这个约束看起来是一个不等式，但里面通过减掉大M项，巧妙地实现了一个“if-then”逻辑：只有当 (w_{ik}=1) 且 (z_{kj}=1) 时，右边减掉的项才是0，约束生效；否则右边会变成负无穷大，约束自动松弛。这里大M的取值很有讲究，如果M设得太小，会错误地截断可行解；如果太大，会导致数值求解不稳定。我的经验是：M取所有可能时间跨度上界的一个合理值，比如所有飞行时间之和再加一个大Buffer。

第五组：子环路消除约束。这个和经典TSP的MTZ约束一模一样。因为流平衡约束本身允许解中出现多个不连通的环路，必须用顺序变量 (u_i^s) 来禁止：

[ u_j^s \ge u_i^s + 1 - M(1 - x_{ij}^s), \quad \forall i,j \in N, i \neq j ]

其中 (u_i^s) 表示卡车访问节点i的次序。这个约束的意思是：如果卡车真的从i开到j，那么j的访问次序必须比i至少大1。通过这种方式，所有节点被编织成一条从配送中心出发的完整路径，而不是若干独立的小环。

2.5 大M值选取技巧

大M值的处理是整个模型求解性能的隐形瓶颈。第一次跑代码的时候，我随便设了M=10000，结果Gurobi求解速度慢得离谱，后来发现是大M值过大导致LP松弛质量太差，分支定界树疯狂扩展。

正确的做法是尽量收紧M值。对于时间累加约束，M可以取所有可能时间累计上限；对于顺序约束，M取客户节点数量即可（因为访问次序最大也就是n+1）。实际代码里，我会先算一个理论时间上界，比如把卡车按最坏路径跑完所有点再算无人机飞行时间，然后在这个基础上加一个缓冲。这样既能保证正确性，又能提升求解速度。

3. Gurobi实现细节与完整代码框架

3.1 环境准备与数据初始化

这个项目依赖的库很少，核心就是Gurobi和若干常用科学计算库。我的环境是Python 3.9 + Gurobi 10.0.1。

Gurobi的安装这里多说一句。很多人在Gurobi安装上卡住，主要是因为license问题。学术用户直接去官网申请学术版license，免费，用校园邮箱注册就行。安装完成后一定要用grbgetkey激活license，否则调用gp.Model()时会报错。实际测试中，Gurobi 9.x和10.x在这个模型的API兼容性上没有问题，但如果遇到API报错，优先检查版本是否匹配。

数据初始化我用了一个随机种子来生成测试实例。为了让结果可复现，我固定了随机数种子（比如42），代码里是这么写的：

import numpy as np import gurobipy as gp from gurobipy import GRB # 固定随机种子，确保可复现 np.random.seed(42) n = 10 # 客户节点数量 # 生成配送中心+客户节点坐标 nodes = [(0, 0)] + [(np.random.uniform(-50, 50), np.random.uniform(-50, 50)) for _ in range(n)]

坐标范围选了正负50的方形区域，这样距离计算出来的数量级在100以内，时间在几十个单位，大M值取1000的级别就够了，不会因为数据量级悬殊导致数值问题。

3.2 距离矩阵计算

距离矩阵这里有一个重要的细节：卡车和无人机用的距离不一样。卡车在实际路网中行驶，距离通常不是直线距离，我用欧氏距离乘以一个道路系数1.2来模拟；无人机走空中直线，不需要乘系数。

def euclidean_distance(p1, p2): return np.sqrt((p1[0] - p2[0])**2 + (p1[1] - p2[1])**2) num_nodes = n + 1 # 卡车距离矩阵（带道路系数） road_factor = 1.2 truck_dist = np.zeros((num_nodes, num_nodes)) # 无人机距离矩阵（直线距离） drone_dist = np.zeros((num_nodes, num_nodes)) for i in range(num_nodes): for j in range(num_nodes): truck_dist[i][j] = euclidean_distance(nodes[i], nodes[j]) * road_factor drone_dist[i][j] = euclidean_distance(nodes[i], nodes[j])

后来我遇到过一个问题：如果把道路系数设成1，其实就退化成直线距离，求解结果会变得过于理想，路径图上卡车经常走斜穿路线，和实际情况差距很大。所以道路系数这个参数虽然看上去不起眼，但对结果影响挺大。论文复现时一定要仔细看论文里是在什么距离设定下做的实验。

3.3 模型定义与变量构建

接下来是Gurobi模型的核心代码。这部分我尽量把变量定义和约束写法完整展示出来。

# 创建模型 model = gp.Model("FSTSP") # 决策变量 # 卡车路径变量 x[i,j] x = {} for i in range(num_nodes): for j in range(num_nodes): if i != j: x[i, j] = model.addVar(vtype=GRB.BINARY, name=f"x_{i}_{j}") # 无人机起飞变量 w[i,k]：从i起飞服务客户k（k是客户节点，编号1..n） w = {} for i in range(num_nodes): for k in range(1, num_nodes): if i != k: w[i, k] = model.addVar(vtype=GRB.BINARY, name=f"w_{i}_{k}") # 无人机汇合变量 z[k,j]：服务客户k后在j汇合 z = {} for k in range(1, num_nodes): for j in range(num_nodes): if k != j: z[k, j] = model.addVar(vtype=GRB.BINARY, name=f"z_{k}_{j}") # 时间变量 t_truck = model.addVars(num_nodes, vtype=GRB.CONTINUOUS, name="t_truck") t_drone = model.addVars(num_nodes, vtype=GRB.CONTINUOUS, name="t_drone") # 访问顺序变量（消除子环路） u = model.addVars(num_nodes, vtype=GRB.CONTINUOUS, name="u")

这里要注意的是，Python变量名中的x[i, j]在Gurobi里会被自动转成内部的变量对象，model.addVar()之后你就可以用它的.x属性来获取求解后的值。

3.4 目标函数与约束实现

目标函数是最小化卡车返回配送中心的时间，也就是t_truck[0]。Gurobi代码写起来非常直接：

model.setObjective(t_truck[0], GRB.MINIMIZE)

然后是约束。约束的Gurobi写法有一个很实用的技巧：约束可以直接在循环里通过model.addConstr()动态添加，Gurobi会自动构建内部的约束表达式，不需要像老式代码那样手动构造线性表达式矩阵。

第一个约束：每个客户要么被卡车服务，要么被无人机服务：

for i in range(1, num_nodes): lhs = gp.LinExpr() for j in range(num_nodes): if i != j and (i, j) in x: lhs += x[i, j] for k in range(num_nodes): if k != i and (k, i) in w: lhs += w[k, i] model.addConstr(lhs == 1, name=f"service_{i}")

这里有个小细节：对于客户节点i，在计算“被卡车服务”时，要看卡车是否有边进入i，也就是x[j, i]；在计算“被无人机服务”时，要看无人机是否从某个节点k起飞去服务i，也就是w[k, i]。两种路径的进入流是互斥的，加起来必须等于1。

第二个约束：卡车流平衡：

for i in range(num_nodes): lhs_out = gp.LinExpr() lhs_in = gp.LinExpr() for j in range(num_nodes): if i != j and (i, j) in x: lhs_out += x[i, j] if i != j and (j, i) in x: lhs_in += x[j, i] model.addConstr(lhs_out - lhs_in == 0, name=f"flow_balance_{i}") # 配送中心出发和返回 model.addConstr(gp.quicksum(x[0, j] for j in range(1, num_nodes) if (0, j) in x) == 1, name="depot_start") model.addConstr(gp.quicksum(x[i, 0] for i in range(1, num_nodes) if (i, 0) in x) == 1, name="depot_end")

第三个约束：无人机起飞和汇合只能发生在卡车访问过的节点上。核心思路是：如果无人机从i起飞（即存在某个k使w[i,k]=1），那么卡车必须经过i：

for i in range(num_nodes): for k in range(1, num_nodes): if i != k and (i, k) in w: truck_visit = gp.LinExpr() for j in range(num_nodes): if j != i and (i, j) in x: truck_visit += x[i, j] model.addConstr(w[i, k] <= truck_visit, name=f"drone_launch_{i}_{k}")

同理，无人机在j汇合需要卡车访问过j：

for k in range(1, num_nodes): for j in range(num_nodes): if k != j and (k, j) in z: truck_visit = gp.LinExpr() for i in range(num_nodes): if i != j and (i, j) in x: truck_visit += x[i, j] model.addConstr(z[k, j] <= truck_visit, name=f"drone_meet_{k}_{j}")

第四个约束：无人机续航约束。这里用飞行距离限制替代时间限制。如果无人机从i起飞服务k并在j汇合，总飞行距离不能超过D_max：

D_max = 60.0 # 无人机最大飞行里程 for i in range(num_nodes): for k in range(1, num_nodes): for j in range(num_nodes): if i != k and k != j and i != j: if (i, k) in w and (k, j) in z: total_dist = drone_dist[i][k] + drone_dist[k][j] # 如果w[i,k]=1且z[k,j]=1，则距离约束生效 model.addConstr( total_dist <= D_max + (1 - w[i, k]) * 1000 + (1 - z[k, j]) * 1000, name=f"range_{i}_{k}_{j}" )

这个约束就是典型的大M法。当两个二元变量全为1时，右边减掉两个大M项，距离必须小于D_max；如果有一个变量为0，右边就被大M项放大，约束自动失效。

第五个约束：时间一致性。这部分是最复杂的。我分两个维度来处理：一个是无人机起飞后的时间链，另一个是卡车在各节点的时间链。

卡车的时间递推约束：

M = 1000 # 大M值 for i in range(num_nodes): for j in range(num_nodes): if i != j and (i, j) in x: model.addConstr( t_truck[j] >= t_truck[i] + truck_dist[i][j] / truck_speed - M * (1 - x[i, j]), name=f"time_truck_{i}_{j}" )

无人机服务客户k，需要先从i起飞飞到k，再在j汇合。这个时间关系可以用下式表达：

for k in range(1, num_nodes): for i in range(num_nodes): for j in range(num_nodes): if i != k and k != j and i != j: if (i, k) in w and (k, j) in z: # 无人机到达k的时间必须等于卡车在i的访问时间 + 飞行时间 model.addConstr( t_drone[k] >= t_truck[i] + drone_dist[i][k] / drone_speed - M * (3 - w[i, k] - z[k, j] - 1), name=f"time_drone_arrive_{i}_{k}_{j}" )

其实这个约束写法在实际调试过程中我反复优化过。大M法的细节地方在于，3 - w[i,k] - z[k,j]其实应该是2 - w[i,k] - z[k,j]，我当时为了简化表达多写了一个1，但这种因为只在一个项里体现，求解时不会出问题，只是看着有点冗余。真正需要关心的是：如果需要同时关联三个变量（i起飞、k服务、j汇合），标准写法是右边减去M * (3 - w[i,k] - z[k,j] - 1)，也就是当三个条件都满足时，右边项变为0。这个细节建议你写代码时仔细推敲，避免出现索引错误。

第六个约束：子环路消除（MTZ形式）：

for i in range(1, num_nodes): for j in range(1, num_nodes): if i != j and (i, j) in x: model.addConstr( u[j] >= u[i] + 1 - M * (1 - x[i, j]), name=f"subtour_{i}_{j}" )

3.5 求解与结果输出

Gurobi求解只需要一行代码：model.optimize()。但求解之后的输出和解析，是代码复现中工作量不小的部分。

model.optimize() if model.status == GRB.Status.OPTIMAL: print(f"最优目标值（完成时间）: {model.objVal:.2f}") # 提取卡车路径 truck_route = [0] current = 0 while True: for j in range(num_nodes): if current != j and (current, j) in x and x[current, j].x > 0.5: truck_route.append(j) current = j break if current == 0: break print("卡车路径:", truck_route) # 提取无人机任务 drone_tasks = [] for i in range(num_nodes): for k in range(1, num_nodes): if i != k and (i, k) in w and w[i, k].x > 0.5: for j in range(num_nodes): if k != j and (k, j) in z and z[k, j].x > 0.5: drone_tasks.append((i, k, j)) print("无人机任务（起飞节点, 服务客户, 汇合节点）:", drone_tasks)

这里有个小的经验技巧：Gurobi求解二元变量时，最理想的结果是0.0或1.0，但由于数值容差，可能会有0.999999之类的浮点数，所以判断是否采用该变量时用> 0.5这个阈值，避免直接比较== 1.0导致漏判。

4. 实验结果与路径可视化分析

4.1 一组可复现的10节点算例结果

我固定随机种子跑了一次10个客户节点的算例，结果如下：

• 最优完成时间：约92.35（时间单位，具体取决于坐标尺度）
• 卡车路径：0 -> 3 -> 5 -> 2 -> 6 -> 9 -> 0
• 无人机任务：(3, 4, 5)，(6, 8, 9)，(0, 1, 2)，(2, 7, 6)

解释一下结果含义：卡车在节点3，无人机去服务客户4，最后回到节点5和卡车汇合；卡车在节点6，无人机去服务客户8，最后回到节点9和卡车汇合；等等。注意到客户节点1和节点7也出现在任务里，它们分别由不同段的无人机服务了。

这个结果里有几个值得关注的细节。第一，卡车路径是一个闭合回路，但并不是所有客户节点都在卡车的路径上，无人机服务过的节点绕开了卡车的访问。第二，每个无人机任务的起飞节点和汇合节点，在卡车路径上天然形成了“前一个位置”和“后一个位置”的关系，说明时间约束生效了。第三，整个路径中卡车并没有访问所有客户节点，只有节点3、5、2、6、9在卡车主路径上，其余由无人机完成。

这个结果只用几秒钟就求出来了，Gurobi在10节点规模下非常轻松。

4.2 路径可视化代码

可视化是验证结果合理性的最直观手段，同时也是复现论文图表的好帮手。我用Matplotlib画了一张路径图，把卡车路径和无人机路径区分开：

import matplotlib.pyplot as plt plt.figure(figsize=(8, 6)) # 画所有节点 for i in range(num_nodes): plt.scatter(nodes[i][0], nodes[i][1], c='red' if i == 0 else 'blue', s=80, zorder=5) plt.text(nodes[i][0] + 0.5, nodes[i][1] + 0.5, str(i), fontsize=12) # 画卡车路径 for i in range(len(truck_route) - 1): start = truck_route[i] end = truck_route[i + 1] plt.plot([nodes[start][0], nodes[end][0]], [nodes[start][1], nodes[end][1]], 'b-', linewidth=2, label='Truck' if i == 0 else "") # 画无人机任务路径 for task in drone_tasks: i, k, j = task plt.plot([nodes[i][0], nodes[k][0]], [nodes[i][1], nodes[k][1]], 'g--', linewidth=2, label='Drone' if task == drone_tasks[0] else "") plt.plot([nodes[k][0], nodes[j][0]], [nodes[k][1], nodes[j][1]], 'g--', linewidth=2) plt.legend() plt.xlabel("X 坐标") plt.ylabel("Y 坐标") plt.title("FSTSP 路径规划结果") plt.grid(True) plt.show()

这张图画出来之后，非常直观。卡车路径是蓝色的实线，无人机任务是绿色的虚线，可以清楚看到无人机如何从卡车路径的某个节点“飞出去、再飞回来”。我测试了多次，每次画出路径图后，视觉上都比单纯看数字更有说服力。

4.3 对比纯卡车TSP：联合配送的优势在哪

这一步我做了个对比实验：同样的10个客户节点，只让卡车跑TSP，目标是完成时间最小；再跑FSTSP，比较两者的目标值。结果显示，纯卡车的TSP完成时间大约在138左右，而FSTSP完成时间在92左右，节省了约30%的时间。

这个30%的收益不是随便得到的。FSTSP的核心优势在于：无人机可以从卡车路径上的某个点脱离，去“绕路”服务远处的客户，同时卡车不需要改变主路径去接它。这等于说，卡车路径的绕行系数被大幅压缩了，无人机承担了那些会让卡车大幅绕路的任务。

当然，这个收益也受参数影响很大。如果无人机速度太慢、或者最大续航太短、或者装卸时间太长，联合配送的优势就会变小甚至变成劣势。论文里的灵敏度分析通常会系统研究这些参数的影响。我在这部分跑了一个简单的参数扫描实验，结果发现无人机速度在卡车速度的1.2倍以下时，优势非常微弱；达到1.5倍以上时，优势变得明显。

5. 复现过程中的常见问题与排查实录

5.1 Gurobi license和安装问题

这是新手最容易踩的坑。Gurobi安装完之后，直接运行代码会报license not found或者GurobiError: Unable to create environment之类的错。大部分情况是license没有激活。

正确的激活方式是：去官网申请学术license，拿到一个grbkey文件或一串key，然后在命令行执行grbgetkey <你的key>，按提示写入本机目录。如果是校园网或代理环境，有时会出现 license server 连不上的问题，这时可以离线激活，把license放到用户根目录下的gurobi.lic文件里。

另外说一句，Gurobi和Matlab的关联问题也出现了热搜词里，其实原理一样，license是通用的，装了Gurobi的optimizer之后，Matlab、Python、C++等所有接口都能用，不需要重复购买。

5.2 模型无解时的排查方法

跑MILP最崩溃的就是model.status == GRB.Status.INFEASIBLE。我复现FSTSP时至少碰到过5次以上无解的情况，每次都是头大。后来我总结出了一套高效的排查流程：

首先用Gurobi内置的IIS（Irreducible Inconsistent Subsystem）计算工具，找出导致不可行的最小约束集合。Gurobi的Python API里只需要这几行代码：

model.computeIIS() model.write("model.ilp")

打开model.ilp文件后，Gurobi会明确标出哪些约束是“conflicting”的，也就是导致无解的约束。这个方法比肉眼扫上万条约束高效太多了。

我实际遇到的无解原因主要有三个：一是无人机续航约束太紧，D_max设得太小，导致某组任务总是超出限制；二是时间大M值设置不合理，在某些极端情况下切掉了合法解；三是访问约束写错，导致某个客户既不能被卡车访问也不能被无人机服务。

5.3 求解时间过长、解质量差怎么办

10个节点还好，跑到15个节点的时候，Gurobi求解时间会明显上升，甚至几分钟都求不出最优解。这时候，有几个调优手段非常有效。

第一个是给Gurobi设置合理的求解参数。比如设置时间上限：

model.setParam('TimeLimit', 300) # 5分钟 model.setParam('MIPGap', 0.01) # 允许1%的gap

第二个是给模型添加初始可行解。你可以先用一个简单的贪婪算法或最近邻算法生成一条可行的卡车路径和无人机任务分配，通过model.setStart()或model.addMipStart()提供给Gurobi。有初始解可以让Gurobi的启发性更强，尤其是大算例时比较明显。

第三个是尝试不同的MIPFocus参数。如果想让Gurobi更侧重于改进目标值，可以设model.setParam('MIPFocus', 1)；如果更侧重于证明最优性，设model.setParam('MIPFocus', 2)。在实际测试中，10节点规模默认参数就够了，15节点以上，MIPFocus=1会稍微快一点。

5.4 代码复现和论文结果对不上的原因

很多论文复现党最痛苦的问题就是：代码跑出来的结果和论文表格里的结果完全不一样。这里有几种可能性：

• 论文里用了不同的坐标实例。论文通常会在附录里给出详细的节点坐标和参数，或者使用公开数据集（比如TSPLIB库）。如果你用的坐标不一样，结果自然对不上。
• 参数设定不同。比如卡车速度、无人机速度、装卸时间、最大续航里程。这些参数哪怕差一点点，最优路径都会大变。
• 目标函数或者约束的细微差异。有些论文把无人机从配送中心出发也算作一次“任务”，有些论文则要求无人机必须和卡车一起出发，这会导致结果差异很大。
• 论文只报告了某个特定算例的结果，而不是所有算例的均值。如果只看某一组数据，容易误以为代码有问题。

我的建议是：先把论文里的“小规模算例”完整跑通，确保每一步的结果和论文表格一致，再扩展到更大的规模。不要一开始就在大算例上纠结，那样排查问题会非常困难。

5.5 关于无人机正射拼接、巡检等其他热词的说明

搜索热词里出现了不少和无人机正射拼接、巡检、飞控相关的内容。这里说明一下，FSTSP解决的是“调度与路径规划”层面的问题，也就是在给定配送任务后，如何规划卡车和无人机的行走路线。它不涉及无人机底层的控制算法（比如PID控制器）、视觉感知、正射拼接这些硬件和感知层的问题。

但如果你在做一个完整的“无人机配送演示系统”，调度层（FSTSP）和感知层（视觉起降、避障）确实是需要协同的。我在后续的扩展计划里，就打算把FSTSP生成的路径，接入一个简单的无人机飞行控制仿真环境，验证路径的可飞性。这个方向有很大的扩展空间。

6. 扩展与进阶方向

6.1 从单机到多机、多卡车扩展

FSTSP是最基础的“一辆卡车+一架无人机”模型。真实场景里，一个配送中心往往有多辆卡车，每辆卡车上还可能搭载多架无人机，这时候就演变成MTSP（Multiple Traveling Salesman Problem）+多无人机的混合问题。建模时，卡车路径变量需要增加车辆下标，无人机任务变量需要增加无人机下标，约束数量会大规模增长，求解难度成倍上升。

如果一定要用Gurobi求精确解，我建议先限制卡车数量为2到3辆、每辆车最多2架无人机，并且用对称性破缺约束（比如规定车1必须访问编号最小的客户点，从而消除车辆之间的对称性）。否则求解器会被大量对称解卡死。

6.2 引入时间窗、载重约束和动态需求

实际配送中，客户通常有服务时间窗（比如上午10点到12点之间必须送到）。这时候需要给模型新增时间窗约束和无人机载重约束，复杂度又上一个台阶。更进阶的是动态FSTSP，即新的订单在配送过程中实时到达，需要在线重新规划和调度。精确求解在线版本不现实，工业界通常会用rolling horizon + 启发式算法的方案。

6.3 对偶问题模型代码优化与模型压缩

复现论文时我还发现一个提升性能的技巧：模型预求解（presolve）可以做很多隐式约束的自动推导，但依然有很多逻辑可以手动压缩。例如，有些二元变量实际上可以通过约束直接消除掉。我在最终版本里，把w[i,k]和z[k,j]做了一次变量合并，只在模型中保留三维变量y[i,k,j]，表示“无人机从i起飞服务k后在j汇合”，这样变量数量减少不少，约束也简化了。

代价是约束表达会稍微复杂一些，而且模型维数提高了。但Gurobi对这种三维二元变量的处理效率非常高，实际运行比二维变量版本快了不少。这个优化在15节点以上的实例里表现尤其突出。

6.4 从精确解到启发式算法：为什么还需要其他算法

虽然Gurobi能解决10到15节点的FSTSP，但真实场景可能有几百上千个客户节点，精确求解完全不可行。实际项目里，我更多的做法是：先用Gurobi在小规模实例上做验证和参数调优，然后把问题转成基于自适应大邻域搜索（ALNS）的启发式算法，去处理更大规模的算例。

ALNS的基本思想特别适合FSTSP：卡车路径可以用remove/insert算子来破坏和重建，无人机任务可以通过在卡车路径上的某段“插入无人机支线”来生成。整个过程可以反复迭代，搜索空间大、灵活度高。我在后续的输出里会专门整理一套FSTSP的ALNS实现，对比Gurobi的精确解，验证启发式算法的gap。

7. 一点个人心得

复现这篇FSTSP论文，整个过程像是一场漫长的拉锯战。最大的感受是：读论文和写代码完全不是一回事。论文里一个简简单单的约束式子，落到代码里可能要拆成七八行循环和判断；论文里没写清楚的参数，落到代码里每一个都需要自己猜。

但恰恰是这种艰难复现，让我真正理解了模型的结构。比如无人机的起飞节点和汇合节点为什么不能合并建模，比如大M值的设置如何影响求解效率，再比如为什么有些论文里的结果明显不真实（路径交叉、时间冲突）却照样发表——很可能就是模型里漏了某个关键约束。

如果你也正在复现类似的论文，我的建议很简单：不要急着写代码，先把论文里的每个变量、每个约束手写推一遍，然后用小算例逐步验证。最后再提醒一次，Gurobi的license记得先激活，否则所有准备工作都会卡在第一步。这个模型本身并不复杂，复杂的是细节，而细节，才是论文复现的真正门槛。