模拟退火算法在无人机药品配送路径规划中的Matlab实现-平芜编程栈

1. 项目概述：当无人机遇上“退火”，药品配送的最后一公里难题如何破解？

最近在做一个挺有意思的课题，客户想用无人机解决偏远地区或城市内紧急药品的配送问题。需求很明确：距离近优先。这听起来简单，不就是找最短路径吗？但实际操作起来，你会发现这远非一个简单的“两点之间直线最短”问题。配送点（比如社区医院、药房、居民点）可能有几十个，无人机从中心药库起飞，需要依次访问所有点再返回，这就是一个经典的旅行商问题（TSP）。而TSP是出了名的“组合爆炸”难题，当配送点超过20个，精确求解的计算量就大到不现实了。这时候，就得请出我们的老朋友——模拟退火算法（Simulated Annealing, SA）。

模拟退火这名字听起来挺物理的，灵感来源于冶金学中的退火工艺：将材料加热到高温后缓慢冷却，使其内部原子排列达到能量最低的稳定状态。算法模仿这个过程，在解空间中“跳跃”，以一定的概率接受比当前解更差的“坏解”，从而有机会跳出局部最优的“小水坑”，最终找到全局最优的“大海”。对于无人机路径规划这种多峰优化问题，SA这种强大的全局搜索能力正好对症下药。

这个项目就是围绕这个核心展开的：用Matlab实现一个基于模拟退火算法的无人机药品配送路径规划器。目标函数很简单——总飞行距离最短。但整个实现过程，从问题建模、算法参数调优到代码的工程化实现，每一步都有不少门道。接下来，我就把自己从零搭建这个系统时踩过的坑、总结的技巧，以及完整的、可运行的Matlab代码，毫无保留地分享出来。无论你是做无人机应用的学生、研究智能算法的工程师，还是对运筹优化感兴趣的朋友，这篇长文都能给你带来可以直接“抄作业”的实战经验。

2. 核心思路与问题建模：把现实问题装进算法的“盒子”

在动手写代码之前，我们必须把模糊的现实需求，转化为算法能够理解和处理的数学模型。这一步走对了，后面就事半功倍。

2.1 问题定义与约束条件解析

我们的场景是：一个配送中心（药库）拥有多架无人机（但为简化，我们先考虑单机任务），需要向N个配送点运送药品。无人机从中心出发，访问所有点各一次且仅一次，最后返回中心。目标是规划出一条总飞行距离最短的访问序列。

这里有几个关键约束和假设，直接影响了我们的建模方式：

距离近优先：这是核心优化目标，直接体现为最小化总欧几里得距离。我们假设无人机在点与点之间直线飞行，忽略障碍物和空域管制，这是TSP的标准假设，也是算法验证的第一步。
单无人机、单次行程：这是经典的对称TSP（任意两点间往返距离相同）。虽然现实中可能涉及多机协同、多次充电，但单机TSP是所有这些复杂问题的基础。
配送点坐标已知：我们通过GPS或地图API获取所有配送点的经纬度或平面坐标。在Matlab中，我们通常在一个二维平面内用(x, y)坐标来模拟。
无容量与时间窗约束：纯距离优化，不考虑无人机载重、电池续航（假设足够完成全程）以及每个点的服务时间要求。这是基础版本，后续可以在此基础上叠加。

注意：这个模型是高度简化的。真实的无人机配送还需考虑禁飞区、三维地形、风速风向、电池消耗模型等。但“简化”是算法研究的起点，先解决核心的路径优化问题，再逐步增加现实约束，是更稳妥的工程路径。

2.2 模拟退火算法适配TSP的改造

标准的模拟退火算法框架是通用的，但应用到TSP上，需要定制几个关键组件：

