SD-MTSP：利用KOA算法优化单仓库多旅行商问题的MATLAB实现，可灵活调整数据集与参数
SD-MTSP：开普勒优化算法KOA求解单仓库多旅行商问题MATLAB（可更改数据集，旅行商的数量和起点）
最近在折腾多旅行商问题的时候，发现了一种挺有意思的优化算法——开普勒优化算法（KOA）。这玩意儿把天体运动规律套用到优化问题上，特别是解决需要同时考虑路径总长度和任务均衡度的单仓库多旅行商问题（SD-MTSP），效果还挺让人惊喜的。今天咱们就手把手实现个MATLAB版本，顺便聊聊代码里那些有意思的设计点。
先简单说下SD-MTSP的应用场景：假设你手头有五个外卖小哥，要从同一个奶茶店出发给全城客户送奶茶，怎么规划路线能让总电瓶车里程最少，同时保证每个小哥送的单量差不多？这就是典型的单仓库多旅行商问题。
直接上核心代码结构。咱们的脚本主要包含这几个部分：
% 基础参数设置 m = 3; % 外卖小哥数量 Citys = rand(50,2)*100; % 50个随机客户坐标 start_point = [50,50]; % 奶茶店位置 % KOA参数 Max_iter = 100; % 迭代次数 N = 20; % 行星数量 dim = length(Citys); % 问题维度=客户数量 % 初始化种群 Population = initialization(N, dim, Citys, m, start_point); % 迭代优化 for iter = 1:Max_iter % 计算适应度（总距离+均衡惩罚项） [fitness, totalDist] = fitnessFunction(Population, Citys, m, start_point); % 行星位置更新（开普勒运动规则） Population = updatePopulation(Population, fitness, iter, Max_iter); % 保留历史最优解 [bestFitness, idx] = min(fitness); BestSol(iter) = struct('Routes',Population(idx).Routes,... 'totalDist',totalDist(idx)); end适应度函数的设计是算法的灵魂。这里我们不仅要考虑总配送距离，还得加上路径均衡度的惩罚项：
function [fitness, totalDist] = fitnessFunction(Population, Citys, m, start_point) % 计算各路径总距离 totalDist = arrayfun(@(x) calculateTotalDistance(x.Routes, Citys, start_point), Population); % 计算均衡惩罚系数（标准差越大惩罚越重） balancePenalty = std(totalDist) * 0.8; % 综合适应度 = 总距离 + 均衡惩罚 fitness = sum(totalDist) + balancePenalty; end这里有个小技巧——用标准差乘以0.8作为均衡惩罚系数。这个0.8是调参经验值，太大容易陷入局部最优，太小又起不到均衡作用。实际使用中可以根据数据集规模调整。
SD-MTSP：开普勒优化算法KOA求解单仓库多旅行商问题MATLAB（可更改数据集，旅行商的数量和起点）
行星位置更新规则是KOA的精髓所在，对应着开普勒定律中的椭圆轨道特性：
function Population = updatePopulation(Pop, fitness, iter, Max_iter) % 按适应度排序确定行星层级（类似引力权重） [~, sortIndex] = sort(fitness); for i = 1:size(Pop,1) % 计算当前轨道离心率（控制搜索范围） eccentricity = 0.5*(1 + cos(pi*iter/Max_iter)); % 按不同轨道层级的行星更新策略 if i <= size(Pop,1)/3 % 内层行星策略：精细搜索 step = eccentricity * randn(size(Pop(i).Routes)); else % 外层行星策略：大范围探索 step = (1 - eccentricity) * rand(size(Pop(i).Routes)); end % 位置更新并施加边界约束 newRoutes = Pop(i).Routes + step; Pop(i).Routes = mod(round(newRoutes), length(Citys)) + 1; end end这里有几个设计亮点：
- 离心率参数eccentricity随迭代次数动态变化，前期鼓励探索，后期偏向开发
- 内层行星（适应度好的解）采用高斯扰动进行精细搜索
- 外层行星使用均匀分布扰动扩大搜索范围
- mod运算确保路径编码始终在有效客户编号范围内
可视化部分最能体现算法效果。咱们用动态图展示迭代过程：
% 实时绘制路径变化 figure; for iter = 1:Max_iter routes = BestSol(iter).Routes; subplot(1,2,1); hold off; for k = 1:m plotRoute(routes{k}, Citys, start_point); hold on; end title(['Iteration: ',num2str(iter),' 总距离: ',num2str(BestSol(iter).totalDist)]); % 绘制收敛曲线 subplot(1,2,2); plot([BestSol(1:iter).totalDist], 'LineWidth',2); xlabel('迭代次数'); ylabel('最优总距离'); drawnow; end实际跑起来会发现，算法前20轮左右总距离下降明显，之后进入微调阶段。比如在50个客户点、3个旅行商的场景下，典型收敛曲线长这样：
!示例收敛曲线
几个实用调试技巧：
- 遇到早熟收敛时，适当增大行星数量N（别超过50，否则计算量爆炸）
- 调整fitnessFunction中的0.8系数可改变均衡权重
- 客户坐标建议做归一化处理，避免维度差异影响距离计算
修改起来也方便。比如想改成北京朝阳区真实坐标：
% 替换Citys变量即可 Citys = [116.4831, 39.9219; % 三里屯 116.4612, 39.9375; % 国贸 116.4993, 39.9107; % 南锣鼓巷 ... ]; % 其他坐标或者测试不同小哥数量对结果的影响：
for m = 2:5 % 重新运行算法... compareResult(m) = BestSol(end).totalDist; end这种基于天体物理启发的算法，在处理多目标优化问题时展现出了独特的优势。相比传统遗传算法，KOA在路径均衡性方面表现更稳定，可能是因为离心率机制有效平衡了探索与开发的关系。不过要注意，当客户点超过200个时，编码方式需要改用更高效的矩阵运算，否则迭代速度会明显下降。