gasystem.h

蚁群算法的C++实现,让读者更清晰的了解蚁群算法的完整思路

#ifndef GATSP_H
#define GATSP_H

#include <vector>
#include <algorithm>
#include "../MemDC.h"
#include "../utils.h"
#include "../Chart.h"

#include "GAConstants.h"
extern CAutoFont gFont;
using namespace std;

namespace GA{
	extern CGAConstants gGAConstants;
	

	//-------------------------------------------------
	// GENOME 
	//-------------------------------------------------
	class CGenome
	{

	public:
		vector<int>		vecCityPath;	//path
		double			  _Fitness;        //gene fitness

		//-------------------------------------------------
		// Default Constructor 
		//-------------------------------------------------
		CGenome():_Fitness(0){}

		//-------------------------------------------------
		// Contructor 
		//-------------------------------------------------
		CGenome(int nc): _Fitness(0)
		{
			vecCityPath = CreateRandomTour(nc);
		}

		//-------------------------------------------------
		// CreateRandomTour 
		//-------------------------------------------------
		vector<int>	CreateRandomTour(int &limit)
		{
			vector<int> vecPerm;

			for (int i=0; i<limit; i++)
			{
				//we use limit-1 because we want ints numbered from zero
				int NextPossibleNumber = RandomInteger(0, limit-1);

				while(TestNumber(vecPerm, NextPossibleNumber))
				{
					NextPossibleNumber = RandomInteger(0, limit-1);
				}
				vecPerm.push_back(NextPossibleNumber);
			}
			return vecPerm;
		}


		//-------------------------------------------------
		// TestNumber 
		//-------------------------------------------------
		bool		TestNumber(const vector<int> &vec, const int &number)
		{
			for (int i=0; (i<(int)vec.size()); ++i)
			{
				if (vec[i] == number)
				{
					return true;
				}
			}
			return false;
		}

	};

	//-------------------------------------------------
	// CLASS CGASystem 
	//-------------------------------------------------
	class CGASystem
	{
	public:
		double	MUTATION_RATE;
		double	CROSSOVER_RATE;	
        int MAX_CITIES;
		int  POPULATION_SIZE;	
		double	_TotalFitness;
		double	_bestSoFar;
		double	_worstSoFar;	
		double		_BestPossible;
		double _AverageWorst;
		double _ElapsedTime;		
		int _ChromosomeLength;	
		int _FittestGene;	
		int _Generation;	
		bool _ProblemSolved;
		vector<CGenome> _Population;
		vector<CPoint>	_Cities;

		char buffer[1000];
		CPen* oldPen;
		CPen blackPen;
		CPen bluePen;

		CChart* _graph;
		//-------------------------------------------------
		// Constructor 
		//-------------------------------------------------
		CGASystem(int NumCities=30,
			double	  mut_rat=0.2,
			double	  cross_rat=0.75,
			int		  pop_size=40					
			):
		MAX_CITIES(NumCities),
				MUTATION_RATE(mut_rat),
			CROSSOVER_RATE(cross_rat),
			POPULATION_SIZE(pop_size),
			_FittestGene(0),
			_Generation(0),
			_bestSoFar(CGlobals::INFINITY),
			_BestPossible(CGlobals::INFINITY),
			_worstSoFar(0),
			_ChromosomeLength(NumCities),
			_ProblemSolved(false),
			_ElapsedTime(0)
		{
			bluePen.CreatePen(PS_SOLID,1,CColor::blue);
			blackPen.CreatePen(PS_SOLID,1,CColor::black);			
		}

		//-------------------------------------------------
		// Destructor 
		//-------------------------------------------------
		~CGASystem()
		{

		}

		void Init(CRect rect)
		{
			//_Cities.clear();
			int width = rect.right - rect.left;
			int height = rect.bottom-rect.top;
			int radius;
			int offset=10;

			if(width < height)
			{
				radius = width-offset;
			}
			else
			{
				radius=height-offset;
			}

			//CPoint origin(radius/ 2, radius / 2);
			CPoint midPoint;
			//midPoint.x = rect.left + width/2 ;
			//midPoint.y = rect.top + height/2 ;
			midPoint.x = rect.left + width/2+offset;
			midPoint.y = rect.top + height/2 ;

			//calculate angle division between adjacent cities.
			double SegmentSize = 2 * (CGlobals::PI) / MAX_CITIES;
			double angle = 0;

			vector<CPoint> vecCities;

			while (angle < 2 * CGlobals::PI)
			{
				CPoint tempCity;
				tempCity.x = static_cast<long>((radius/2 )* sin(angle) + midPoint.x);
				tempCity.y = static_cast<long>((radius/2) * cos(angle) + midPoint.y);
				_Cities.push_back(tempCity);

				angle += SegmentSize;
			}
			CalculateBestPossibleRoute();
			_AverageWorst=DistanceA_to_B(_Cities[0],_Cities[MAX_CITIES/2])*MAX_CITIES;
			CreateStartingPopulation();
		}

		void Resize(CRect rect)
		{
			_Cities.clear();
			Init(rect);
			//CalculateBestPossibleRoute();
			//CreateStartingPopulation();
		}