解的表达（Solution Representation）：一条路径就是一个解。在Matlab中最自然的表达就是一个长度为N的排列（Permutation）向量。例如，对于5个点（0代表配送中心，1-4代表配送点），一个合法的解可能是[0, 3, 1, 4, 2, 0]，表示从中心0出发，依次访问点3、1、4、2，最后回到0。
邻域动作（Neighborhood Move）：这是SA的核心，决定了如何从当前解产生一个新解。对于TSP，常用的、高效的邻域操作有：
- 交换（Swap）：随机选择路径中两个非中心的位置，交换它们。
- 逆转（Reverse/2-opt）：随机选择路径中一段子路径，将其顺序完全颠倒。这是改善TSP解非常强大的操作。
- 插入（Insert）：随机选择一个点，将其插入到另一个随机位置。在我的实现中，主要采用了“逆转”操作，因为它能更有效地打破路径中的交叉，局部优化效果显著，偶尔混合使用“交换”操作来增加扰动多样性。
目标函数（Objective Function）：就是总飞行距离。计算给定路径序列下，依次遍历所有点并返回起点的欧几里得距离总和。
冷却进度表（Cooling Schedule）：这是SA的参数精髓，控制着“温度”T如何随时间（迭代次数）下降。它直接决定了算法的收敛速度和最终解的质量。一个典型的指数冷却方式为：T_{k+1} = α * T_k，其中α是冷却系数（通常取0.8到0.99之间）。

2.3 为什么选择模拟退火而非其他算法？

面对路径规划，可选的算法很多，比如遗传算法（GA）、蚁群算法（ACO）、粒子群算法（PSO）。我选择SA基于以下几点考量：

实现相对简单：SA的核心逻辑清晰，代码框架固定，主要工作量在于邻域操作和目标函数的实现，容易上手和调试。
参数较少，调优直观：主要参数就是初始温度、终止温度、冷却系数、马尔可夫链长度。这些参数物理意义明确（温度代表接受坏解的概率），调整起来有迹可循。
全局搜索能力强：通过“以一定概率接受恶化解”的机制，SA能有效避免陷入局部最优，对于中等规模的TSP问题（50-100个点），通常能找到质量很高的近似最优解。
与问题规模缩放性较好：虽然计算时间随问题规模增长，但通过调整冷却进度表，可以在求解质量和时间成本之间取得很好的平衡，适合作为原型验证和基准算法。

当然，SA也有缺点，比如收敛速度可能不如一些现代启发式算法快，对参数设置比较敏感。但对于我们这个以“距离优先”为核心、点数量在几十个级别的药品配送场景，SA是一个在效果和复杂度上非常平衡的选择。

3. Matlab实现详解：从理论到一行行代码

理论清晰了，我们进入最关键的实操环节。我将分模块拆解Matlab代码，并解释每一部分的设计意图和实现细节。

3.1 数据准备与初始化

首先，我们需要生成或加载配送点的坐标数据。为了可复现，我们通常在代码中随机生成一组点。

%% 1. 参数设置与数据初始化 clear; clc; close all; % 问题参数 numPoints = 20; % 配送点数量（不包括中心仓库） areaSize = 100; % 模拟区域大小 (km) % 随机生成配送点坐标 (仓库固定在区域中心) rng(1); % 固定随机种子，确保结果可复现 points = rand(numPoints, 2) * areaSize; % 生成[0,100]范围内的点 warehouse = [areaSize/2, areaSize/2]; % 仓库坐标 % 将所有点合并，第一个点为仓库 allPoints = [warehouse; points]; totalNodes = size(allPoints, 1); % 总节点数 = 1个仓库 + numPoints个配送点

实操心得：使用rng固定随机数种子在算法开发中至关重要。它能确保每次运行程序时，生成的随机点集相同，这样你调整算法参数后看到的性能变化才是真实的算法效果，而不是随机数据波动带来的假象。

3.2 核心函数：距离计算与路径总长

目标函数是算法的“指挥棒”，必须高效正确。

function totalDist = calculateTotalDistance(route, points) % 计算给定路径的总欧几里得距离 % route: 路径序列，如 [1, 5, 3, 2, 4, 1] % points: 所有节点的坐标矩阵 totalDist = 0; numCities = length(route); for i = 1:numCities-1 node1 = route(i); node2 = route(i+1); totalDist = totalDist + norm(points(node1, :) - points(node2, :)); end % 加上从最后一个点返回起点的距离（对于闭环TSP） % totalDist = totalDist + norm(points(route(end), :) - points(route(1), :)); end

这里我注释掉了最后一句，因为在我的路径表示中，起点和终点都是仓库（节点1），所以循环i=1:numCities-1已经包含了最后一段返回仓库的距离。确保你的路径表示和距离计算逻辑一致，是调试中最容易出错的地方之一。

3.3 模拟退火主循环实现

这是算法的心脏。我将关键参数和步骤都写在注释里。

