~swaits/differential-evolution

6e2d4bf3cab5304fdc5e1cee23a291cb149151e5 — swaits 19 years ago 70650ba
adjust eol style
7 files changed, 838 insertions(+), 838 deletions(-)

M de.h
M de.inl
M derand.cpp
M derand.h
M derand.inl
M devec.h
M devec.inl
M de.h => de.h +80 -80
@@ 1,80 1,80 @@
#ifndef __de_h__
#define __de_h__


#include "devec.h"
#include "derand.h"


namespace DE
{

	
	template <unsigned int DIM, unsigned int POP=(DIM*10)>
	class Engine
	{
	public:

		Engine();
		virtual ~Engine();

		void Reset();

		const double* GetBest() const;

		void SetRange(double minimum, double maximum);
		void SetRange(const Vector<DIM>& minimum, const Vector<DIM>& maximum);
		
		virtual double CalculateError(const double testcase[DIM], bool& stop) = 0;

		bool Solve(unsigned int maxgenerations = 1000000000);

	private:

		void Select(unsigned int candidate, unsigned int *r1 = 0, unsigned int *r2 = 0, unsigned int *r3 = 0, unsigned int *r4 = 0, unsigned int *r5 = 0);

		void Dump(const double& fitness, Vector<DIM>& individual) const;

