Source Discussion

The source code for the genetic algorithm, as well as the virtual computer that evaluates the chromosomes, is very simple. Let's first look at the virtual computer implementation ( otherwise known as a virtual machine).

Virtual Machine Implementation

The virtual machine (or VM) is a very simple architecture that supports a single integer stack and six instructions that use operands on the stack. Listing 6.1 shows the VM symbolics (from  common.h ) as well as the virtual machine (stm.c).

Listing 6.1: Virtual Machine Implementation.

 #define DUP              0x00 #define SWAP             0x01 #define MUL              0x02 #define ADD              0x03 #define OVER             0x04 #define NOP              0x05 #define MAX_INSTRUCTION          (NOP+1) #define NONE                     0 #define STACK_VIOLATION          1 #define MATH_VIOLATION           2 #define STACK_DEPTH       25 int stack[STACK_DEPTH]; int stackPointer; #define ASSERT_STACK_ELEMENTS(x) \      if (stackPointer < x) { error = STACK_VIOLATION ; break; } #define ASSERT_STACK_NOT_FULL \      if (stackPointer == STACK_DEPTH) { error = STACK_VIOLATION ;        break; } #define SPUSH(x) (stack[stackPointer++] = x) #define SPOP     (stack[--stackPointer]) #define SPEEK    (stack[stackPointer-1]) /*  * interpretSTM  *  *   program    - sequence of simple instructions  *   progLength - length of program  *   args       - arguments on the virtual stack  *   argsLength - number of arguments on the virtual stack  *  */ int interpretSTM(const int *program, int progLength,                    const int *args, int argsLength) {   int pc = 0;   int i, error = NONE;   int a, b;   stackPointer = 0;   /* Load the arguments onto the stack */   for (i = argsLength-1 ; i >= 0 ; i--) {     SPUSH(args[i]);   }   /* Execute the program */   while ((error == NONE) && (pc < progLength)) {     switch(program[pc++]) {       case DUP:         ASSERT_STACK_ELEMENTS(1);         ASSERT_STACK_NOT_FULL;         SPUSH(SPEEK);         break;       case SWAP:         ASSERT_STACK_ELEMENTS(2);         a = stack[stackPointer-1];         stack[stackPointer-1] = stack[stackPointer-2];         stack[stackPointer-2] = a;         break;       case MUL:         ASSERT_STACK_ELEMENTS(2);         a = SPOP; b = SPOP;         SPUSH(a * b);         break;       case ADD:         ASSERT_STACK_ELEMENTS(2);         a = SPOP; b = SPOP;         SPUSH(a + b);         break;       case OVER:         ASSERT_STACK_ELEMENTS(2);         SPUSH(stack[stackPointer-2]);         break;     } /* Switch opcode */   } /* Loop */   return(error); } 

 

