simple-gp.c

其中提到遺傳學的程式碼與應用提供給次淚相向的研究者參考下載

/*   simple-gp.c
// Simple Genetic Program
//
// W. B. Langdon, cwi.nl, 21 February 1999
//
// Modifications:
// Author:  WBL
// Date  :  21 Feb 1999
// Change:  Increase portability by replacing drand48() by park-miller, 
//          and renaming functions to avoid name clashes

// Author:  WBL
// Date  :  18 Feb 1994, 21 April 1994
// Change:  Replace rand() with drand48(),
//          reduce POPULATION_SIZE from 500 to 16
*/
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <time.h>
#include <float.h>
#include <assert.h>

/* SUNos v4 */
#define EXIT_FAILURE 1
#define FILENAME_MAX 80

#define POPULATION_SIZE  16
#define PCROSS            0.6
#define PMUT               .01
#define MAX_GEN          50
#define MAX_PROG_SIZE    25
#define STACK_SIZE       25
#define LINE_WIDTH      200
#define THRESHOLD          0.01

#define MAX_TEST_POINTS   10

float x[MAX_TEST_POINTS];
float y[MAX_TEST_POINTS];
int   num_x;
char output_file[FILENAME_MAX];

struct prog
{
  char string [MAX_PROG_SIZE+1];
  float fitness;
  int  hits;
};

struct prog     population [POPULATION_SIZE];
struct prog new_population [POPULATION_SIZE];

/*Data for display purposes only */
#define COPY_MUM -1
#define MUTATE   -1
struct par
{
  int mum, mumcross;
  int dad, dadcross;
};
struct par parents [POPULATION_SIZE];
int num_tests  =  0;

int best_prog, max_hits, required_hits;

float new_sumfitness;
float total_fitness;

int generation;


/* start code */

int    seed;

read_inputs()
/*read pairs of numbers from stdin
//
//outputs:    required_hits
//            x, y, num_x
*/
{
int i;
time_t now;

  printf ("Enter up to %d x y pairs of numbers> ", MAX_TEST_POINTS );
  fflush(stdout);

  for (i = 0; i < MAX_TEST_POINTS; i++ )
    { if ( scanf("%e %e", &x[i], &y[i]) != 2 ) break; }

  if (i <= 0) {printf ("\nNo data pairs! Stopping!\n"); exit (EXIT_FAILURE);}

  num_x = i;
  required_hits = num_x;
  strcpy (output_file, "gp.c");


  printf ("\nGP regression on %d points\n", num_x);
  for( i=0; i < num_x; i++ )
  {
    printf("(%e,%e) ",x[i],y[i]);
    if ( (i%2) == 1 ) printf("\n");
  }

  time (&now);
  seed = now;

  printf ("\nPop = %d, PCross %.5f, PMutate %.5f, Max gens %d, seed %d\n",
          POPULATION_SIZE, PCROSS, PMUT, MAX_GEN, seed );

  printf ("\nGP starting\n");

} /*end read_inputs() */


/*park-miller.cc   Park-Miller random number generator
 W. Langdon cs.ucl.ac.uk 5 May 1994*/
int intrnd () /* 1<=seed<=m */
{
double const a    = 16807;      /* ie 7**5 */
double const m    = 2147483647; /* ie 2**31-1 */

        double temp = seed * a;
        seed = (int) (temp - m * floor ( temp / m ));
        return seed;
}
float rand_0to1() { return (float)(intrnd()) / 2147483647.0; }

float random_number()  /* single precision ok */
{
float r;

for (; (r = rand_0to1() ) >= 1.0; )
{printf("random number >=1.0 %e\n",r);} /* ensure we dont return 1.0 */

return ( r );

}/* end random_number()*/

int random_point ( char string[] )
{
/*return a random point in the string, ie from pos zero to length of string -1
*/

  return  ((int) (random_number() * (float) (strlen(string) )  ));

}/*end random_point ()*/


int matching_point ( char string[], int subtree )
/*Return: Point in string immediately after subtree starting at subtree*/
{
int i;
int dangling_limbs;

if ( string[subtree] == '\0' )
  dangling_limbs = 0;         /*already at end of tree */
else
  dangling_limbs = 1;

for (i = subtree; dangling_limbs > 0; i++ )
  {
      switch (string [i])
      {
        case '+':
        case '-':
        case '*':
        case '/':         dangling_limbs++; /*function (all have two limbs)*/
                          break;
      
        case '\0':        assert (1 == 0);  /*oh dear*/
                          break;

        default :         dangling_limbs--; /*terminal*/
                          break;
      }
  }/*end for*/

return ( i );

}/*end matching_point()*/



