annealing.cpp

模拟退火算法求解TSP问题,求解TSP问题的模拟退火算法

//模拟退火算法求解旅行商问题（TSP）
//问题描述:设有n个城市和距离矩阵D=[dij]，其中dij表示城市i到城市j的距离，
//i、j = 1、2 … n，则问题是要找出遍访每个城市恰好一次的一条回路，并使
//其路径长度为最短。

#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <time.h>

#define CITY_NUM  10            //城市个数，即问题规模
#define INITIAL_T 100           //控制参数初值
#define LENGTH_K  10000          //Mapkob链长

int distance[CITY_NUM][CITY_NUM]={0,8,14,22,30,30,25,20,15,12,  //距离矩阵
                                  8,0,19,10,15,24,20,15,10,10,
								  14,19,0,10,10,25,17,20,18,18,
                                  22,10,10,0,15,22,20,18,17,20,
                                  30,15,10,15,0,20,15,30,23,24, 
                                  30,24,25,22,20,0,16,21,19,30,
								  25,20,17,20,15,16,0,18,16,26,
								  20,15,20,18,30,21,18,0,17,22,
								  15,10,18,17,23,19,16,17,0,15,
								  12,10,18,20,24,30,26,22,15,0};
int pathlength=0; //路径长度
//二变换法
void GENERATE(int cnum,int *cha)
{
    int i=0;
	cha[0]=cha[1]=0;
	while((cha[0]==cha[1])||(cha[0]>cha[1])||cha[0]>=CITY_NUM||cha[1]>=CITY_NUM)         
	{
		cha[0]=(int)(rand()/(1024*3.2*(10/CITY_NUM)));  //产生随机数
		cha[1]=(int)(rand()/(1024*3.2*(10/CITY_NUM)));
    }
	cha[2]=(cha[0]-1+cnum)%cnum;
	cha[3]=(1+cha[1])%cnum;
}

//计算路径差
void CACULATE(int *c,int *pa,int &dff)
{
	int c0[4];
	for(int i=0;i<4;i++)
	      c0[i]=pa[c[i]];
	dff=0;
	if(c[1]-c[0]==CITY_NUM-1)  dff=0;    
	else  dff = distance[c0[0]][c0[3]]+distance[c0[1]][c0[2]]
		       -distance[c0[0]][c0[2]]-distance[c0[1]][c0[3]];
}

//判断是否接受新解
bool ACCEPT(float t,int dl)
{   
	if((dl<=0)||((dl/t<88)&&(dl>0)&&(exp(-dl/t)>(float)rand()/32768)))
	   	return true;
	else return false;
}

//接受新解，生成新的路径
void TWOCHAIN(int cnum,int *c,int *path_next)
{
     int i;
	 int u,v,w;
	 i=(1+(c[1]-c[0]+cnum)%cnum)/2;
	 for(int j=0;j<i;j++)
	 {
         u=(c[0]+j)%cnum;
		 v=(c[1]-j+cnum)%cnum;
         w=path_next[u]; path_next[u]=path_next[v]; path_next[v]=w;
	 }
}

//退火过程
void ANNEALING(int cnum,int lini,int sflag,float tini,float dt,int *path)
{
	int schange=0;                         //相同路径数 
    int flag,ntime=0;  
	int c[4];                              //二变换法中逆转元素的下标值
	int df;                                //路径差
	int bondtime=0;
//	FILE *fp;
	while(schange!=sflag)                  //若相等则算法停止
	{
		flag=0;
        for(int k=0;k<lini;k++)        
		{   
			GENERATE(cnum,c);              //采用二变换法产生一个新解
            CACULATE(c,path,df);           //计算路径长度差
			if(ACCEPT(tini,df))            //判断是否接受新解
			{
				TWOCHAIN(cnum,c,path);     //生成新的路径
				if(df==0) bondtime+=1;
				  else bondtime=0;
				pathlength=pathlength+df;  //计算新的路径长度
				flag=1;
			}
		
		}
	//	fp=fopen("result.txt","a");
	//	fprintf(fp,"第%d次迭代后 路径长度:%d 温度:%.10f  bondtime:%d\n",ntime++,pathlength,tini,bondtime);
	//  fclose(fp);
		tini=tini*dt;                       //降温，控制参数改变
		if(flag==0)  schange+=1;
		  else       schange=0;
		if(bondtime>=LENGTH_K*100)  schange=2;   //时间约束
	}
}

void main()
{
  int Cnum,Lini,Sflag;            
  float Tini,dt=(float)0.9;                      //衰减量
  //初始化过程
  Cnum=CITY_NUM;                                 //城市个数
  int path[CITY_NUM]={0,1,2,3,4,5,6,7,8,9};      //初始解
  clock_t start,finish;
  double duration;
  for(int i=0;i<CITY_NUM-1;i++)
  {
     pathlength+=distance[i][i+1];               //初始路径长度
  }
     pathlength+=distance[CITY_NUM-1][0];
  
  Lini=LENGTH_K;                                 //解的变换次数
  Tini=INITIAL_T;                                //初始控制参数
  Sflag=2;                 //停止准则，连续两个Mapkob链路无变化
  
  //调用退火过程开始退火
  start=clock();
  ANNEALING(Cnum,Lini,Sflag,Tini,dt,path);
  
  finish=clock();
  duration=(double)(finish-start)/CLOCKS_PER_SEC;
    
  //输出最后结果
  printf("Result :\n\t");
  for(int j=0;j<CITY_NUM;j++)
      printf("%d->",path[j]);
  printf("%d\n",path[0]);
  printf("The length of the path : %d\n",pathlength);
  printf("Using time: %.3f seconds\n", duration );

}