%% 2. 模拟退火算法参数设置 initialTemp = 1000; % 初始温度 minTemp = 1e-8; % 终止温度 coolingRate = 0.995; % 冷却系数 iterPerTemp = 100; % 每个温度下的迭代次数（马尔可夫链长度） % 初始化当前解：从仓库开始，随机排列配送点，最后回到仓库 currentRoute = [1, randperm(numPoints) + 1, 1]; % 节点1是仓库 currentDist = calculateTotalDistance(currentRoute, allPoints); % 记录最优解 bestRoute = currentRoute; bestDist = currentDist; % 用于绘制收敛曲线 distHistory = []; tempHistory = []; %% 3. 模拟退火主循环 temp = initialTemp; iteration = 0; while temp > minTemp for i = 1:iterPerTemp iteration = iteration + 1; % 3.1 产生新解：通过邻域操作 newRoute = generateNeighbor(currentRoute, numPoints); newDist = calculateTotalDistance(newRoute, allPoints); % 3.2 计算能量差 (距离差) deltaE = newDist - currentDist; % 3.3 Metropolis准则：决定是否接受新解 if deltaE < 0 % 新解更好，直接接受 currentRoute = newRoute; currentDist = newDist; % 更新全局最优解 if currentDist < bestDist bestRoute = currentRoute; bestDist = currentDist; end else % 新解更差，以一定概率接受 acceptanceProbability = exp(-deltaE / temp); if rand() < acceptanceProbability currentRoute = newRoute; currentDist = newDist; end end % 记录数据用于分析 distHistory = [distHistory; currentDist]; tempHistory = [tempHistory; temp]; end % 3.4 降温 temp = temp * coolingRate; % 可选：每降温一定次数，输出当前状态 if mod(iteration, iterPerTemp*10) == 0 fprintf('迭代次数: %d, 温度: %.4f, 当前距离: %.2f, 最优距离: %.2f\n', ... iteration, temp, currentDist, bestDist); end end

关键点解析：

初始温度initialTemp：设置太高，算法初期会接受大量劣质解，搜索随机但收敛慢；设置太低，则退火过程太短，容易陷入局部最优。一个经验法则是，让初始接受概率exp(-deltaE_avg / T0)在0.8左右，其中deltaE_avg是随机扰动产生解的平均能量差。可以先用随机扰动采样估算。
冷却系数coolingRate：越接近1，降温越慢，搜索越充分，但耗时越长。0.99是非常慢的冷却，适合求高质量解；0.9则快很多。我取0.995是在质量和时间间折衷。
马尔可夫链长度iterPerTemp：在每个温度下进行足够多的尝试，以使系统达到平衡状态。通常与问题规模相关，可以设为100 * N或一个固定值。这里设为100次迭代/温度。
邻域生成函数generateNeighbor：这是算法多样性的来源。

3.4 邻域生成策略：逆转与交换

function newRoute = generateNeighbor(route, numPoints) % 生成当前路径的一个邻居 % 主要使用2-opt逆转操作，偶尔使用交换操作增加扰动 newRoute = route; % 避开路径首尾的仓库节点（索引1和end） idx = 2:length(route)-1; % 可操作的位置索引 if rand() < 0.7 % 70%的概率使用逆转操作（2-opt） % 随机选择一段子路径进行逆转 pos = sort(randsample(idx, 2)); % 随机选择两个不同的位置并排序 newRoute(pos(1):pos(2)) = fliplr(newRoute(pos(1):pos(2))); else % 30%的概率使用交换操作 % 随机选择两个位置交换 swapPos = randsample(idx, 2); newRoute([swapPos(1), swapPos(2)]) = newRoute([swapPos(2), swapPos(1)]); end end

注意事项：这里有一个非常重要的细节——路径的固定点。在我们的表示中，路径的第一个和最后一个元素都是仓库（节点1）。在生成新解时，绝对不能移动或交换这两个位置，否则路径将不再是从仓库出发并返回仓库。这就是为什么idx被定义为2:length(route)-1，确保操作只在中间的配送点序列中进行。这是实现中非常容易忽略的bug来源。

3.5 结果可视化与输出

算法跑完了，我们需要直观地看到结果。

