旅行推销员问题:基于求解器

打开直播脚本

此示例显示如何使用二进制整数程序来解决经典旅行推销员问题。这个问题涉及通过一组停止（城市）找到最短的闭幕之旅（路径）。在这种情况下，有200个停止，但你可以轻松改变nStops变量得到不同的问题大小。您将解决初始问题，并看到解决方案具有子图。这意味着找到的最佳解决方案不会给出贯穿所有点的连续路径，而是有几个不相连的循环。然后，您将使用一个迭代过程来确定子团，添加约束，并重新运行优化，直到子团被消除。

对于基于问题的方法，请参阅旅行推销员问题:基于问题．

问题公式化

整数线性规划的旅行商问题如下:

生成所有可能的旅行，这意味着所有不同的站点都有。
计算每次旅行的距离。
要最小化的成本函数是旅游中每一趟的旅行距离之和。
决策变量是二进制的，并与每个行程相关联，其中每个1表示行程中存在的行程，每个0表示行程中不存在的行程。
为了确保行程包括每一站，要考虑到每一站正好是两次行程的线性约束。这意味着到站一次，离站一次。

生成停止

在美国大陆大陆的原油多边形表示内产生随机停止。

负载(“usborder.mat”，“x”，'是'，'xx'，'yy'）;rng (3“旋风”）％在缅因州和佛罗里达州的停止发出剧情，并且是可重复的nstops = 200;％您可以使用任何数字，但问题大小尺度为n ^ 2stopslon = zeros（nstops，1）;分配nStops的x坐标stopsLat = stopsLon;%分配坐标n = 1;尽管(n <= nStops) xp = rand*1.5;yp =兰德;如果inpolygon (xp, yp, x, y)%测试是否在边界内stopslon（n）= xp;stopslat（n）= yp;n = n + 1;结束结束

计算点之间的距离

因为有200个停止，有19,900个旅行，意思是19,900个二进制变量（＃变量= 200选择2）。

生成所有的行程，即所有对的站点。

idx = nchoosek (1: nStops, 2);

计算所有行程距离，假设地球是平坦的，以便使用毕达哥拉斯族规则。

dist = armot（stopslat（Idxs（：，1）） -  stopslat（Idxs（:,2）），．..stopslon（IDXS（：，1）） -  stopslon（IDXS（:,2）））;lendist =长度（dist）;

的这个定义dist向量，行程的长度是

dist'* x_tsp

在哪里X_TSP.是二进制解决方案矢量。这是游览您尝试最小化的距离。

创建图形和绘制地图

将问题表示为图形。创建一个图表，停止是节点，并且跳闸是边缘。

g =图表（IDXS（：，1），IDXS（：，2））;

使用曲线图显示停止。绘制没有图边的节点。

图= plot（g，“XData”，stopslon，“YData”stopsLat,'linestyle'，“没有”，'nodelabel'，{}）;持有上绘制外边框绘图（x，y，的r -） 抓住离开

约束

创建每个站点有两个相关行程的线性约束，因为每个站点必须有一个行程和从每个站点出发的行程。

Aeq = spalloc (nStops、长度(idx) nStops * (nStops-1));%分配稀疏矩阵为II = 1：nstops哪个idxs =（idxs == ii）;％找到包括停止II的旅行whichIdxs =稀疏(sum (whichIdxs, 2));％包括跳闸，其中II在任何一端Aeq (ii):) = whichIdxs ';％包括在约束矩阵中结束说真的= 2 * 1 (nStops 1);

二进制界限

所有决策变量都是二进制的。现在，设置intcon参数为决策变量的数量，将每个变量的下界设为0，上界设为1。

intcon = 1: lendist;1磅= 0 (lendist);乌兰巴托= 1 (lendist, 1);

使用Intlinprog优化

这个问题随时可以解决。要抑制迭代输出，请关闭默认显示。

选择= optimoptions ('intlinprog'，'展示'，'离开'）;[x_tsp, costopt exitflag、输出]= intlinprog(经销、intcon [], [], Aeq,说真的,磅,乌兰巴托,选择);

使用解决方案TRIPS创建一个新图形作为边缘。为此，请在某些值不完全整数的情况下绕过解决方案，并将结果值转换为逻辑．

x_tsp =逻辑(圆(x_tsp));Gsol =图(idx (x_tsp, 1), idx (x_tsp 2));

可视化解决方案

持有上突出(hGraph Gsol,'linestyle'，“- - -”） 标题（“解决方案与Subtours”）

可以在地图上看到，解决方案有几个子流。指定到目前为止指定的约束不会阻止这些子房子发生。为了防止任何可能的子学会发生，您需要一个令人难以置信的大量不等式约束。

Subtour约束

因为您无法添加所有子资金的约束，所以采取迭代方法。检测当前解决方案中的子流，然后添加不等式约束，以防止这些特定的子流程发生。通过这样做，您可以在几个迭代中找到适当的巡演。

消除带有不等式约束的子团。举个例子，如果你在一个subtour中有5个点，那么你有5条线连接这些点来创建这个subtour。通过实现一个不等式约束来消除这个subtour，即在这五个点之间必须有小于或等于四条线。

甚至更多，在这五点之间找到所有线条，并限制解决方案不超过四条线的存在。这是一个正确的约束，因为如果在解决方案中存在五个或更多行，则解决方案将有一个子流（图表 $n$ 节点和 $n$ 边总是包含一个循环)。

通过识别中的连接组件来检测子图GSOL.，使用当前解决方案中的边缘构建的图形。Conncomp.返回带有每个边缘所属的子流量的向量的向量．

tourIdxs = conncomp (Gsol);numtours = max (tourIdxs);%次参观的数量流('子房子数量：％d \ n', numtours);

子房子数量：27

包括线性不等式约束来消除子遍历，并反复调用求解器，直到只剩下一个子遍历。

lendist = spalloc (0, 0);%分配一个稀疏线性不等式约束矩阵b = [];尽管numtours > 1重复，直到只剩下一个小巡视%添加subtour约束b = [b; 0 (numtours, 1)];%分配bA = [; spalloc (numtours、lendist nStops)];％猜测有多少欧塞斯分配为II = 1：numtours rowidx = size（a，1）+ 1;索引的％计数器subTourIdx = find(tourIdxs == ii);%提取当前子巡回％下一行找到与之关联的所有变量％特定的子，然后添加不等式约束禁止百分比和使用这些停止的所有子流。变型= nchoosek（1：长度（亚尺寸），2）;为JJ = 1：长度（变型）哪个var =（sum（idxs == subtouridx（变体（Jj，1）），2）））和．..（SUM（IDXS == SubtouridX（变体（JJ，2）），2））;a（rowidx，whenvar）= 1;结束b(rowIdx) = length(subTourIdx) - 1;％少的行程比子房子停止结束％尝试再次优化[x_tsp, costopt exitflag、输出]= intlinprog(经销、intcon A、b Aeq,说真的,磅,乌兰巴托,选择);x_tsp =逻辑(圆(x_tsp));Gsol =图(idx (x_tsp, 1), idx (x_tsp 2));％可视化结果hGraph。线型=“没有”；%删除前面突出显示的路径突出(hGraph Gsol,'linestyle'，“- - -”）绘制这次有多少次短途旅行?tourIdxs = conncomp (Gsol);numtours = max (tourIdxs);%次参观的数量流('子房子数量：％d \ n'numtours)结束

子房子数量：20

子房子数量：7

子房子数量：9

子房子数量：9

副数量：3

副数量：2

子房子数量：7

副数量：2

子房子数量：1

标题（“消除子图的解决方案”）;持有离开

解决方案质量

因为它是一个单闭环，所以解决方案代表了一个可行的行程。但是这是一个最低费用的旅行吗?一种方法是检查输出结构。

disp (output.absolutegap)

绝对间隔的小意味着解决方案要么是最优的，要么总长度接近最优。