The symbolics section of Listing 6.1 shows the creation of the instructions and their values that will be used by the algorithm. A number of other constants are defined to understand how the program exited (whether successfully NONE , or by an error). We create a simple integer stack with a defined depth and a number of macros that simplify the development of the VM. The assert_ macros are used to identify violations within the program (such as blowing out the stack or trying to use an element that's not there). The macros SPEEK , SPUSH , and SPOP simplify stack-based operations.

Function interpretSTM (or interpret Stack Machine) is the virtual machine. This function accepts a program ( program ) and the length of the passed program ( progLength ) as well as a set of arguments ( args ), and the number passed ( argsLength ). The first operation of the VM is to store the passed arguments onto the stack. We walk through the list backwards , so whatever element the user specified as element zero is at the top of the stack.

Finally, we walk through the program that was passed in, executing each of the instructions. When we reach the end of the user program, we return with the error status. If an error was encountered along the way, the error variable is set and we return automatically.

Genetic Algorithm Implementation

The genetic algorithm that will be implemented here follows the basic flow discussed earlier in this chapter (see Figure 6.2). We'll discuss the implementation from the top down, starting with the main program, initialization, evaluation, selection, and recombination.

C main Function

The main function, shown in Listing 6.2, illustrates the basic flow of the algorithm.

Listing 6.2: C Main Routine for the Genetic Algorithm.

 int main() {   int generation = 0, i;   FILE *fp;   extern float minFitness, maxFitness, avgFitness;   extern int curCrossovers, curMutations;   extern int curPop;   void printProgram( int, int );   /* Seed the random number generator */   srand(time(NULL));   curPop = 0;   fp = fopen("stats.txt", "w");   if (fp == NULL) exit(-1);   /* Initialize the initial population and check each element's    * fitness.    */   initPopulation();   performFitnessCheck( fp );   /* Loop for the maximum number of allowable generations */   while (generation < MAX_GENERATIONS) {     curCrossovers = curMutations = 0;     /* Select two parents and recombine to create two children */     performSelection();     /* Switch the populations */     curPop = (curPop == 0) ? 1 : 0;     /* Calculate the fitness of the new population */     performFitnessCheck( fp );     /* Emit statistics every 100 generations */     if ((generation++ % 100) == 0) {       printf("Generation %d\n", generation-1);       printf("\tmaxFitness = %f (%g)\n", maxFitness, MAX_FIT);       printf("\tavgFitness = %f\n", avgFitness);       printf("\tminFitness = %f\n", minFitness);       printf("\tCrossovers = %d\n", curCrossovers);       printf("\tMutation  = %d\n", curMutations);       printf("\tpercentage = %f\n", avgFitness / maxFitness);     }     /* Check the diversity of the population after _ of the      * generations has completed. If the population has prematurely      * converged, exit out to allow the user to restart.      */     if ( generation > (MAX_GENERATIONS * 0.25) ) {       if ((avgFitness / maxFitness) > 0.98) {         printf("converged\n");         break;       }     }     if (maxFitness == MAX_FIT) {       printf("found solution\n");       break;     }   }   /* Emit final statistics */   printf("Generation %d\n", generation-1);   printf("\tmaxFitness = %f (%g)\n", maxFitness, MAX_FIT);   printf("\tavgFitness = %f\n", avgFitness);   printf("\tminFitness = %f\n", minFitness);   printf("\tCrossovers = %d\n", curCrossovers);   printf("\tMutation = %d\n", curMutations);   printf("\tpercentage = %f\n", avgFitness / maxFitness);   /* Emit the highest fit chromosome from the population */   for (i = 0 ; i < MAX_CHROMS ; i++) {     if (populations[curPop][i].fitness == maxFitness) {       int index;       printf("Program %3d : ", i);       for (index = 0 ; index < populations[curPop][i].progSize ;            index++) {         printf("%02d ", populations[curPop][i].program[index]);       }       printf("\n");       printf("Fitness %f\n", populations[curPop][i].fitness);       printf("ProgSize %d\n\n", populations[curPop][i].progSize);       printProgram(i, curPop);       break;     }   }   return 0; } 

 

A summarized flow of the main function is as follows. We begin by seeding the random number generator with the srand call. We also open a file to store the fitness information (for which we'll see a plot in the next section). The population is then initialized with random chromosomes ( initpopulation ), and the fitness check is performed on the population ( performFitnessCheck ). This fitness check is done now because the fitness of each chromosome is necessary for selection ( performSelection ). This is because we're using probabilistic selection based upon the fitness value. Recombination is not shown specifically in this listing, as it is performed with the selection process.

Once selection has been performed, we switch the variable called curPop. This variable defines the current population. The population of chromosomes is defined as a 2-dimensional array that defines two populations (symbolically, old and new). See Listing 6.3 for this declaration.

Listing 6.3: Declaration of Populations.

 typedef struct population {   float fitness;   int   progSize;   int   program[MAX_PROGRAM]; } POPULATION_TYPE; POPULATION_TYPE populations[2][MAX_CHROMS]; int curPop; 

 

A population is made up of a number of chromosomes (as defined by MAX_CHROMS ), where each chromosome is a program , program size ( progSize ), and fitness . The actual populations declaration shows two populations, as discussed before representing the old and new. The curPop variable defines which population is currently under view; the alternate is the target for any changes. When recombination occurs, we select parents from the current population (as defined by curPop ) and introduce the children into the alternate population (as defined by ! curPop ).

Returning to the main function, the main loop (while loop of generations) simply performs selection, checks the fitness, and continues (swapping populations each iteration). The remainder of the loop emits debugging information and checks for a couple of exit conditions. If the average fitness is 98% of the maximum fitness, the population has converged and we exit. If the maximum fitness is the actual maximum allowable by the objective function, we also exit since we've found a solution to the problem at hand. Once a solution has been found, we emit more information to the user and display the highest fit program from the pool.

Initialization

Initializing the population is a simple process. We simply walk through the available chromosomes of the population, and assign each of the genes with a random instruction. We also initialize the fitness of each chromosome to zero and set the program size to the maximum (see Listing 6.4).

Listing 6.4: Initializing the Population.

 void initMember( pop, index ) {   int progIndex;   populations[pop][index].fitness = 0.0;   populations[pop][index].progSize = MAX_PROGRAM-1;   /* Randomly create a new program */   progIndex = 0;   while (progIndex < MAX_PROGRAM) {     populations[pop][index].program[progIndex++] =                                      getRand(MAX_INSTRUCTION);   } } void initPopulation( void ) {   int index;   /* Initialize each member of the population */   for (index = 0 ; index < MAX_CHROMS ; index++) {     initMember(curPop, index);   } } 

 

The initPopulation function is called from our main function, which in turn calls the initMember function to initialize each chromosome within the population.

Evaluating the Fitness

Fitness evaluation is performed through the performFitnessCheck function (shown in Listing 6.5).

Listing 6.5: Evaluating the Fitness of the Population.

 float maxFitness; float avgFitness; float minFitness; extern int stackPointer; extern int stack[]; static int x = 0; float totFitness; int performFitnessCheck( FILE *outP ) {   int chrom, result, i;   int args[10], answer;   maxFitness = 0.0;   avgFitness = 0.0;   minFitness = 1000.0;   for ( chrom = 0 ; chrom < MAX_CHROMS ; chrom++ ) {   populations[curPop][chrom].fitness = 0.0;   for ( i = 0 ; i < COUNT ; i++ ) {     args[0] = (rand() & 0x1f) + 1;     args[1] = (rand() & 0x1f) + 1;     args[2] = (rand() & 0x1f) + 1;     /* Sample Problem: x^3 + y^2 + z */     answer = (args[0] * args[0] * args[0]) +              (args[1] * args[1]) + args[2];     /* Call the virtual stack machine to check the program */     result = interpretSTM(populations[curPop][chrom].program,                           populations[curPop][chrom].progSize,                           args, 3);     /* If no error occurred, add this to the fitness value */     if (result == NONE) {       populations[curPop][chrom].fitness += TIER1;     }     /* If only one element is on the stack, add this to the       fitness      * value.      */     if (stackPointer == 1) {       populations[curPop][chrom].fitness += TIER2;     }     /* If the stack contains the correct answer, add this to the      * fitness value.      */     if (stack[0] == answer) {       populations[curPop][chrom].fitness += TIER3;     }   }   /* If this chromosome exceeds our last highest fitness, update    * the statistics for this new record.    */   if (populations[curPop][chrom].fitness > maxFitness) {     maxFitness = populations[curPop][chrom].fitness;     } else if (populations[curPop][chrom].fitness < minFitness) {       minFitness = populations[curPop][chrom].fitness;     }     /* Update the total fitness value (sum of all fitnesses) */     totFitness += populations[curPop][chrom].fitness;   }   /* Calculate our average fitness */   avgFitness = totFitness / (float)MAX_CHROMS;'   if (outP) {     /* Emit statistics if we have an output file pointer */     fprintf(outP, "%d %6.4f %6.4f %6.4f\n",              x++, minFitness, avgFitness, maxFitness);   }   return 0; } 

 

Function performFitnessCheck walks through every chromosome in the current population and evaluates its fitness. The fitness is then stored with the chromosome in the populations structure.

The algorithm begins by clearing out the fitness measures in preparation for evaluating the current population. A loop is then defined to walk through each chromosome. Once the current fitness is cleared, another loop walks through the tests for the chromosome. To avoid a chromosome providing the correct answer, but working in only one case, we test the chromosome a number of times (in this case 10 times). The args array contains the arguments for the problem. This array is loaded onto the stack upon a call to the evaluation function. In each case, the args array is loaded with a random set of parameters. This is to ensure that it actually solves the problem at hand.

Once the arguments for the problem have been defined, we calculate the value that should result and load it into the answer variable. We then call the function interpretSTM (interpret stack machine) to evaluate the fitness (see Listing 6.1). From the result of the evaluation function, we calculate the fitness for the chromosome. The fitness is based upon a three-tiered measure. If the program exited successfully (no math error or program error), then we give it a TIER1 value. If only one value was left on the stack, we add in TIER2. Finally, if the top of the stack is the correct value (expected value of answer ) we add in TIER3 . The tier values are defined such that TIER3 is greater than (more important than) TIER2 , and so on for TIER1 . This measure gives the genetic algorithm information to work with during the selection process. The genetic algorithm works best with gradual improvement of the calculated to expected response of the chromosome, this is not quite the case with this problem but it still works.

Once we finish checking the fitness for the current chromosome, we check to see if it is either the highest or lowest fit chromosome. This indication is stored within the maxFitness and minFitness variables accordingly . We then calculate the totFitness and once the loop to check all of the chromosomes is complete, we calculate the average fitness of the population (stored in avgFitness ).

Finally, we store the current generation's information to our output file. This can later be used to plot the fitnesses of the populations to monitor their progress.

Selection and Recombination

During selection and recombination, we select parents from the current population to recombine as the children in the next. This is performed very succinctly in Listing 6.6.

Listing 6.6: Performing Selection.

 int performSelection( void ) {   int par1, par2;   int child1, child2;   int chrom;   /* Walk through the chromosomes, two at a time */   for (chrom = 0 ; chrom < MAX_CHROMS ; chrom+=2) {     /* Select two parents, randomly */     par1 = selectParent();     par2 = selectParent();     /* The children are loaded at the current index points */     child1 = chrom;     child2 = chrom+1;     /* Recombine the parents to the children */     performReproduction( par1, par2, child1, child2 );]   }   return 0; } 

 

This function operates from the indices into the population tables. The variables par1 and par2 are indexes into the current population, while child1 and child2 are indexes into the new population. The algorithm simply selects child indexes starting at zero and incrementing by two each count. The parent indexes are selected using the selectParent function. Once we have four index values, the performReproduction function is called to perform the actual recombination. The selectParent call is shown in Listing 6.7.

Listing 6.7: Selecting a Parent.

 int selectParent( void ) {   static int chrom = 0;   int ret = -1;   float retFitness = 0.0;   /* Roulette-wheel selection process */   do {     /* Select the target fitness value */     retFitness = (populations[curPop][chrom].fitness /                   maxFitness);     if (chrom == MAX_CHROMS) chrom = 0;     /* If we've walked through the population and have reached our      * target fitness value, select this member.      */     if (populations[curPop][chrom].fitness > minFitness) {       if (getSRand() < retFitness) {         ret = chrom++;         retFitness = populations[curPop][chrom].fitness;         break;       }     }     chrom++;   } while (1);   return ret; } 

 

Parent selection operates on the principal that a chromosome's chances of being selected are proportional to that chromosome's fitness compared to the overall population. We calculate this probability (stored within retFitness ) at the start of the do loop. We then check to ensure that the current chromosome is greater than the minimum fitness of the population. In other words, we don't select from the least fit of the population (elitist in some manner). We then calculate a random number (from 0 to 1) and compare this to our fitness value. If the random number is less than our fitness value, we select the parent and allow the function to return. Otherwise, we continue to the top of the loop and calculate the fitness ratio for the next chromosome.

Once we've selected both parents, we perform recombination. This process is shown in Listing 6.8.

Listing 6.8: Recombining Parent Chromosomes to Create Two New Children.

 int performReproduction( int parentA, int parentB,                          int childA, int childB ) {   int crossPoint, i;   int nextPop = (curPop == 0) ? 1 : 0;   int mutate( int );   /* If we meet the crossover probability, perform crossover of the    * two parents by selecting the crossover point.    */   if (getSRand() > XPROB) {     crossPoint =          getRand(MAX(populations[curPop][parentA].progSize-2,                      populations[curPop][parentB].progSize-2))+1;     curCrossovers++;   } else {     crossPoint = MAX_PROGRAM;   }   /* Perform the actual crossover, in addition to random mutation     */   for (i = 0 ; i < crossPoint ; i++) {     populations[nextPop][childA].program[i] =         mutate(populations[curPop][parentA].program[i]);     populations[nextPop][childB].program[i] =         mutate(populations[curPop][parentB].program[i]);   }   for ( ; i < MAX_PROGRAM ; i++) {     populations[nextPop][childA].program[i] =         mutate(populations[curPop][parentB].program[i]);     populations[nextPop][childB].program[i] =         mutate(populations[curPop][parentA].program[i]);   }   /* Update the program sizes for the children (based upon the    * parents).    */   populations[nextPop][childA].progSize =     populations[curPop][parentA]   progSize; populations[nextPop][childB].progSize =     populations[curPop][parentB].progSize;   return 0; } 

 

We first check to see if we are to perform the crossover operator (random number greater than XPROB ). If so, we calculate the crossover point ( crossPoint ) based upon the maximum length of the chromosome. This is based upon the larger of the two parents. The chromosomes may be different sizes (though in this simulation, they are always the maximum). The crossover point is ensured never to be the first or last gene of the chromosome (since then no crossover is actually occurring). If no crossover is to be performed, the crossover point is defined as the size of the program.

The next step is to perform the crossover (or tail swapping) of the two chromosomes. We copy parent A to child A and parent B to child B, up until the crossover point. At this point (the next for loop), we copy parent A to child B and parent B to child A. Note that if the crossover point was the size of the program, the first for loop would have copied the entire chromosome and the second for loop would have nothing to copy (the desired result in this case).

Finally, the children in the new population inherit the program sizes of their parents. In this case, the program size is constant, but this can easily be changed during the init process.

Note that we create a variable called nextPop , which is based upon curPop . Variable nextPop defines the population that the children will be created in; curPop defines where the parents are read. This provides a ping-pong population buffer from old population to new.

The final function is mutate , shown in Listing 6.9. This function simply redefines a gene in the current chromosome to a new instruction, based upon the mutation probability.

Listing 6.9: Mutation Operator.

 int mutate(int gene) {   float temp = getSRand();   /* If we've met the mutation probability, randomly mutate this    * gene to a new instruction.    */   if (temp > MPROB) {       gene = getRand(MAX_INSTRUCTION);       curMutations++;   }   return gene; }