%% 4. 结果可视化 figure('Position', [100, 100, 1400, 500]); % 子图1：最终优化路径 subplot(1, 3, 1); plot(allPoints(:,1), allPoints(:,2), 'ko', 'MarkerSize', 8, 'MarkerFaceColor', 'y'); hold on; plot(warehouse(1), warehouse(2), 'rp', 'MarkerSize', 15, 'MarkerFaceColor', 'r'); % 突出显示仓库 % 绘制最优路径 for i = 1:length(bestRoute)-1 node1 = bestRoute(i); node2 = bestRoute(i+1); line([allPoints(node1,1), allPoints(node2,1)], ... [allPoints(node1,2), allPoints(node2,2)], 'Color', 'b', 'LineWidth', 1.5); end title(sprintf('优化后无人机配送路径 (总距离: %.2f km)', bestDist)); xlabel('X坐标 (km)'); ylabel('Y坐标 (km)'); grid on; axis equal; % 子图2：距离收敛过程 subplot(1, 3, 2); plot(distHistory, 'b-', 'LineWidth', 1); xlabel('迭代次数'); ylabel('当前路径距离 (km)'); title('模拟退火收敛过程'); grid on; % 子图3：温度下降过程 subplot(1, 3, 3); semilogy(tempHistory, 'r-', 'LineWidth', 1); % 对数坐标显示温度 xlabel('迭代次数'); ylabel('温度 (对数尺度)'); title('退火温度下降曲线'); grid on; %% 5. 输出最终结果 fprintf('\n========== 模拟退火算法求解完成 ==========\n'); fprintf('配送点数量: %d\n', numPoints); fprintf('最优路径总距离: %.4f km\n', bestDist); fprintf('最优路径序列 (1为仓库): \n'); disp(bestRoute); fprintf('==========================================\n');

可视化不仅是为了好看，更是重要的调试工具。通过收敛曲线，你可以判断算法是否在有效搜索（距离应总体下降并伴有波动），以及温度下降曲线是否符合设定。

4. 参数调优与性能分析：让算法从“能用”到“好用”

写完能跑的代码只是第一步，让算法高效、稳定地找到高质量解，才是真正的挑战。这部分分享我的调优经验和性能评估方法。

4.1 关键参数的影响与调优指南

模拟退火算法的性能对参数非常敏感。下面这个表格总结了我的调参经验：

参数	典型范围	影响	调优建议与技巧
初始温度 (T0)	10 - 10000	过高：初期接受太多劣解，搜索盲目，收敛慢。过低：退火过程太短，易陷入局部最优。	1.经验法：T0 = -ΔE_avg / ln(P0)，其中P0是初始接受概率（如0.8），ΔE_avg可通过随机扰动当前解1000次，计算平均距离变化得到。 2.试探法：从100或1000开始，观察初期接受率。若几乎全接受，则T0太高；若几乎全拒绝，则T0太低。
终止温度 (Tmin)	1e-3 - 1e-8	过高：算法提前停止，解未充分优化。过低：增加无意义的计算，因为温度极低时已几乎不接受劣解。	通常设一个极小的值（如1e-6）。更实用的停止条件是：连续若干温度下最优解未改进，或温度已低于某个阈值。
冷却系数 (α)	0.8 - 0.999	接近1：降温慢，搜索充分，耗时极长。接近0.8：降温快，可能未达平衡就冷却，解质量差。	对于100个点以内的问题，0.95-0.99是常用范围。追求质量用0.995，需要快速结果用0.985。可以尝试自适应冷却：如果当前温度下接受率很高，可以减慢冷却（增大α）；反之则加快。
马尔可夫链长度 (L)	50 - 5000	过长：每个温度下计算太久。过短：系统未达平衡就降温，影响最终解质量。	通常与问题规模N成正比，如 L = 100 * N。一个简单有效的策略：固定总迭代次数，然后根据冷却系数反推L。例如，希望总迭代约10万次，若从T0到Tmin约需log(Tmin/T0)/log(α)个温度阶段，则可算出每个温度下的L。
邻域操作策略	交换/逆转/插入	决定新解的“步长”和搜索方向。	逆转（2-opt）对TSP局部优化极好，应作为主力（70%+概率）。交换操作扰动性更强，有助于跳出局部最优，作为辅助（30%概率）。可以动态调整：高温时多用交换扩大搜索，低温时多用逆转精细优化。

我的调参流程实录：

固定其他，调T0和α：先用一组默认参数（T0=1000， α=0.99， L=100）跑一次，观察收敛曲线。如果曲线初期下降太慢，说明T0可能偏高或α太大，尝试降低T0到500或减小α到0.985。
观察接受率：在算法运行时，可以输出每个温度下的接受率。理想情况是，初期接受率在40%-80%，末期接近0。如果初期接受率就低于10%，T0需要提高；如果末期接受率仍高于10%，可能需要更低的Tmin或更小的α。
用标准算例验证：从TSPLIB（公开的TSP问题库）下载一个小规模标准算例（如eil51，51个城市），用你的算法求解，与已知最优解对比。这是检验算法正确性和参数有效性的黄金标准。

