利用遗传算法求解TSP问题

转载: http://blog.csdn.net/lalor/article/details/7704011

1.遗传算法

        遗传算法是受大自然的启发,模拟生物在自然环境中的遗传和进化过程而形成的一种自适应、具有全局优化能力的随机搜索算法。

自然界的进化包括3个原则:

(1)适者生存原则,这意味着适应能力强的物种,会在残酷的竞争中生存下来,而适应能力差的物种会逐渐地消亡。

(2) 两性繁殖。这意味着种群中性别不同的个体,生活在一起,产生新的个体。

(3) 变异。 由于环境的变化,新物种的出现,以及不同物种的交互都会引起种群的变异。

        遗传算法的思路是通过从给定一个初始群体出发,利用选择算子、杂交算子以及变异算子来模拟自然进化的三种原则,逐步改进种群,一步步逼近最优解,以达到求解最优华问题的目的。


GA算法的计算步骤:




        记住遗传算法的过程很重要,首先是初始化一群解,然后再在这些解中选择较优的一部分,将选择的这一部分解进行交叉,且以一定概率变异,(交叉一般能得到比当前解更好的解,而变异则很可能让结果变差,所以变异的概率一般不是很大,但是这样有助于我们跳出局部最优)。交叉变异以后进行群体更新,对于TSP问题,群体更新时保存这一次迭代产生的最好的解,然后继续进行下一次迭代,直到满足终结条件为止。


GA的算法过程:


        初始化t,t代表while循环已经迭代了多少次。其中f(pop(t))是指这个解的适应度,对于TSP问题,适应度就是它的代价,第8行是按一定的概率选择较有的解。第9行以Pc概率进行交叉,第10行以Pm概率进行变异。

2. 问题建模

        遗传算法其实很简单,就是初始化一群解,然后选择这一群里面较优的解,在较优的解里面,让其中的个体交叉,使得交叉后得到更好的解,再按一定概率进行变异,希望变异能跳出局部最优。对于遗传算法求解TSP问题,最难的地方在于问题建模,刚开始根本不知道如何用遗传算法来求解旅行商问题,如何交叉,如何变异。

        首先初始化一群解,可以通过C++提供的库函数来产生一个城市的随机排列,每一个排列代表一个解random_shuffle(temp.path, temp.path + nCities)。然后以一定概率选择较优的解,选择的方法有很多,我们不一定非要按照上面伪代码的方式来选择,比如我们希望每次保存当前这群解中的前60%,则我们可以按解的适应度排序,然后取前60%的解,对于后40%的解,我们可以用前40%的解去覆盖它,则前40%的解就有2个副本,只要我们交叉的时候不要让相同的两个副本交叉就行了,因为相同的两个解交叉,不会让结果变得更好。

        变异也很简单,只需要在一个解中随机的选择两个城市,然后交换它们即可。注意,变异的概率不宜太大。

        最难的部分是交叉,我们要如何用两个解得到一个更好的解?这就是交叉,让一代比一代强,我们才可能慢慢接近最优解。交叉的方法有很多,可以参考http://blog.csdn.net/xuyuanfan77/article/details/6726477


源码中采用类似于三交换启发交叉(THGA),我把它改成了二交叉的。

三交换启发交叉方法的基本思想如下:

选3个参加交配的染色体作为父代,以8个城市为例来说明这一过程,其中dij由前面的表1给出,父代染色体为

A = 3 2 1 4 8 7 6 5

B = 2 4 6 8 1 3 5 7

C = 8 7 5 6 4 3 2 1

SUM1=42,SUM2=40,SUM3=46(SUM1,SUM2,SUM3分别为这3种排法所走的距离总和数).

随机选出初始城市j=1,Sj=3右转动,使3成为3父代的第1位置.

A = 3 2 1 4 8 7 6 5

B = 3 5 7 2 4 6 8 1

C = 3 2 1 8 7 5 6 4

由于d(3,2)>d(3,5),所以有:

A = × 5 2 1 4 8 7 6

B = × 5 7 2 4 6 8 1

C = × 5 6 4 2 1 8 7

由此规则计算可得:

O = 3 5 7 6 8 4 2 1

