sa.cpp

利用模拟退火算法解决TSP问题.很经典的一个算法.使用中国144个城市测试,在很短的时间内找到的比较好的解.

#include "stdafx.h"

#include "sa.h"
#include <math.h>
#include <time.h>
#include <stdio.h>
#include <string.h>
#include <stdlib.h>

CSaTsp::CSaTsp(char *datafile)
{
	FILE *file=fopen(datafile, "r");
	fscanf(file, "%d", &m_nCity);

	m_curPath = new Path;
	m_curPath->Citys = new int[m_nCity];

	m_newPath = new Path;
	m_newPath->Citys = new int[m_nCity];

	m_globalPath = new Path;
	m_globalPath->Citys = new int[m_nCity];

	m_xyPosition = new Position[m_nCity];

	for(int i=0; i<m_nCity; i++)
	{
		fscanf(file, "%d%d", &(m_xyPosition[i].x), &(m_xyPosition[i].y));
	}
	fclose(file);
	file = NULL;

	srand( (unsigned)time(NULL) );
}
CSaTsp::~CSaTsp()
{
	delete[] m_curPath->Citys;
	m_curPath->Citys = NULL;
	delete m_curPath;

	delete[] m_newPath->Citys;
	m_newPath->Citys = NULL;
	delete m_newPath;

	delete[] m_globalPath->Citys;
	m_globalPath->Citys = NULL;
	delete m_globalPath;

	delete[] m_xyPosition;
	m_xyPosition = NULL;
}

//////////////////////////////////////////////////////////////////////////
void CSaTsp::Init()
{
	//初始化温度
	m_Temperature = 150.0;	
	//初始化路径
	for(int i=0; i<m_nCity; i++)
	{
		m_curPath->Citys[i] = i;
	}

	//随机排列
	for(int i=0; i< m_nCity; i++)		
	{
		int index, temp;

		//生成i到(n-1)的随机数
		index = i + rand()%( m_nCity - i );	

		//交换
		temp = m_curPath->Citys[i];
		m_curPath->Citys[i] = m_curPath->Citys[index];
		m_curPath->Citys[index] = temp;
	}

	//计算路径长度
	m_curPath->Length = 0.0;
	double x = m_xyPosition[m_curPath->Citys[0]].x-m_xyPosition[m_curPath->Citys[m_nCity-1]].x;
	double y = m_xyPosition[m_curPath->Citys[0]].y-m_xyPosition[m_curPath->Citys[m_nCity-1]].y;
	m_curPath->Length += sqrt(x*x +y*y);
	for( int i=1; i<m_nCity; i++)
	{
		x = m_xyPosition[m_curPath->Citys[i]].x-m_xyPosition[m_curPath->Citys[i-1]].x;
		y = m_xyPosition[m_curPath->Citys[i]].y-m_xyPosition[m_curPath->Citys[i-1]].y;
		m_curPath->Length += sqrt(x*x +y*y);
	}
	m_globalPath->Length = m_curPath->Length;
	memcpy(m_globalPath->Citys, m_curPath->Citys, m_nCity*sizeof(int));
	
}