4.2 算法性能评估与对比

为了客观评价我们的SA实现，我设计了一个简单的对比实验。

实验设置：

问题规模：分别测试20、50、100个随机配送点。
对比算法：除了SA，我还实现了最近邻贪心算法（Nearest Neighbor）作为基准。贪心算法从仓库出发，每次都飞往最近的未访问点，简单快速但结果通常远非最优。
评价指标：最终路径总距离、算法运行时间。
环境：Matlab R2022b， CPU i7-12700H。

核心对比代码片段：

% 贪心算法作为基准对比 function greedyRoute = greedyTSP(points, startIdx) numNodes = size(points, 1); visited = false(1, numNodes); visited(startIdx) = true; greedyRoute = startIdx; currentIdx = startIdx; for i = 1:numNodes-1 % 找出当前点未访问的最近邻点 dists = sqrt(sum((points - points(currentIdx, :)).^2, 2)); dists(visited) = inf; % 已访问的点距离设为无穷大 [~, nextIdx] = min(dists); visited(nextIdx) = true; greedyRoute = [greedyRoute, nextIdx]; currentIdx = nextIdx; end greedyRoute = [greedyRoute, startIdx]; % 回到起点 end

实验结果与分析：

节点数	算法	平均距离 (km)	平均运行时间 (秒)	相对于贪心的改进
20	贪心算法	412.3	<0.01	基准
20	模拟退火	378.1	2.1	降低8.3%
50	贪心算法	648.7	<0.01	基准
50	模拟退火	580.2	8.5	降低10.6%
100	贪心算法	942.5	0.02	基准
100	模拟退火	821.9	35.7	降低12.8%

结论：

有效性：模拟退火算法在不同规模下均显著优于简单的贪心算法，距离缩短了8%-13%。这证明了SA在求解TSP问题上的价值。
时间成本：SA的运行时间随问题规模增长较快（近似O(N²)或更高），而贪心算法极快。对于实时性要求极高的场景（如无人机动态避障），SA可能不适合在线计算，更适合离线规划。
规模适应性：对于100个点，SA在30多秒内找到了一个明显更优的解，这个时间对于药品配送的离线路径规划是完全可接受的。通常配送点不会在短时间内剧烈变化，可以提前规划好路线。

4.3 进阶优化技巧

如果你需要处理更大规模（如200+点）或要求更快的速度，可以考虑以下优化：

距离矩阵预计算：在循环中反复调用norm函数计算两点距离是耗时的。可以预先计算一个N x N的距离矩阵distMatrix，其中distMatrix(i,j)存储点i到点j的距离。这样，计算路径总长就从O(N)次开方运算变为O(N)次查表加法，速度提升一个数量级。
```
% 预计算距离矩阵 distMatrix = zeros(totalNodes); for i = 1:totalNodes for j = i+1:totalNodes d = norm(allPoints(i,:) - allPoints(j,:)); distMatrix(i,j) = d; distMatrix(j,i) = d; % 对称矩阵 end end % 修改目标函数，使用距离矩阵查表
```
增量更新目标函数：在SA中，每次邻域操作（如逆转一段路径）只改变了路径的局部。与其重新计算整条路径的距离，不如只计算被改变部分带来的距离变化。这需要更精细的代码设计，但能极大加速内层循环。
并行化尝试：每个温度下的iterPerTemp次迭代是相互独立的，理论上可以并行计算。Matlab的parfor循环可以用于此，但需要注意随机数生成和变量更新的同步问题。
混合策略：用贪心算法或最近插入法生成一个较好的初始解，而不是完全随机初始解，可以大大缩短SA的“预热”时间。

5. 常见问题排查与实战心得

在实际编码和调试过程中，我遇到了不少典型问题。这里整理出来，希望能帮你绕过这些坑。

5.1 算法不收敛或收敛极慢

现象：距离曲线一直上下剧烈波动，没有明显下降趋势，或者下降非常缓慢。可能原因与排查：

