tspsa.cpp

使用模拟退火算法解决TSP问题

//author: yinhui
//Date  : 2007.06.15

//TSPSA.cpp

#include "StdAfx.h"
#include "TSPSA.h"
#include <math.h>
#include <algorithm>

#define MAX_SAME 40

void Annealing::Init()
{
	srand((unsigned)time(NULL));
	ComputeTotalDistace();
	CreateInitSolution(m_vecRouter);
}

void Annealing::ObtainParams(float startT, float endT, float decRate, 
							 int malkovLen, VecCity vCitys, int nCityNum)
{
	m_fStaTemper = startT;
	m_fEndTemper = endT;
	m_fDecRatio = decRate;
	m_nMalkovLength = malkovLen;
	m_vecCitys = vCitys;
	m_nCityNum = nCityNum;
}

void Annealing::CreateInitSolution(vecCityRouter& initRouter)
{
	for (int i=1; i<=m_nCityNum; i++)
	{
		initRouter.push_back(i);
		//m_vecRouter.push_back(i);
	}
	std::random_shuffle(initRouter.begin(), initRouter.end());
}

vecCityRouter Annealing::NextSolution()
{
	int fSwap, eSwap;
	vecCityRouter tempRouter;
	tempRouter = m_vecRouter;

	while (1)
	{
		fSwap = (rand()%m_nCityNum)-1;
		if (fSwap>=0)
			break;
	}

	while (1)
	{
		eSwap = (rand()%m_nCityNum)-1;
		if ((eSwap!=fSwap) && (eSwap>=0))
		{
			int temp = tempRouter[fSwap];
			int xx = tempRouter[eSwap];
			tempRouter[fSwap] = tempRouter[eSwap];
			tempRouter[eSwap] = temp;
			break;
		}
	}
	return tempRouter;
}

void Annealing::ComputeTotalDistace()
{
	std::vector<City>::iterator fCityIter;
	std::vector<City>::iterator eCityIter;

	Distance tempDist;

	for (fCityIter=m_vecCitys.begin(); fCityIter!=m_vecCitys.end(); fCityIter++)
	{
		for (eCityIter=m_vecCitys.begin(); eCityIter!=m_vecCitys.end(); eCityIter++)
		{
			ComputeDistance(*fCityIter, *eCityIter, tempDist);
			m_vecDistance.push_back(tempDist);
		}
	}
}

double Annealing::ComputeDistance(City& fCity, City& eCity, Distance& cityDist)
{
	double tempDist;

	float fCityX = fCity.m_corCity.m_fX;
	float fCityY = fCity.m_corCity.m_fY;
	float eCityX = eCity.m_corCity.m_fX;
	float eCityY = eCity.m_corCity.m_fY;

	tempDist = sqrt((fCityX-eCityX)*(fCityX-eCityX) +
		            (fCityY-eCityY)*(fCityY-eCityY));

	cityDist.m_nFromCity = fCity.m_nIndex;
	cityDist.m_nEndCity = eCity.m_nIndex;
	cityDist.m_fDist = (float)tempDist;

	return tempDist;
}

double Annealing::FindDistance(int fIndex, int eIndex)
{
	VecDistance::iterator tempIter;

	for (tempIter=m_vecDistance.begin(); tempIter!=m_vecDistance.end(); tempIter++)
	{
		if ((tempIter->m_nFromCity==fIndex) && (tempIter->m_nEndCity==eIndex))
		{
			return tempIter->m_fDist;
		}
	}
	return 0.0;
}

double Annealing::ComputeRouterDistance(vecCityRouter& curCityRou)
{
	vecCityRouter::iterator tempRouterIter;
	double allDist=0.0;

	for (tempRouterIter=curCityRou.begin(); tempRouterIter!=(curCityRou.end()-1); tempRouterIter++)
	{
		allDist += FindDistance(*tempRouterIter, *(tempRouterIter+1));
	}

	allDist += FindDistance(*(curCityRou.end()), *(curCityRou.begin()));
	return allDist;
}

float Annealing::GetNowDistance()
{
	float dist;
	dist = (float)ComputeRouterDistance(m_vecRouter);
	return dist;
}

float Annealing::DoOneSearch(float t)
{
	double currDist = ComputeRouterDistance(m_vecRouter);
	vecCityRouter nextRouter = NextSolution();
	double nextDist = ComputeRouterDistance(nextRouter);

	double random = (double)rand()/(double)RAND_MAX;
	if ((nextDist-currDist) < 1.0e-6)
	{
		m_vecRouter = nextRouter;
	}
	/*else if (exp((currDist-nextDist)/t) > random)
	{
		m_vecRouter = nextRouter;
	}*/
		
	return (float)nextDist;
}

void Annealing::TSPAnneal(CListBox *pLB)
{
	pLB->ResetContent();
	CString strDisplay, strValue;
	strDisplay = "城市个数：";
	strValue.Format("%d\n", m_nCityNum);
	strDisplay += strValue;
	pLB->InsertString(0,strDisplay);

	int count = 0;    //统计相同解出现的次数
	float currT = m_fStaTemper;    //当前的温度
	double nextDist;
	double currDist;
	double bestDist;
	double countDist;    //记录最优解出现次数用

	currDist = ComputeRouterDistance(m_vecRouter);
	bestDist= currDist;
	countDist = currDist;

	strDisplay = "初始距离：";
	strValue.Format("%f\n", bestDist);
	strDisplay += strValue;
	pLB->InsertString(0, strDisplay);

	vecCityRouter nextRouter;

	int innerSame = 0;    //内循环出现不改进情况的计数器

	while ((currT>m_fEndTemper) && (count<MAX_SAME))
	{
		for (int i=0; i<m_nMalkovLength; i++)
		{
	        nextRouter = NextSolution();
	        nextDist = ComputeRouterDistance(nextRouter);

			double random = (double)rand()/(double)(1000*RAND_MAX);
			if (nextDist<=currDist)
			{
				m_vecRouter = nextRouter;
				currDist = nextDist;
			}else if (exp((currDist-nextDist)/currT) > random)
			{
				m_vecRouter = nextRouter;
				currDist = nextDist;
			}
				
			if (currDist>bestDist)
			{
				innerSame++;
			}
			else
			{
				innerSame = 0;
				bestDist = currDist;
			}
			

			//if (innerSame > (int)(0.95*m_nMalkovLength))
		//	{
		//		break;
		//	}
			
		}

		strDisplay = "当前距离：";
	    strValue.Format("%f\n", bestDist);
	    strDisplay += strValue;
	    pLB->InsertString(0, strDisplay);

		strDisplay = "当前温度：";
	    strValue.Format("%f\n", currT);
	    strDisplay += strValue;
	    pLB->InsertString(0, strDisplay);
		
		if (abs(bestDist-countDist) < 1.0e-6)
		{
			count++;
		}
		else
		{
			countDist = bestDist;
			count = 0;
		}
		
		currT *= m_fDecRatio;
	}
}