基于遗传算法解决TSP问题

1、基于遗传算法解决TSP问题摘要 题目要求给出环游全国全部省会的最短路径方案，是传统的TSP问题，本文将图表数据数字化后，将其转变成为线性规划问题，进而采取遗传算法用Matlab求解出理论上的最短路径与路线图。通过第一问求出的路线顺序结合实际情况求解出实际情况下的最短路径与最短时间。针对第一问，首先建立基本TSP模型，求出其线性规划方程组，用Matlab对地图做出基本处理，求出其像素坐标的矩阵。将省会城市初始化为种群数据，用遗传算法求解出模型最优解，即最短路径大小与旅游城市顺序。针对第二问，由于遗传算法求出的是近似最优解，以及实际道路情况不可能是直线距离，所以理论数据与实际有一定差别。将旅行顺序

2、求解出后，需要根据实际道路情况重新求解出最短路径大小，并根据题目所给条件求解出最短时间。根据实际情况，求得最后的最短路径长度为20402.9公里，时间为45天。题目中给出城市转化为图集，一共有33个顶点，数据较大，用Dijkstra算法和Floyd算法虽然答案可能更精确，但数据处理量大，时间复杂度高。采用遗传算法可以得到近似最优解，并且精简了时间复杂度。关键词： 最短路径 TSP问题 线性规划 遗传算法问题重述 如果从杭州出发，要想开车走遍全国所有的省会城市，而且在到了每个省会城市以后都必须住一晚，第二天早上才能出发，安全起见每天开车时间最多在8小时左右，车速视实际路况而定，一般高速公路可以取

3、平均车速100公里/小时左右，不是高速公路平均车速取60公里/小时左右，从杭州出发要走遍所有省会城市以后回到杭州，开车里程最少需要多少公里？最少需要多少时间？（暂不考虑台湾省）问题分析题目要求求出最短路径与时间是典型的TSP问题，即要用最短的总路径走遍所有城市，此处需要设立一个二元组，并且求出一个点的邻接矩阵。根据经典的TSP模型，得出一般的线性规划方程组。解TSP模型的一般算法有Dijkstra算法和Floyd算法，由于TSP问题是典型的NP难题，其可能路径数目与城市总数目n是呈指数型增长，并不适合用以上两种算法。另外常用来解决TSP问题的退火算法受概率和退火过程影响，可能运算起来十分缓慢。

4、为了精简时间复杂度，故选取了遗传算法将图形数据初始为种群数据进行计算。并且随着种群规模的增大，其结果越逼近最优值。对于第二问，最短时间的计算，按照第一问求出的最优路线行走，并且结合实际高速公路分布与走向来计算最短路径的真实值，并且统计出高速公路与普通公路的实际数值。题目中给出了车辆在高速公路上行车速度为100km/h，普通公路上为60km/h，统计出实际数值后，计算出总行车时间即可得出最短时间。基本假设题目给出的速度符合实际。无论第一天何时到达，第二天均会离开该城市。每天开车时间为8小时左右，对于9到10小时内可以开到的城市均算为8小时左右，不需要中途休息。开车行程在10小时以后的，每天开车八

5、小时即停下并有地方休息。第一问用TSP模型求解时，所有城市之间全部算为直线相连。符号说明Z： 走遍所有省会城市的最短路径长度C： 路网有向图矩阵； ：从节点i到节点j的弧，C；：路段的长度； M： 种群个数N： 染色基因个数（即城市个数）c： 迭代次数Pc： 交叉概率Pm： 变异概率len： 某组解的总距离minlen：该次迭代中的最短距离maxlen：该次迭代中的最长距离m： 适应度归一化淘汰加速指数模型建立与求解5.1线性规划的约束条件与目标函数： 已设为两城市之间的距离，根据TSP问题的一般约束条件可知，记为赋权图G=(V,E)，V为顶点集，E为边集，各顶点间的距离已知。设目标函数与约束

6、条件： 5.2遗传算法前的基本处理 将全国城市的经纬度制成表格，用excel将其制成散点图，方便接下来的数据处理与运算。城市东经北纬城市东经北纬城市东经北纬北京116.4639.92沈阳123.3841.8长沙11328.21天津117.239.13石家庄114.4838.03武汉114.3130.52上海121.4831.22太原112.5337.87广州113.2323.16重庆106.5429.59西宁101.7436.56海口110.3520.02拉萨91.1129.97济南11736.65兰州103.7336.03乌鲁木齐87.6843.77郑州113.6534.76西安108.95