我们本来是3个不同的解,现在得到了一个比三个解都优的解,总不能让原来的三个解都等于现在的这个局部最优解吧,这样不利于下次交叉,我们可以用如下的方法改变另外两个解的路径:

  1. Rotate(q.path, nCities, rand() % nCities);  
        

        上行代码执行以后,它的代价还是和原来一样的,路径也是一样,只是起点变了,这样有什么好处呢?有利于下次交叉的时候,原来的两个相同代价,不同路径的解能和其他解交叉出不同的结果,这样有利于找到更好的解。


3. 代码实现

  1. /* 
  2.  * * 
  3.  * * 
  4.  * * Copyright(c) Computer Science Department of XiaMen University  
  5.  * * 
  6.  * * 
  7.  * * Authored by lalor on: 2012年 06月 29日 星期五 23:49:57 CST 
  8.  * * 
  9.  * * 
  10.  * * Email: mingxinglai(at)gmail.com 
  11.  * * 
  12.  * * 
  13.  * * @desc: 
  14.  * * 
  15.  * * 
  16.  * * @history 
  17.  * * 
  18.  * * 
  19.  * * 说明:本程序使用的测试数据来自权威的benchmark,其最优解是1211.数据保存在source.txt 
  20.  * * 本例的测试数据来自http://www.iwr.uni-heidelberg.de/groups/comopt/software/TSPLIB95/tsp/ 
  21.  * * rat99.tsp.gz 
  22.  * * 数据如下 
  23.  * ﹡格式:(城市编号,横坐标,纵坐标) 
  24.  1  6  4 
  25.   2 15 15 
  26.   3 24 18 
  27.   4 33 12 
  28.   5 48 12 
  29.   6 57 14 
  30.   7 67 10 
  31.   8 77 10 
  32.   9 86 15 
  33.  10  6 21 
  34.  11 17 26 
  35.  12 23 25 
  36.  13 32 35 
  37.  14 43 23 
  38.  15 55 35 
  39.  16 65 36 
  40.  17 78 39 
  41.  18 87 35 
  42.  19  3 53 
  43.  20 12 44 
  44.  21 28 53 
  45.  22 33 49 
  46.  23 47 46 
  47.  24 55 52 
  48.  25 64 50 
  49.  26 71 57 
  50.  27 87 57 
  51.  28  4 72 
  52.  29 15 78 
  53.  30 22 70 
  54.  31 34 71 
  55.  32 42 79 
  56.  33 54 77 
  57.  34 66 79 
  58.  35 78 67 
  59.  36 87 73 
  60.  37  7 81 
  61.  38 17 95 
  62.  39 26 98 
  63.  40 32 97 
  64.  41 43 88 
  65.  42 57 89 
  66.  43 64 85 
  67.  44 78 83 
  68.  45 83 98 
  69.  46  5 109 
  70.  47 13 111 
  71.  48 25 102 
  72.  49 38 119 
  73.  50 46 107 
  74.  51 58 110 
  75.  52 67 110 
  76.  53 74 113 
  77.  54 88 110 
  78.  55  2 124 
  79.  56 17 134 
  80.  57 23 129 
  81.  58 36 131 
  82.  59 42 137 
  83.  60 53 123 
  84.  61 63 135 
  85.  62 72 134 
  86.  63 87 129 
  87.  64  2 146 
  88.  65 16 147 
  89.  66 25 153 
  90.  67 38 155 
  91.  68 42 158 
  92.  69 57 154 
  93.  70 66 151 
  94.  71 73 151 
  95.  72 86 149 
  96.  73  5 177 
  97.  74 13 162 
  98.  75 25 169 
  99.  76 35 177 
  100.  77 46 172 
  101.  78 54 166 
  102.  79 65 174 
  103.  80 73 161 
  104.  81 86 162 
  105.  82  2 195 
  106.  83 14 196 
  107.  84 28 189 
  108.  85 38 187 
  109.  86 46 195 
  110.  87 57 194 
  111.  88 63 188 
  112.  89 77 193 
  113.  90 85 194 
  114.  91  8 211 
  115.  92 12 217 
  116.  93 22 210 
  117.  94 34 216 
  118.  95 47 203 
  119.  96 58 213 
  120.  97 66 206 
  121.  98 78 210 
  122.  99 85 204 
  123.  ﹡* 
  124.  * * 
  125.  * */  
  126.   
  127. #include <iostream>  
  128. #include <string.h>  
  129. #include <fstream>  
  130. #include <iterator>  
  131. #include <algorithm>  
  132. #include <limits.h>  
  133. #include <math.h>  
  134. #include <stdlib.h>  
  135.   
  136. using namespace std;  
  137.   
  138. const int nCities = 99; //No. of node  
  139. //const double PC = 0.9; //交叉概率  
  140. double PM = 0.1; //变异概率  
  141. double PS = 0.8;//保留概率  
  142. int GEN_MAX = 50; //最大代数  
  143. const int UNIT_NUM = 5000; //群体规模为50  
  144. double length_table[nCities][nCities];//distance  
  145.   
  146.   
  147.   
  148. //城市  
  149. struct node  
  150. {  
  151.     int num;//城市的编号  
  152.     int x;//横坐标  
  153.     int y;//纵坐标  
  154.   
  155. }nodes[nCities];  
  156.   
  157.   
  158.   
  159. struct unit  
  160. {  
  161.     double length;//代价,总长度  
  162.     int path[nCities];//路径  
  163.   
  164.     bool operator < ( const struct unit &other) const //用于群体的排序  
  165.     {  
  166.         return length < other.length;  
  167.     }  
  168.   
  169. };  
  170.   
  171.   
  172. //群体规模(群体规模是指有 UNIT_NUM 个不同的解,而bestone 用于保存最好的一个解)  
  173. struct unit group[UNIT_NUM];  
  174.   
  175.   
  176. //保存最好的一个解  
  177. unit bestone = {INT_MAX, {0} };  
  178.   
  179.   
  180. // create matrix to storage the Distance each city  
  181. void init_dis();   
  182.   
  183. //计算 unit 中的length, 也就是群体的一个个体(一个解)的长度  
  184. void CalCulate_length(unit &p);  
  185.   
  186. //查找id (代表城市) 在当前解中的位置,用于两个解的交叉  
  187. int search_son(unit &p, int id);  
  188.   
  189. //打印一个解  
  190. void print( unit &p);  
  191.   
  192. //初始化群体,由C++ 中的 random_shuff 产生一个随机排列  
  193. void Initial_group( unit group[]);  
  194.   
  195. //开始进化,在本函数中执行群体中个体的交叉和变异  
  196. void Evolution_group(unit group[]);  
  197.   
  198. //变异,随机的选择一个群体,然后随机选择两个点,交换它们的位置  
  199. void Varation_group(unit group[]);  
  200.   
  201. //交叉  
  202. void Cross_group( unit &p, unit &q);  
  203.   
  204. int main(int argc, char* argv[])  
  205. {  
  206.     srand(time(NULL));  
  207.   
  208.     init_dis();  
  209.   
  210.     //初始化种群  
  211.     Initial_group( group );  
  212.   
  213.   
  214.     //种群进化:选择,交叉,变异  
  215.     Evolution_group( group );  
  216.   
  217.     cout << "变异概率PM = " << PM << endl;  
  218.     cout << "保留概率PS = " << PS << endl;  
  219.     cout << "最大代数 = " << GEN_MAX << endl;  
  220.     cout << "群体规模 = " << UNIT_NUM  << endl;  
  221.     cout << "代价是: = " << bestone.length << endl;  
  222.   
  223.     print(bestone);  
  224. }  
  225.   
  226.   
  227.   
  228.   
  229.   
  230. void init_dis() // create matrix to storage the Distance each city  
  231. {  
  232.     int i, j;  
  233.   
  234.     ifstream in("source.txt");  
  235.   
  236.     for (i = 0; i < nCities; i++)  
  237.     {  
  238.         in >> nodes[i].num >> nodes[i].x >> nodes[i].y;  
  239.     }  
  240.   
  241.     for (i = 0; i < nCities; i++)  
  242.     {  
  243.         length_table[i][i] = (double)INT_MAX;  
  244.         for (j = i + 1; j < nCities; j++)  
  245.         {  
  246.             length_table [i][j] = length_table[j][i] =sqrt(   
  247.                     (nodes[i].x - nodes[j].x) * (nodes[i].x - nodes[j].x) +  
  248.                     (nodes[i].y - nodes[j].y) * (nodes[i].y - nodes[j].y) );  
  249.         }  
  250.   
  251.     }  
  252.   
  253. }  
  254.   
  255.   
  256. void CalCulate_length(unit &p)  
  257. {  
  258.     int j = 0;  
  259.   
  260.     p.length = 0;  
  261.   
  262.     for (j = 1; j < nCities; j++)   
  263.     {  
  264.         p.length += length_table[ p.path[j-1] ][ p.path[j] ];  
  265.     }  
  266.     p.length += length_table[ p.path[nCities - 1] ][ p.path[0] ];  
  267.   
  268. }  
  269.   
  270.   
  271.   
  272. void print( unit &p)  
  273. {  
  274.     int i;  
  275.     cout << "代价是:" << p.length << endl << "路径是:";  
  276. //  for (i = 0; i < nCities; i++)   
  277. //  {  
  278. //      cout << p.path[i] << " ";  
  279. //  }  
  280.   
  281.   
  282.     copy(p.path, p.path + nCities, ostream_iterator<int>(cout, " -> "));  
  283.     cout << p.path[0] << endl;  
  284. }  
  285.   
  286.   
  287. //函数对象,给generate 调用  
  288. class GenByOne  
  289. {  
  290.     public:  
  291.         GenByOne (int _seed = -1): seed(_seed)  
  292.         {  
  293.           
  294.         }  
  295.   
  296.         int operator() ()  
  297.         {  
  298.             return seed += 1;  
  299.         }  
  300.   
  301.   
  302.     private:  
  303.         int seed;  
  304. };  
  305.   
  306.   
  307. //随机产生 UNIT_NUM 个解空间  
  308. void Initial_group( unit group[])  
  309. {  
  310.     int i, j;  
  311.     unit temp;  
  312.   
  313.   
  314.     //1, 2, 3, 4 ...... nCities  
  315.     generate(temp.path, temp.path + nCities, GenByOne(0));  
  316.   
  317.     // 产生 UNIT_NUM 个解,也就是群体  
  318.     for (i = 0; i < UNIT_NUM; i++)   
  319.     {  
  320.         //产生一个随机排列,也就是初始化一个解  
  321.         random_shuffle(temp.path, temp.path + nCities);  
  322.         memcpy(&group[i], &temp, sizeof(temp));  
  323.         CalCulate_length(group[i]);  
  324.     }  
  325. }  
  326.   
  327. void Evolution_group(unit group[])  
  328. {  
  329.     int i, j;  
  330.     int n = GEN_MAX;  
  331.     int num1, num2;  
  332.   
  333.     //以PS 的概率选择前 num2 个解,抛弃其后的num1 个解。  
  334.     num1 = UNIT_NUM * ( 1 - PS);  
  335.     num2 = UNIT_NUM * PS;  
  336.   
  337.     //迭代几次,即繁衍多少代  
  338.     while (n-- ) //循环GEN-MAX次  
  339.     {  
  340.         //选择部分优秀的种群  
  341.         sort(group, group + UNIT_NUM);  
  342.         if (group[0].length < bestone.length)   
  343.         {  
  344.             memcpy(&bestone, &group[0], sizeof(unit));  
  345.         }  
  346.   
  347.   
  348.         for (j = 0; j <=  num1 - 1; j++)   
  349.         {  
  350.             memcpy(&group[ num2 + j], &group[j], sizeof(unit));  
  351.         }  
  352.   
  353.         //交叉  
  354.         for (j = 0; j < UNIT_NUM / 2; j+= 1)   
  355.         {  
  356.             Cross_group(group[j], group[ UNIT_NUM - j -1]);  
  357.         }  
  358.   
  359.         //变异  
  360.         Varation_group(group);  
  361.     }  
  362.   
  363.     //保存已找最好的解  
  364.     sort(group, group + UNIT_NUM);  
  365.     if (group[0].length < bestone.length)   
  366.     {  
  367.         memcpy(&bestone, &group[0], sizeof(unit));  
  368.     }  
  369. }  
  370.   
  371. void Varation_group(unit group[])  
  372. {  
  373.     int i, j, k;  
  374.     double temp;  
  375.     //变异的数量,即,群体中的个体以PM的概率变异,变异概率不宜太大  
  376.     int num = UNIT_NUM * PM;  
  377.   
  378.     while (num--)   
  379.     {  
  380.   
  381.         //确定发生变异的个体  
  382.         k = rand() % UNIT_NUM;  
  383.   
  384.         //确定发生变异的位  
  385.         i = rand() % nCities;  
  386.         j = rand() % nCities;  
  387.   
  388.         //exchange  
  389.         temp  = group[k].path[i];  
  390.         group[k].path[i] = group[k].path[j];   
  391.         group[k].path[j] = temp;  
  392.   
  393.         CalCulate_length(group[k]);  
  394.     }  
  395. }  
  396.   
  397. int Search_son( int path[], int len, int city)  
  398. {  
  399.     if (city <= 0 || city > nCities)   
  400.     {  
  401.         cout << "city outfiled, city = " << city << endl;  
  402.         return -1;  
  403.     }  
  404.   
  405.     int i = 0;  
  406.     for (i = 0; i < len; i++)   
  407.     {  
  408.         if (path[i] == city)   
  409.         {  
  410.             return i;  
  411.         }  
  412.     }  
  413.     return -1;  
  414. }  
  415.   
  416. //reverse a array  
  417. //it's a auxiliary function for Rotate()   
  418. void Reverse(int path[], int b, int e)  
  419. {  
  420.     int temp;  
  421.   
  422.     while (b < e)   
  423.     {  
  424.         temp = path[b];   
  425.         path[b] = path[e];  
  426.         path[e] = temp;  
  427.   
  428.         b++;  
  429.         e--;  
  430.     }  
  431. }  
  432.   
  433.   
  434. //旋转 m 位  
  435. void Rotate(int path[],int len, int m)  
  436. {  
  437.     if( m < 0 )  
  438.     {  
  439.         return;  
  440.     }  
  441.     if (m > len)   
  442.     {  
  443.         m %= len;  
  444.     }  
  445.   
  446.     Reverse(path, 0, m -1);  
  447.     Reverse(path, m, len -1);  
  448.     Reverse(path, 0, len -1);  
  449. }  
  450.   
  451.   
  452. void Cross_group( unit &p, unit &q)  
  453. {  
  454.     int i = 0, j ,k;  
  455.     int pos1, pos2;  
  456.     int len = nCities;  
  457.     int first;  
  458.   
  459.     double len1 = length_table[p.path[0] ][ p.path[1] ];  
  460.     double len2 = length_table[q.path[0] ][ q.path[1] ];  
  461.   
  462.     if (len1 <= len2)   
  463.     {  
  464.         first = p.path[0];  
  465.     }  
  466.     else  
  467.     {  
  468.         first = q.path[0];  
  469.     }  
  470.   
  471.   
  472.   
  473.     pos1 = Search_son( p.path + i, len, first);  
  474.     pos2 = Search_son( q.path + i, len, first);  
  475.   
  476.     Rotate(p.path + i, len, pos1);  
  477.     Rotate(q.path + i, len, pos2);  
  478.   
  479.     while ( --len > 1)   
  480.     {  
  481.         i++;  
  482.         double span1  = length_table[ p.path[i - 1] ][ p.path[i] ];  
  483.         double span2  = length_table[ q.path[i - 1] ][ q.path[i] ];  
  484.   
  485.         if ( span1 <= span2 )  
  486.         {  
  487.             pos2 = Search_son( q.path + i, len, p.path[i]);  
  488.             Rotate(q.path + i, len, pos2);  
  489.         }  
  490.         else  
  491.         {  
  492.             pos1 = Search_son( p.path + i, len, q.path[i]);  
  493.             Rotate(p.path + i, len, pos1);  
  494.         }  
  495.     }  
  496.   
  497.     Rotate(q.path, nCities, rand() % nCities);  
  498.   
  499.     CalCulate_length(p);  
  500.     CalCulate_length(q);  
  501. }  

4.参考资料:

[1] http://blog.csdn.net/xuyuanfan77/article/details/6726477

[2] 算法设计与分析(高级教程),国防工业出版社

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值