灾情巡视路线最优化方案(刘).doc

约束最优路线的模拟退火解法说明：以98年全国大学生数模竞赛中的B题（即“灾情巡视路线”）为例，介绍能解一类较广泛的约束最优路线问题的方法模拟退火法1。该法对“灾情巡视路线”这类有约束以及“（一般）旅行推销员”、“中国邮递员”等无约束组合优化问题均能求得较好的近似解，具有适用范围广和可拓展的优点。一、问题描述对于最短路、最大流、中国邮递员、旅行推销员等最优路线问题，常采用各自不同的方法求解。若在这些问题中再加入一些约束条件，则原方法往往不再有效，如98年大学生数模竞赛中的B题就是如此。我们设计的方法较好地解决了这一问题。现以98年B题为例，介绍该法及其实现。下面为该题文字部分，并称其四问分别为问题1至问题4：下图（略）为某县的乡（镇）、村公路网示意图，公路边的数字为该路段的公里数。今年夏天该县遭受水灾。为考察灾情、组织自救，县领导决定，带领有关部门负责人到全县各乡（镇）、村巡视。巡视路线指从县政府所在地出发，走遍各乡（镇）、村，又回到县政府所在地的路线。若分三组（路）巡视，试设计总路程最短且各组尽可能均衡的巡视路线。假定巡视人员在各乡（镇）停留时间T=2小时，在各村停留时间t=1小时，汽车行驶速度V=35公里小时。要在24小时内完成巡视，至少应分几组；给出这种分组下你认为最佳的巡视路线。在上述关于T，t和V的假定下，如果巡视人员足够多，完成巡视的最短时间是多少；给出在这种最短时间完成巡视的要求下，你认为最佳的巡视路线。若巡视组数已定（如三组），要求尽快完成巡视，讨论T，t和V改变对最佳巡视路线的影响。二、问题分析及模型的建立因为是分组巡视（不妨设分N 组），要直接确定一个组巡视哪些地点是困难的。由于将各组巡视的路线连接起来可看成一条N 次相继从县城出发又回到县城的路线，这样，多组巡视就化成了单组巡视。经分析，我们认为前3问及第4问计算部分都是组合规划中的约束优化问题，均属以模型 （I）为基础的约束最优路线模型。下面根据各问的要求，分别对4个问题进行具体讨论。对于问题1，如果选取总路程最短的所有巡视路线中最均衡的，一般这一路线仍会很不均衡。故除了要总路程短，另需“均衡”提出一定的要求，即组间巡视路线的长度差不大于某给定值L。还有路线能够分成3次从县城O出发再回到O、各组经过地点的并集为所有顶点的集合只之约束。模型如下： (II)其中F为巡视总路程，N为要求的分组数（本问N=3），n是优化过程中路线的实际分组数，fmax和fmin分别为n组中最长和最短组的巡视路程，Pi 为第i组巡视地点的集合，A是所有顶点的集合。约束条件(fmax-fmin)L用来保证各组路程基本均衡，目标函数中加入l(fmax-fmin)是为了使各组的路程尽可能均衡，这是一个约束多目标规划。权重l的取值应远小于目标函数中F的权系数1，否则会因离散问题函数值的跳跃现象而导致优化困难。这里，我们取l=1/L。 对问题2和3，因在原图中不是任意两顶点间均有边，故在多组巡视时可能存在二组以上经过的公共点，从而使顶点赋权遇到如下困难：若先赋点权，则这些公共点的权会被重复计算；而若在优化过程中赋点权，则又有公共点究竟应该由哪次（组）巡视的问题。对此，我们采用Dijkstra法算出任意两点间的最短路程并除以速度V转化为时间，并用其作为两顶点间边的权构造一新图。第三个约束条件保证了每点至少经过一次，而这时若路线中含除县城外的重复地点将至少增加1个小时（除县城为0外，其它点权为1或2小时），因而能保证优化结果中不出现县城以外的重复点。经分析，我们认为这两问同属分组数和路线双重优化问题，具体模型如下： (III)其中N、n 、A、Pi与模型（II）相同，F为各组巡视时间之和， fmax为n组中时间最长组的巡视时间，TMAX是巡视允许的最长时间。对问题2，TMAX=24；对问题3，因人员足够多，故每个地点可由一组单独巡视，因此完成巡视的最短时间是这些单点巡视的最佳路线中花时间最长的组对应的巡视时间，其值TMAX由程序计算确定。对最少分组数的优化，模型中没有体现，我们是通过主程序与辅助控制函数协同工作来实现该目的的，即在优化过程中，若路线调整发现最短与次短组的时间之和不大于TMAX，且满足约束，则记录该路线后返回，并由主程序将分组数减1后重新优化，如果重新优化找不到满足约束的可行路线，则认为找到了最优解；而若优化函数没有找到满足约束的可行路线，主程序则将分组数加1后进行下一轮优化。我们用m和M分别表示求解时尝试的最小和最大分组数，对问题2，取m=4，M=8，而对问题3，m=20，M=35。对于问题4，其最佳路线是指在分组数已确定的情况下，巡视时间最长的组所对应的时间尽量短的路线。下面给出其模型，各参数的含义同模型（III），关于T、t和V对路线的影响，我们将在结果与分析部分给出。 (IV)三、模型的求解 (一）算法选择 由以上分析可知，本题的所有问题都可视为“在一定约束下，从一点出发再回到该点的最优路线寻求问题”。如问题1可看成有约束的“一般旅行推销员”问题，但即使是无约束的“旅行推销员”问题，也尚未找到解这类问题最优解的多项式算法3，对本题n=53这种情形，无论是“分枝与界限法”4还是动态规划5等方法，想得到最优解一般是不可能的，只能用近似算法。选择模拟退火法，是因为它不同于一般的单纯下降法，它具有一定条件下的随机上升过程，有可能从一极小区域跃过波峰到过另一极小区域，因而具有解更接近全局最优解的优点。 (二)计算方法与步骤 1方法概述 下面介绍模拟退火法函数，主函数部分略。模拟退火法是一种较新的算法，尚无固定的计算步骤1。我们设计的基本操作包括路线调整、约束条件判断、调整许可预测、内部优化、与已有最佳路线的比较和最佳路线及相关数据的记录等。路线调整在给定的工作点集序列上进行，其前部分为当前路线，后部分为剩余点集。如此设计是因为在许多问题中某些点在路线中可能出现多次，如“田”字型图的“一般旅行推销员”问题，此时可将每一可能重复的点或所有点以最大可能重复出现的次数置于工作点集中。路线的调整有5种（所有操作对象都是随机选取的）：交换路线内两路段（含点）互换；反转路线中选取一路段反向；删除路线内一路段移至路线外；插入路线外一路段插入到路线中；替换路线内外各取一路段互换。约束条件判断由用户函数完成；调整允许预测由内部函数完成，它决定是否接受一上升过程；内部优化是仅有交换和反转操作的优化过程，与上层模拟退火法函数基本相同。2基本步骤 (1) 由主程序传送工作点集、初始路线终点、记录器数目、温度初值t0、温度下降率FT、温度循环上限maxtk、每个温度的最大循环次数maxk和最多成功次操作次数max_suc、连续温度无改进时结束优化参数K等，并设置寻找可行路线的不成功次数n_suc:=0、最优路线条数rec_n:=0、优化值fx0:=maxfx、温度循环变量tk:=1，优化连续无改进温度次数same_tn:=0等内部变量初值；(2) 计算初始路线的目标函数值fx、优化控制函数值fx1。若fx1maxfx，则记录路线并置fx0:=fx+fx1、rec_n:=0；(3) 如果tk=maxtk或same_tn=K，转步(15)；否则，置循环变量kk:=1，成功操作次数计数器sucn:=0等；(4) 如果kk=maxk或sucn=max_suc，转步(14)；否则，计算当前路线的优化控制函数值fx1、在复制路线上调整路线并计算目标函数值fx和优化控制函数值fx2，得到两者的改变量dfx和dfx1，置dfx12:=dfx+dfx1，kk:=kk+1；(5) 若fx20或rec_n=0且fx2fx1，转步(4)；(6) 执行100次内部优化并更改路线，若某次内部优化后fx230时，转步(15)；其它情况转步(4)；(7)计算路线的辅助控制函数值fx3，若fx3fx0，即比已有的最佳路线差，转步(4)；否则检查是否与记录的路线重复，若重复则转步(4)；(10)若fx+fx2fx0，则记录器中原最佳路线后移，最佳路线数rec_n:=1，更新优化值fx0:=fx+fx2，路线及相关数据记入第一个位置，转步(4)；(11) rec_n:=rec_n+1，并计算记录器中最后一条最优路线的辅助控制函数值fx4；(12)若fx3fx4，与记录器中各路线的辅助控制函数值比较并据此确定记录器指针，当前路线插入相应位置，转步(4)；(13)若记录器没满，则路线及相关参数记在已有最优路线之后；转步(4)；(14)若tk0，则找到了最优路线，可从记录器中获取最优和次优路线及相关数据。3编程提示 对分组的处理及模型中的变量n、Pi、fmax、fmin、F和下面用到的次短组巡视时间fmin1均由分组处理函数完成，因此它是解决分组问题的关键，分组方法是：当自然分组数不大于N时，取自然分组；否则将最短两组合并，直到组数等于N为止。与总的设计思想相对应，程序设计应考虑多功能与通用性，以适应起点与终点相同和不同、单组和多组等各种要求。由于有上升过程，故需一组变量存放当前最佳路线、最优值和路线终点位置等。记录器参数组可满足人们想获取多条最优路线的要求。参数maxfx之值为“无穷大”，也是判断约束是否满足的依据，即满足约束时优化控制函数值小于maxfx。三个用户函数分别为目标函数、优化控制函数和辅助控制函数。优化控制函数在未找到可行路线时为罚函数（值大于等于maxfx），否则为0或多目标规划中减去实际目标函数值的剩余部分。辅助控制函数实现两种控制，一是优先记录某种（不影响优化过程的）特性的路线，如优先记录均衡路线等；其二是有些问题仅要求实现一固定目标，如本题最少分组数的寻找，只要两最短组的时间之和不超过要求的时间，它们即应合并成一组，此时继续优化是多余的，应将分组数减1后重新优化。函数返回值大于或小于等于-maxfx可区分上述两种情况。另外，定义了一称为优化值的变量，它是目标函数值与优化函数值之和。程序在没有找到和已有可行路线时的处理不同，前者主要是寻找可行路线，只关心优化控制函数值是否下降，而已有可行路线时的关注点在优化值的改变。使用优化值而不是目标函数值判断路线的优劣，是因为在许多多目标规中常常只要一个目标的结果，且把其余部分放在优化控制函数中可以仅计算路线改变部分的目标函数值，因而在处理多目标规划时可有更好的性能。程序设有奇偶两种处理方式，其中偶方式为标准多目规划或只计算路线改变部分得不到目标函数值（如第4问）的情形而设。设计奇偶方式而非固定参数，可以使程序具有对多个任务的处理能力。（三）求解 对本题的求解，我们根据各问的要求和约束条件，设计相应的目标函数、优化控制函数和辅助控制函数，用Borland C+ 3.1编程，在多台486、586及PII机上运行通过，解见结果与分析。对所有问题，编制同一组（三个）用户函数，根据处理方式参数fxmode之值判断处理：对问题1, fxmode1；对问题2和3， fxmode3；对问题4, fxmode2。目标函数是明显的，现按单问题形式将其它两函数介绍如下：问题1 优化控制函数为：f = maxfx (r1+ r2)+ r3 maxfx+1000( fmax - fmin -L )；辅助控制函数为：f = fmax - fmin ，起优先记录均衡路线的作用。其中r1是路线中缺少地点的数目，控制必须巡视每个地点；r2 =N-n为要求与实际分组的差值，控制必须分成N个组，在下面的各问题中r1、r2的意义和作用完全相同；而则用于实现均衡的基本要求。问题2和3 优化控制函数为：maxfx(r1+ r2) + r3 maxfx+1000 ( fmax -TMAX)；辅助控制函数为：（1-r4）( fmax - fmin) - r4 maxfx；这里，控制每组巡视时间不超过要求；主要用于优化最少分组数。辅助控制函数在本模型中具双重功能，一是在优化过程中遇到巡视时间最短的两个组的时间之和不超过要求的时间时返回-maxfx，使优化函数记录该路线并返回，实现优化分组数的目的；其次是优先记录均衡路线。问题4 优化控制函数为：maxfx；辅助控制函数为：F ，使最佳路线中总时间短的路线优先记录（F 为总路程）。四、结果与分析为使记录的路线不重复，且结果输出时能按原图输出等，我们编制了原图与新图路程数据切换、路线是否相同的检查、相邻两点用其间的最短路线替换，以及将编号转化为原地点名等多个辅助函数。为防止内存溢出，路程数据只存储矩阵的上三角部分。为使解尽量好，还采用了多轮算法，即先以任给初始路线进行若干次优化，再把结果中最好的路线作为初始路线进一步优化，如果找到了更好的路线，则再以该路线为初始路线进行优化，直到没有改进时为止。所有计算的记录器均设为10组，一般都有多条最佳路线。（一）运行结果 问题1 为实现均衡要求，我们令L=70，即与完全均衡时相比，路程差不超过1小时的路程。运行结果为：总路程568.7公里，平均每组巡视189.57公里，组间最大路程差值58.00公里，与平均值的最大差值为31.67公里。选1条如下：第1组(157.9公里)：OR29Q282728Q303231333534AB1O；第2组(194.9公里)：OM2520L19J18I15I161722K212324N26PO；第3组(215.9公里)：OC3D48E11G1314H12F10F9E7652O. 问题2 至少需要4组巡视，最短时间23小时，4组合计时间87.066小时，组间最大时间差3.45小时，列1条如下：第1组(21.5小时)：O2567E9F10F9E84D3CO；第2组(19.6小时)：O1BA343533313230Q29RO；第3组(23.0小时)：O2567E11G12H1413J19L2025N26PO；第4组(22.9小时)：OP28272423221716I15I18K2125MO.其中用括号“”括起来的部分为经过而不巡视的路段，问题3亦是如此。 问题3 巡视最短时间TMAX=6.429小时，至少需要22组巡视，22组合计时间130.231小时，组间最大时间差值有两种。各选一条且第二条只列出与第一条有差别的第8组和第12组： (1)最长时间为6.429小时（第16组），最短时间是4.777小时（第4组），最大组间时间差值为1.651小时： 第 1组（6.383小时）：O2567E7652O； 第 2组（5.583小时）：O2567E11G11E7652O； 第 3组（6.309小时）：OM2521K1716I15I18K2125MO； 第 4组（4.777小时）：OP26N26PO； 第 5组（5.717小时）：OP26N23242726PO； 第 6组（6.120小时）：OP26N232217K2125MO； 第 7组（5.391小时）：OM52O； 第 8组（5.354小时）：O1A33A1O； 第 9组（6.111小时）：OR29Q30Q28PO； 第10组（5.994小时）：O23D4D32O； 第11组（6.366小时）：O256L2025MO； 第12组（4.983小时）：O1BCO； 第13组（6.154小时）：O256L19J131413J19L652O； 第14组（6.317小时）：OP29RO； 第15组（6.103小时）：O256L19J19L652O； 第16组（6.429小时）：O2567E9F12H12F9E7652O； 第17组（6.317小时）：O1A34353231RO； 第18组（6.149小时）：O2567E9F9E7652O； 第19组（5.491小时）：OM2521K18I18K2125MO； 第20组（6.023小时）：OM2521K18K2125MO； 第21组（5.937小时）：O2567E9F12G11E7652O. 第22组（6.223小时）：O2567E9F10F9E8E7652O. (2)最长时间为6.429小时（第16组），最短时间是4.354小时（第8组），最大组间时间差值为2.074小时：第 8组(4.354小时)：O1A33A1O；第12组(5.983小时)：O1BCO. 问题4 受篇幅限制，仅给出结论：无论是T、t 还是V，都有增加时巡视总路程增加，而减少时巡视总路程变短的趋势。运行结果符合分析结论，且以县城为参照点，路线大体上分成左上、左下和右（含左上几个点）三部分。 （二）结果分析 由于我们手中没有其它程序和结果可供比较，故仅能以多次运行的情况进行分析。 对问题1，若将其理解为“寻找最短路线中的最均衡的路线”，则计算结果为：总路程值533公里，第一组12公里，第二组23.8公里，第三组497.2公里。走法为：第1组：O1O；第2组：O1B1O；第3组：OP29R3133A34353230Q282726N2423221716I15I18K2125 20L19J11G1314H12F10F9E84D76M 52 3CO.这显然与“尽量均衡”的要求相去甚远，与我们在“问题分析及模型的建立”中的预见相符。如果用多目标规划而不加差值约束，当总路