灾情巡视路线
灾情巡视路线
摘要
本文讨论了由某县的赋权交通网络图确定分组进行灾情巡视的最佳路线规划问题。
对于问题一,对所给的赋权网络图要求分三组(路)巡视,得到总路程最短且各组尽可能均衡的巡视路线,可转化为三个售货员的最佳旅行售货员问题。本文先用MATLAB 软件编程计算得到加权网络图的最小生成树,按每块近似有相等总路程的标准将最小生成树分成三块,每一块都转化为一个最佳旅行售货员问题。再确定总路程最短且满足各组尽可能均衡的路线的目标函数,然后对目标函数进行调整改进,比较得到最终的双目标最优化模型,最后确定的最佳行驶路线
目要求。因此首先考虑分四组的情况。在分三组的基础上根据分组原则将图大致划分为4个子图。同样以巡视路程最短和时间均衡度最小为目标函数,各巡视时间小于24小时作为约束条件建立多目标规划模型。然后采用类似于问题一的分
间为最短时间,根据最短路径树,从远到近,依次巡视各村镇,在所用时间不大于最短时间的前提下,各组尽可能多的巡视几个村镇,进行分组确立巡视点,并对已巡视过的点进行逐个剔除。通过分组优化对比等,得出最佳分组情况及巡视
响的只是整个巡视时间。我们利用MATLAB 编程画图得到T 、t 、V 与巡视时间的关系曲线。观察曲线发现:(1)当停留时间t 不变时,随着行驶速度的增加,巡视所用的时间越来越短。在保证安全的行驶过程中,适当加快行驶速度,可以减少花费在路途中的时间。(2)停留时间t 与巡视时间呈线性关系,在保持V 不变时,巡视最短的时间T min 与停留时间t 成正比关系,即无论停留时间t 取何值,对T min 的影响都较大。
关键词:最小生成树 多目标规划 最佳路线
一、问题重述
如图为某县的乡(镇)、村公路网示意图。该县遭受水灾,为考察灾情、组织自救,县领导决定带领有关部门负责人到全县各乡(镇)、村巡视。巡视路线指从县政府所在地出发,走遍各乡(镇)、村,又回到县政府所在地的路线。现有如下问题:
县公路网络示意图
1. 若分三组(路)巡视,试设计总路程最短且各组尽可能均衡的巡视路线。2. 假定巡视人员在各乡(镇)停留时间 T=2 小时,在各村停留时间 t=1小时,汽车行驶速度 V=35公里/小时。要在24小时内完成巡视,至少应分几组;给出这种分组下你认为最佳的巡视路线。
3. 在上述关于 T,t和 V的假定下,如果巡视人员足够多,完成巡视的最短时间是多少;给出在这种最短时间完成巡视的要求下,你认为最佳的巡视路线。4. 若巡视组数已定(比如三组),要求尽快完成巡视,讨论T ,t 和 V改变对最佳巡视路线的影响。
二、模型假设
1、若汽车只是路过乡村,则不考虑汽车穿越乡村的时间。
2、假设汽车行驶速度稳定不变。 3、假设无其它意外情况
三、符号说明
四、问题的分析
本文是领导巡视受灾县, 并求最佳巡视路线的数学建模问题,题中已给出该县公路的网络图,要求在不同的题目要求下,得到灾情巡视的最佳分组方案和路线。
若将每个乡(镇)或村看作一个图的顶点,各乡镇、村之间的公路看作此图对应顶点间的边,各条公路的长度(或行驶时间)看作对应边上的权,所给公路网就转化为加权网络图,问题就转化图论中旅行售货员的一类问题,即在给定的加权网络图中寻找从给定点O 出发,行遍所有顶点至少一次再回到点O ,使得总权(路程或时间)最小. 本题所求的分组巡视的最佳路线,也就是m 条经过同一点并覆盖所有其他顶点又可使边上的权之和达到最小的闭路,即最佳旅行售货员问题。
对于问题一, 要求分三组(路)巡视,得到总路程最短且各组尽可能均衡的巡视路线,可转化为三个售货员的最佳旅行售货员问题。先用MATLAB 软件编程计算得到加权网络图的最小生成树,按每块近似有相等总路程的标准将最小生成树分成三块,每一块都转化为一个最佳旅行售货员问题。再确定总路程最短且满足各组尽可能均衡的路线的目标函数,最后对目标函数进行调整改进,比较得到最终的双目标最优化模型。
对于问题二,由于T =2小时,t =1小时,V =35公里/小时, 需访问的乡镇共有17个,村共有35个。计算出在乡(镇) 及村的总停留时间为17⨯2+35=69小时,要在24小时内完成巡回,若不考虑行走时间,有:69/i
对于问题三,考虑在人员足够多的情况下,求出最短的巡视时间。假设一个小组只巡视一个点的情况下,则去巡视离O 点最远的点所花时间最长。我们以巡视小组中所耗时间最长的小组所用时间作为这次整个巡视的最短时间。要使这次巡视时间最短,则要求去巡视离点O 最远的点所花时间最小,由题目所给图可知,离O 点最远的点为H ,所以就以巡视H 所花时间作为T min 。当此小组只巡视时,T min 最小。在不超过T min 的情况下,根据其他小组的剩余时间确定沿途是否巡视其他点。其中巡视原则为:
(1)当一组人员巡视完规定点后时,在剩余时间允许的情况下,优先考虑原巡视点附近而距离O 较远的点;
(2)最大限度使用剩余时间,主要考虑原则(1)。按照此原则,逐个巡视,直至巡视完所有点。
对于问题四,若巡视组数已定,则每个小组的最短路径就已确定,T 、t 、V 改变只影响的是整个的巡视时间。要求尽快完成巡视,即巡视时间要尽可能小。巡视时间包括到巡视点的行驶时间和在巡视点的停留时间。行驶时间主要取决于速度V ,而停留时间由T 、t 决定。所以此问题讨论的是当T 、t 、V 改变时对巡视时间的影响,即对T 、t 、V 的影响分析。
五、模型的建立与求解
由题目可知,共有53个乡(镇)和村,即在赋权网络图G (V ,E ,W )中共有53个点. 其中V ={v 1, v 2, v n }(n =53) 表示图G 的53个节点,E ={e ij }表示相关联的两点i , j 的边集,W ={w ij }表示相关联两点i , j 间的权值。定义决策变量r ij :
⎧1
r ij =⎨
⎩0
表示i , j 两点相关联
表示i , j 两点不相关
其中相关联表示i , j 两点之间有权值,不相关表示i , j 之间没有权值。将这些
点和权值生成图G 的矩阵,对于不相关的两点权值作为无穷大处理。用MATLAB 编写程序,得出该网络图G 的最小生成树。画出的图形如下图所示。
最小生成树图
5.1问题一模型的建立与求解 5.1.1模型的建立
通过问题一的分析,建立多目标规划模型。 (1)三组巡视的总路线最短:
min
(2)巡视路线尽量均衡:
∑c
i =1
3
i
max(c i ) -min(c i )
min α=
max(c i )
上式两个目标函数相互矛盾,不可能同时达到最好。而α越小,表示均衡度越好。设当α
通过增环、扩环、换枝等方法,对子图内部进行适当调整优化,找出多种方案进行对比,寻找出最佳巡视回路,得到结果如下表:
203.5-200.4
⨯100%=1.52%.
203.5
本题的分组路程均衡度为1.52%,说明本题分组路程均衡性很好,满足题目要求的均衡路线,总路程为607.4公里。各小组的巡视路线如下图. 其中绿色表示第一组,红色表示第二组,蓝色表示第三组。
根据上表数据,得到分组路程均衡度为α=
各小组巡视路线图
5.2问题二模型的建立与求解 5.2.1模型的建立
由前面对问题二的分析可知,要满足条件,最少要分成4组进行灾情巡视。先按4组建立模型并求解,若分成4组仍无法满足,则继续考虑增加巡视组数。要求巡视时间与总路程长度均达到最小. 得出目标函数:
4
⎧
⎪min ∑c i
⎨i =1
⎪min α⎩
t i =mT +nt +外加巡视时间不超过24小时的约束条件:
c i
5.2.2模型的求解
依据该县地图的最小生成树图首先将图形大致的分为四组, 然后采用类似于问题一的分组规则进行优化调整. 显然,该问要求在时间尽量短的条件下,使总的行驶路程最短。同样可以得到多个满足条件的解, 现在仅将多组巡视方案中时
定义最小生成树的分支上未分组的点中到O 最远的点为V N ,次最远点为
V N -1,依此类推,直至离O 最近的点V 1。O 点到点N 的时间记为T n ,停留时间记为t N 。将V N 看成第一组,再根据约束条件判断V 1, V 2,... V N 分组情况,判别方法如下:
(1)若T N +t N ≤T N ≤T N +t N +t N -1,则点V N 应单独分在一组;若
T min >T N +t N ≥T N +t N +t N -1,则V N 和V N -1分为一组。同理,可依此类推判断V N -2, V N -3,... 直至待判点不能和前面已判点分在一组。
(2)再在分支上未分组的点中找到离O 最远的点,作为下一组,用上述同样的方法判断,直至所有节点分组完。按以上求解过程,逐步求解可得以下共分
5.4问题四模型的建立与求解 完成所有巡视所用最短时间:
T min
S min =+nT +mt V
以问题一中的第一组为例,其中S 1=200. 4km , m =10, n =7。为了方便讨论,将T ,t 统一作为时间变量,设T =bt (b =2) ,则上式可以改写为:
T =
S min
+(m +bn ) t V
(1)当停留时间t 不变时,利用matlab 作图得到巡视所用最短时间T min 与行驶速度V 的关系图如下:
由T min -V 图可以知道,随着行驶速度的增加,巡视所用的最短时间越来越短。在保证安全的行驶过程中,适当加快行驶速度,可以减少花费在路途中的时间。而速度过小时,对T min 的影响较大,这个时候可以考虑重新分组。 (2)当保持V =35km /h 不变时,利用matlab 作图得到的巡视所用最短时间T min 与停留时间t 的关系图如下:
T min -t 图
由T min -t 图可以知道,在保持V 不变时,巡视最短的时间T min 与停留时间t 成正比关系,即无论停留时间t 取何值,对T min 的影响都较大,当停留时间小范围波
动时,也需重新考虑分组。
六、模型的评价与推广
6.1优点:
1)本模型运用最小生成树的方法,将区域大致分为三个部分,使求解简化; 2)采用均衡度的表达式,使模型检验合理方便; 3)模型运用网络图的表示方法给出了灾情最佳巡视路线,不仅缩短了时间,还考虑了各个组巡视的均衡度,对考察灾情、组织自救起到了很大的作用。 6.2缺点:
1)运用MATLAB 软件求得最小生成树图后,寻找各个组的巡视路线得最优解工作量较大。
2)问题四中,在对T 、t 、V 的灵敏度分析中,只对其变化作了定性分析,而未作定量分析。
七、参考文献
[1]姜启源,数学建模(第三版)[M],北京:高等教育出版社,2003 [2]韩中庚,数学建模方法及其应用[M],北京:高等教育出版社,2005 [3]王庆波,最佳灾情巡视路线[J],哈尔滨工业大学学报,第31卷第2期,1999年4月
八、附录
程序一
最小生成树的程序: clc,clear
a=zeros(53);
a(50,1 )=6.0;a(50,53)=12.9;a(50,38)=11.5;a(50,2)=9.2;a(50,48)=19.8;a(50,51)=10.1; a(1,36)=10.3; a(1,37)=5.9;a(1,38)=11.2; a(2,3)=4.8;a(2,5)=8.3;
a(3,38)=7.9;a(3,39)=8.2;a(4,39)=12.7;a(4,8)=20.4; a(5,48)=11.4;a(5,39)=11.3;a(5,6)=9.7; a(6,48)=9.5;a(6,7)=7.3;a(6,47)=11.8;
a(7,39)=15.1;a(7,40)=7.2;a(7,47)=14.5;a(8,40)=8.0; a(9,40)=7.8;a(9,41)=5.6;a(10,41)=10.8;
a(11,45)=13.2;a(11,40)=14.2;a(11,42)=6.8; a(12,42)=7.8;a(12,41)=12.2;a(12,43)=10.2;
a(13,44)=16.4;a(13,45)=9.8;a(13,42)=8.6;a(13,14)=8.6;
a(14,15)=15;a(14,43)=9.9;a(15,44)=8.8;
a(16,17)=6.8;a(16,44)=11.8;a(17,22)=6.7;a(17,46)=9.8; a(18,46)=9.2;a(18,45)=8.2;a(18,44)=8.2; a(19,20)=93;a(19,47)=7.2;a(19,45)=8.1; a(20,21)=7.9;a(20,25)=6.5;a(20,47)=5.5; a(21,23)=9.1;a(21,25)=6.5;a(21,46)=4.1;
a(22,23)=10.0;a(22,46)=10.1;a(23,24)=8.9;a(23,49)=7.9; a(24,27)=18.8;a(24,49)=13.2;a(25,49)=8.8;a(25,48)=12.0; a(26,27)=7.8;a(26,51)=10.5;a(26,49)=10.5;a(27,28)=7.9; a(28,52)=8.3;a(28,51)=12.1;
a(29,52)=7.2;a(29,53)=7.9;a(29,51)=15.2; a(30,32)=10.3;a(30,52)=7.7;
a(31,32)=8.1;a(31,33)=7.3; a(31,53)=9.2; a(32,33)=19;a(32,35)=14.9;a(33,36)=7.4; a(34,35)=8.2;a(34,36)=11.5;a(34,13)=17.6;
a(37,38)=12.2;a(36,53)=8.8;a(37,38)=11.0;a(44,45)=15.8;a(48,49)=14.2; a=a+a';
a(find(a==0))=inf; result=[]; p=1;
tb=2:length(a);
while length(result)~=length(a)-1 temp=a(p,tb);temp=temp(:); d=min(temp);
[jb,kb]=find(a(p,tb)==d); j=p(jb(1));k=tb(kb(1)); result=[result,[j;k;d]]; p=[p,k];
tb(find(tb==k))=[]; end result
程序二
题二. 县城到O 到各点的距离程序及结果 model: sets:
cities/A35,A34,A33,A32,A31,A30,A29,A28,A27,A26,A25,A24,A23,A22,A21,A20,A19,A18,A17,A16,A15,A14,A13,A12,A11,A10,A9,A8,A7,A6,A5,A4,A3,A2,A1,R,Q,P,N,M,L,K,J,I,H,G,F,E, D,C,B,A,O/:fl;
roads(cities,cities)/A35,A34 A35,A33 A35,A32 A34,B A34,A A32,A33 A31,A33 A33,A A32,A31 A30,A32 A31,R A30,Q A29,R Q,A29 A29,P A27,A28 Q,A28 A28,P A27,A26 A24,A27 A26,P N,A26 A21,A25 A20,A25 N,A25 A25,M A24,A23 A24,N A22,A23 A23,A21 A23,N A17,A22 A22,K A21,A20 K,A21 A19,A20 A20,L A19,L J,A19 A18,K A18,J I,A18 A16,A17
A17,K A16,I A15,A14 A15,I A14,A13 H,A14 A13,J A13,I A13,G H,A12 A12,G A12,F J,A11 G ,A11 A11,E A10,F F,A9 A9,E A8,A4 A8,E A7,A6 L,A7 E,A7 A7,D A6,A5 M,A6 L,A6 A5,A2 M,A5 D,A5 A4,D A3,A2 D,A3 A3,C A2,O C,A1 B,A1 A,A1 A1,O A,R R,O P,O N,M M,O I,J C,B C,O/:w,p; endsets data:
w=8.2 20.3 14.9 17.6 11.5 19 7.3 7.4 8.1 10.3 9.2 7.7 7.9 7.2 15.2 7.9 8.3 12.1 7.8 18.8 10.5 10.5 7.8 6.5 8.8 12 8.9 13.2 10 9.1 7.9 6.7 10.1 7.9 4.1 9.3 5.5 7.2 8.1 9.2 8.2 8.2 6.8 9.8 11.8 15 8.8 8.6 9.9 9.8 16.4 8.6 10.2 7.8 12.2 13.2 6.8 14.2 10.8 5.6 7.8 20.4 8 7.3 14.5 7.2 15.1 9.7 9.5 11.8 8.3 11.4 11.3 12.7 4.8 8.2 7.9 9.2 11.2 5.9 10.3 6 8.8 12.9 10.1 14.2 19.8 15.8 11 11.5; enddata
n=@size(cities); fl(n)=0;
@for(cities(i)|i#lt#n:fl(i)=@min(roads(i,j):w(i,j)+fl(j))); @for(roads(i,j):p(i,j)=@if(fl(i)#eq#w(i,j)+fl(j),1,0)); end
Feasible solution found.
Total solver iterations: 0
V ariable Value N 53.00000 FL( A35) 36.00000 FL( A34) 27.80000 FL( A33) 23.70000 FL( A32) 30.20000 FL( A31) 22.10000 FL( A30) 35.70000 FL( A29) 20.80000 FL( A28) 22.20000 FL( A27) 28.40000 FL( A26) 20.60000 FL( A25) 31.80000 FL( A24) 44.30000 FL( A23) 39.00000 FL( A22) 49.00000 FL( A21) 39.60000 FL( A20) 38.30000 FL( A19) 46.20000 FL( A18) 52.90000 FL( A17) 53.50000 FL( A16) 60.30000 FL( A15) 69.90000 FL( A14) 72.70000
FL( A13) 64.10000 FL( A12) 67.30000 FL( A11) 55.90000 FL( A10) 65.90000 FL( A9) 49.50000 FL( A8) 49.70000 FL( A7) 34.50000 FL( A6) 27.20000 FL( A5) 17.50000 程序三
问题四编程(以第一问第一组为例) t1=1; s=200.4; m=7; n=10; b=2;
v=15:0.5:100;
T=s./v+(a*m+n)*t1; plot(v,t); xlabel('V'); ylabel('T');
FL( A4) FL( A3) FL( A2) FL( A1) FL( R) FL( Q) FL( P) FL( N) FL( M) FL( L) FL( K) FL( J) FL( I) FL( H) FL( G) FL( F) FL( E) FL( D) FL( C) FL( B) FL( A) FL( O) 34.90000 14.00000 9.200000 6.000000 12.90000 28.00000 10.10000 31.10000 19.80000 39.00000 43.70000 54.30000 61.10000 77.50000 62.70000 55.10000 41.70000 22.20000 11.50000 11.90000 16.30000 0.000000
v=35;
s=200.4;m=7; n=10; b=2;
t1=0:0.2:2;
T=s./v+(a*m+n)*t1; plot(t1,T); xlabel('t1'); ylabel('Tˊ')