基于边缘重组的遗传算法求解TSP问题

0.引言

旅行商问题是一个经典的优化组合问题，它可以扩展到很多问题，如电路布线、输油管路铺设等，但是，由于TSP问题的可行解数目与城市数目N是成指数型增长的，是一个NP难问题，因而一般只能近似求解，遗传算法（GA）是求解该问题的较有效的方法之一，当然还有如粒子群算法，蚁群算法，神经网络算法等优化算法也可以进行求解。遗传算法是美国学者Holland根据自然界“物竞天择，适者生存”现象而提出的一种随机搜索算法，本文采用C语言来实现遗传算法解决TSP问题。

1.旅行商问题

旅行商问题可以具体描述为：已知n个城市之间的相互距离，现有一个推销员从某一个城市出发，必须遍访这n个城市，并且每个城市只能访问一次，最后又必须返回到出发城市，如何安排他对这些城市的访问次序，可使其旅行路线的总长度最短。用图论术语来表示，就是有一个图g=(v,e),其中v是点，e是边集，设d=(dij)是有顶点i和顶点j之间的距离所组成的距离矩阵，旅行商问题就是求出一条通过所有顶点且每个顶点只通过一次的最短距离的回路。若对与城市v={v1,v2,v3…vn}的一个访问顺序为t=(t1,t2,t3…,tn)，其中ti∈v(i=1,2,..n),且记tn+1=t1

2.遗传算法

2.1遗传算法

遗传算法是一种借鉴生物界自然选择和自然遗传机制的随机化搜索算法，由美国J.Holland教授提出，其主要特点是群体搜索策略和群体中个体之间的信息交换，搜索不依赖于梯度信息。它是一种全局化搜索算法，尤其适用于传统搜索算法难于解决的复杂和非线性问题。
选择(selection)算子、交叉(crossover)算子和变异(mutation)算子是遗传算法的3个主要操作算子。遗传算法中包含了如下5个基本要素:
(1)对参数进行编码；

(2)设定初始种群大小；
(3)适应度函数的设计；
(4)遗传操作设计；
(5)控制参数设定(包括种群大小、最大进化代数、交叉率、变异率等)。

2.2遗传算法的过程

遗传算法的基本过程是：
初始化群体。
计算群体上每个个体的适应度值
由个体适应度值所决定的某个规则选择将进入下一代个体。
按概率Pc进行交叉操作。
按概率Pm进行变异操作。
没有满足某种停止条件，则转第2步，否则进入第7步。
输出种群中适应度值最优的染色体作为问题的满意解或最优界。
停止条件有两种：一是完成了预先给定的进化代数则停止；二是种群中的最优个体在连续若干代没有改进或平均适应度在连续若干代基本没有改进时停止。
遗传算法过程图如图1：



图1：遗传算法过程框图

3.遗传算法C代码实现

遗传算法中控制参数如下：

#defineNUMCITIES 51 /*城市数量*/
#definePOPULATION_SIZE 10000 /*初始种群数量*/
#defineEVOLVING_POPULATION_SIZE 500 /*进化种群数量*/
#defineELITISM_PCT 0.1 /*精英种群比例*/
#defineNUMBER_OF_GENERATIONS 400 /*迭代次数*/
#defineMUTATION_RATE 0.05 /*变异概率*/
#defineCROSSOVER_RATE 0.9 /*交叉概率*/
#defineTOURNAMENT_SIZE 10 /*竞赛次数*/

3.1开始准备

先读入文件，读入51个城市坐标读入矩阵，并计算城市距离。详见附件中源代码函数：
LoadData与CreateMartix