int random_function()
{
#define N_FUNCTIONS 4
int function [ N_FUNCTIONS ] = { '+', '-', '*', '/'};

   return ( function [ (int) (random_number() * N_FUNCTIONS)] );

}/*end random_function()*/


int random_terminal()
{
   return ( 'x' );

}/*end random_terminal()*/


int random_tree (char buff[], int size)
/*generate random program, place in buff
//
//Inputs: size of buff
//
//Returns: length of tree
*/
{
int dangling_limbs = 1;
int i;

  for ( i = 0; (dangling_limbs > 0) && (i <= size) ; i++ )
  {  /*chose function or terminal at random*/

     if( random_number() > (float)(dangling_limbs*dangling_limbs+1)/(float)(size-i) )
       {
         buff[i] = random_function();
         dangling_limbs++;                 /*all operators have two limbs*/
       }
     else
       {
         buff[i] = random_terminal();
         dangling_limbs--;
       }

  }/*end for*/

  assert ( dangling_limbs == 0 ); /*Opps...*/

  return (i);

}/*end random_tree()*/

int splice  ( char output[], int outsize,
              char buff1[],  int end1,
              char buff2[],  int end2,
              char buff3[],  int end3 )
/*returns: 0 if ok*/
{
  if (( end1 + end2 + end3 ) >= outsize ) return ( 1 );

  memcpy ( output,             buff1, end1 );
  memcpy ( &output[end1],      buff2, end2 );
  memcpy ( &output[end1+end2], buff3, end3 );

  output [end1+end2+end3] = '\0'; /*terminate string*/

  return ( 0 ); /*ie ok*/

}/*end splice()*/


mutate(int mum, int child)
/*make a single random change to mum from old population and store
//child in new_population */
{
int mum_end1;
int mum_start2, mum_end2;
int mut_length;

int sanity = 0;

char buffer [MAX_PROG_SIZE];

do {
  assert (sanity++ < 1000);

  mum_end1   = random_point ( population[mum].string );
  mum_start2 = matching_point ( population[mum].string, mum_end1 ); 
  mum_end2   = strlen ( population[mum].string );

  mut_length = random_tree ( buffer, MAX_PROG_SIZE );

} while (splice ( 
         new_population[child].string, MAX_PROG_SIZE+1,
         population[mum].string, mum_end1,
         buffer, mut_length,
         &population[mum].string[mum_start2], mum_end2 - mum_start2 ) != 0 );

/*display data */

parents[child].mum       = mum;
parents[child].mumcross  = mum_end1;
parents[child].dad       = MUTATE;
parents[child].dadcross  = mut_length;


}/*end mutate()*/


int   stack;
float  value_stack [STACK_SIZE];
char      op_stack [STACK_SIZE];
char  output_stack [STACK_SIZE][LINE_WIDTH];


push_operator ( char op )
{
  op_stack [stack++]      = op;
  value_stack [stack]     = 0.0;   /* keep it clean */
  output_stack [stack][0] = '\0';  /* and tidy      */

}/*end push_operator()*/


push_value ( float value )
{
float a;
float result;

if ( (stack <= 0) || (op_stack[stack-1] != 0) )  /* operator on top of stack */
  {                                              /* or stack empty           */
    op_stack [ stack ]    = 0; /* operand on top of stack */
    value_stack [stack++] = value;
  }
else
  {
    a = value_stack[--stack];  /* pop first operand */
    switch (op_stack[--stack])   /* pop operator */
      {
        case '+':  result = a + value;
                   break;

        case '-':  result = a - value;
                   break;

        case '*':  result = a * value;
                   break;

        case '/':  if (value == 0.0)
                     result = 1.0;
                   else
                     result = a / value;
                   break;

        default:
                   assert (10 == 0 ); /* opps */
                   break;
      }

    push_value ( result );
  }

}/*end push_value()*/

float prog_value (char string[], float input )
/* return value of program child in new_population with input*/
{
int i;

stack = 0;

for ( i = 0; string[i] != 0; i++ )
  {
      switch (string [i])
      {
        case '+':
        case '-':
        case '*':
        case '/':         push_operator(string[i]);
                          break;
      
        default:          push_value ( input );
                          break;
      }
  }/*end for*/

assert (stack == 1);

/*printf("prog_value: %s = %e, input = %e\n", string, value_stack[0], input );
*/

return (value_stack[--stack]); /*pop value*/

} /*end prog_value()*/



update_hits_etc ( float fit, int n_hits, int child )
/*
//store fitness and hits of child in new_population
//update best program statistics.
//