7、34.27银川106.2738.47南京118.7832.04成都104.0630.67呼和浩特111.6540.82合肥117.2731.86贵阳106.7126.57南宁108.3322.84杭州120.1930.26昆明102.7325.04哈尔滨126.6345.75福州119.326.08香港114.122.2长春125.3543.88南昌115.8928.68澳门113.3322.13 用基于Java的经纬度坐标与地图容器像素坐标转换器，求出33个省会城市的像素坐标如下：1.西藏 2.云南 3.四川 4.青海 5.宁夏 6.甘肃 7.内蒙古 8.黑龙江 9.吉林 10.辽宁 11.

8、北京 12 天津 13.河北 14.山东 15.河南 16.山西 17. 陕西 18.安徽 19.江苏 20.上海 21.浙江 22.江西 23.湖北 24.湖南 25,贵州 26. 广西 27.广东 28.福建 29.海南 30.澳门 31.香港 32.重庆 33.新疆 像素坐标如下： Columns 1 through 11 100 187 201 187 221 202 258 352 346 336 290 211 265 214 158 142 165 121 66 85 106 127 Columns 12 through 22 297 278 296 274 265 239 30

9、2 316 334 325 293 135 147 158 177 148 182 203 199 206 215 233 Columns 23 through 33 280 271 221 233 275 322 250 277 286 220 104216 238 253 287 285 254 315 293 290 226 775.3遗传算法求解最优值1. 程序初始化 程序首先读入 33 个省会城市坐标，计算任意两个城市的距离； 设置遗传算法控制参数。 clear all;clc;clf; load(testdata.mat); nlen=length(x1); xy=x1;y1; n