int LoadData(char const *path)
{
int i = 0;
int j;
char line[100];
char buffTemp[5][10];
char split[] = " ";
char *token;
FILE *fp;
fp = fopen(path, "r");
if (fp == NULL)
{
perror(path);
exit(EXIT_FAILURE);
}
while(!feof(fp) && i<NUMCITIES && fgets(line,sizeof(line),fp))
{
j = 0;
token = (char *)strtok(line, split);
while(token != NULL)
{
strcpy(buffTemp[j], token);
++j;
token = (char *)strtok(NULL, split);
}
strcpy(data[i].name, buffTemp[0]); //printf("%s\n", data[i].name);
data[i].id = i; //printf("%d\n", i);
data[i].x = strtod(buffTemp[1], NULL); //printf("%f\n", data[i].x);
data[i].y = strtod(buffTemp[2], NULL); //printf("%f\n", data[i].y);
++i;
}
return i;
}

void CreateMartix(void)
{
int i,j;
double distance;
City city1,city2;// = (City *)malloc(sizeof(City));
//City *city2 = (City *)malloc(sizeof(City));
for (i = 0; i < NUMCITIES - 1; ++i)
{
martix[i][i] = 0.0;
city1 = data[i];
for (j = i + 1; j < NUMCITIES; ++j)
{
city2 = data[j];
distance = sqrt(pow((city1.x - city2.x), 2) + pow((city1.y - city2.y), 2));
martix[i][j] = distance;
martix[j][i] = distance;
}
}
martix[NUMCITIES-1][NUMCITIES-1] = 0.0;

}

3.2初始化种群

遗传算法对求解问题本身是一无所知的，这里初始化种群分由两部分组成：采用随机（ShuffleArray函数）生成初始化种群和通过最近距离函数（NearestNeighborTour）生成部分初始化种群，核心代码如下：

City *ShuffleArray(City *dest,const City *src, int n)//, unsigned int seed
{
assert(dest != NULL && src !=NULL);
int i,j, t=NUMCITIES;
City temp;
City *c = dest;
CitynCpy(dest, src, n);
//srand(time(NULL)+seed);
while(t--)
{
do{
i=rand()%n;//printf("i%d\n", i);//i = genrand_int32()%n;
j=rand()%n;//printf("j%d\n", j);//j = genrand_int32()%n;
}while(i==j);
temp = *(dest+i);
*(dest+i) = *(dest+j);
*(dest+j) = temp;
}
return c;
}

/*
**as a city and calc the nearest Neighbor tour
**return a Route
*/
Route *NearestNeighborTour(City const *c)
{
int i,j;
int id;
int index;
double minDist, val;
City cityList[NUMCITIES];
int  citiesVisited[NUMCITIES];
Route *r;
//City *temp;// = (City *)malloc(sizeof(City));
//temp = c;
City temp;
temp = *c;
for (i = 0; i < NUMCITIES; ++i)
{
cityList[i] = temp;
id = temp.id;
//citiesVisited = (int *)malloc(sizeof(int));
citiesVisited[i] = id;
minDist = DBL_MAX; //DBL_MAX defined in float.h
index = 0;
for (j = 0; j < NUMCITIES; ++j)
{
val = martix[id][j];
if (!Search(citiesVisited,j,i+1) && val < minDist)
{
minDist = val;
index = j;
}
}
temp = data[index];//data + index;
}
//free(temp);
r = (Route *)malloc(sizeof(Route));
if (r == NULL)
{
printf("malloc error\n");
exit(EXIT_FAILURE);
}
CitynCpy(r->cities, cityList, NUMCITIES);
r->length = CalcRouteLength(cityList);
return r;
}

static void CreatePopulation()
{
int i;
unsigned int randint;
Route *r;
Route *init_routes;//[POPULATION_SIZE];
City *cities;
init_routes = (Route *)malloc((POPULATION_SIZE+NUMCITIES)*sizeof(Route));
cities = (City *)malloc(sizeof(City)*NUMCITIES);
if (init_routes == NULL || cities == NULL)
{
printf("malloc error\n");
exit(EXIT_FAILURE);
}
CitynCpy(cities, data, NUMCITIES);
srand((unsigned)time(NULL));
for (i = 0; i < POPULATION_SIZE; ++i)
{
//randint = rand();//genrand_int32();
ShuffleArray((init_routes+i)->cities, cities, NUMCITIES);//, randint
(init_routes + i)->length = CalcRouteLength((init_routes+i)->cities);
}
for (i = POPULATION_SIZE; i < POPULATION_SIZE + NUMCITIES; ++i)
{
*(init_routes + i) = *(NearestNeighborTour(&data[i-POPULATION_SIZE]));
}
qsort(init_routes, POPULATION_SIZE+NUMCITIES, sizeof(Route), CmpArrayASC);
routes = (Route *)malloc(sizeof(Route)*EVOLVING_POPULATION_SIZE);
if (routes == NULL)
{
printf("Out of memory\n");
exit(EXIT_FAILURE);
}
RoutenCpy(routes, init_routes, EVOLVING_POPULATION_SIZE);
//printf("%f\n", routes->length);
free(init_routes);
free(cities);
}