		void CreateStartingPopulation()
		{
			//make sure the vector of genomes is empty before we
			//start
			_Population.clear();

			//create a new population of random genomes
			for (int i=0; i<POPULATION_SIZE; ++i)
			{
				_Population.push_back(CGenome(_ChromosomeLength));
			}

			//make sure variables are initialised correctly
			_Generation	 = 0;
			_bestSoFar = CGlobals::INFINITY;
			_worstSoFar=0;
			CalculateBestPossibleRoute();
			_FittestGene = 0;
			_ProblemSolved			 = false;

		}

		CGenome& RouletteWheelSelection()
		{
			double fSlice	= RandomFloat() * _TotalFitness;
			double cfTotal	= 0.0;
			int	SelectedGenome = 0;
			for (int i=0; i<POPULATION_SIZE; ++i)
			{
				cfTotal += _Population[i]._Fitness;
				if (cfTotal > fSlice) 	
				{
					SelectedGenome = i;
					break;
				}
			}

			return _Population[SelectedGenome];
		}




		void Mutate(vector<int> &chromo)
		{

			//return dependent upon mutation rate
			if (RandomFloat() > 	MUTATION_RATE) return;

			//choose first gene
			int pos1 = RandomInteger(0, static_cast<int>(chromo.size()-1));

			//choose second
			int pos2 = pos1;

			while (pos1 == pos2)
			{
				pos2 = RandomInteger(0, static_cast<int>(chromo.size()-1));
			}

			//swap their positions
			swap(chromo[pos1], chromo[pos2]);
		}

		//-------------------------------------------------
		// Crossover
		//-------------------------------------------------
		void	Crossover(const vector<int>	&mother, 	const vector<int>	&father, 
								vector<int>	      &baby1, vector<int>       &baby2)
		{
			baby1 = mother;
			baby2 = father;

			//just return dependent on the crossover rate or if the
			//chromosomes are the same.
			if ( (RandomFloat() > CROSSOVER_RATE) || (mother == father)) 
			{
				return;
			}

			//first we choose a section of the chromosome
			int beg = RandomInteger(0, static_cast<int>(mother.size()-2));

			int end = beg;

			//find an end
			while (end <= beg)
			{
				end = RandomInteger(0, static_cast<int>(mother.size()-1));
			}

			//now we iterate through the matched pairs of genes from beg
			//to end swapping the places in each child
			vector<int>::iterator posGene1, posGene2;

			for (int pos = beg; pos < end+1; ++pos)
			{
				//these are the genes we want to swap
				int gene1 = mother[pos];
				int gene2 = father[pos];

				if (gene1 != gene2)
				{
					//find and swap them in baby1
					posGene1 = find(baby1.begin(), baby1.end(), gene1);
					posGene2 = find(baby1.begin(), baby1.end(), gene2);

					swap(*posGene1, *posGene2);

					//and in baby2
					posGene1 = find(baby2.begin(), baby2.end(), gene1);
					posGene2 = find(baby2.begin(), baby2.end(), gene2);

					swap(*posGene1, *posGene2);
				}

			}//next pair
		}	

		void CalculatePopulationsFitness()
		{
			//for each chromo
			for (int i=0; i<POPULATION_SIZE; ++i)
			{

				//calculate the tourlength for each chromosome
				double TourLength = 	GetTourLength(_Population[i].vecCityPath);
				_Population[i]._Fitness = TourLength;

				//keep a track of the shortest route found each generation
				if (TourLength <_bestSoFar)
				{
					_bestSoFar = TourLength;
					_FittestGene = i;
				}

				//keep a track of the worst tour each generation
				if (TourLength > _worstSoFar)
				{
					_worstSoFar = TourLength;
				}

			}//next chromo

			//Now we have calculated all the tour lengths we can assign
			//the fitness scores
			for (int i=0; i<POPULATION_SIZE; ++i)
			{
				_Population[i]._Fitness = _worstSoFar - _Population[i]._Fitness;
			}

		}

		double	DistanceA_to_B(const CPoint &city1, const CPoint &city2)
		{
			double xDist = city1.x - city2.x;
			double yDist = city1.y - city2.y;

			return sqrt(xDist*xDist + yDist*yDist);
		}

		double GetTourLength(const vector<int> &route)
		{
			double TotalDistance = 0;

			for (int i=0; i<static_cast<int>(route.size()-1); ++i)
			{
				int city1 = route[i];
				int city2 = route[i+1];

				TotalDistance += DistanceA_to_B(_Cities[city1], _Cities[city2]);
			}

			//don't forget this is a closed loop so we need to add in the distance 
			//from the last city visited back to the first
			int last  = route[route.size()-1];
			int first = route[0];

			TotalDistance += DistanceA_to_B(_Cities[last], _Cities[first]);

			return TotalDistance;
		}


		void CalculateBestPossibleRoute()
		{
			_BestPossible = 0;

			for (int city=0; city<static_cast<int>(_Cities.size()-1); ++city)
			{
				_BestPossible += DistanceA_to_B(_Cities[city], _Cities[city+1]);

				//add in a small amount to cover any precision errors we may have made
				//m_dBestPossibleRoute += gGAConstants.EPSILON;
			}

			//add in the distance from the last to the first
			_BestPossible += DistanceA_to_B(_Cities[_Cities.size()-1], _Cities[0]);
		}




