1、数学建模暑期培训,旅行商问题(TSP),主要内容,基本概念,算法简介,TSP模型的应用,最佳灾情巡视路线的模型的建立与求解,引 例,:,引 例,1. 98年全国大学生数学建模竞赛B题“最佳灾,今年(1998年)夏天某县遭受水灾. 为考察灾情、,组织自救,县领导决定,带领有关部门负责人到,全县各乡(镇)、村巡视. 巡视路线指从县政府,所在地出发,走遍各乡(镇)、村,又回到县政,府所在地的路线.,情巡视路线”中的前两个问题:,:,引 例,1)若分三组(路)巡视,试设计总路程最,短且各组尽可能均衡的巡视路线.,2)假定巡视人员在各乡(镇)停留时间T=2,小时,在各村停留时间t=1小时,汽车行驶速度V
2、,=35公里/小时. 要在24小时内完成巡视,至少应分,几组;给出这种分组下最佳的巡视路线.,引 例,公路边的数字为该路段的公里数.,引 例,2. 问题分析,本题给出了某县的公路网络图,要求在不同的条件下,,灾情巡视的最佳分组方案和路线.,将每个乡(镇)或村看作一个图的顶点,各乡镇、村之,间的公路看作此图对应顶点间的边,各条公路的长度(或,行驶时间)看作对应边上的权,所给公路网就转化为加权,图,问题就转化为图论中一类称之为旅行推销员问题,,即在给定的加权网络图中寻找从给定点O出发,行遍所有,顶点至少一次再回到点O,使得总权(路程或时间)最小.,旅行商问题(TSP,traveling sales
3、man problem)一个商人欲到n个城市推销商品,每两个城市i和j之间的距离为dij,如何选择一条道路使得商人每个城市正好走一遍后回到起点且所走路线最短。,:,原始问题,图论模型 构造一个图G=(V,E),顶点表示城市,边表示连接两城市的路,边上的权W(e)表示距离(或时间或费用)。于是旅行推销员问题就成为在加权图中寻找一条经过每个顶点正好一次的最短圈的问题,即求最佳Hamilton 圈的问题。,:,原始问题,基本概念,1) 哈米尔顿路径(H路径): 经过图G每个顶点正好一次的路径; 2) 哈米尔顿圈(H圈);经过G的每个顶点正好一次的圈; 3) 哈米尔顿图(H图): 含H圈的图。 4)
4、最佳H圈: 在加权图G=(V,E)中,权最小的H圈; 5) 最佳推销员回路: 经过每个 顶点一次的权最小闭通路; 6) TSP问题: 在完备加权图中求最佳H圈的问题。,:,TSP问题举例,工件排序设有n个工件等待在一台机床上加工,加工完i,接着加工j,这中间机器需要花费一定的准备时间tij,问如何安排加工顺序使总调整时间最短?此问题可用TSP的方法求解,n个工件对应n个顶点,tij表示(i,j) 上的权,因此需求图中权最小的H路径。 计算机布线 一个计算机接口含几个组件。每个组件上都置有若干管脚。这些管脚需用导线连接。考虑到以后改变方便和管脚的细小。要求每个管脚最多连两条线。为避免信号干扰以及
5、布线的简洁,要求导线总长度尽可能小。,TSP问题举例,计算机布线(续)问题容易转化为TSP问题,每个管脚对应于图的顶点,d(x,y)代表两管脚x与y的距离,原问题即为在图中寻求最小权H路径。 电路板钻孔 MetelcoSA是希腊的一个印刷电路板(PCCB)制造商。在板子上对应管脚的地方必须钻孔,以便以后电子元件焊在这板上。典型的电路板可能有500个管脚位置,大多数钻孔都由程序化的钻孔机完成,求最佳钻孔顺序。此问题其实就是求500个顶点的完备加权图的最佳H圈的问题,即TSP问题。用求解出的H圈来指导生产,使Metclo的钻孔时间缩短了30%,提高了生产效率。,算法简介,TSP问题是NP-hard
6、问题,即不存在多项式时间算法.,也就是说,对于大型网络(赋权图),目前还没有一个精确求解,TSP问题的有效算法,因此只能找能求出相当好(不一定最优)的解的算法.,算法简介,近似算法或启发式算法 一般是以构造型算法得到一个初始解,然后再用改进型算法逐步改进。 常见的构造型算法有两种:Christofides最小权匹配算法 ,对角线完全算法。 常见的改进型算法有两种:二次逐次修正法,Feiring矩阵逐次改进法。,分枝定界法,算法简介,二次逐次修正法,(1)任取初始H圈C0=v1,v2,vi,vj,vm,v1 (2)对所有的i,j,1i+1jm,若w(vi,vj)+w(vi+1,vj+1) w(v
7、i,vi+1)+w(vj,vj+1) 则在C0中删去边(vi,vi+1)和(vj,vj+1)而加入边(vi,vj)和(vi+1,vj+1),形成新的H圈C,即C=v1,v2,vi,vj,vj-1,vi+1,vj+1,vi,v1 (3)C0C,重复步骤(2),直到条件不满足为止,最后得到的C即为所求。,例对下图的K6,用二边逐次修正法求较优H圈.,较优H圈:,其权为W(C3)=192,分析: 这个解的近似程度可用最优H圈的权的下界与,其比较而得出.即利用最小生成树可得最优H圈的一个下界.,设C是G的一个最优H圈,则对G的任一顶点v, C-v是,G-v的生成树.如果T是G-v的最小生成树,且e1是
8、e2与v关联,的边中权最小的两条边,则w(T)+w(e1)+w(e2)将是w(C)的一个下界.,取v=v3,得G-v3的一最小生成树(实线),其权,w(T)=122,与v3关联的权最,小的两条边为v1v3和v2v3,w(T)+w(v1v3)+w(v2v3),=178. 故最优H圈的权,故w(C),应满足178 w(C) 192.,对角线完全算法,结论: 若能在nn距离矩阵中找出n个不同行也不同列的元素,使它们的和为最小值。若这n个元素构成一条哈米尔顿圈时,此圈便是最佳H圈。,矩阵的简化 :将矩阵的每一行的各元素减去本行的最小元素称为对行简化,从第i行减去的数称为第i行的约数,记为R(i)。将矩
9、阵的每一列的各元素减去本列的最小元素称为对列简化,从第j列减去的数称为第j列的约数,记为R(j)。各行各列的约数之和称为矩阵的约数,记为R。矩阵经行的简化和列的简化后所得矩阵称为该矩阵的简化矩阵。,对角线完全算法,例 求下列距离矩阵D的简化矩阵及各行,各列的约数。其中,D中各行的约数 R(1)=25,R(2)=5,R(3)=1,R(4)=6,R(5)=7,对角线完全算法,D经行简化得,求出D各列的约数 R(1)=0,R(2)=0,R(3)=0,R(4)=3,R(5)=0,对角线完全算法,D经列简化得,由于简化矩阵同一行各元素减同一数,同列各元素也是减同一数,因而每个H圈的权都减少同一数即R.,
10、对角线完全算法,定理 设D是图G的距离矩阵D的简化矩阵,则D对应的图G的最佳(有向)H圈也是原图G的最佳(有向)H圈。G只是边权与G不同,去掉权之后完全一样。因此当简化矩阵中的零元素构成H圈时,该H圈也是原问题的最佳H圈。,罚数: 在图G的距离矩阵的简化矩阵D中,第i行的最小元素与次小元素之差称为第i行的罚数,记为P(i)。第j列的最小元素与次小元素之差称为第j列的罚数,记为P(j),某行(或列)的罚数即是若H圈不选择该行(或列)的最小元素会使其权增加的最小值。,对角线完全算法,算法步骤 输入:图的距离矩阵D (1)求D的简化矩阵D以及各行各列的约数R(i),R(j),罚数P(i),P(j)
11、(2)计算在简化矩阵中零元素所在行与列的罚数和,即P(i,j)=P(i)+P(j)。将P(i,j)由大到小排列后,依次选取可作为可行部分路的边(i,j)。这些边对应的零元素记为0*。用这些选择出来的边构成可行部分路。,定义 一个H圈的一些不相交路径称为可行部分路。,对角线完全算法,(3)构造新的距离矩阵称为重构距离矩阵:按上述可行部分路的顶点序重新排列简化距离D的行,列也按使上述所有“0*”位于对角线上的次序重新排列。 (4)产生D的子阵:设重构矩阵对角线上m个非零元素对应的边为(i1,j1),(i2,j2),(im,jm),则从D中取出相应的m行,m列构成一个mm子阵D1。为保证选出的边与原
12、来的可行部分路不形成子循环,有m条边不能选择,将其对应的元素置为。并将列作适当调整使对角线元素为。 (5)对D1重复(1)(4)步,直到重构矩阵对角线上的元素全为0为止,这时便可得到一个H圈。,对角线完全算法,例 用对角线完全法求加权图K10的较佳H圈,对角线完全算法,对角线完全算法,2 求出以上第一级简化阵中零元素对应的罚数, 如P(1,2)=P(1,2)=P(1)+P(2)=116+52=168, 并将这些罚数按由大到小的次序排列如下: P(1,2)=168, P(2,1)=168, P(8,6)=124 P(6,8)=103, P(4,5)=67, P(7,3)=52 P(10,9)=4
13、5, P(9,10)=34, P(5,4)=30 P(2,7)=6, P(3,4)=6, P(2,3)=0 P(6,4)=0, P(9,5)=0,对角线完全算法,现依次从上述各边选择能形成可行部分路的边,先选第一边(1,2), 之后(2,1),(2,1)不行,因(1,2),(2,1)形成子循环;接着(8,6),但(6,8)不行,(6,8)会与(8,6)构成子循环;再下是(4,5),(7,3),(10,9),但(9,10)不能入选;(5,4)也不能入选,因会和前面选中的(4,5)构成子循环;接着是(2,7),(3,4),但(2,3)不能入选,否则与(2,7)会使2的出次大于1。同理(6,4),(
14、9,5)也不能入选。综上得到可形成可行部分路的边为:(1,2), (8,6), (4,5), (7,3),(10,9),(2,7)和(3,4)。它们对应的零元素为0*。 可行部分路:1-2-7-3-4-5,8-6和10-9,三条不相交路径。,对角线完全算法,3. 以1,2,7,3,4,5,8,6,10,9的顺序重新排列D的行,为使所有0*位于对角线上,D的列按2,7,3,4,5,8,6,10,9,1的顺序排列,这样便得到第一级重构距离矩阵,对角线完全算法,4. 对角线上有三个非零元素281,477和644,其对应的边为: (5,8),(6,10)和(8,1),相应的行数为5,6,9; 列数为8
15、,10,1。从原始权矩阵D中取出这三行,三列构成 一个33型的D的子阵:,对角线完全算法,5. 重复步骤(1)-(4),得到第二级简化矩阵及相应的约数,罚数,计算各零元素的罚数并按由大到小排列如下: P(6,1)=243 P(9,8)=243 P(5,8)=54 P(6,10)=54 依次选择可行边:(6,1)和(9,8)。它们对应的零元素记为0*。 与原来已经选出的边一起形成可行部分路如下:1-2-7-3-4-5;8-6-1;10-9-8, 或更清楚地应为:10-9-8-6-1-2-7-3-4-5。,对角线完全算法,上述顶点序得到第二级重构距离矩阵,最后还剩一个元素(5,10)不为零,没有选
16、择余地,第三迭代必定是选(5,10)与前面得到的可行部分将一起构成H圈10-9-8-6-1-2-7-3-4-5-10.,对角线完全算法,因此第三级重构距离矩阵只需在第二级距离矩阵中335换成0*即可。第三级重构距离矩阵为:,故所求H圈为:10-9-8-6-1-2-7-3-4-5-10,其权:W=3022。,旅行商问题的数学规划模型,旅行商问题的数学规划模型,或,旅行商问题的数学规划模型,例 用LINGO软件求解,MODEL:1sets:2 cities/110/:level; !level(i)= the level of city;3 link(cities, cities):4 dista
17、nce, !The distance matrix;5 x; ! x(i,j)=1 if we use link i,j;6endsets7data: !Distance matrix, it need not be symmetirc;8 distance = 0 8 5 9 12 14 12 16 17 22,旅行商问题的数学规划模型,9 8 0 9 15 16 8 11 18 14 2210 5 9 0 7 9 11 7 12 12 1711 9 15 7 0 3 17 10 7 15 1512 12 16 9 3 0 8 10 6 15 1513 14 8 11 17 8 0 9 14
18、 8 1614 12 11 7 10 10 9 0 8 6 1115 16 18 12 7 6 14 8 0 11 1116 17 14 12 15 15 8 6 11 0 1017 22 22 17 15 15 16 11 11 10 0;18enddata19n=size(cities); !The model size;20! Minimize total distance of the links;,旅行商问题的数学规划模型,21min=sum(link(i,j)|i #ne# j: distance(i,j)*x(i,j);22!For city i;23for(cities(i)
19、:24! It must be entered;25 sum(cities(j)| j #ne# i: x(j,i)=1;26! It must be departed;27 sum(cities(j)| j #ne# i: x(i,j)=1;28! level(j)=levle(i)+1, if we link j and i;29 for(cities(j)| j #gt# 1 #and# j #ne# i :30 level(j) = level(i) + x(i,j)31 - (n-2)*(1-x(i,j) + (n-3)*x(j,i);32 );33);,旅行商问题的数学规划模型,3
20、4! Make the xs 0/1;35for(link : bin(x);36! For the first and last stop;37for(cities(i) | i #gt# 1 :38 level(i)=1+(n-2)*x(i,1);40);END,水平变量(level)仍然是用来保证所选的边除第1 点外不构成圈的.计算结果如下:,旅行商问题的数学规划模型,Global optimal solution found at iteration: 90 Objective value: 73.00000Variable Value Reduced CostX( 1, 2) 1.0
21、00000 8.000000X( 2, 6) 1.000000 8.000000X( 3, 1) 1.000000 5.000000X( 4, 3) 1.000000 7.000000X( 5, 4) 1.000000 3.000000X( 6, 9) 1.000000 8.000000X( 7, 10) 1.000000 11.00000X( 8, 5) 1.000000 6.000000X( 9, 7) 1.000000 6.000000X( 10, 8) 1.000000 11.00000,旅行商问题的数学规划模型,旅行商经过10 个城镇的最短距离为73公里,其连接情况如图所示.,最佳灾
22、情巡视路线的模型的建立与求解,问题转化为在,给定的加权网,络图中寻找从,给定点O出发,行遍所有顶点,至少一次再回,回到点O ,使得,总权(路程或时,时间)最小,即,最佳旅行推销,员问题.,最佳旅行推销员问题是NP完全问题,采用一种,近似算法求其一个近似最优解,来代替最优解.,算法一 求加权图的最佳旅行推销员回路近似算法:,1) 用图论软件包求出G中任意两个顶点间的最短路,构造出完全图,2) 输入图 的一个初始H圈;,3) 用对角线完全算法产生一个初始圈;,4) 随机搜索出 中若干个H圈,例如2000个;,5) 对第2),3),4)步所得的每个H圈,用二边逐次,修正法进行优化,得到近似最优H圈;
23、,6) 在第5)步求出的所有H圈中,找出权最小的一个,此即要找的最优H圈的近似解.,因二边逐次修正法的结果与初始圈有关,故本算法,第2),3),4)步分别用三种方法产生初始圈,以保,证能得到较优的计算结果.,问题一 若分为三组巡视,设计总路程最短且各,组尽可能均衡的巡视路线.,此问题是多个推销员的最佳旅行推销员问题.,4),此问题包含两方面:a)对顶点分组, b)在每组中,求(单个推销员)最佳旅行推销员回路.,因单个推销员的最佳旅行推销员回路问题不存,存在多项式时间内的精确算法.,而图中节点数较多,为53个,我们只能去寻求,一种较合理的划分准则,对图1进行初步划分后,求,出各部分的近似最佳旅行
24、推销员回路的权,再进一,步进行调整,使得各部分满足均衡性条件3).,从O点出发去其它点,要使路程较小应尽量走,O点到该点的最短路.,故用软件包求出O点到其余顶点的最短路.,这些最短路构成一棵O为树根的树.,将从O点出发的树枝称为干枝.,从O点出发到其它点共有6条干枝,它们的名称,分别为,.,在分组时应遵从准则:,准则1 尽量使同一干枝上及其分枝上的点分在同一组.,准则2 应将相邻的干枝上的点分在同一组;,准则3 尽量将长的干枝与短的干枝分在同一组.,由上述分组准则,我们找到两种分组形式如下:,分组1:(,),(,),(,),分组2:(,),(,),(,),分组1极不均衡,故考虑分组2.,分组2
25、:(,),(,),(,),对分组2中每组顶点的生成子图,用算法一求出近似最优解及相应的巡视路线.,在每个子图所构造的完全图中,取一个尽量包含上图中树上的边的H圈作为其第2)步输入的初始圈.,分组2的近似解,因为该分组的均衡度,. 所以此分法的均衡性很差.,为改善均衡性,将第组中的顶点C,2,3,D,4,分给第组(顶点2为这两组的公共点),重新分,分组后的近似最优解见表2.,因该分组的均衡度,所以这种分法的均衡性较好.,问题二 当巡视人员在各乡(镇)、村的停留,停留时间一定,汽车的行驶速度一定,要在24小时内,完成巡视,至少要分几组及最佳旅行的巡视路线.,由于T=2小时,t=1小时,V=35公里
26、/小时,需访问,的乡镇共有17个,村共有35个. 计算出在乡(镇)及,村的总停留时间为17 2+35=69小时,要在24小时,内完成巡回,若不考虑行走时间,有: 69/i24(i为,分的组数). 得i最小为4,故至少要分4组.,由于该网络的乡(镇)、村分布较为均匀,故有可,能找出停留时间尽量均衡的分组,当分4组时各组停,停留时间大约为69/4=17.25小时,则每组分配在路,路途上的时间大约为24-17.25=6.75小时.而前面讨,论过,分三组时有个总路程599.8公里的巡视路线,,分4组时的总路程不会比599.8公里大太多,不妨以,599.8公里来计算.路上约用599.8/35=17小时,
27、若平,均分配给4个组,每个组约需17/4=4.25小时6.75小,小时,故分成4组是可能办到的.,现在尝试将顶点分为4组.分组的原则:除遵从,前面准则1、2、3外,还应遵从以下准则:,准则4 尽量使各组的停留时间相等.,用上述原则在下图上将图分为4组,同时计算,各组的停留时间,然后用算法一算出各组的近似最,最佳旅行推销员巡回,得出路线长度及行走时间,从而得出完成巡视的近似最佳时间. 用算法一计,计算时,初始圈的输入与分三组时同样处理.,这4组的近似最优解见表3.,表3符号说明:加有底纹的表示前面经过并停留过,此次只经过不停留;加框的表示此点只经过不停留.,可看出,表3分组的均衡度很好,且完全满足24,小时完成巡视的要求.,该分组实际均衡度 4.62%,