3.3选择操作

选择的目的是为了从当前群体中选出优良的个体，本文使用竞赛机制进行选择，使它们有机会作为父代产生后代个体，核心代码（详见函数evolve）：

possibleParents = (Route *)malloc(sizeof(Route)*TOURNAMENT_SIZE);
for (i = 0; i < TOURNAMENT_SIZE; ++i)
{
index = rand() % EVOLVING_POPULATION_SIZE;//genrand_int32()%(EVOLVING_POPULATION_SIZE);
//printf("%d\n",index);
*(possibleParents+i) = *(routes+index);
}
//printf("%f\n", possibleParents->length);
qsort(possibleParents, TOURNAMENT_SIZE, sizeof(Route), CmpArrayASC);
possibleParents =(Route *)realloc(possibleParents,sizeof(Route)*2);

3.4交叉操作

许多生物的繁衍是通过染色体的交叉完成的，在遗传算法中使用这一概念，并把交叉作为遗传算法的一个操作算子，其实现过程是对选中用于繁殖下一代的个体，随机地选择两个个体的相同位置，按交叉概率CROSSOVER_RATE，在选择的位置实行交换，目的在于产生新的基因组合，即新的个体。本文的算法设计了两种：
一种是部分映射交叉算法但由于顺序编码交叉会产生不合法的编码，因此本文采取修复策略：部分映射交叉（PMX）。部分映射交叉(PMX)算子，整个交叉操作由两步来完成:首先对个体编码串进行常规的双点交叉操作。然后根据交叉区域内各个基因值之间的映射关系来修改交叉区域之外的各个基因座的基因值。交叉对象一般是符号编码表示的个体。（源代码太多，详见附件中函数PMXCrossover）。
另一种是改进的边缘重组算法，详见ERCrossover

