//geneticEngine.mel by Martin Hemberg 2004

//This is a genetic engine, ie a simple genetic algorithm for
//evolutionary search implemented in MEL. In order to make it useful,
//one must define a fitness function.

//Some global variables. Unfortunately, matrix sizes can not be set
//dynamically in MEL. Moreover, one can't declare constants. Thus I'm
//using this method, I think it's marginally better than not having
//them. Thus, we are stuck with a fixed size, if you want to change pop
//size or length, change these numbers and do a search and replace on
//the rest of the file 

global int $gGenomeLength = 100; 
global int $gPopulationSize = 50; 
global int $gMaxGeneValue = 42; //For a limit when generating random numbers, you may need to adjust this depending on your problem.

global proc evolve(int $generations, int $elites, int $tournamentSize, float $mutationRate)
{
	//Call this function in order to do an evolutionary run for a
	//fixed number of generations. Unfortunately, there is no way to
	//continue running with the same population. However, that should not be
	//to hard to fix. It's just a matter of saving the variables in the
	//scene adn writing a function for reading them again at re-start.
	matrix $genomes[50][100] = initializePopulation(); //[popsize][geneomelength] Each row represents an individual
	float $fitness[]; //An array containing the fitness values
	int $i;
	//Iterate through the generations
	for($i=0; $i<$generations; $i++){
		$fitness = evaluate($genomes); //Find out who's best
		$genomes = breed($genomes, $fitness, $elites, $tournamentSize, $mutationRate); //Create the population for the next generation
	}
}

//This function creates a population of individuals with random gene values.
proc matrix initializePopulation()
{
	global int $gPopulationSize;
	global int $gGenomeLength;
	global int $gMaxGeneValue;
	matrix $genomes[50][100];
	int $i, $j;
	//Set each gene to a random value
	for($i=0; $i<$gPopulationSize; $i++){
		for($j=0; $j<$gGenomeLength; $j++)
			$genomes[$i][$j] = rand ($gMaxGeneValue);
	}
	return $genomes;
}

//This function is used to evaluate the fitness of the population. It
//contains a large gap - you'll have to insert your own fitness function
//that maps the array of doubles (the genome) to a scalar value (the
//fitness). The fitness should be set so that a low value indicates a fit
//individual.
proc float[] evaluate(matrix $genomes)
{
	global int $gPopulationSize;
	int $i;
	float $fitness[50];
	for($i=0; $i<$gPopulationSize; $i++){
		//Here each individual should be evaluated
	}
	return $fitness;
}

//This function is used to create the population for the next generation
proc matrix breed(matrix $genomes, float $fitness[], int $elites, int $tournamentSize, float $mutationRate)
{
	global int $gPopulationSize;
	global int $gGenomeLength;
	matrix $newGenomes[50][100]; //The next generation
	int $lowest = -1;
	//Copy the elites to the new generation, they will automatically be the first individuals in the new generation.
	for($i=0; $i<$elites; $i++){
		$lowest = findLowestFitness($fitness, $lowest);
		for($j=0; $j<$gGenomeLength; $j++)
			$newGenomes[$i][$j] = $genomes[$lowest][$j];
	}
	//Create the rest through tournament selection
	for($i=$elites; $i<$gPopulationSize; $i++){
		int $father = tournamentSelect($fitness, $tournamentSize);
		int $mother = tournamentSelect($fitness, $tournamentSize);
		$newGenomes = crossover($genomes, $newGenomes, $father, $mother, $i);
	}
	return mutate($newGenomes, $elites, $mutationRate);
}

//Returns the index of the individual with the lowest fitness above the one indicated by the argument $lowest.
//If $lowest==-1, return the global minimum
proc int findLowestFitness(float $fitness[], int $lowest)
{
	int $i, $index;
	float $indexFitness = 100000; //Very high number
	float $lowestFitness = $fitness[$lowest];
	for($i=0; $i<size($fitness); $i++){
		if($fitness[$i]<$indexFitness && $fitness[$i]>=$lowestFitness && $i>$lowest){
			$index = $i;
			$indexFitness = $fitness[$i];
		}
	}
	return $index;
}

//Choose a parent using tournament selection. $tournamentSize
//individuals are randomly chosen from the population. The best one of
//these is then used as one of the parents for a member of the next
//generation. A high value of $tournamentSize (compared to the pop size)
//means that there will be less variation since the fittest individuals
//are more likely to be picked in every tournament.

proc int tournamentSelect(float $fitness[], int $tournamentSize) {
	global int $gPopulationSize;
	int $parent = (int)rand ($gPopulationSize), $i, $tmp;
	for($i=1; $i<$tournamentSize; $i++){
		$tmp = (int)rand ($gPopulationSize);
		if($fitness[$tmp]<$fitness[$parent]) //Fitness minimization
			$parent = $tmp;
	}
	return $parent;
}

//Combine to individuals to produce a new one for the next generation.
proc matrix crossover(matrix $genomes, matrix $newGenomes, int $father, int $mother, int $child)
{
	global int $gGenomeLength;
	int $i, $crossoverPoint = (int)rand ($gGenomeLength);
	for($i=0; $i<$crossoverPoint; $i++)
		$newGenomes[$child][$i] = $genomes[$father][$i];
	for($i=$crossoverPoint; $i<$gGenomeLength; $i++)
		$newGenomes[$child][$i] = $genomes[$mother][$i];
	return $newGenomes;
}

//Random mutation
proc matrix mutate(matrix $newGenomes, int $elites, float $mutationRate)
{
	global int $gPopulationSize;
	global int $gGenomeLength;
	int $i, $j;
	for($i=$elites; $i<$gPopulationSize; $i++){ //Don't mutate the elites
		for($j=0; $j<$gGenomeLength; $j++){
			if($mutationRate>`rand 100.0`){ //Mutate
				if(`rand 2.0`>1.0)
					$newGenomes[$i][$j] = $newGenomes[$i][$j] + 1;
				else
					$newGenomes[$i][$j] = $newGenomes[$i][$j] - 1;
			}
		}
	}
	return $newGenomes;
}