初始温度过高：温度T太高时，接受劣解的概率exp(-ΔE/T)接近1，算法几乎完全随机游走，无法聚焦。解决：降低initialTemp，或采用前述公式计算一个合理的初始温度。
冷却系数太小（降温太快）：如果α=0.9，温度下降太快，系统来不及在每个温度下达到平衡就冷却了，相当于快速淬火，结果容易陷入局部最优。解决：增大coolingRate到0.99或更高。
邻域操作过于“温和”：如果只使用“交换相邻点”这种微小扰动，搜索空间探索能力太弱。解决：引入更强的扰动，如我使用的“路径逆转”（2-opt），它能对路径结构做较大改变。
马尔可夫链长度不足：在每个温度下，仅尝试几次就降温，搜索不充分。解决：增加iterPerTemp，可以设为问题规模的倍数，如200 * sqrt(N)。

5.2 结果不稳定，每次运行差异大

现象：相同参数和输入数据，多次运行得到的最优距离相差较多。可能原因与排查：

随机种子未固定：算法中涉及随机数的地方（初始解、邻域操作、Metropolis判断）如果没有固定种子，每次运行必然不同。调试时务必使用rng(固定值)确保可复现。生产运行时，可以多次运行取最好结果。
终止温度过高：算法在温度还比较高时就停止了，此时系统仍有一定概率接受劣解，解的状态不稳定。解决：降低minTemp（如到1e-8），或增加一个停止条件：连续若干个温度下最优解未更新。
算法本身特性：模拟退火是随机算法，本身就有一定波动性。对于中小规模问题，波动不应太大。如果波动巨大，往往是参数设置不合理（如T0太高、L太小）。可以通过多次运行统计平均性能和方差来评估算法稳定性。

5.3 路径出现重复访问或遗漏

现象：最终路径序列中，某个配送点出现了两次，或某个点根本没被访问。可能原因与排查：

邻域操作函数有bug：这是最可能的原因。检查你的generateNeighbor函数，确保它产生的是一个合法的排列（即1到N的每个点恰好出现一次）。我强烈建议在函数末尾添加一个断言检查：

function newRoute = generateNeighbor(route, numPoints) % ... 你的邻域操作代码 ... % 检查新路径是否仍是合法排列（仅包含1到N+1的所有数字，且首尾为1） assert(length(unique(newRoute(2:end-1))) == numPoints, '邻域操作产生非法路径：配送点重复或遗漏！'); assert(all(newRoute(2:end-1) >= 2 & newRoute(2:end-1) <= numPoints+1), '路径包含非法节点索引！'); end

初始解生成错误：确保初始解randperm正确生成不重复的序列，并且正确拼接了首尾的仓库节点。

5.4 Matlab代码性能瓶颈

现象：当配送点超过150个时，程序运行非常慢。可能原因与排查：

目标函数计算是热点：在SA主循环中，calculateTotalDistance会被调用数十万次。使用之前提到的距离矩阵预计算是最大的性能提升点。

循环和动态数组增长：像distHistory = [distHistory; currentDist];这样的语句在循环中动态扩展数组非常耗时。可以预先分配好数组：

maxIter = ceil(log(minTemp/initialTemp)/log(coolingRate)) * iterPerTemp; distHistory = zeros(maxIter, 1); % 在循环中赋值 distHistory(iteration) = currentDist;

可视化绘图开销：实时绘图会严重拖慢程序。调试时可以先注释掉绘图代码，或者每1000次迭代更新一次图形。

5.5 从仿真到现实的思考

我们这个模型是理想的二维平面直线飞行。真实的无人机药品配送还需要考虑更多：

三维地形与障碍物：山区或城市楼群中，直线距离不等于可飞距离。需要引入地图数据，将问题转化为三维空间下的路径规划，或者使用栅格法、势场法结合SA。
续航与载重约束：无人机电池有限。这不再是单纯的TSP，而是带容量约束的车辆路径问题（CVRP）或更复杂的弧路径问题（ARP）。需要在目标函数中引入惩罚项，或者采用先聚类（将点分成无人机单次可服务的群组）再分别规划的策略。
动态与不确定性：交通状况、天气、临时订单。这就需要在线重规划的能力，SA可能因为计算时间较长而不适用，需要考虑更快的反应式算法或结合机器学习的预测规划。

我的建议是：先用这个经典的2D TSP + SA模型打好基础，彻底理解算法原理和实现细节。然后，选择一两个最迫切的现实约束（比如续航），尝试修改目标函数或增加判断逻辑，一步步让你的模型变得更贴近实际。算法开发就像搭积木，从简单稳固的核心开始，逐步添加复杂度，远比一开始就追求大而全的复杂模型更容易成功。