		void MakeTrial_best1exp(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_rand1exp(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_randtobest1exp(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_best2exp(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_rand2exp(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_best1bin(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_rand1bin(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_randtobest1bin(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_best2bin(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_rand2bin(unsigned int candidate, Vector<DIM>& trial);

	private:

		static const double DEFAULTSCALE;
		static const double DEFAULTCROSSOVER;
		static const double BIGDOUBLE;


	private:

		Random      prng;

		Vector<DIM> best;
		Vector<DIM> minimum;
		Vector<DIM> maximum;
		Vector<DIM> population[POP];
		double      fitness[POP];

		double      scale;
		double      crossover;
		double      bestfitness;

	};


};


#include "de.inl"


#endif // __de_h__

#ifndef __de_h__
#define __de_h__


#include "devec.h"
#include "derand.h"


namespace DE
{

	
	template <unsigned int DIM, unsigned int POP=(DIM*10)>
	class Engine
	{
	public:

		Engine();
		virtual ~Engine();

		void Reset();

		const double* GetBest() const;

		void SetRange(double minimum, double maximum);
		void SetRange(const Vector<DIM>& minimum, const Vector<DIM>& maximum);
		
		virtual double CalculateError(const double testcase[DIM], bool& stop) = 0;

		bool Solve(unsigned int maxgenerations = 1000000000);

	private:

		void Select(unsigned int candidate, unsigned int *r1 = 0, unsigned int *r2 = 0, unsigned int *r3 = 0, unsigned int *r4 = 0, unsigned int *r5 = 0);

		void Dump(const double& fitness, Vector<DIM>& individual) const;

		void MakeTrial_best1exp(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_rand1exp(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_randtobest1exp(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_best2exp(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_rand2exp(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_best1bin(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_rand1bin(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_randtobest1bin(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_best2bin(unsigned int candidate, Vector<DIM>& trial);
		void MakeTrial_rand2bin(unsigned int candidate, Vector<DIM>& trial);

	private:

		static const double DEFAULTSCALE;
		static const double DEFAULTCROSSOVER;
		static const double BIGDOUBLE;


	private:

		Random      prng;

		Vector<DIM> best;
		Vector<DIM> minimum;
		Vector<DIM> maximum;
		Vector<DIM> population[POP];
		double      fitness[POP];

		double      scale;
		double      crossover;
		double      bestfitness;

	};


};


#include "de.inl"


#endif // __de_h__


M de.inl => de.inl +407 -407
@@ 1,407 1,407 @@
#ifndef __de_inl__
#define __de_inl__


#include <cstdio>


#if !(defined(__de_h__))
#	 error "de.inl should only be included by de.h"
#endif


template <unsigned int DIM, unsigned int POP>
const double DE::Engine<DIM,POP>::DEFAULTSCALE = 0.85;

template <unsigned int DIM, unsigned int POP>
const double DE::Engine<DIM,POP>::DEFAULTCROSSOVER = 1.0;

template <unsigned int DIM, unsigned int POP>
const double DE::Engine<DIM,POP>::BIGDOUBLE = 1.0e30;

template <unsigned int DIM, unsigned int POP>
inline DE::Engine<DIM,POP>::Engine()
{
	scale     = DEFAULTSCALE;
	crossover = DEFAULTCROSSOVER;

	for ( unsigned int i=0;i<DIM;i++ )
	{
		this->minimum[i] = -BIGDOUBLE;
		this->maximum[i] = -(this->minimum[i]);
	}

	Reset();

	return;
}

template <unsigned int DIM, unsigned int POP>
inline DE::Engine<DIM,POP>::~Engine()
{
	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::Reset()
{
	// reset our best
	bestfitness = BIGDOUBLE;
	for ( unsigned int i=0;i<DIM;i++ )
	{
		best[i] = 0.0;
	}

	// reset all vectors
	for ( unsigned int i=0;i<POP;i++ )
	{
		// energy tied for worst
		fitness[i] = bestfitness;

		// make a new individual
		for ( unsigned int j=0;j<DIM;j++ )
		{
			population[i][j] = prng.RandDouble(minimum[j],maximum[j]);
		}
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::SetRange(double minimum, double maximum)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		this->minimum[i] = minimum;
		this->maximum[i] = maximum;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::SetRange(const Vector<DIM>& minimum, const Vector<DIM>& maximum)
{
	this->minimum = minimum;
	this->maximum = maximum;

	return;
}

template <unsigned int DIM, unsigned int POP>
inline const double* DE::Engine<DIM,POP>::GetBest() const
{
	return best;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::Select(unsigned int candidate, unsigned int *r1, unsigned int *r2, unsigned int *r3, unsigned int *r4, unsigned int *r5)
{
	if ( r1 )
	{
		do
		{
			*r1 = prng.RandU32(0,POP-1);
		}
		while ( *r1 == candidate );
	}

	if ( r2 )
	{
		do
		{
			*r2 = prng.RandU32(0,POP-1);
		}
		while ( (*r2 == candidate) || (*r2 == *r1) );
	}

	if ( r3 )
	{
		do
		{
			*r3 = prng.RandU32(0,POP-1);
		}
		while ( (*r3 == candidate) || (*r3 == *r2) || (*r3 == *r1) );
	}

	if ( r4 )
	{
		do
		{
			*r4 = prng.RandU32(0,POP-1);
		}
		while ( (*r4 == candidate) || (*r4 == *r3) || (*r4 == *r2) || (*r4 == *r1) );
	}

	if ( r5 )
	{
		do
		{
			*r5 = prng.RandU32(0,POP-1);
		}
		while ( (*r5 == candidate) || (*r5 == *r4) || (*r5 == *r3) || (*r5 == *r2) || (*r5 == *r1) );
	}

	return;
}

// TODO: change this to run one generation at a time???
template <unsigned int DIM, unsigned int POP>
inline bool DE::Engine<DIM,POP>::Solve(unsigned int maxgenerations)
{
	bool success = false;

	Reset();

	for (unsigned int generation=0;generation<maxgenerations && !success;generation++)
	{
		printf("====== GENERATION %07d =============================================\n\n",generation);

		for (int candidate=0;candidate<POP && !success;candidate++)
		{
			Vector<DIM> trial;
			MakeTrial_randtobest1bin(candidate,trial);

			double trialfitness = CalculateError(trial,success);

			if ( trialfitness < fitness[candidate] )
			{
				fitness[candidate]    = trialfitness;
				population[candidate] = trial;

				if ( trialfitness < bestfitness )
				{
					bestfitness = trialfitness;
					best        = trial;
					
					Dump(bestfitness,best);
				}
			}
		}
	}

	printf("====== RUN COMPLETE       =============================================\n\n");
	Dump(bestfitness,best);

	return success;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::Dump(const double& fitness, DE::Vector<DIM>& individual) const
{
	printf(" fitness = %60.30f\n",fitness);
	for (unsigned int i=0;i<DIM;i++)
	{
		printf("   [%03d] = %60.30f\n",i,individual[i]);
	}
	printf("\n");
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_best1exp(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2;
	unsigned int n;

	Select(candidate, &r1, &r2);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; (prng.RandDouble(0.0,1.0) < crossover) && (i < DIM); i++) 
	{
		trial[n] = best[n] + scale * (population[r1][n] - population[r2][n]);
		n        = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_rand1exp(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; (prng.RandDouble(0.0,1.0) <= crossover) && (i < DIM); i++) 
	{
		trial[n] = population[r1][n] + scale * (population[r2][n] - population[r3][n]);
		n        = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_randtobest1exp(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2;
	unsigned int n;

	Select(candidate, &r1, &r2);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; (prng.RandDouble(0.0,1.0) <= crossover) && (i < DIM); i++) 
	{
		trial[n] += scale * (best[n] - trial[n]) + scale * (population[r1][n] - population[r2][n]);
		n         = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_best2exp(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3, r4;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3, &r4);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; (prng.RandDouble(0.0,1.0) <= crossover) && (i < DIM); i++) 
	{
		trial[n] = best[n] + scale * (population[r1][n] + population[r2][n] - population[r3][n] - population[r4][n]);
		n        = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_rand2exp(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3, r4, r5;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3 ,&r4, &r5);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; (prng.RandDouble(0.0,1.0) <= crossover) && (i < DIM); i++) 
	{
		trial[n] = population[r1][n] + scale * (population[r2][n] + population[r3][n] - population[r4][n] - population[r5][n]);
		n        = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_best1bin(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2;
	unsigned int n;

	Select(candidate, &r1, &r2);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; i < DIM; i++) 
	{
		if ( (prng.RandDouble(0.0,1.0) <= crossover) || (i == (DIM - 1)) )
		{
			trial[n] = best[n] + scale * (population[r1][n] - population[r2][n]);
		}
		n = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_rand1bin(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; i < DIM; i++) 
	{
		if ( (prng.RandDouble(0.0,1.0) <= crossover) || (i == (DIM - 1)) )
		{
			trial[n] = population[r1][n] + scale * (population[r2][n] - population[r3][n]);
		}
		n = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_randtobest1bin(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2;
	unsigned int n;

	Select(candidate, &r1, &r2);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; i < DIM; i++) 
	{
		if ( (prng.RandDouble(0.0,1.0) <= crossover) || (i == (DIM - 1)) )
		{
			trial[n] += scale * (best[n] - trial[n]) + scale * (population[r1][n] - population[r2][n]);
		}
		n = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_best2bin(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3, r4;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3, &r4);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; i < DIM; i++) 
	{
		if ( (prng.RandDouble(0.0,1.0) <= crossover) || (i == (DIM - 1)) )
		{
			trial[n] = best[n] + scale * (population[r1][n] + population[r2][n] - population[r3][n] - population[r4][n]);
		}
		n = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_rand2bin(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3, r4, r5;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3, &r4, &r5);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; i < DIM; i++) 
	{
		if ( (prng.RandDouble(0.0,1.0) <= crossover) || (i == (DIM - 1)) )
		{
			trial[n] = population[r1][n] + scale * (population[r2][n] + population[r3][n] - population[r4][n] - population[r5][n]);
		}
		n = (n + 1) % DIM;
	}

	return;
}

#endif // __de_inl__

#ifndef __de_inl__
#define __de_inl__


#include <cstdio>


#if !(defined(__de_h__))
#	 error "de.inl should only be included by de.h"
#endif


template <unsigned int DIM, unsigned int POP>
const double DE::Engine<DIM,POP>::DEFAULTSCALE = 0.85;

template <unsigned int DIM, unsigned int POP>
const double DE::Engine<DIM,POP>::DEFAULTCROSSOVER = 1.0;

template <unsigned int DIM, unsigned int POP>
const double DE::Engine<DIM,POP>::BIGDOUBLE = 1.0e30;

template <unsigned int DIM, unsigned int POP>
inline DE::Engine<DIM,POP>::Engine()
{
	scale     = DEFAULTSCALE;
	crossover = DEFAULTCROSSOVER;

	for ( unsigned int i=0;i<DIM;i++ )
	{
		this->minimum[i] = -BIGDOUBLE;
		this->maximum[i] = -(this->minimum[i]);
	}

	Reset();

	return;
}

template <unsigned int DIM, unsigned int POP>
inline DE::Engine<DIM,POP>::~Engine()
{
	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::Reset()
{
	// reset our best
	bestfitness = BIGDOUBLE;
	for ( unsigned int i=0;i<DIM;i++ )
	{
		best[i] = 0.0;
	}

	// reset all vectors
	for ( unsigned int i=0;i<POP;i++ )
	{
		// energy tied for worst
		fitness[i] = bestfitness;

		// make a new individual
		for ( unsigned int j=0;j<DIM;j++ )
		{
			population[i][j] = prng.RandDouble(minimum[j],maximum[j]);
		}
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::SetRange(double minimum, double maximum)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		this->minimum[i] = minimum;
		this->maximum[i] = maximum;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::SetRange(const Vector<DIM>& minimum, const Vector<DIM>& maximum)
{
	this->minimum = minimum;
	this->maximum = maximum;

	return;
}

template <unsigned int DIM, unsigned int POP>
inline const double* DE::Engine<DIM,POP>::GetBest() const
{
	return best;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::Select(unsigned int candidate, unsigned int *r1, unsigned int *r2, unsigned int *r3, unsigned int *r4, unsigned int *r5)
{
	if ( r1 )
	{
		do
		{
			*r1 = prng.RandU32(0,POP-1);
		}
		while ( *r1 == candidate );
	}

	if ( r2 )
	{
		do
		{
			*r2 = prng.RandU32(0,POP-1);
		}
		while ( (*r2 == candidate) || (*r2 == *r1) );
	}

	if ( r3 )
	{
		do
		{
			*r3 = prng.RandU32(0,POP-1);
		}
		while ( (*r3 == candidate) || (*r3 == *r2) || (*r3 == *r1) );
	}

	if ( r4 )
	{
		do
		{
			*r4 = prng.RandU32(0,POP-1);
		}
		while ( (*r4 == candidate) || (*r4 == *r3) || (*r4 == *r2) || (*r4 == *r1) );
	}

	if ( r5 )
	{
		do
		{
			*r5 = prng.RandU32(0,POP-1);
		}
		while ( (*r5 == candidate) || (*r5 == *r4) || (*r5 == *r3) || (*r5 == *r2) || (*r5 == *r1) );
	}

	return;
}

// TODO: change this to run one generation at a time???
template <unsigned int DIM, unsigned int POP>
inline bool DE::Engine<DIM,POP>::Solve(unsigned int maxgenerations)
{
	bool success = false;

	Reset();

	for (unsigned int generation=0;generation<maxgenerations && !success;generation++)
	{
		printf("====== GENERATION %07d =============================================\n\n",generation);

		for (int candidate=0;candidate<POP && !success;candidate++)
		{
			Vector<DIM> trial;
			MakeTrial_randtobest1bin(candidate,trial);

			double trialfitness = CalculateError(trial,success);

			if ( trialfitness < fitness[candidate] )
			{
				fitness[candidate]    = trialfitness;
				population[candidate] = trial;

				if ( trialfitness < bestfitness )
				{
					bestfitness = trialfitness;
					best        = trial;
					
					Dump(bestfitness,best);
				}
			}
		}
	}

	printf("====== RUN COMPLETE       =============================================\n\n");
	Dump(bestfitness,best);

	return success;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::Dump(const double& fitness, DE::Vector<DIM>& individual) const
{
	printf(" fitness = %60.30f\n",fitness);
	for (unsigned int i=0;i<DIM;i++)
	{
		printf("   [%03d] = %60.30f\n",i,individual[i]);
	}
	printf("\n");
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_best1exp(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2;
	unsigned int n;

	Select(candidate, &r1, &r2);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; (prng.RandDouble(0.0,1.0) < crossover) && (i < DIM); i++) 
	{
		trial[n] = best[n] + scale * (population[r1][n] - population[r2][n]);
		n        = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_rand1exp(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; (prng.RandDouble(0.0,1.0) <= crossover) && (i < DIM); i++) 
	{
		trial[n] = population[r1][n] + scale * (population[r2][n] - population[r3][n]);
		n        = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_randtobest1exp(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2;
	unsigned int n;

	Select(candidate, &r1, &r2);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; (prng.RandDouble(0.0,1.0) <= crossover) && (i < DIM); i++) 
	{
		trial[n] += scale * (best[n] - trial[n]) + scale * (population[r1][n] - population[r2][n]);
		n         = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_best2exp(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3, r4;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3, &r4);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; (prng.RandDouble(0.0,1.0) <= crossover) && (i < DIM); i++) 
	{
		trial[n] = best[n] + scale * (population[r1][n] + population[r2][n] - population[r3][n] - population[r4][n]);
		n        = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_rand2exp(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3, r4, r5;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3 ,&r4, &r5);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; (prng.RandDouble(0.0,1.0) <= crossover) && (i < DIM); i++) 
	{
		trial[n] = population[r1][n] + scale * (population[r2][n] + population[r3][n] - population[r4][n] - population[r5][n]);
		n        = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_best1bin(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2;
	unsigned int n;

	Select(candidate, &r1, &r2);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; i < DIM; i++) 
	{
		if ( (prng.RandDouble(0.0,1.0) <= crossover) || (i == (DIM - 1)) )
		{
			trial[n] = best[n] + scale * (population[r1][n] - population[r2][n]);
		}
		n = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_rand1bin(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; i < DIM; i++) 
	{
		if ( (prng.RandDouble(0.0,1.0) <= crossover) || (i == (DIM - 1)) )
		{
			trial[n] = population[r1][n] + scale * (population[r2][n] - population[r3][n]);
		}
		n = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_randtobest1bin(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2;
	unsigned int n;

	Select(candidate, &r1, &r2);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; i < DIM; i++) 
	{
		if ( (prng.RandDouble(0.0,1.0) <= crossover) || (i == (DIM - 1)) )
		{
			trial[n] += scale * (best[n] - trial[n]) + scale * (population[r1][n] - population[r2][n]);
		}
		n = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_best2bin(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3, r4;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3, &r4);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; i < DIM; i++) 
	{
		if ( (prng.RandDouble(0.0,1.0) <= crossover) || (i == (DIM - 1)) )
		{
			trial[n] = best[n] + scale * (population[r1][n] + population[r2][n] - population[r3][n] - population[r4][n]);
		}
		n = (n + 1) % DIM;
	}

	return;
}

template <unsigned int DIM, unsigned int POP>
inline void DE::Engine<DIM,POP>::MakeTrial_rand2bin(unsigned int candidate, DE::Vector<DIM>& trial)
{
	unsigned int r1, r2, r3, r4, r5;
	unsigned int n;

	Select(candidate, &r1, &r2, &r3, &r4, &r5);
	n = prng.RandU32(0,DIM-1);

	trial = population[candidate];
	for (unsigned int i=0; i < DIM; i++) 
	{
		if ( (prng.RandDouble(0.0,1.0) <= crossover) || (i == (DIM - 1)) )
		{
			trial[n] = population[r1][n] + scale * (population[r2][n] + population[r3][n] - population[r4][n] - population[r5][n]);
		}
		n = (n + 1) % DIM;
	}

	return;
}

#endif // __de_inl__


M derand.cpp => derand.cpp +1 -1
@@ 1,1 1,1 @@



M derand.h => derand.h +76 -76
@@ 1,76 1,76 @@
#ifndef __derand_h__
#define __derand_h__

namespace DE
{


	typedef unsigned int u32;
	typedef signed   int s32;
	
	
	class Random
	{
	public:
	
		void Seed(u32 seed);
		
		u32 Rand();
		u32 operator()();
		
		float RandFloat();
		float RandFloat(float min, float max);
	
		double RandDouble();
		double RandDouble(double min, double max);
	
		u32 RandU32();
		u32 RandU32(u32 min, u32 max);
	
		s32 RandS32();
		s32 RandS32(s32 min, s32 max);
	
	

	private:
	
		static const u32 U32_MAX = (u32)0xffffffff;

		/**
		 * LCG multiplier by M. Lavaux and F. Janssens
		 */
		static const u32 MULTIPLIER = 1664525UL;
	
		/**
		 * LCG increment by M. Lavaux and F. Janssens
		 */
		static const u32 INCREMENT = 1013904223UL;
	
		/**
		 * LCG modulus set to 2^32, but not used due to MODULUS_ASSUMED
		 * above (since we're using a 32 bit word here)
		 *
		 * Set this to 0 if the modulus falls out automatically
		 * due to your data type (i.e. "addition modulo".  For
		 * instance, if you're using a 32 bit integer, and your
		 * modulus is 2^32, then there's no need to perform the
		 * modulo since   a == a % (2^32). (so the modulo is "free")
		 * 
		 */
		static const u32 MODULUS = 0xffffffff;
	
		/**
		 * LCG state
		 */
		u32 seed;
	};


};


#include "derand.inl"


#endif // __derand_h__

#ifndef __derand_h__
#define __derand_h__

namespace DE
{


	typedef unsigned int u32;
	typedef signed   int s32;
	
	
	class Random
	{
	public:
	
		void Seed(u32 seed);
		
		u32 Rand();
		u32 operator()();
		
		float RandFloat();
		float RandFloat(float min, float max);
	
		double RandDouble();
		double RandDouble(double min, double max);
	
		u32 RandU32();
		u32 RandU32(u32 min, u32 max);
	
		s32 RandS32();
		s32 RandS32(s32 min, s32 max);
	
	

	private:
	
		static const u32 U32_MAX = (u32)0xffffffff;

		/**
		 * LCG multiplier by M. Lavaux and F. Janssens
		 */
		static const u32 MULTIPLIER = 1664525UL;
	
		/**
		 * LCG increment by M. Lavaux and F. Janssens
		 */
		static const u32 INCREMENT = 1013904223UL;
	
		/**
		 * LCG modulus set to 2^32, but not used due to MODULUS_ASSUMED
		 * above (since we're using a 32 bit word here)
		 *
		 * Set this to 0 if the modulus falls out automatically
		 * due to your data type (i.e. "addition modulo".  For
		 * instance, if you're using a 32 bit integer, and your
		 * modulus is 2^32, then there's no need to perform the
		 * modulo since   a == a % (2^32). (so the modulo is "free")
		 * 
		 */
		static const u32 MODULUS = 0xffffffff;
	
		/**
		 * LCG state
		 */
		u32 seed;
	};


};


#include "derand.inl"


#endif // __derand_h__


M derand.inl => derand.inl +69 -69
@@ 1,69 1,69 @@
#ifndef __derand_inl__
#define __derand_inl__


#if !(defined(__derand_h__))
#	 error "derand.inl should only be included by derand.h"
#endif


inline void DE::Random::Seed(u32 seed)
{
	this->seed = seed;
}

inline DE::u32 DE::Random::Rand()
{
	seed = ( (seed * DE::Random::MULTIPLIER) + DE::Random::INCREMENT );
	seed %= DE::Random::MODULUS;
	
	return seed;
}

inline DE::u32 DE::Random::operator()()
{
	return this->Rand();
}

inline float DE::Random::RandFloat()
{
	return (float)Rand() / (float)DE::Random::U32_MAX;
}

inline float DE::Random::RandFloat(float min, float max)
{
	return ( RandFloat() * (max - min) + min );
}

inline double DE::Random::RandDouble()
{
	return (double)Rand() / (double)DE::Random::U32_MAX;
}

inline double DE::Random::RandDouble(double min, double max)
{
	return ( RandDouble() * (max - min) + min );
}

inline DE::u32 DE::Random::RandU32()
{
	return Rand();
}

inline DE::u32 DE::Random::RandU32(DE::u32 min, DE::u32 max)
{
	return ( ( Rand() % (max - min) ) + min );
}

inline DE::s32 DE::Random::RandS32()
{
	return (DE::s32)Rand();
}

inline DE::s32 DE::Random::RandS32(DE::s32 min, DE::s32 max)
{
	return (DE::s32)( ( Rand() % (max - min) ) + min );
}

#endif // __derand_inl__

#ifndef __derand_inl__
#define __derand_inl__


#if !(defined(__derand_h__))
#	 error "derand.inl should only be included by derand.h"
#endif


inline void DE::Random::Seed(u32 seed)
{
	this->seed = seed;
}

inline DE::u32 DE::Random::Rand()
{
	seed = ( (seed * DE::Random::MULTIPLIER) + DE::Random::INCREMENT );
	seed %= DE::Random::MODULUS;
	
	return seed;
}

inline DE::u32 DE::Random::operator()()
{
	return this->Rand();
}

inline float DE::Random::RandFloat()
{
	return (float)Rand() / (float)DE::Random::U32_MAX;
}

inline float DE::Random::RandFloat(float min, float max)
{
	return ( RandFloat() * (max - min) + min );
}

inline double DE::Random::RandDouble()
{
	return (double)Rand() / (double)DE::Random::U32_MAX;
}

inline double DE::Random::RandDouble(double min, double max)
{
	return ( RandDouble() * (max - min) + min );
}

inline DE::u32 DE::Random::RandU32()
{
	return Rand();
}

inline DE::u32 DE::Random::RandU32(DE::u32 min, DE::u32 max)
{
	return ( ( Rand() % (max - min) ) + min );
}

inline DE::s32 DE::Random::RandS32()
{
	return (DE::s32)Rand();
}

inline DE::s32 DE::Random::RandS32(DE::s32 min, DE::s32 max)
{
	return (DE::s32)( ( Rand() % (max - min) ) + min );
}

#endif // __derand_inl__


M devec.h => devec.h +63 -63
@@ 1,63 1,63 @@
#ifndef __devec_h__
#define __devec_h__

namespace DE
{
	
	
	template <unsigned int DIM>
	class Vector
	{
	public:

		// default constructor
		Vector();

		// destructor
		~Vector();

		// copy constructor
		Vector(const Vector& other);

		// assignment
		Vector& operator= (const Vector& other);

		// addition
		Vector& operator+= (const Vector& other);
		Vector  operator+  (const Vector& other);

		// subtraction
		Vector& operator-= (const Vector& other);
		Vector  operator-  (const Vector& other);

		// multiply by scalar
		Vector& operator*= (const double& multiplier);
		Vector  operator*  (const double& multiplier);

		// divide by scalar
		Vector& operator/= (const double& divisor);
		Vector  operator/  (const double& divisor);

		// element access
		const double& operator[] (unsigned int i) const;
		      double& operator[] (unsigned int i);

		// data access
		operator const double* () const;
		operator       double* ();

	private:

		// data lives here
		double tuple[DIM];
	};


};


#include "devec.inl"


#endif // __devec_h__

#ifndef __devec_h__
#define __devec_h__

namespace DE
{
	
	
	template <unsigned int DIM>
	class Vector
	{
	public:

		// default constructor
		Vector();

		// destructor
		~Vector();

		// copy constructor
		Vector(const Vector& other);

		// assignment
		Vector& operator= (const Vector& other);

		// addition
		Vector& operator+= (const Vector& other);
		Vector  operator+  (const Vector& other);

		// subtraction
		Vector& operator-= (const Vector& other);
		Vector  operator-  (const Vector& other);

		// multiply by scalar
		Vector& operator*= (const double& multiplier);
		Vector  operator*  (const double& multiplier);

		// divide by scalar
		Vector& operator/= (const double& divisor);
		Vector  operator/  (const double& divisor);

		// element access
		const double& operator[] (unsigned int i) const;
		      double& operator[] (unsigned int i);

		// data access
		operator const double* () const;
		operator       double* ();

	private:

		// data lives here
		double tuple[DIM];
	};


};


#include "devec.inl"


#endif // __devec_h__


M devec.inl => devec.inl +142 -142
@@ 1,142 1,142 @@
#ifndef __devec_inl__
#define __devec_inl__


#if !(defined(__devec_h__))
#	 error "devec.inl should only be included by devec.h"
#endif


template <unsigned int DIM>
inline DE::Vector<DIM>::Vector()
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] = 0.0;
	}
	return;
}

template <unsigned int DIM>
inline DE::Vector<DIM>::~Vector()
{
	// empty
	return;
}

template <unsigned int DIM>
inline DE::Vector<DIM>::Vector(const Vector& other)
{
	*this = other;
	return;
}

template <unsigned int DIM>
inline DE::Vector<DIM>& DE::Vector<DIM>::operator= (const Vector& other)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] = other.tuple[i];
	}
	return *this;
}

template <unsigned int DIM>
inline DE::Vector<DIM>& DE::Vector<DIM>::operator+= (const Vector& other)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] += other.tuple[i];
	}
	return *this;
}

template <unsigned int DIM>
inline DE::Vector<DIM> DE::Vector<DIM>::operator+ (const Vector& other)
{
	Vector tmp(*this);
	tmp += other;
	return tmp;
}

template <unsigned int DIM>
inline DE::Vector<DIM>& DE::Vector<DIM>::operator-= (const Vector& other)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] -= other.tuple[i];
	}
	return *this;
}

template <unsigned int DIM>
inline DE::Vector<DIM> DE::Vector<DIM>::operator- (const Vector& other)
{
	Vector tmp(*this);
	tmp -= other;
	return tmp;
}

template <unsigned int DIM>
inline DE::Vector<DIM>& DE::Vector<DIM>::operator*= (const double& multiplier)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] *= multiplier;
	}
	return *this;
}

template <unsigned int DIM>
inline DE::Vector<DIM> DE::Vector<DIM>::operator* (const double& multiplier)
{
	Vector tmp(*this);
	tmp *= multiplier;
	return tmp;
}

template <unsigned int DIM>
inline DE::Vector<DIM>& DE::Vector<DIM>::operator/= (const double& divisor)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] /= divisor;
	}
	return *this;
}

template <unsigned int DIM>
inline DE::Vector<DIM> DE::Vector<DIM>::operator/ (const double& divisor)
{
	Vector tmp(*this);
	tmp /= divisor;
	return tmp;
}

template <unsigned int DIM>
inline const double& DE::Vector<DIM>::operator[] (unsigned int i) const
{
	return tuple[i];
}

template <unsigned int DIM>
inline double& DE::Vector<DIM>::operator[] (unsigned int i)
{
	return tuple[i];
}

template <unsigned int DIM>
inline DE::Vector<DIM>::operator const double* () const
{
	return tuple;
}

template <unsigned int DIM>
inline DE::Vector<DIM>::operator double* ()
{
	return tuple;
}


#endif // __devec_inl__

#ifndef __devec_inl__
#define __devec_inl__


#if !(defined(__devec_h__))
#	 error "devec.inl should only be included by devec.h"
#endif


template <unsigned int DIM>
inline DE::Vector<DIM>::Vector()
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] = 0.0;
	}
	return;
}

template <unsigned int DIM>
inline DE::Vector<DIM>::~Vector()
{
	// empty
	return;
}

template <unsigned int DIM>
inline DE::Vector<DIM>::Vector(const Vector& other)
{
	*this = other;
	return;
}

template <unsigned int DIM>
inline DE::Vector<DIM>& DE::Vector<DIM>::operator= (const Vector& other)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] = other.tuple[i];
	}
	return *this;
}

template <unsigned int DIM>
inline DE::Vector<DIM>& DE::Vector<DIM>::operator+= (const Vector& other)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] += other.tuple[i];
	}
	return *this;
}

template <unsigned int DIM>
inline DE::Vector<DIM> DE::Vector<DIM>::operator+ (const Vector& other)
{
	Vector tmp(*this);
	tmp += other;
	return tmp;
}

template <unsigned int DIM>
inline DE::Vector<DIM>& DE::Vector<DIM>::operator-= (const Vector& other)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] -= other.tuple[i];
	}
	return *this;
}

template <unsigned int DIM>
inline DE::Vector<DIM> DE::Vector<DIM>::operator- (const Vector& other)
{
	Vector tmp(*this);
	tmp -= other;
	return tmp;
}

template <unsigned int DIM>
inline DE::Vector<DIM>& DE::Vector<DIM>::operator*= (const double& multiplier)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] *= multiplier;
	}
	return *this;
}

template <unsigned int DIM>
inline DE::Vector<DIM> DE::Vector<DIM>::operator* (const double& multiplier)
{
	Vector tmp(*this);
	tmp *= multiplier;
	return tmp;
}

template <unsigned int DIM>
inline DE::Vector<DIM>& DE::Vector<DIM>::operator/= (const double& divisor)
{
	for (unsigned int i=0;i<DIM;i++)
	{
		tuple[i] /= divisor;
	}
	return *this;
}

template <unsigned int DIM>
inline DE::Vector<DIM> DE::Vector<DIM>::operator/ (const double& divisor)
{
	Vector tmp(*this);
	tmp /= divisor;
	return tmp;
}

template <unsigned int DIM>
inline const double& DE::Vector<DIM>::operator[] (unsigned int i) const
{
	return tuple[i];
}

template <unsigned int DIM>
inline double& DE::Vector<DIM>::operator[] (unsigned int i)
{
	return tuple[i];
}

template <unsigned int DIM>
inline DE::Vector<DIM>::operator const double* () const
{
	return tuple;
}

template <unsigned int DIM>
inline DE::Vector<DIM>::operator double* ()
{
	return tuple;
}


#endif // __devec_inl__