10、 = 500; %种群数目 C = 5000; %进化迭代次数 m=2; %适应度归一化淘汰加速指数，取值不宜太大 alpha=0.8; %淘汰保护指数，范围01，为1时关闭保护 a = meshgrid(1:nlen); %生成 n x n矩阵 dmat = reshape(sqrt(sum(xy(a,:)-xy(a,:).2,2),nlen,nlen); %计算城市距离矩阵 遗传算法对求解问题本身是一无所知的，这里采用随机生成初始化种群，如 下： N,NN=size(dmat); farm=zeros(n,N); %用于存储种群 for i=1:n farm(i,:)=randperm(N

11、); %随机生成初始化种群 end 2 计算适应度 本程序目标函数为经过 33 省会的总距离，适应度与目标函数的正相关，取 值范围 01，适应度计算公式为： 3.4 变异操作 采用变异操作解 决同一问题。算法核心思路是，在每次迭代中，解的个体随机按 4 个为 1 组，每组中只保留最优解，然后对此最优解进行左右翻转、交换、向前移位三种变异操作，生成三个新个体，再参与下次迭代。整个算法不需要计算归一化适应度。在此算法中，每次迭代淘汰率固定为 75%，三种变异操作覆盖面比较广，直 接以最短距离为适应函数，省去每次适应度的计算。初始值种群数目设为 500， 进化迭代次数为 5000，寻得最优解为 12

12、95.72。5.4算法求解结果图一 散点图最优表示图二 仿真图表示5.4根据实际情况的求解答案路线：杭州上海南京合肥武汉南昌长沙重庆成都西安郑州济南天津沈阳长春哈尔滨北京石家庄太原呼和浩特银川兰州西宁乌鲁木齐拉萨昆明贵阳南宁海口澳门广州香港福州杭州实际数据图起点终点总距离/公里公路/公里时间/小时高速/公里时间/小时总时间/小时天数杭州上海177.600177.61.7761.7761上海南京299.200299.22.9922.9921南京合肥175.600175.61.7561.7561合肥武汉378003783.783.781武汉南昌355.400355.43.5543.5541南昌长沙

13、345.300345.33.4533.4531长沙重庆892.400892.48.9248.9241重庆成都311.700311.73.1173.1171成都西安712.900712.97.1297.1291西安郑州482004824.824.821郑州济南436.600436.64.3664.3661济南天津325.400325.43.2543.2541天津沈阳664.500664.56.6456.6451沈阳长春380.400380.43.8043.8041长春哈尔滨271.400271.42.7142.7141哈尔滨北京1241.6001241.612.41612.4162北京石家庄29

14、2.200292.22.9222.9221石家庄太原240.800240.82.4082.4081太原呼和浩特468.667.51.125401.14.0115.1361呼和浩特银川718.200718.27.1827.1821银川兰州440.700440.74.4074.4071兰州西宁225.429.60.493333333195.81.9582.4513333331西宁乌鲁木齐1740.8339.35.6551401.514.01519.673乌鲁木齐拉萨2667.12667.144.451666670044.451666676拉萨昆明2223.82223.837.06333333003

15、7.063333335昆明贵阳507.500507.55.0755.0751贵阳南宁577.800577.85.7785.7781南宁海口474004744.744.741海口澳门583.600583.65.8365.8361澳门广州136.700136.71.3671.3671广州香港172.700172.71.7271.7271香港福州864.500864.58.6458.6451福州杭州618.500618.56.1856.1851总计20402.95327.388.7883333315075.6150.756150.75645 最终的最优解，最短路径长度为20402.9公里，时间为45

16、天。模型评价本文采用MATLAB实现遗传算法求解TSP问题，对结果进行了分析。遗传算法是一种智能优化算法。该模型最大的缺点是，求出的不是最优解，而是近似最优解。同时，在做该模型时，是通过经纬度计算，算出的距离接近于两地直线距离。而真实情况中，大部分公路都不是笔直的，而是弯弯曲曲的，因此该模型求解答案与真实值有差距。同时，建最短时间模型过程中，忽略了很多客观的不定因素，比如天气问题等，会使正常计划受到一定的影响。还有驶入高速前的一段公路上的速度也被近似看成高速车速。但是本模型用区域性划分和动态分析法解决游遍全国的省会城市、直辖市、香港、澳门，相对于传统的动态规划解法，达到了省时、简便的效果，大大

17、降低了计算的复杂性。我们还需再靠近实际情况，再精确该模型，减少误差。参考文献1铁菊红，彭辉：一种改进的基于高斯分布拟合的提取标志点像素坐标方法。Computer and Modernization，2008（4）。2解晨，韦雄亦：模拟退火算法和遗传算法的比较与思考。电脑知识与技术，2013,9（19）。3王勇：用遗传算法求解中国旅行商问题。哈尔滨商业大学学报：自然科学版，2005（4）。4刘英：遗传算法中适应度函数的研究。兰州工业高等专科学校学报，2006.9，vol13（3）。5周凯，宋全军，邬学军：数学建模竞赛与提高，图与网络模型。6 王剑文，戴光明，谢柏桥，张全元。求解 TSP 问题算法

18、综述J。计算机 工程与科学.2008(02)。 八、附录基于Java的经纬度坐标与地图容器像素坐标转换器经纬度坐标与地图容器像素坐标相互转换bodymargin:0;height:100%;width:100%;position:absolute;#mapContainerposition: absolute;top:0;left: 0;right:0;bottom:150px;#btnDivposition:absolute;bottom: 0;left:20px;height:150px;font-size:12px;#btnDiv inputheight:23px;outline:non

19、e;border:1px solid #ddd;padding-left:5px;border-radius:3px;#btnDiv inputtype=buttonheight:28px;line-height:28px;outline:none;text-align:center;padding-left:5px;padding-right:5px;color:#FFF;background-color:#0D9BF2;border:0;border-radius: 3px;margin-top:5px;margin-left:5px;cursor:pointer;margin-right

20、:10px;地图经纬度坐标：(鼠标左键在地图上单击获取经纬度坐标)X：Y：地图容器像素坐标：X：Y： /初始化地图对象，加载地图var map = new AMap.Map(mapContainer, resizeEnable: true);AMap.event.addListener(map,click,getLnglat);/鼠标在地图上点击，获取经纬度坐标function getLnglat(e) document.getElementById(lngX).value = e.lnglat.getLng(); document.getElementById(latY).value = e

21、.lnglat.getLat(); /经纬度坐标转换为容器像素坐标document.getElementById(lng2x).onclick = function()var px = document.getElementById(lngX).value; var py = document.getElementById(latY).value; var pixel = map.lnglatTocontainer(new AMap.LngLat(px,py); document.getElementById(pixelx).value = pixel.getX(); document.get

22、ElementById(pixely).value = pixel.getY(); ;/容器像素坐标转换为经纬度坐标document.getElementById(x2lng).onclick = function()var lnglatX = parseInt(document.getElementById(pixelx).value); var lnglatY = parseInt(document.getElementById(pixely).value); var ll = map.containTolnglat(new AMap.Pixel(lnglatX,lnglatY); doc

23、ument.getElementById(lngX).value = ll.getLng(); document.getElementById(latY).value = ll.getLat(); ;变异法核心代码for p = 4:4:pop_size rtes = pop(rand_pair(p-3:p),:); dists = total_dist(rand_pair(p-3:p); ignore,idx = min(dists); best_of_4_rte = rtes(idx,:); ins_pts = sort(ceil(n*rand(1,2); %生成 1x2 每一列元素 %按

24、照升序排列矩阵 I = ins_pts(1); J = ins_pts(2); for k = 1:4 %保留最佳个体，繁殖三个新个体 tmp_pop(k,:) = best_of_4_rte; switch k case 2 %左右翻转 tmp_pop(k,I:J) = fliplr(tmp_pop(k,I:J); case 3 %交换 tmp_pop(k,I J) = tmp_pop(k,J I); case 4 %向前移动一位 tmp_pop(k,I:J) = tmp_pop(k,I+1:J I); otherwise end end new_pop(p-3:p,:) = tmp_pop

25、; end pop = new_pop; 附录资料:不需要的可以自行删除c语言典型问题处理方法小结循环问题（1）、数论问题1、求素数 for(i=2;i1，如果它仅有平凡约数1和a，则我们称a为素数(或质数)。整数 1 被称为基数，它既不是质数也不是合数。整数 0 和所有负整数既不是素数，也不是合数。 2、求最大公约数和最小公倍数a、 if(ab) t=a; a=b; b=t; for(i=a;i=1;i-) if(a%i=0&b%i=0) break; printf(largest common divisor:%dn,i); printf(least common multiple:%dn

26、,(a*b)/is);b、辗转相除法求解 a1=a; b1=b; while(a%b!=0) t=a%b; a=b; b=t; printf(largest common divisor:%dnleast common multiple:%d,b,a1*b1/b);3、求完数 一个数如果恰好等于它的因子之和，这个数就称为“完数”。 例如：6的因子为1、2、3，而6123，因此6是“完数”。for(a=1;a=1000;a+) s=0; for(i=1;i=a) break; if(s=a) printf(%dt,a);注意S=0所放的位置 4、分解质因数 将一个整数写成几个质因数的连乘积，如：

27、 输入36，则程序输出36=2*2*3*3 。解一、看似简单，但要自己完整地写出来还真不容易！竟然还动用了goto语句，正好可以熟悉一下goto语句的用法！main() int a,z,i; clrscr(); scanf(%d,&a);判断下一个数开始有要重新从2开始了。所以用loop语句回到for语句，这是for语句仍从2初始化。从2开始的原则不变，变的是a的值。 loop: for(z=2;z=a;z+)判断是否为质数for(i=2;i=z;i+) if(z%i=0) break;判断是否为a的质因数 if(z=i) if(a%z=0) k+; if(k=1) printf(%d=%d,

28、a1,z);用计数器来解决每行输入不同的问题。 else printf(*%d,z); a/=z; goto loop; 解二：main() int n, k=2, isfirst=1; printf(Input n=); scanf(%d,&n); while(k=n) if(n%k=0) if(isfirst) printf(%d=%d, n, k); isfirst=0; else printf(*%d,k); n/=k; else k+; printf(n);5、从键盘输入两个整数，输出这两个整数的商的小数点后所有1000位整数 for(i=1;i=2;i-) if(fm%i=0&fz

29、%i=0) fz/=i; fm/=i; z=fz/fm; fzx=fz%fm; if(fzx=0) printf(%d%d/%d-%d%d/%d=%dn,z1,fz1,fm1,z2,fz2,fm2,z); else if(z=0) printf(%d%d/%d-%d%d/%d=%d/%dn,z1,fz1,fm1,z2,fz2,fm2,fzx,fm); else printf(%d%d/%d-%d%d/%d=%d%d/%dn,z1,fz1,fm1,z2,fz2,fm2,z,fzx,fm);（2）近似问题1、书P122习题4-6注意千万不要忘记添加#include “math.h”#include

30、 math.hmain() float x,j=1,k,s,so; int n; scanf(%f,&x); s=x; so=x+1; for(n=1;fabs(s-so)1e-6;n+) for(k=1;k1e-6) x=(x1+x2)/2; f=x*x*x+4*x*x-10;可以用/*if(f*f10) x2=x; else x1=x; printf(%fn,x);（3）枚举法（4）数列问题二、数组问题（1）排序问题1、从小到大排序main() int a10,i,j,t; for(i=0;i10;i+) scanf(%d,&ai); for(i=1;i10;i+) for(j=0;jaj

31、+1) t=aj+1;aj+1=aj;aj=t; for(i=0;i10;i+) printf(%d ,ai); printf(n);注意排序问题：1、须迅速，熟练，无差错经常插入在程序中间2、现使用最大数下沉冒泡法还可以使用最小数上浮冒泡法3、j控制前面一个数和后面一个数一一比较。由于是最大数下沉，i+1后j仍要从0开始。4、i控制这样的操作一共要做多少次5、注意i j的控制次数2、从大到小排序main()现使用最大数上浮冒泡法还可使用最小数下沉冒泡法 int a10,i,j,t; for(i=0;i10;i+) scanf(%d,&ai); for(i=1;i=i;j-) if(ajaj-

32、1) t=aj-1; aj-1=aj; aj=t; for(i=0;i10;i+) printf(%d ,ai);（2）二维数组三、字符或字符串输入输出问题（1）字符打印1、打印*此类题的溯源为书P122 4.11（1），其他题都是它的拓展 for (i=1;i=n;i+) 一共要输出的行数 for(j=1;j=i;j+) 每行要打印的*数 printf(*); printf(n); a、*解题要点：此类题关键在于找到每行要打印的个数和行数的关系。此题j=i j=n-i+1b、* for(i=1;i=n;i+) 一共要输出的行数 for(j=1;j=n-i;j+) 控制空格数 printf(

33、); for(k=1;k=i;k+) 每行要打印的*数 printf(*); printf(n); c、 * * *解题要点：在出现空格的时候，在找到每行要打印的*个数和行数的关系后，还应找到空格和行数的关系，分不同的参数进行循环。此题k=i j=n-i j=i-1k=n-i+1d、* * for(i=1;i=n;i+) for(j=1;j=n-i;j+) printf( ); for(k=1;k=2*i-1;k+) printf(*); printf(n); * *e、 * * * for(i=1;i=n-1;i+) for(j=1;j=i;j+) printf( ); for(k=1;k=

34、2*(n-1-i)+1;k+) printf(*); printf(n); for(i=1;i=n;i+) for(j=1;j=n-i;j+) printf( ); for(k=1;k=2*i-1;k+) printf(*); printf(n); * * * * * *2、打印9*9乘法表解题要点：注意寻找行与列的规律。i*ji代表列j代表行for(i=1;i=9;i+) for(j=1;j=9;j+) printf(%-3d ,i*j); 注意输出格式的控制 printf(n); 3、九九乘法表1 2 3 4 5 6 7 8 92 4 6 8 10 12 14 16 183 6 9 12

35、15 18 21 24 27 9 18 27 36 45 54 63 72 814、杨晖三角形11 11 2 11 3 3 11 4 6 4 11 5 10 10 5 11 6 15 20 15 6 1（2）字符串打印问题for(i=1;i=7;i+) ai1=1; aii=1; for(i=3;i=7;i+) for(j=2;j=i-1;j+) aij=ai-1j-1+ai-1j; gets(a); puts(a); for(i=1;i0;j-) aj=aj-1; a0=t; for(k=0;k=a&ai=z) ai-=32; puts (a);3、逆序输出gets (a); c=strle

36、n(a); for(i=0;i=0;i-) 蓝色部分可以简写为绿色部分coutai-1;4、如输入：ab1 3,;z 输出：ab1注意点：1、= =2、while语句的使用处体会3、全面考虑问题 3,;zgets(a); while(a0= ) for(i=0;ai!=0;i+) ai=ai+1; for(i=0;ai!=0;i+)if(ai= &ai+1!= ) printf(n); else if(ai= &ai+1= ) for(k=i;ak!=0;k+) ak+1=ak+2; i-; elseprintf(%c,ai);5、输入3个字符串，按从小到大排序输出这3个字符串。 使用一个两维

37、数组贮存多个字符串char a8181;注意：如何使用一个两维数组贮存多个字符串 int i,j; for(i=0;i3;i+) gets(ai); for(i=0;i3;i+) puts(ai);注意：1、scanf(%d%s,&n,str) 其中%s为字符串格式2、逐个给字符串赋值的方法见书140页。 不可for（i=0；ai!=0;i+）3、stri=stri-A+10;4、pow函数5、任何进制转为十进制的方法6、输入一个整数n和一个字符串str，计算并输出n进制数str的值。 如输入：7 16则输出：13(16)7=(13)10如输入：16 3A则输出：58(3A)16=(58)10

38、#include stdio.h#include math.hmain() char str81; int n,i,s=0,t; clrscr(); scanf(%d%s,&n,str); for(i=0;stri!=0;i+) if(stri=A) stri=stri-A+10; else stri=stri-0; t=strlen(str); for(i=0;stri!=0;i+) s+=strt-i-1*pow(n,i); printf(%d,s);编写程序，将一个十进制正整数转换成十六进制数。 注意类比#include main()char a20;int x,i=0,j;clrscr

39、();scanf(%d,&x);while(x) if(x%16=10&x%16=0;j-)printf(%c,aj);printf(n);7、输入一个字符串，将其中的缩写形式展开，并输出展开后的该字符串。所谓展开缩写形式就是将其中由大小写字母或数字构成的形如a-f、U-Z、3-8 的形式展开成为 abcdef 、UVWXYZ 、345678，若出现f-a、A-7、9-5等形式则不予理睬。例如： 输入：qwe246e-hA-d$-%4-7A-Dz-xp-R4-0输出：qwr246efghA-d$-%4567ABCDz-xp-R4-0main() char a81; int i,c,s,k,t;

40、 gets(a); for(i=0;ai!=0;i+) if(ai=-) if(ai-1=A&ai+1=a&ai+1=0&ai+1i;k-)ak+c-2=ak;as-1+c-2+1=0; for(;i=t;i+) ai=ai-1+1; puts(a);补充：循环：求：a+aa+aaa+.的值#includevoid main()int a,n,i=1,sn=0,tn=0;coutinput a and nan;while(i=n)tn=tn+a;sn+=tn;a*=10;i+;coutthe answer is snendl;两个乒乓球队进行比赛，各出3人。甲队为A,B,C;已对是X,Y,Z；

41、已经抽签决定比赛名单。有人向队员大厅比赛的名单。A说他不和X比，C说他不和X,Z比。请编程序找出3对赛手的名单。#includevoid main()char i,j,k;for(i=X;i=Z;i+)for(j=X;j=Z;j+)if(i!=j)for (k=X;k=Z;k+)if(i!=k&j!=k)if(i!=X&k!=X&k!=Z)coutA-i B-j C-kendl;枚举口袋中有红，黄，蓝，白，黑5种颜色的球若干。每次从口袋中任意取出3歌，问得到3种不同颜色球的可能取法，输出每种排列的情况。#include#include /在C语言中不用加这句void main()enum co

42、lorred ,yellow ,blue,white, black;color pri;int i,j,k,n=0,loop;for(i=red;i=black;i+)for(j=red;j=black;j+)if(i!=j)for (k=red;k=black;k+)if(k!=i)&(k!=j)n+;coutsetw(3)n; /setw是输出格式的限定for(loop=1;loop=3;loop+)switch(loop)case 1:pri=color(i);break;case 2:pri=color(j);break;case 3:pri=color(k);break;default:break;switch(pri)case red:coutsetw(8)red;break;case yellow:coutsetw(8)yellow;break;case blue:coutsetw(8)blue;break;case white:coutsetw(8)white;break;case black:coutsetw(8)black;break;defa