home *** CD-ROM | disk | FTP | other *** search
/ NetNews Usenet Archive 1992 #19 / NN_1992_19.iso / spool / comp / lang / cplus / 13066 < prev    next >
Encoding:
Text File  |  1992-08-30  |  10.2 KB  |  536 lines

  1. Path: sparky!uunet!charon.amdahl.com!pacbell.com!sgiblab!darwin.sura.net!jvnc.net!rutgers!sun-barr!olivea!hal.com!decwrl!access.usask.ca!kakwa.ucs.ualberta.ca!alberta!ubc-cs!news.UVic.CA!sol.UVic.CA!schumach
  2. From: schumach@sol.UVic.CA (Ian Schumacher)
  3. Newsgroups: comp.lang.c++
  4. Subject: c++ genetic algorithm example
  5. Summary: example genetic algorithm program in c++
  6. Keywords: gnetic
  7. Message-ID: <schumach.715194212@sol.UVic.CA>
  8. Date: 30 Aug 92 17:03:32 GMT
  9. Sender: news@sol.UVic.CA
  10. Distribution: all
  11. Organization: University of Victoria
  12. Lines: 521
  13. Nntp-Posting-Host: sol.uvic.ca
  14.  
  15.  
  16. I have seen a lot of interest in Gene algorithms these days.  I was
  17. pretty interested myself, so I threw together a Gene and Gene_Pool
  18. Class in C++.  Figured I send it out, so that anyone that was interested
  19. could have a example to work with, improve upon, whatever.
  20. Bugs?  You tell me.  It hasn't been very thoroughly tested, but there's
  21. not that much that could go wrong.  It's in three parts gene.h gene.m and
  22. imp.m.  Enjoy.
  23.  
  24. --------------------------------------------------------------------
  25.  
  26. // class gene.h: gene structure
  27. // Programmer : Ian Schumacher
  28. // Comments welcome
  29. // email: schumach@sol.UVic.CA
  30. // Genie : I.SCHUMACHER
  31.  
  32. class Gene
  33. {
  34.     unsigned short *gene;
  35.  
  36. protected:
  37.  
  38.     int length;
  39.     double mutate_rate;
  40.     double value;
  41.  
  42. public:
  43.     Gene(int);
  44.     Gene();
  45.     ~Gene();
  46.     Gene& operator+(Gene &);
  47.     unsigned short& operator[](int);
  48.     Gene& operator=(const Gene &);
  49.     Gene& operator=(const double);
  50.     Gene& mutate(void);
  51.     void set_rate(double);
  52.     int get_length(void);
  53.     double get_value(void) {return value;}
  54. };
  55.  
  56. class Gene_Pool
  57. {
  58.     int size;
  59.     int plength;
  60.     Gene *gene;
  61.     int *index;
  62.     
  63. public:
  64.     Gene_Pool(int number=1,int a=1);
  65.     ~Gene_Pool() {delete gene;delete index;}
  66.     Gene_Pool& operator+(Gene_Pool &);
  67.     Gene& operator[](int);
  68.     Gene_Pool& operator=(const Gene_Pool &);
  69.     Gene_Pool& mutate(void);
  70.     int get_size(void) { return size;}
  71.     int get_length(void) {return plength;}
  72.     Gene_Pool& evolve(void);
  73.     Gene_Pool& sort(void);
  74.     Gene_Pool& kill(int);
  75.     Gene_Pool& set_rate(double);
  76. };
  77.  
  78. ---------------------------------------------------------------------------
  79.  
  80. // gene.m : Definition of Gene and Gene_Pool functions
  81. // Programmer : Ian Schumacher
  82. // Comments welcome
  83. // email: schumach@sol.UVic.CA
  84. // Genie : I.SCHUMACHER
  85.  
  86. #include "gene.h"
  87. #include <stdio.h>
  88. #include <stdlib.h>
  89. #include <math.h>
  90.  
  91. #define RATE 0.01            // default mutation rate
  92. #define RAND_MATE 1        // mate randomly or in order switch
  93.  
  94. void indexx(int ,double *,int *);
  95. int irbit2(void);
  96.  
  97. int pop_length=1;        // gene population value set by Gene_Pool, used by Gene
  98. int iseed=-666;               // random generator seed
  99.  
  100. Gene::Gene(int n)
  101. {
  102.     int i;
  103.     
  104.     if(n>0)
  105.     {
  106.         gene=new unsigned short[n];
  107.         mutate_rate=RATE;
  108.         length=n;
  109.     }
  110.     else
  111.     {
  112.         gene=new unsigned short[1];
  113.         mutate_rate=RATE;
  114.         length=1;
  115.     }
  116.     
  117. // insert random values into gene
  118.     for(i=0;i<length;i++)
  119.         gene[i]=(unsigned short) irbit2();  // don't trust rand() that much, use "NR in C" function
  120. }
  121.  
  122. Gene::Gene()
  123. {
  124. // uses global variable pop_length set by Gene_Pool to determine gene size
  125.     int i;
  126.     
  127.     if(pop_length<1)
  128.         pop_length=1;
  129.     length=pop_length;
  130.     mutate_rate=0.01;
  131.     gene=new unsigned short[pop_length];
  132. // insert random values into gene
  133.     for(i=0;i<pop_length;i++)
  134.         gene[i]=(unsigned short) (irbit2());
  135. }
  136.  
  137. Gene::~Gene()
  138. {
  139.     delete gene;
  140. }
  141.  
  142. Gene& Gene::operator+(Gene &a)
  143. {
  144.     int cp,i;
  145.     int tmp;
  146.  
  147.     static Gene b(length);    //static because its adress is returned
  148.  
  149.     if(this==&a)
  150.         return *this;
  151.     if(length!=a.length)
  152.     {
  153.         printf("Currently cannot mate genes of different lengths. Genes unchanged\n");
  154.         return *this;
  155.     }
  156.     // mating routine
  157.     cp=int(rand()%length);   // Random cross point. I thought this would produce better result
  158.                         // constant half way cross point.
  159.     for(i=0;i<cp;i++)
  160.         b.gene[i]=this->gene[i];
  161.     for(i=cp;i<length;i++)
  162.         b.gene[i]=a.gene[i];
  163.     return b;
  164. }
  165.  
  166. unsigned short& Gene::operator[](int elem)
  167. // Returning address allows gene[j] to be on either side of equals sign
  168. {
  169.     static unsigned short tmp=0;
  170.     
  171.     if(elem>=0&&elem<length)
  172.         return gene[elem];
  173.     else
  174.     {
  175.         printf("Gene access error\n");
  176.         return tmp;
  177.     }
  178. }
  179.  
  180. Gene& Gene::operator=(const Gene& a)
  181. // set one gene equal to another
  182. {
  183.     int i;
  184.  
  185.     if(this==&a)
  186.         return *this;
  187.     length=a.length;
  188.     mutate_rate=a.mutate_rate;
  189.     value=a.value;
  190.     delete gene;
  191.     gene=new unsigned short[length];
  192.     for(i=0;i<length;i++)
  193.         gene[i]=a.gene[i];
  194.     return *this;
  195. }
  196.  
  197. Gene& Gene::operator=(const double x)
  198. // set gene value
  199. {
  200.     value=x;
  201.     return *this;
  202. }
  203.  
  204. Gene& Gene::mutate(void)
  205. {
  206.     int i;
  207.     
  208.     for(i=0;i<length;i++)
  209.     {
  210.         if(rand()/(RAND_MAX+1.0)<mutate_rate)
  211.         {
  212.             gene[i]=gene[i]<1 ? 1:0;
  213.         }
  214.     }
  215.     return *this;
  216. }
  217.  
  218. void Gene::set_rate(double r)
  219. {
  220.     if(r>0.0)
  221.         mutate_rate=r;
  222.     return;
  223. }
  224.  
  225. int Gene::get_length(void)
  226. {
  227.     return length;
  228. }
  229.  
  230. // Gene_Pool definitions
  231.  
  232. Gene_Pool::Gene_Pool(int number,int a)
  233. // uses global variable pop_length to set gene size
  234. {
  235.     int i;
  236.     
  237.     pop_length=a;
  238.     plength=a;
  239.     size=number;
  240.     index=new int[size];
  241.     for(i=0;i<size;i++)
  242.         index[i]=i;
  243.     gene=new Gene[size];
  244. }
  245.  
  246. Gene_Pool& Gene_Pool::operator+(Gene_Pool &a)
  247. {
  248.     int i;
  249.     
  250.     static Gene_Pool b(size+a.size,plength);
  251.     
  252.     if(this==&a)
  253.     {
  254.         delete &b;
  255.         return *this;
  256.     }
  257.     if(plength!=a.plength)
  258.     {
  259.         printf("incompatible genes\n");
  260.         delete &b;
  261.         return *this;
  262.     }
  263.     for(i=0;i<a.size;i++)
  264.         b.gene[i]=gene[i];
  265.     for(i=size;i<b.size;i++)
  266.         b.gene[i]=a.gene[i-size];
  267.     return b;
  268. }
  269.  
  270. Gene& Gene_Pool::operator[](int elem)
  271. {
  272.     if(elem>=0&&elem<size)
  273.         return gene[index[size-elem-1]];
  274.     else
  275.     {
  276.         printf("Gene access error\n");
  277.         return gene[0];
  278.     }
  279. }
  280.  
  281. Gene_Pool& Gene_Pool::operator=(const Gene_Pool& a)
  282. {
  283.     int i;
  284.  
  285.     if(this==&a)
  286.         return *this;
  287.     size=a.size;
  288.     pop_length=a.plength;
  289.     plength=a.plength;
  290.     delete gene;
  291.     delete index;
  292.     gene=new Gene[size];
  293.     index=new int[size];
  294.     for(i=0;i<pop_length;i++)
  295.     {
  296.         gene[i]=a.gene[i];
  297.         index[i]=a.index[i];
  298.     }
  299.     return *this;
  300. }
  301.  
  302. Gene_Pool& Gene_Pool::mutate(void)
  303. {
  304.     int i;
  305.  
  306.     for(i=0;i<size;i++)
  307.         gene[i].mutate();
  308.     return *this;
  309. }
  310.  
  311. Gene_Pool& Gene_Pool::set_rate(double r)
  312. {
  313.     int i;
  314.  
  315.     for(i=0;i<size;i++)
  316.         gene[i].set_rate(r);
  317.     return *this;
  318. }
  319.  
  320. Gene_Pool& Gene_Pool::sort(void)
  321. {
  322.     int i,j;
  323.     double *values;
  324.  
  325.     values=new double[size];
  326.  
  327.     for(i=0;i<size;i++)
  328.         values[i]=gene[i].get_value();
  329.     indexx(size,values,index);    // "NR in C" heap sort function.  Sorts lowest to highest
  330.     delete values;
  331.     return *this;
  332. }
  333.  
  334. Gene_Pool& Gene_Pool::evolve(void)
  335. // allows for random mating of top half or ordered mating (you decide which is better)
  336. {
  337.     int i;
  338.  
  339.     this->sort();
  340.     if(RAND_MATE)
  341.     {
  342.         // randomly mate top half of population and store in bottom half
  343.         for(i=0;i<size/2;i++)
  344.         {
  345.             gene[index[i]]=gene[index[size-rand()%size-1]]+gene[index[size-rand()%size-1]];
  346.         }
  347.         this->mutate();
  348.     }
  349.     else
  350.     {
  351.         // or mate top half of population together in order
  352.         for(i=0;i<size/2;i++)
  353.         {
  354.             gene[index[i]]=gene[index[size-i-1]]+gene[index[size-i-2]];
  355.         }
  356.     }
  357.     return *this;
  358. }
  359.  
  360. Gene_Pool& Gene_Pool::kill(int n)
  361. // if after evolve for a while you want to kill off the weak to speed things up
  362. {
  363.     int i;
  364.     Gene *tmp;
  365.  
  366.     tmp=new Gene[size-n];
  367.     if(n>=size)
  368.         n=size-1;
  369.     size=size-n;
  370.     for(i=0;i<size;i++)
  371.         tmp[i]=gene[index[size+n-i]-1];
  372.     delete gene;
  373.     delete index;
  374.     gene=new Gene[size];
  375.     index=new int[size];
  376.     for(i=0;i<size;i++)
  377.         gene[i]=tmp[i];
  378.     this->sort();
  379.     delete tmp;
  380.     return *this;
  381. }
  382.  
  383. // "Numerical Recipes is C" heapsort of arrin with corresponding reordering of indx
  384. void indexx(int n,double *arrin,int *indx)
  385. {
  386.     int l,j,ir,indxt,i;
  387.     float q;
  388.  
  389.     indx--;
  390.     for (j=1;j<=n;j++) indx[j]=j-1;
  391.     l=(n >> 1) + 1;
  392.     ir=n;
  393.     for (;;)
  394.     {
  395.         if (l > 1)
  396.             q=arrin[(indxt=indx[--l])];
  397.         else
  398.         {
  399.             q=arrin[(indxt=indx[ir])];
  400.             indx[ir]=indx[1];
  401.             if (--ir == 1)
  402.             {
  403.                 indx[1]=indxt;
  404.                 indx++;
  405.                 return;
  406.             }
  407.         }
  408.         i=l;
  409.         j=l << 1;
  410.         while (j <= ir)
  411.         {
  412.             if (j < ir && arrin[indx[j]] < arrin[indx[j+1]]) j++;
  413.             if (q < arrin[indx[j]])
  414.             {
  415.                 indx[i]=indx[j];
  416.                 j += (i=j);
  417.             }
  418.             else j=ir+1;
  419.         }
  420.         indx[i]=indxt;
  421.     }
  422. }
  423.  
  424. #define IB1 1
  425. #define IB2 2
  426. #define IB5 16
  427. #define IB18 131072
  428. #define MASK IB1+IB2+IB5
  429.  
  430. // "Numerical Recipes in C" random bit function (uses periodic maximal length sequences)
  431. int irbit2(void)
  432. {
  433.     if (
  434.     iseed & IB18) {
  435.         iseed=((iseed ^ MASK) << 1) | IB1;
  436.         return 1;
  437.     } else {
  438.         iseed <<= 1;
  439.         return 0;
  440.     }
  441. }
  442.  
  443. #undef MASK
  444. #undef IB18
  445. #undef IB5
  446. #undef IB2
  447. #undef IB1
  448.  
  449. --------------------------------------------------------------------------
  450.  
  451. // imp.m : Example implementation of Gene_Pool class
  452. // Programmer : Ian Schumacher
  453. // Comments welcome 
  454. // email : schumach@sol.UVic.CA
  455. // Genie : I.SCHUMACHER
  456.  
  457. // Pathetically unrealistic problem to be solved by using gene algorithms,
  458. // but illustrates principles.
  459.  
  460. #include "gene.h"
  461. #include <stdio.h>
  462. #include <stdlib.h>
  463. #include <math.h>
  464.  
  465. #define size 500        // Gene_Pool population size
  466. #define length 30        // Gene length
  467. #define steps 100        // maximum number of evolution steps
  468. #define sqr(x) x*x        // simple squaring function
  469.  
  470. int main()
  471. {
  472.     int i,j,k,t;
  473.     Gene_Pool a(size,length);
  474.     double value1,value2,value3;
  475.     double tmp1,tmp2;
  476.     int tmp;
  477.  
  478.     a.set_rate(0.05/length);        // set mutation rate, (in this case to 5% of population)
  479.  
  480.     for(i=0;i<steps;i++)        // number of evolution steps
  481.     {
  482.         for(j=0;j<a.get_size();j++)
  483.         {
  484.             // this is where the value of gene is determined
  485.             value1=0;
  486.             value2=0;
  487.             value3=0;
  488.             tmp=a[j].get_length()/3;
  489.             for(t=0;t<tmp;t++)
  490.             {
  491.                 // interperate gene as three long integers
  492.                 value1+=(1L<<t)*a[j][t];
  493.                 value2+=(1L<<t)*a[j][t+tmp];
  494.                 value3+=(1L<<t)*a[j][t+2*tmp];
  495.             }
  496.             // normalize values to 1
  497.             value1/=1L<<tmp;
  498.             value2/=1L<<tmp;
  499.             value3/=1L<<tmp;
  500.             // maximize simple nonlinear equation
  501.             a[j]=exp(-sqr((value1-0.33))-sqr((value2-0.5))-sqr((value3-0.7)));  // set gene[j] value
  502.         }
  503. // break out of evolution early if genes starting to look the same
  504.         if((a[0].get_value()==a[1].get_value())&&a[1].get_value()==a[2].get_value()&&tmp1==tmp2)
  505.             break;
  506.         tmp1=a[0].get_value();
  507.         tmp2=a[1].get_value();
  508.         a.evolve();
  509.         // display sample of every fifth generation
  510.         if(!(i%5))
  511.         {
  512.             for(j=0;j<5;j++)
  513.             {
  514.                 for(k=0;k<length;k++)
  515.                 {
  516.                     printf("%d ",a[j][k]);
  517.                 }
  518.                 printf("%f\n",a[j].get_value());
  519.             }
  520.             printf("\n");
  521.         }
  522.     }
  523.     // display final result
  524.     for(t=0;t<tmp;t++)
  525.     {
  526.         value1+=(1L<<t)*a[0][t];
  527.         value2+=(1L<<t)*a[0][t+tmp];
  528.         value3+=(1L<<t)*a[0][t+2*tmp];
  529.     }
  530.     value1/=1L<<tmp;
  531.     value2/=1L<<tmp;
  532.     value3/=1L<<tmp;
  533.     printf("value1 = %f value2= %f value3= %f\n",value1,value2,value3);
  534.     return 1;
  535. }
  536.