//返回f(s)-f(s')
double CSaTsp::NewPath()
{
	memcpy(m_newPath->Citys, m_curPath->Citys, m_nCity*sizeof(int));

	int c1,c2 = rand()%m_nCity;
	do 
	{
		c1 = rand()%m_nCity;
	} while(c1==c2);

	int temp;
	if(c1>c2)
	{
		temp = c1;
		c1 = c2;
		c2 = temp;
	}
	for(int i=0; i<=c2-c1; i++)
	{
		m_newPath->Citys[c1+i]=m_curPath->Citys[c2-i];
	}

	int cp,cn;
	double delta = 0.0, x, y;

	cp = (c1+m_nCity-1)%m_nCity;
	cn = (c2+1)%m_nCity;

	x = m_xyPosition[m_curPath->Citys[c1]].x-m_xyPosition[m_curPath->Citys[cp]].x;
	y = m_xyPosition[m_curPath->Citys[c1]].y-m_xyPosition[m_curPath->Citys[cp]].y;
	delta += sqrt(x*x+y*y);
	x = m_xyPosition[m_newPath->Citys[c1]].x-m_xyPosition[m_newPath->Citys[cp]].x;
	y = m_xyPosition[m_newPath->Citys[c1]].y-m_xyPosition[m_newPath->Citys[cp]].y;
	delta -= sqrt(x*x+y*y);
	x = m_xyPosition[m_curPath->Citys[c2]].x-m_xyPosition[m_curPath->Citys[cn]].x;
	y = m_xyPosition[m_curPath->Citys[c2]].y-m_xyPosition[m_curPath->Citys[cn]].y;
	delta += sqrt(x*x+y*y);
	x = m_xyPosition[m_newPath->Citys[c2]].x-m_xyPosition[m_newPath->Citys[cn]].x;
	y = m_xyPosition[m_newPath->Citys[c2]].y-m_xyPosition[m_newPath->Citys[cn]].y;
	delta -= sqrt(x*x+y*y);

	m_newPath->Length = m_curPath->Length-delta;
	return delta;
}
/*
double CSaTsp::NewPath()
{
	memcpy(m_newPath->Citys, m_curPath->Citys, m_nCity*sizeof(int));

	int c1,c2 = rand()%m_nCity;
	do 
	{
		c1 = rand()%m_nCity;
	} while(c1==c2);
	int temp=m_newPath->Citys[c1];
	m_newPath->Citys[c1] = m_newPath->Citys[c2];
	m_newPath->Citys[c2] = temp;

	int cp,cn;
	double delta = 0.0, x, y;
	
	cp = (c1+m_nCity-1)%m_nCity;
	cn = (c1+1)%m_nCity;
	x = m_xyPosition[m_curPath->Citys[c1]].x-m_xyPosition[m_curPath->Citys[cp]].x;
	y = m_xyPosition[m_curPath->Citys[c1]].y-m_xyPosition[m_curPath->Citys[cp]].y;
	delta += sqrt(x*x+y*y);
	x = m_xyPosition[m_curPath->Citys[c1]].x-m_xyPosition[m_curPath->Citys[cn]].x;
	y = m_xyPosition[m_curPath->Citys[c1]].y-m_xyPosition[m_curPath->Citys[cn]].y;
	delta += sqrt(x*x+y*y);

	x = m_xyPosition[m_newPath->Citys[c1]].x-m_xyPosition[m_newPath->Citys[cp]].x;
	y = m_xyPosition[m_newPath->Citys[c1]].y-m_xyPosition[m_newPath->Citys[cp]].y;
	delta -= sqrt(x*x+y*y);
	x = m_xyPosition[m_newPath->Citys[c1]].x-m_xyPosition[m_newPath->Citys[cn]].x;
	y = m_xyPosition[m_newPath->Citys[c1]].y-m_xyPosition[m_newPath->Citys[cn]].y;
	delta -= sqrt(x*x+y*y);

	//////////////////////////////////////////////////////////////////////////
	cp = (c2+m_nCity-1)%m_nCity;
	cn = (c2+1)%m_nCity;
	x = m_xyPosition[m_curPath->Citys[c2]].x-m_xyPosition[m_curPath->Citys[cp]].x;
	y = m_xyPosition[m_curPath->Citys[c2]].y-m_xyPosition[m_curPath->Citys[cp]].y;
	delta += sqrt(x*x+y*y);
	x = m_xyPosition[m_curPath->Citys[c2]].x-m_xyPosition[m_curPath->Citys[cn]].x;
	y = m_xyPosition[m_curPath->Citys[c2]].y-m_xyPosition[m_curPath->Citys[cn]].y;
	delta += sqrt(x*x+y*y);

	x = m_xyPosition[m_newPath->Citys[c2]].x-m_xyPosition[m_newPath->Citys[cp]].x;
	y = m_xyPosition[m_newPath->Citys[c2]].y-m_xyPosition[m_newPath->Citys[cp]].y;
	delta -= sqrt(x*x+y*y);
	x = m_xyPosition[m_newPath->Citys[c2]].x-m_xyPosition[m_newPath->Citys[cn]].x;
	y = m_xyPosition[m_newPath->Citys[c2]].y-m_xyPosition[m_newPath->Citys[cn]].y;
	delta -= sqrt(x*x+y*y);

	m_newPath->Length = m_curPath->Length-delta;
	return delta;
}
*/
int CSaTsp::GetNumCitys()
{
	return m_nCity;
}
void CSaTsp::UpdateBestPath()
{
	Path *temp	= m_curPath;
	m_curPath	= m_newPath;
	m_newPath	= temp;

	if( m_curPath->Length < m_globalPath->Length)
	{
		m_globalPath->Length = m_curPath->Length;
		memcpy(m_globalPath->Citys, m_curPath->Citys, m_nCity*sizeof(int));
	}
}
Position* CSaTsp::GetPosition()
{
	return m_xyPosition;
}
void CSaTsp::GetBestPath(Path *path)
{
	int Threshold = 0;
	double Length=0.0;

	Init();
	while( Threshold<4 )	//比较简单的条件
	{
		for(int i=0; i<20000; i++)
		{
			double delta = NewPath();
			if( delta>=0 )	//原来路径的长度不小于新路径的长度
				UpdateBestPath();
			else
			{
				double p = (rand()%100001)/100000.0;
				double e = exp( delta/m_Temperature );
				if( e>p )	//按概率接受
					UpdateBestPath();
			}
		}
		if ( m_curPath->Length != Length )
		{
			Length = m_curPath->Length;
			Threshold = 0;
		}
		else
		{
			Threshold++;
		}
		m_Temperature *= 0.9;
	}
	path->Length = m_globalPath->Length;
	memcpy(path->Citys, m_globalPath->Citys, m_nCity*sizeof(int));
}