/*
**Change Crossover to ERCrossOver
**Enhanced Edge Recombination (ER) Algorithm.
*/
void ERCrossover(Route *child, Route *parents)
{
assert(validRoute(parents,2)==TRUE);
assert(child != NULL);
int i, j, dadIndex, momIndex, cityId;
int tmprand;
int dadInitialCityConnctions, momInitialCityConnctions;
int nextCityId;
//City childCityArray[NUMCITIES];
//City *dadArray, *momArray;
City c;
City currentCity;
List unvisitedCityIds;
Map *edgeMap;
edgeMap =(Map *)calloc(NUMCITIES, sizeof(Map));
Route *dad;
Route *mom;
dad	= parents;
mom = parents + 1;
// Create the edge map.
// <CityId : List of neighboring CityId's in mom and dad>.
for (dadIndex = 0; dadIndex < NUMCITIES; ++dadIndex)
{
c = dad->cities[dadIndex];//dadArray[dadIndex];
cityId = c.id;
// Get location of current city in dad--in mom.
momIndex = 0;
for (j = 0; j < NUMCITIES; ++j)
{
if (cityId == mom->cities[j].id)
{
momIndex = j;
break;
}
}
(edgeMap+dadIndex)->key = cityId;
(edgeMap+dadIndex)->edge = getEdges(dad->cities, mom->cities, dadIndex, momIndex);
}
//City[] childCityArray = new City[NUMCITIES];
//ArrayList<Integer> unvisitedCityIds = new ArrayList<Integer>(NUMCITIES);
unvisitedCityIds = InitEdges();
for (i = 0; i < NUMCITIES; ++i)
EdgesAddTail(&unvisitedCityIds, i);

// Pick start city...
dadInitialCityConnctions = GetMap(edgeMap, NUMCITIES, dad->cities[0].id)->size;
momInitialCityConnctions = GetMap(edgeMap, NUMCITIES, mom->cities[0].id)->size;
if (dadInitialCityConnctions >= momInitialCityConnctions)
currentCity = dad->cities[0];
else
currentCity = mom->cities[0];

child->cities[0] = currentCity;//childCityArray[0] = currentCity;

i = 1;
EdgesRemove(&unvisitedCityIds, currentCity.id);
ListRemove(edgeMap, NUMCITIES, currentCity.id);

while (unvisitedCityIds->size > 0)
{
if (GetMap(edgeMap, NUMCITIES, currentCity.id)->size > 0) // Step 4.
currentCity = PickNextCity(edgeMap, currentCity.id);
else
{ // Step 5.
if ( unvisitedCityIds->size == 1 )
nextCityId = unvisitedCityIds->next->data;
else
{
tmprand = rand() % unvisitedCityIds->size;
///printf("%d\n", tmprand);
nextCityId = EdgesnNode(unvisitedCityIds, tmprand)->data;
}
currentCity = data[nextCityId];
}
child->cities[i] = currentCity;//childCityArray[i] = currentCity;
++i;

//printf("%d\n", unvisitedCityIds->size);
assert(unvisitedCityIds->size > 0);

EdgesRemove(&unvisitedCityIds, currentCity.id);
ListRemove(edgeMap, NUMCITIES, currentCity.id);
}

//CitynCpy(child->cities, childCityArray, NUMCITIES);
child->length = CalcRouteLength(child->cities);
//printf("%f\n", child->length);
assert(validRoute(child,1)==TRUE);
free(edgeMap);
//return child;
}

3.5变异操作

是按一定概率随机改变某个个体的基因值，以变异概率MUTATION_RATE，随机地在染色体上选取两个位置，交换基因的位置如代码：

/*
**Randomly swap two citites in the route
*/
void Mutate(Route *chromosome)
{
int c1, c2;
City tmp;
//srand((unsigned)time(NULL));
do{
c1 = rand()%NUMCITIES;//genrand_int32()%NUMCITIES;
c2 = rand()%NUMCITIES;//genrand_int32()%NUMCITIES;
}while(c1==c2);
//printf("%d&&%d\n",c1,c2);
tmp = chromosome->cities[c1];
chromosome->cities[c1] = chromosome->cities[c2];
chromosome->cities[c2] = tmp;
chromosome->length = CalcRouteLength(chromosome->cities);
assert(validRoute(chromosome,1)==TRUE);
}

4实验分析

实验硬件环境为普通笔记本电脑，内存2GB，CPU频率2.0GHz。软件环境为WindowsXP
Professional SP3操作系统，TinyCC。实验对象是51个城市的旅行商问题。
得出实验结果为：
PMX：


图二：PMX实验结果路径
ER：



图三：ER实验结果路径

5总结

本文采用ANSIC实现遗传算法求解TSP问题，并根据实验结果进行了分析，遗传算法是一种智能优化算法，它的实现有些关键点，一是串的编码方式，本质就是问题编码，串长度及编码形式对算法收敛影响极大；二是适应函数的确定，这是选择的基础；三是自身参数的设定，其中重要的是群体大小，最大迭代次数，交叉概率和变异概率，通过实验我们可以看到最大迭代次数对问题求解的精度有影响，交叉概率和变异概率的设定对问题的收敛速度和求解精度都有极大的影响，目前很多研究都是根据具体的领域问题，改进交叉算子，变异算子，寻找最优的参数设定来提高算法收敛速度和保证最优解的得到，对算子的改进和参数值的设定这是将来的研究工作。

参考文献：

[1]http://en.wikipedia.org/wiki/Edge_recombination_operator
[2]http://en.wikipedia.org/wiki/Category:Genetic_algorithms