无人机配送路径规划：K-means聚类与遗传算法Matlab实战详解

简介本资源面向物流优化、智能算法与无人机系统研究领域的高校师生及工程技术人员聚焦多无人机协同包裹递送中的区域划分与路径规划两大核心问题。通过融合K-means聚类实现中心辐射式配送区域自动划分结合遗传算法求解全局最优配送路径在载荷、续航与地理约束下提升整体配送效率并降低成本。压缩包仅含1个HTML文件10KB完整嵌入算法原理说明、流程图解、Matlab代码实现及仿真结果可视化代码支持地理坐标输入、动态任务调整与多场景路径重优化可直接运行验证算法有效性。目前已有23人学习下载适用于算法教学演示、课程设计实践及科研原型快速验证为无人配送系统建模与仿真提供即开即用的技术支撑。 最近带学生做无人机物流仿真项目时把以前那套“中心辐射模型 K-means聚类 遗传算法优化”的Matlab代码翻出来重写了一遍。这三个词组合在一起几乎是无人机配送路径规划领域的“老三样”论文里一抓一大把但真到了自己上手从零写到能跑、能调、能可视化我才发现里面全是细节K-means的K到底怎么定遗传算法的染色体怎么编码才能表达“多架无人机同时配送”适应度函数里的续航惩罚项该给多大权重这些问题在大多数参考代码里都不会直说只有踩过坑才明白。这篇文章就把我重构三版之后稳定运行的完整思路整理出来。从需求建模、聚类选址细节、GA路径规划的编码与算子设计到Matlab调试过程里的高频报错和避坑方法一次性讲清楚。代码结构也会拆开说明适合正在做物流无人机、智能优化算法课程设计或者想用K-means 遗传算法解决组合优化问题的读者对照参考。1. 中心辐射为什么无人机送货要先“聚”再“送”中心辐射模型不是无人机行业发明的它最早大规模用在航空客运网络里几个大型枢纽机场之间用主干航线连接乘客再从枢纽坐支线飞机去中小城市。这套思路放在快递物流里就是经典的hub-and-spoke结构包裹先集中到区域分拣中心再分拨到末端网点最后送到收件人手里。无人机配送之所以天然适配这种结构核心原因只有一个电池续航扛不住点对点长距离飞行。1.1 无人机配送的续航悖论假设一个城市里有100个客户点仓库在城郊如果让一架无人机从仓库出发挨个配送路径总长度很可能超过100公里而市面上多数多旋翼物流无人机的有效航程也就15到30公里。哪怕任务分配给多架无人机每架都从仓库直飞客户点再返回仍然会出现大量长距离冗余飞行。实际项目中我们通常会在配送区域内设置若干个起降站点无人机从站点出发只负责服务站点周边半径3到5公里内的客户。这样虽然多了一次“中转”但每架无人机的飞行距离大幅缩短换电、装载、充电都可以在站点完成整体效率反而更高。1.2 中心辐射网络的分层逻辑这套模型在我们的Matlab代码里是这样分层的第一层大型分拨中心。相当于“根节点”负责接收所有包裹并向各站点补货。第二层配送站点。通过K-means聚类生成每个站点服务一片客户簇选址目标是最小化簇内客户到站点的平均距离。第三层末端客户点。即需要交付包裹的具体坐标它们被聚到某个站点名下再由站点派出无人机逐一访问。那么整个路径规划问题就拆解为两个子问题一是站点选址与客户分配K-means负责二是每个站点内多无人机的访问顺序遗传算法负责。实际代码中我们把这两个问题串成一个Pipeline先用kmeans函数把客户点聚成K类类中心作为初始站点然后把每类客户看作一个TSP旅行商问题实例用GA规划无人机访问顺序。1.3 模型假设与输入数据设计为了让问题可解我采用了以下假设条件这些也是大多数同类仿真论文的默认设置所有客户点坐标已知无人机飞行路径为直线等效为欧氏距离。每架无人机运载能力相同单次最多携带一定数量的包裹后必须返回站点补货。无人机飞行速度恒定续航里程有上限超过上限的路线会被惩罚。站点位置可自由选择暂不考虑空域限制、楼顶平台、信号覆盖等真实条件。这些假设在第一章先定下来后面写遗传算法适应度函数时就不至于反复推翻重来。实际代码中我定义了一个全局结构体param把无人机数量、最大航程、载重上限、客户坐标等放进去后续所有脚本都从param读取避免硬编码。2. K-means选址聚类中心不等于起降点很多初学K-means的人以为把客户坐标丢进kmeans函数输出聚类中心就可以直接当无人机起降点了。实际使用中完全不是这么回事这里有两个很容易被忽略的点一是K值怎么确定二是聚类中心是几何中心不代表现实中可以落地的位置。2.1 K-means在这一步到底在算什么我们的目标是把N个客户点划分到K个站点使得每个客户到它所属站点的距离之和最小。用数学语言描述就是[ \min \sum_{i1}^{K} \sum_{x \in C_i} |x - \mu_i|^2 ]其中(C_i)是第i个簇(\mu_i)是第i个簇的中心。这个目标函数对应的是“簇内平方和”也就是Matlab里kmeans默认优化的目标。在无人机场景里站点是终点也是起点这一步的核心价值是把“服务范围”圈定出来让每个站点的客户数量大致均衡从而避免某架无人机累死、另一架闲死。2.2 K值选取肘部法则里的数学K值不是越大越好。K值太小站点覆盖范围太大无人机单程飞不到K值太大站点数量多基础设施成本过高而且每个簇里客户太少遗传算法容易退化。项目里我一般用肘部法则辅助判断。Matlab中可以用以下代码绘制簇内平方和随K值变化的曲线rng(42); X [customers.x, customers.y]; K_list 1:10; sumd_list zeros(size(K_list)); for i 1:length(K_list) [~, ~, sumd] kmeans(X, K_list(i), Replicates, 10, MaxIter, 300); sumd_list(i) sum(sumd); end plot(K_list, sumd_list, o-); xlabel(K); ylabel(Sum of within-cluster distances);这里我刻意使用了Replicates参数。kmeans默认只跑一次初始中心完全随机结果方差很大。设置Replicates10表示从10组不同的初始中心出发返回目标函数最小的一次结果这样结果更稳定。肘部法则的“肘”出现在曲线下降趋势明显变缓的位置。比如我测试某个城市数据时K从1到4下降非常陡峭从4到10逐渐平缓那个拐点K4就是合理的选择。当然把K取自变量还要结合实际业务约束比如手里只有5架无人机那K就不能超过5否则站点多但无人机调度不过来。2.3 Matlab聚类实现与初始中心敏感性K-means对初始中心敏感这是老话题了。Matlab的kmeans虽然内置了k-means策略Start参数默认为plus但依然可能陷入局部最优。常见的优化方式是多组初始中心并行跑然后取目标函数最小的解这正是前面代码里Replicates10的作用。另一个容易被忽略的细节是坐标尺度对聚类结果的影响。如果客户点分布在经度、纬度这种量纲不一致的坐标上直接聚会产生严重偏差。我的做法是把经纬度投影成平面直角坐标常见选择是Web Mercator或者高斯-克吕格投影然后再做标准化处理X [customers.x, customers.y]; X zscore(X); % 标准化注意zscore会把坐标中心化、缩放成单位方差这样聚出来的簇结构更稳健但聚类中心返回的是标准化后的坐标使用时必须反向还原才能得到真实坐标。2.4 从簇心到真实起降点的修正这一步是我在实际项目中最想强调的经验K-means给出的聚类中心是一个“几何质心”但在真实城市场景里这个位置可能位于马路中央、湖面上或者建筑群正中间无人机根本无法起降。所以我会在聚类之后增加一个“最近可用点匹配”步骤提前准备一份“候选起降点”列表例如城市中允许无人机起降的屋顶平台、空地和站场然后对每个聚类中心在候选点列表里找距离最近且未被占用的点作为实际站点。这个逻辑在Matlab里用dsearchn就能实现candidate_sites [sites.x, sites.y]; centroids C; selected_sites zeros(K,2); for i 1:K nearest_idx dsearchn(candidate_sites, centroids(i,:)); selected_sites(i,:) candidate_sites(nearest_idx,:); end这一步虽然属于工程细节但在论文里写“站点选址采用K-means聚类与最近可用点匹配相结合”会让评审觉得你考虑到了真实约束而不只是跑了个数学模型。3. 遗传算法跑路线编码、适应度、算子三件套站点分好了接下来每个站点内的客户点顺序就交给遗传算法。为什么不用精确算法因为这是个NP难问题站点内客户一旦超过15个暴力枚举就变得不现实。遗传算法虽然没有理论上的最优解保证但胜在通用、容易实现而且和Matlab的矩阵运算天然契合。3.1 为什么不用精确算法而选GA对于单个站点10个客户TSP的路径数量是10! 3,628,800条计算机还可以枚举但站点内客户到20个时20!已经超出任何民用机器的计算能力。精确算法里像分支定界法能在小规模问题下找到最优解但问题规模稍微变大就会指数爆炸。遗传算法不一样它不追求枚举全部路径而是从一组随机初始解出发通过选择、交叉、变异一步步逼近较优解虽然无法保证找到全局最优但在工程上“够用”且计算时间可接受。3.2 多无人机任务分配的双染色体编码单站点只有一架无人机时染色体用一条整数排列表示客户访问顺序即可。但实际场景往往是多个站点、每个站点有多架无人机怎么设计编码就各有门道了。我采用了一种双染色体编码方案染色体1任务分配染色体长度等于客户数每个基因位存储该客户被分配给哪个站点。染色体2访问顺序染色体长度等于客户数存储一个全局客户序号的排列。解码时先根据染色体1把客户按站点分组再按染色体2中各客户的相对顺序进行排序从而得到每个站点内各无人机的访问序列。这个方案的好处是交叉和变异可以分别作用在两条染色体上既保证了站点分配的变化又保证了访问顺序的多样性。举个具体例子假设有6个客户2个站点。染色体1为[1,2,1,2,2,1]染色体2为[4,6,2,1,5,3]。那么客户1、3、6属于站点1客户2、4、5属于站点2。按染色体2的顺序排序站点1的访问顺序是6→1→3站点2的访问顺序是2→4→5。3.3 适应度函数距离只是底线约束才是关键适应度函数直接决定GA搜索方向我只把“总飞行距离短”作为主要目标还不够必须把续航约束、载重约束一起纳入否则算法会一厢情愿地规划出“理论距离最短但无人机根本飞不到”的路线。设某条完整路径的总距离为(D)超出最大续航的距离为(E_{exceed})超出载重量的包裹数为(P_{exceed})那么适应度函数定义为[ fitness D \lambda_1 \cdot E_{exceed} \lambda_2 \cdot P_{exceed} ]其中(\lambda_1)和(\lambda_2)是惩罚系数。代码中我会这样处理function fit compute_fitness(chrom1, chrom2, param) % 解码 routes decode(chrom1, chrom2, param); total_distance 0; penalty 0; for i 1:length(routes) route routes{i}; if isempty(route) continue; end dist route_distance(route, param.distance_matrix); total_distance total_distance dist; if dist param.max_range penalty penalty (dist - param.max_range); % 超出部分线性惩罚 end % 载重约束 package_num param.package_demand(route); if sum(package_num) param.max_payload penalty penalty (sum(package_num) - param.max_payload) * 10; end end fit total_distance penalty; end这里有个关键点惩罚系数的量纲要和距离保持兼容。param.max_payload的单位是“个包裹”如果罚款系数设得太大比如1000那么算法会完全无视距离而只想着别超载如果设得太小比如0.1那么超载的路线会被当作好路线保留下来。经验上我会先统计初始种群中的距离均值再把超载惩罚系数设置成距离均值的1/10左右比如距离均值是1000那超载1个包裹罚100。实际上在测试中我还会对超出续航的部分采取“硬惩罚”和“软惩罚”结合超出1%以内只加一个固定惩罚值超出超过1%则返回一个极大值直接淘汰该个体。这样处理的好处是避免算法在续航边界反复试探。具体阈值要根据无人机型号参数来设我这里用的是10%。3.4 交叉与变异针对顺序编码的算子设计遗传算法的算子必须和编码方式匹配。对顺序编码的染色体2不能直接使用单点交叉否则很容易产生重复客户点、丢失部分客户点的非法解。我常用的顺序交叉Ordered CrossoverOX思路是随机选择两个交叉点将父代1两交叉点之间的基因段保留到子代1对应位置子代1剩余位置按照父代2中非重复基因的相对顺序依次填充父代2同理生成子代2。这种操作保证了每一个客户点只出现一次适合路径顺序类编码。MATLAB里我写了一个简单的OX交叉函数核心代码如下function child ox_crossover(p1, p2) n length(p1); point1 randi(n - 1); point2 point1 randi(n - point1); child zeros(1, n); child(point1:point2) p1(point1:point2); % 从p2中提取非重复基因 remaining p2(~ismember(p2, child)); idx 1; for i 1:n if child(i) 0 child(i) remaining(idx); idx idx 1; end end end变异算子我用的是“两点交换片段倒置”混合变异以一定概率随机交换两个基因位或者随机选一段基因倒置。两类变异以概率各半触发这样既能带来扰动又能保持路径的局部结构。对染色体1也就是站点分配基因我用的是多点变异随机挑若干基因位把其值改成另一个站点的编号。这个变异力度不能太大否则每个站点客户量严重失衡导致优化过程剧烈震荡。3.5 参数表与调参经验遗传算法参数我常用以下初始值实际调优建议在此基础上微调参数取值说明种群规模100客户多时可增加到200最大迭代次数200看收敛曲线决定交叉概率0.85不宜超过0.9避免随机搜索变异概率0.15太大容易破坏优秀解精英保留数2每代直接复制到下一代惩罚系数λ110按续航超出距离比例调惩罚系数λ2100按平均距离的1/10设定调参的经验是先固定变异概率和交叉概率只调惩罚系数等路线不再出现超续航的个体之后再调种群规模和迭代次数。不要一上来就同时改所有参数否则你根本不知道是哪个参数让结果变好的。4. 调试与避坑Matlab实现中的高频问题这部分是我实际调试代码过程中遇到的最多、最典型的坑。每一个我都在代码注释里标了醒目的警示符号这里专门整理出来。4.1 坐标格式带来的距离失真第一版代码里我直接用经纬度做欧氏距离运算导致结果完全失真。经度1度在赤道附近大约是111公里但在高纬度地区可能只有几十公里纬度1度则基本恒定在111公里左右。如果直接把经纬度丢进sqrt((x1-x2)^2(y1-y2)^2)在哈尔滨和在三亚算出来的同度数距离会差很多。解决办法是先用deg2km函数或者等距投影转成平面坐标。Matlab中可以用lldistkm函数需要Mapping Toolbox或者自定义投影函数% 简单等距投影适用于城市尺度 lat0 mean(customers.lat); lon0 mean(customers.lon); x (customers.lon - lon0) * 111.32 * cosd(lat0); y (customers.lat - lat0) * 110.94;如果你的客户坐标只是仿真生成的二维坐标这一步可以省略但代码里我还是保留了投影函数入口为的是以后接入真实试验数据时不用改大框架。4.2 K-means空簇问题kmeans在某些随机初始中心下可能产生空簇也就是某个聚类中心没有分配到任何客户。这种情况下遗传算法里对应站点的客户列表为空计算适应度时会出现NaN或索引越界报错。我在实际代码中的处理方式是在kmeans之后加一个检查如果任何簇为空重新运行聚类并递增Replicates。还有一个更实用的思路先给每个簇强制分配距离该簇中心最近的客户然后再跑一轮迭代修正簇心。这有点类似带种子约束的K-means变体在Matlab中实现也不复杂。4.3 GA早熟如何判断和干预遗传算法最常见的现象是迭代到50代左右最优适应度就不再下降了但解明显还不是最优解。这就是早熟收敛。我通常用两个指标判断最优个体适应度连续50代没有变化种群平均适应度和最优适应度差值小于某个阈值说明种群多样性丧失。遇到早熟我不会盲目加大变异概率而是采用“自适应变异”当检测到最优值停滞时把变异概率临时提升到0.3维持若干代后再降回来。另外还可以引入“移民策略”每20代随机生成一批新个体加入种群替代最差的一些个体。代码里我用一个if mod(gen, 20) 0的触发器来做简单有效。4.4 惩罚项权重或无穷大警告适应度值出现Inf或者NaN时大概率是惩罚系数设置不合理。在较早期版本中我直接用fit 1/(total_distance penalty)作为适应度当total_distance penalty为0时就会爆炸。后来我改成直接使用total_distance penalty作为目标函数遗传算法内以“最小化”为方向只在选择操作时把目标值越小视为适应度越高。如果你坚持用“适应度越大越好”的形式一个更稳妥的转换是fit 1 / (1 total_distance penalty);这样能保证分母恒大于0不会因为惩罚过大出现除零问题。4.5 结果可复现性Matlab的rand和randi每次运行结果不同导致两次跑出的路径差异很大写报告时很尴尬。我在所有随机操作前统一加了一行rng(2025);这样每次跑出的初始种群、聚类中心、交叉变异位置都完全一致方便对比不同参数对结果的影响。注意rng(2025)的位置要放在所有随机操作之前最好放在脚本第一行因为kmeans内部也用了随机数生成器如果你的rng放在kmeans之后才设置控制不了聚类的随机性。可视化部分我也做了个简单处理用不同颜色区分不同站点服务的客户簇用带箭头的线条表示无人机飞行轨迹。Matlab里可以用quiver画箭头也可以用plot加annotation实现。实跑效果来看一张图能把站点、客户点、航线全部展示出来比单纯输出数字直观得多。5. 从代码到天空这套方案的边界与延伸方向虽然这套仿真跑起来效果不错但我心里清楚它离真实可部署的无人机配送系统还有一段距离。最明显的问题是第一章节我就提到的“直线距离假设”——城市里无人机虽然不受地面道路限制但受空域管制、高楼遮挡、信号干扰、气象条件等多重因素影响实际飞行轨迹永远不可能是标准直线。5.1 当前模型的简化假设除了直线距离还有几个容易忽视的简化电池能耗没有和载重、风速关联真实无人机带2kg包裹和带5kg包裹的续航完全不同客户时间窗没有考虑所有订单都要求“尽快送达”没有优先级站点容量没有约束假设所有站点都能无限停放无人机、无限补货无人机之间的冲突规避没有考虑多机同时飞行时可能存在航线交叉。这些假设决定了这套代码适合做算法验证和技术预研但不适合直接拿去做真实空域规划。如果真想落地至少要把电池能耗模型和空域栅格化建模加进去。5.2 从静态到动态当前流程是“先聚类再规划”属于静态一次性优化。实际配送场景是订单不断到达的所以更合理的方式是多阶段滚动优化每当一批新订单到达重新执行聚类和路径规划保留已经起飞的无人机任务不变只重规划未分配订单。这个方向上可以用Matlab的Timer或实时脚本实现但计算压力会明显增加。5.3 可能的扩展这套代码框架稍微改一改就能扩展到很多方向。比如把目标函数从“总距离最短”改为“总能耗最低”就需要给每个客户点增加“重量需求”属性并在适应度函数中引入载荷相关的能耗模型又比如加入时间窗约束就需要在适应度函数里加入时间惩罚项同时把遗传算法的解码过程改成带时间轴的任务插空排序。我后续打算把多站点之间的补货路径也纳入优化让无人机既送末端客户也承担站点间调拨任务这样整个中心辐射模型就更完整了。如果只是想应付课程设计或论文仿真这套方案已经完全够用如果想进一步挑战自己可以顺着“多目标优化NSGA-II 无人机调度”的方向继续深挖。毕竟K-means和遗传算法都只是工具箱解决实际问题时真正值钱的是对约束的建模能力和对算法行为的理解。本文还有配套的精品资源点击获取