这个例子展示了如何使用二进制整数规划来解决经典的旅行推销员问题。这个问题涉及到通过一组站点(城市)找到最短的封闭旅行(路径)。在这种情况下,有200个停止,但您可以轻松地更改nStops
变量来获得不同的问题大小。您将解决最初的问题,并看到解决方案有子图。这意味着找到的最优解不是给出一条贯穿所有点的连续路径,而是有几个不相连的循环。然后,您将使用一个迭代过程来确定子导览、添加约束并重新运行优化,直到消除子导览。
有关此问题的基于求解器的方法,请参见旅行推销员问题:基于求解器.
对整数线性规划的旅行商问题,公式如下:
生成所有可能的行程,即所有不同的站点对。
计算每次行程的距离。
要最小化的代价函数是行程中每次行程的行程距离之和。
决策变量是二进制的,并且与每次行程相关,其中每个1表示旅行中存在的行程,每个0表示旅行中不存在的行程。
为了确保旅行包含每个站点,需要包含每个站点恰好位于两次行程中的线性约束。这意味着一站到站和一站出站。
在美国大陆的粗略多边形表示中生成随机止损
负载(“usborder.mat”,“x”,“y”,“xx”,“yy”);rng (3“旋风”)在缅因州和佛罗里达州停留,并且是可重复的nStops = 200;您可以使用任何数字,但问题规模的规模为N^2stopsLon = 0 (nStops,1);分配nstop的x坐标stopsLat = stopsLon;分配y坐标N = 1;而(n <= nStops) xp = rand*1.5;Yp =兰特;如果inpolygon (xp, yp, x, y)% tTest如果在边界内stopsLon(n) = xp;stopsLat(n) = yp;N = N +1;结束结束
因为有200个站点,所以有19,900次行程,这意味着19,900个二进制变量(# variables = 200选2)。
生成所有的行程,即所有对站点。
idxs = nchoosek(1:nStops,2);
计算所有的旅行距离,假设地球是平的,以便使用毕达哥拉斯法则。
dist = hypot(stopsLat(idxs(:,1)) - stopsLat(idxs(:,2)),...stopsLon(idxs(:,1)) - stopsLon(idxs(:,2)));Lendist =长度(dist);
用这个定义经销
矢量,旅行的长度是
dist”*旅行
在哪里旅行
表示解决方案所经过的行程的二进制向量。这是你要尽量减少的旅行距离。
用图表表示问题。创建一个图,其中站点是节点,行程是边。
G = graph(idxs(:,1),idxs(:,2));
使用图表显示止损点。绘制没有图边的节点。
figure hGraph = plot(G,“XData”stopsLon,“YData”stopsLat,“线型”,“没有”,“NodeLabel”, {});持有在绘制外部边界情节(x, y,的r -)举行从
创建一个用二进制优化变量表示潜在行程的优化问题。
TSP =优化问题;行程= optimvar(“旅行”lendist 1“类型”,“整数”,下界的0,“UpperBound”1);
在问题中包含目标函数。
tsp。目标= dist'*行程;
创建线性约束,即每个站点有两个相关联的行程,因为每个站点必须有一个行程,并且每个站点必须有一个行程。
使用图形表示法通过查找连接到某个站点的所有边来识别从该站点开始或结束的所有行程。对于每个站点,创建该站点的行程总和等于2的约束。
constr2trips = optimconstr(nStops,1);为stop = 1:nStops which hidxs = outedges(G,stop);识别与该站相关的行程constr2trips(stop) = sum(trips(which hidxs)) == 2;结束tsp.Constraints。Constr2trips = Constr2trips;
问题已经准备好解决了。要抑制迭代输出,请关闭默认显示。
Opts = optimoptions(“intlinprog”,“显示”,“关闭”);Tspsol =解(tsp,“选项”选择)
tspsol =带字段的结构:旅游:[19900×1 double]
创建一个以解行程为边的新图。为此,请四舍五入解决方案(以防某些值不完全是整数),并将结果值转换为逻辑
.
tspsol。旅行= logical(round(tspsol.trips)); Gsol = graph(idxs(tspsol.trips,1),idxs(tspsol.trips,2),[],numnodes(G));% Gsol = graph(idxs(tspsol.trips,1),idxs(tspsol.trips,2));%在大多数情况下也适用
将新图形覆盖在现有图形上,并突出显示其边缘。
持有在突出(hGraph Gsol,“线型”,“- - -”)标题(“分团解决方案”)
从地图上可以看到,解决方案有几个子行程。到目前为止指定的约束并不能阻止这些子巡回的发生。为了防止任何可能发生的子循环,你将需要大量的不等式约束。
因为您不能添加所有的子巡演约束,所以请采用迭代方法。检测当前解决方案中的子巡视,然后添加不等式约束以防止那些特定的子巡视发生。通过这样做,您可以在几次迭代中找到合适的旅行。
消除带有不等式约束的子游历。这是如何工作的一个例子,如果你在一个子漫游中有五个点,那么你有五条线连接这些点来创建子漫游。通过实现一个不等式约束来消除这个子行程,即这五个点之间必须小于或等于四条线。
更重要的是,找到这五个点之间的所有直线,并约束解不超过四条这些直线。这是一个正确的约束,因为如果一个解决方案中存在五条或五条以上的线,那么该解决方案将有一个子巡回(与 节点和 Edges总是包含一个循环)。
中的连接组件来检测子游历Gsol
,即用当前解中的边构建的图。conncomp
返回一个向量,其中包含每条边所属的子游历的编号.
tourIdxs = conncomp(Gsol);numtours = max(tourIdxs);%子行程数流('子旅程#:%d\n', numtours);
子行程数:27
包含线性不等式约束以消除子巡回,并反复调用求解器,直到只剩下一个子巡回。
子行程添加约束的索引K = 1;而Numtours > 1%重复此操作,直到只有一个子巡回添加子行程约束为ii = 1:numtours inSubTour = (tourIdxs == ii);%当前子行程中的边a = all(inSubTour(idxs),2);在子巡回中完成图形索引constrname =“subtourconstr”+ num2str (k);tsp.Constraints.(constrname) = sum(trips(a)) <= (nnz(inSubTour) - 1);K = K + 1;结束再次尝试优化[tspsol,fval,exitflag,output] = solve(tsp,“选项”、选择);tspsol。旅行= logical(round(tspsol.trips)); Gsol = graph(idxs(tspsol.trips,1),idxs(tspsol.trips,2),[],numnodes(G));% Gsol = graph(idxs(tspsol.trips,1),idxs(tspsol.trips,2));%在大多数情况下也适用%绘制新方案hGraph。线型=“没有”;移除之前高亮显示的路径突出(hGraph Gsol,“线型”,“- - -”) drawnow这次有多少次分旅行团?tourIdxs = conncomp(Gsol);numtours = max(tourIdxs);%子行程数流('子旅程#:%d\n'numtours)结束
子行程数:20
子行程:7
子行程:9
子行程:9
子行程:3
子行程:2
子行程:7
子行程:2
子行程#:1
标题(“取消分团的解决方案”);持有从
该解表示可行行程,因为它是一个单闭环。但这是一个费用低廉的旅行吗?一种方法是检查输出结构。
disp (output.absolutegap)
0
绝对间隙的小意味着解决方案要么是最优的,要么有一个接近最优的总长度。