		//-------------------------------------------------
		// RESET 
		//-------------------------------------------------
		void	Reset()
		{
			//just make the shortest route some excessively large distance
			_bestSoFar	= CGlobals::INFINITY;
			_worstSoFar		= 0;
			_TotalFitness		= 0;
		}

		//-------------------------------------------------
		// EPOCH 
		//-------------------------------------------------
		void 	Epoch()
		{

			//first reset variables and calculate the fitness of each genome
			Reset();
			CalculatePopulationsFitness();

			//if a solution is found exit
			//if (_bestSoFar <= _BestPossible)
			if((_bestSoFar-_BestPossible)<CGlobals::EPSILON)
			{
				_ProblemSolved = true;
				return;
			}

			//create a temporary vector for the new population
			vector<CGenome> vecNewPop;

			//First add NUM_BEST_TO_ADD number of the last generation's
			//fittest genome(elitism)
			for (int i=0; i<2; ++i)
			{
				vecNewPop.push_back(_Population[_FittestGene]);
			}


			//now create the remainder of the population
			while (vecNewPop.size() != POPULATION_SIZE)
			{

				//grab two parents
				CGenome mother = RouletteWheelSelection();
				CGenome father = RouletteWheelSelection();

				//create 2 children
				CGenome baby1, baby2;

				//Breed them
				Crossover(mother.vecCityPath,
					father.vecCityPath,
					baby1.vecCityPath,
					baby2.vecCityPath);

				//and mutate them
				Mutate(baby1.vecCityPath);
				Mutate(baby2.vecCityPath);

				//add them to new population
				vecNewPop.push_back(baby1);
				vecNewPop.push_back(baby2);
			}

			//copy into next generation
			_Population = vecNewPop;

			//increment generation counter
			++_Generation;
		}


		//-------------------------------------------------
		// Attach Graph 
		//-------------------------------------------------
		void AttachGraph(CChart* pGraph)
		{
			this->_graph=pGraph;
		}
		//-------------------------------------------------
		//  
		//-------------------------------------------------
		void Render(CMemDC* surface)
		{
			int CITY_SIZE=gGAConstants.CITY_SIZE;
			int NUM_CITIES=gGAConstants.NUM_CITES;
			int offset=100;

			oldPen = (CPen*)surface->SelectObject(&blackPen);
			//draw all the cities
			for (int city_num = 0; city_num <static_cast<int>(_Cities.size()); ++city_num)
			{
				int x = (int)_Cities[city_num].x;
				int y = (int)_Cities[city_num].y;				
				surface->Ellipse(x-CITY_SIZE, y-CITY_SIZE, x+CITY_SIZE, y+CITY_SIZE);
			}

			//draw the fittest route so far
			vector<int> route = _Population[_FittestGene].vecCityPath;

			//only display the routes if we are in a run
			if (_Generation != 0)
			{

				oldPen=surface->SelectObject(&bluePen);
				surface	->MoveTo( (int)_Cities[route[0]].x, (int)_Cities[route[0]].y);

				for (int i=1; i<static_cast<int>(route.size()); ++i)
				{
					int CityX = (int)_Cities[route[i]].x;
					int CityY = (int)_Cities[route[i]].y;
					surface->LineTo(CityX,CityY);
				}

				//now draw the line from the last city visited back to the starting point
				surface->LineTo( (int)_Cities[route[0]].x, (int)_Cities[route[0]].y);
			}	


        
           gFont.SetBold(true);
		   gFont.SetUnderline(true);
           surface->SelectObject(gFont);
			sprintf(buffer,"Epoch = %d",_Generation);
			surface->TextOut(52,10,buffer,(int)strlen(buffer));		
			

			gFont.SetBold(false);
			gFont.SetUnderline(false);				
			surface->SelectObject(&gFont);
			

			if(_bestSoFar == CGlobals::INFINITY)
			{
               strcpy(buffer,"Best So Far = N/A");
			}
			else
			{
				sprintf(buffer,"Best So Far = %.2f",_bestSoFar);				
			}
			surface->SetTextColor(CColor::blue);
			surface->TextOut(52,30,buffer,(int)strlen(buffer));
			
			sprintf(buffer,"Best Possible = %.2f",_BestPossible);
			surface->SetTextColor(CColor::green);
			surface->TextOut(52,50,buffer,(int)strlen(buffer));


			sprintf(buffer,"Elapsed Time = %.2f(ms)",_ElapsedTime);
			surface->SetTextColor(CColor::black);
			surface->TextOut(52,70,buffer,(int)strlen(buffer));


			if(this->_graph!=NULL )
			{								
					_graph->SetRange(0,_Generation+50,0,_AverageWorst);
					_graph->SetXYValue(_Generation,_bestSoFar,_Generation,0);
					_graph->SetXYValue(_Generation,_BestPossible,_Generation,1);								
			}
		}

	};//class
}//namespace
#endif