基于遗传算法的TSP问题求解(C)
TSP问题：
TSP(Travelling salesman problem): 译作“旅行商问题”， 一个商人由于业务的需要，要到n个城市，每个城市之间都有一条路径和其他所有的城市相连。现在要求从一个城市出发，穿越所有其他所有的城市，再回到出发的城市。 出于成本的考虑，要求商人走的路径的长短最短。问能否找到这样的一条路径？
这是个经典的NP-complete问题。 时间复杂度为θ(n!)。 随着城市的数量规模增大，在有限的时间内得不到问题的最优解。 我们只能退而求其次，在有限的时间内(或可接受的时间内)找到一个近似最优解，而且这个近似最优解和最优解的误差我们也是可以接受的。 出于这样的考虑，为求解这类问题的启发式，元启发式算法，演化算法营运而生。
随着研究的深入，TSP问题已经演化成好多版本。本文的C程序对于对称和非对称的版本都适用。
遗传算法：
遗传算法(Genetic Algorithm)，也称进化算法，是依据生物进化的过程而提出的一种启发式算法，由美国的J.Holland于1975年首次提出。其主要特点是依据生物进化中的“优胜劣汰”或者“竞争”法作为问题的解的选择依据。 直接对结构对象进行操作，不存在求导和函数连续性的限定；具有内在的隐并行性和更好的全局寻优能力；采用概率化的寻优方法，能自动获取和指导优化的搜索空间，自适应地调整搜索方向，不需要确定的规则。遗传算法的这些性质，已被人们广泛地应用于组合优化、机器学习、信号处理、自适应控制和人工生命等领域。它是现代有关智能计算中的关键技术之一
遗传算法的基本原理：
首先根据问题产生多个初始可行解()，然后从初始解中选择诺干优异的个体(问题的解)进行各种操作(交叉，变异)用以产生新的后代。再从新的后代中选择优异的个体进行相应的操作产生它们的后代，如此不断循环，直到迭代的次数达到了预先设定的值或者多次迭代以后产生的最优异的个体(最优解)的质量并没有明显的提高，就可以停止迭代，这时的最优个体(最优解)最为问题的解。
从上面的原理我们可以知晓该算法主要涉及的步骤如图 1所示：
- 编码
在解实现初始化之前，如何表示一个解即编码是一个很重要的问题。 一般的编码方式有:
- 基于二进制编码: 作为遗传算法传统的编码方式。对于TSP问题，经过交叉后可能得不到可行解，需要修补操作。
- 基于矩阵的编码: 需要更多的额外空间，而且随着种群规模的增加，存储空间剧增和前期处理工作任务很庞大。后续操作也比较复杂。
- 基于邻近的编码: 采用链表的方式存储,但是变异,交叉操作复杂
- 基于格雷码方法：传统二进制编码的一种改进，容易实现交叉，变异操作，但是对于该问题不是最优的
- 基于符号编码：对于TSP问题，我们直接使用路径来表示一个染色体。即使用一个整数数组，数组长度为TSP城市的数量，数组中存储各个城市编号，从下标为0开始逐个遍历数组元素形成的一个序列即为路径(对于要回到原点的要求，为了方便表示，在这里不作考虑，只是在计算路径长度的时候会添加相应的长度)。
- 形成初始解
采用随机的方式产生诺干个(种群规模)的序列，即产生符合城市编号的随机数存储在数组中，数组中的元素包含所有的城市编号，而且没有重复的编码。数组的个数为种群规模。
- 选择
在形成一定数量的可行解（染色体）后，需要对这些父代的染色体进行筛选。根据生物遗传的基因的优胜劣汰原则，在筛选染色体的我们也会秉承精英保留原则，使得优异的基因有更多的机会得到遗传。
适应度函数： 这里我们选择路径长度的倒数来作为解的适应度
在这个问题中，我们选择“轮盘赌”算法来筛选染色体。
基本原理：计算每个染色体(路径)的长度的倒数，再得到所有路径倒数之和，用每条路径的倒数除以所有所有路径倒数之和即为这条路径被选中的概率(路径越短，概率越高)。
- 交叉
这里我们选择交替位置交叉法（Alternating Position Crossover，APX）来对一对染色体进行交叉操作，其基本原理如下图所示
左边为父代的两个染色体，右边为子代染色体。 将左上的数组第一个元素放入右上数组的第一位置中，再转移到左下数组第一个元素，查看右上数组是否已经包含了该元素，如果未包含将其插入右上数组中，否则插入右下数组中。接着从左上数组的第二个元素开始，到左下第二个元素，和前次同意的判断操作。如此类推直到右边两个数组都被填满了为止。
交叉概率：交叉概率对于解的收敛速度有着重要的影响。一般选择0.6-1。
- 变异
生物的进化，除了遗传父母的基因，还有自身基因有一定的概率突变。基于这个原理，变异操作在一定的概率上是作用于染色体自身的。
变异概率：一定的概率师兄自身基因的改变
在这个问题中，我们选择位置倒换法，即染色体上随机的产生两个位置上数值互换。
- 终止条件
一般有两种方式停止交叉，变异的操作。一，预先设定迭代次数。二，多次跌代后，解的质量得不到一定要求的提高，或者解达到要求的质量，这时都可以停止迭代。这个问题上我们选择第一种。
基于TSP问题的遗传算法代码如下：
1 #include <stdio.h> 2 #include <tchar.h> 3 #include <math.h> 4 #include <stdlib.h> 5 #include <time.h> 6 7 8 int scale; /* 种群规模 */ 9 int cityNum; /* 城市数量 */ 10 int *pathlength; /* 存储种群每个个体路径长度 */ 11 double *cumPropa; /* 存储个体累积概率 */ 12 int bestlength; /* 最佳路径长度 */ 13 int *bestpath; /* 最佳路径 */ 14 float pc; /* 交叉概率 */ 15 float pm; /* 变异概率 */ 16 int count; /* 变异次数 */ 17 int MAX_Gene; /* 迭代次数 */ 18 19 /* 函数声明 */ 20 void APCrossover(int **,int ,int ,int ); 21 void copy(int *,int **,int ,int ); 22 void cumDistance(int **,int **, int *,int ,int ); 23 void cumulatePro(int **,double *,int *,int ); 24 void copytoBestPath(int *,int **,int ,int ); 25 void copy(int *,int **,int ,int ); 26 int *creatArray1(int ); 27 double *creatArraydoub(int ); 28 int **creatArray2(int , int ); 29 void CrossAndMutation(int **,int ); 30 void cumDistance(int **,int **, int *,int ,int); 31 void mutation(int **,int); 32 void Initialize(int **, int , int ); 33 void readfile(int **, int , int ); 34 void rouletteAlgo(int **,double *, int **,int ); 35 void NcopyO(int **,int **); 36 37 //创建一个整形的二维整数数组 38 int **creatArray2(int scale, int cityNum) 39 { 40 int i; 41 42 int **ptemp; 43 44 ptemp=(int **)malloc(sizeof(int *)*scale); 45 for(i=0;i<scale;i++) 46 ptemp[i]=(int *)malloc(sizeof(int)*cityNum); 47 48 return ptemp; 49 } 50 //创建一个整形的一维数组 51 int *creatArray1(int scale) 52 { 53 int *tempcreate; 54 tempcreate=(int *)malloc(sizeof(int)*scale); 55 return tempcreate; 56 } 57 //创建一个双精度型的一维数组 58 double *creatArraydoub(int scale) 59 { 60 double *tempcreate; 61 tempcreate=(double *)malloc(sizeof(double)*scale); 62 return tempcreate; 63 } 64 65 // 二维数组拷贝 66 void NcopyO(int **newgeneration,int **oldgeneratoin) 67 { 68 int i,j; 69 for(i=0;i<scale;i++) 70 for(j=0;j<cityNum;j++) 71 oldgeneratoin[i][j]=newgeneration[i][j]; 72 73 } 74 75 76 77 //计算一次迭代各可行解的路径长度 78 void cumDistance(int **oldgeneration,int **cityDistance, int *pathlength,int cityNum,int scale) 79 { 80 int i,j,k; 81 int length=0; 82 int min=88888; 83 int minp; 84 85 for(i=0;i<scale;i++) 86 { 87 length=0; 88 for(k=0;k<cityNum-1;k++) 89 { 90 j=k+1; 91 length+=cityDistance[oldgeneration[i][k]][oldgeneration[i][j]]; 92 } 93 //题目要求回到原点，所以加上回到原点的距离 94 pathlength[i]=length+cityDistance[oldgeneration[i][cityNum-1]][oldgeneration[i][0]]; 95 96 if(pathlength[i]<min) 97 { 98 min=pathlength[i]; 99 minp=i; 100 } 101 } 102 if(min<bestlength) 103 { 104 bestlength=min; 105 //将每代最优解直接加入下一代中，即“精英保留”原则 106 copytoBestPath(bestpath,oldgeneration,minp,cityNum); 107 } 108 109 } 110 111 112 113 114 void copytoBestPath(int *bestpath,int **oldGeneration,int position,int cityNum) 115 { 116 int i; 117 for(i=0;i<cityNum;i++) 118 bestpath[i]=oldGeneration[position][i]; 119 } 120 121 122 123 //计算一次迭代中各个可行解作为交叉操作的累积概率 124 125 void cumulatePro(int **generation,double *cumPropa,int *pathlength,int scale) 126 { 127 int i; 128 double sumlength=0; 129 130 for(i=0;i<scale;i++) 131 sumlength+=1/(double)pathlength[i]; 132 133 cumPropa[0]=(1/(double)pathlength[0])/sumlength; 134 135 for(i=1;i<scale;i++) 136 { 137 cumPropa[i]=(1/(double)pathlength[i])/sumlength+cumPropa[i-1]; 138 } 139 } 140 141 142 //初始化话生产相应规模的可行解 143 void Initialize(int **geneSolution, int scale, int cityNum) 144 { 145 int i,j,k; 146 int randnum; 147 bool exist; 148 149 srand(time(NULL)); 150 151 for(i=0;i<scale;i++) 152 for(j=0;j<cityNum;j++) 153 { 154 exist=true; 155 while(exist==true){ 156 //注意：rand()/Rand_MAX 结果只能是0， 应该先进行类型转换 157 randnum=(int)(((double)rand()/RAND_MAX)*cityNum); 158 for(k=0;k<j;k++) 159 if(geneSolution[i][k]==randnum) 160 break; 161 if(k==j) 162 exist=false; 163 } 164 geneSolution[i][j]=randnum; 165 } 166 167 } 168 169 //读取文件中的各相邻点的距离信息 170 void readfile(int **cityDistance, int scale, int cityNum) 171 { 172 int i,j; 173 FILE *fp; 174 errno_t err; 175 176 err=fopen_s(&fp,"data.txt","r"); 177 178 if(err!=0) 179 { 180 printf("The file can not be found!\n"); 181 } 182 else 183 { 184 185 for(i=0;i<scale;i++) 186 { 187 for(j=0;j<cityNum;j++) 188 { 189 fscanf_s(fp,"%d ",&cityDistance[i][j]); 190 } 191 } 192 193 194 195 fclose(fp); 196 } 197 198 } 199 200 //拷贝一条路径 201 void copy(int *oldGeneration,int **newGeneration,int position,int cityNum) 202 { 203 int i; 204 for(i=0;i<cityNum;i++) 205 newGeneration[position][i]=oldGeneration[i]; 206 } 207 208 //使用轮盘赌算法选择交叉的对象 209 void rouletteAlgo(int **oldGeneration,double *cumPropa, int **newGeneration,int scale) 210 { 211 int i,j; 212 double randNum; 213 214 for(i=0;i<scale-1;i++) 215 { 216 randNum=(double)rand()/RAND_MAX; 217 if(cumPropa[0]>=randNum) 218 copy(oldGeneration[0],newGeneration,i,scale); 219 else 220 { 221 for(j=0;j<scale;j++) 222 if(randNum>cumPropa[j] && randNum<=cumPropa[j]) 223 break; 224 copy(oldGeneration[i],newGeneration,i,scale); 225 } 226 } 227 228 copy(bestpath,newGeneration,scale-1,scale); 229 230 } 231 232 //交叉操作：交替位置交叉法（Alternating Position Crossover，APX） 233 void APCrossover(int **newgeneration,int p1,int p2,int cityNum) 234 { 235 int i; 236 int s1=1; 237 int s2=0; 238 239 int *tempArray1; 240 int *tempArray2; 241 242 tempArray1=creatArray1(cityNum); 243 tempArray2=creatArray1(cityNum); 244 245 for(i=0;i<cityNum;i++) 246 { 247 tempArray1[i]=newgeneration[p1][i]; 248 tempArray2[i]=newgeneration[p2][i]; 249 } 250 251 int m,n; 252 m=1; 253 n=0; 254 255 while(s1<10 || s2<10) 256 { 257 for(i=0;i<s1;i++) 258 if(tempArray1[m]==newgeneration[p1][i]) 259 break; 260 if(i==s1) 261 { 262 newgeneration[p1][s1]=tempArray1[m]; 263 m++; 264 s1++; 265 } 266 else{ 267 newgeneration[p2][s2]=tempArray1[m]; 268 m++; 269 s2++; 270 } 271 272 273 for(i=0;i<s1;i++) 274 if(tempArray2[n]==newgeneration[p2][i]) 275 break; 276 if(i==s1) 277 { 278 newgeneration[p1][s1]=tempArray2[n]; 279 n++; 280 s1++; 281 } 282 else{ 283 newgeneration[p2][s2]=tempArray2[n]; 284 n++; 285 s2++; 286 } 287 } 288 289 } 290 291 //变异操作 292 void mutation(int **newgeneration,int p1) 293 { 294 int rand1,rand2; 295 int temp; 296 int i; 297 298 srand(time(NULL)); 299 300 for(i=0;i<count;i++) 301 { 302 303 rand1=(int)(((double)rand()/RAND_MAX)*cityNum); 304 rand2=(int)(((double)rand()/RAND_MAX)*cityNum); 305 while(rand1==rand2) 306 { 307 rand2=(int)(((double)rand()/RAND_MAX)*cityNum); 308 } 309 310 temp=newgeneration[p1][rand1]; 311 newgeneration[p1][rand1]=newgeneration[p1][rand2]; 312 newgeneration[p1][rand2]=temp; 313 } 314 315 } 316 317 //对一代群体进行交叉变异操作 318 void CrossAndMutation(int **newgeneration,int cityNum) 319 { 320 float rand1,rand2; 321 int k; 322 for(k=0;k<cityNum;k=k+2) 323 { 324 srand(time(NULL)); 325 rand1=(float)rand()/RAND_MAX; 326 327 if(rand1>pc) 328 { 329 APCrossover(newgeneration,k,k+1,cityNum); 330 } 331 else 332 { 333 rand2=(float)rand()/RAND_MAX; 334 if(rand2>pm) 335 { 336 mutation(newgeneration,k); 337 } 338 rand2=(float)rand()/RAND_MAX; 339 if(rand2>pm) 340 { 341 mutation(newgeneration,k+1); 342 } 343 344 } 345 } 346 347 } 348 349 int main() 350 { 351 int **a; 352 int **oldGeneration; 353 int **newGeneration; 354 int i,j; 355 356 MAX_Gene=20; 357 cityNum=10; 358 scale=10; 359 bestlength=888888; 360 pc=0.6; 361 pm=0.5; 362 count=4; 363 int m; 364 365 366 a=creatArray2(scale,cityNum); 367 pathlength=creatArray1(scale); 368 cumPropa=creatArraydoub(scale); 369 bestpath=creatArray1(cityNum); 370 371 readfile(a,scale,cityNum); 372 373 oldGeneration=creatArray2(scale,cityNum); 374 newGeneration=creatArray2(scale,cityNum); 375 376 Initialize(newGeneration,scale,cityNum); 377 378 for(m=0;m<MAX_Gene;m++) 379 { 380 NcopyO(newGeneration,oldGeneration); 381 cumDistance(oldGeneration,a,pathlength,scale,cityNum); 382 cumulatePro(oldGeneration,cumPropa,pathlength,scale); 383 rouletteAlgo(oldGeneration,cumPropa,newGeneration,scale); 384 CrossAndMutation(newGeneration,cityNum); 385 386 } 387 388 printf("The best path is :\n"); 389 390 391 for(i=0;i<cityNum;i++) 392 { 393 printf("%d ",bestpath[i]); 394 } 395 printf("\n"); 396 397 printf("The minmum length is %d\n",bestlength); 398 399 return 0; 400 }