ga.cpp

该源码是一个用ｃ＋＋编写的经典的遗传算法

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <math.h>
#include <string.h>
#include <assert.h>
#include <iostream>
#include <ctype.h>
#include <cstdlib>

using namespace std;

//用于改变随机二进制位种子iseed的常量
#define IB1 1
#define IB2 2
#define IB5 16
#define IB18 131072L
#define MASK (IB1+IB2+IB5)

//以下常量用于随机数发生器,IA1,IQ1,IR1有三组数可供选择，取决于FPGA易实现性和资源最小化
#define IM1 2147483647L
#define AM1 (1.0/IM1)
#define IA1 16807L
#define IQ1 127773L
#define IR1 2836L

//#define IA1 48271L
//#define IQ1 44488L
//#define IR1 3399L

//#define IA1 69621L
//#define IQ1 30845L
//#define IR1 23902L


//idum为浮点随机数种子,iseed为产生随机二进制位的种子，seed为初始种子
static long idum=0;
static unsigned long iseed;
static unsigned long seed=0;



struct SGenome
{
unsigned char **data;          //data指向一个二维数组，data[种群数][每个基因的染色体数]
double *score, *psum;          //score为每个基因的评估值，psum为每个基因的加权和
unsigned int size, width;      //size是种群数大小，width为每个基因的染色体数
double max_score, min_score, ave_score;  //分别是当前种群评估值的最大、最小和平均值
};


//产生浮点随机数的函数，该函数能保证产生的随机数是均匀的
double garan2()
{
long k;
double temp;

  k=(idum)/IQ1;
  idum=IA1*(idum-k*IQ1)-k*IR1;
  if (idum < 0) idum += IM1;

  temp = AM1*idum;
  return temp;
}

double GARandomFloat()
{
double val;	
	val = garan2();
	return val;
}


//在low和high之间产生一个整型随机数的函数
int GARandomInt(int low, int high)
{ 
double val;

  val = (double)high-(double)low+(double)1; 
  val*=garan2(); 
  return ((int)val+low);
}


//随机产生一个0或1的函数
int GARandomBit()
{
  if (iseed & IB18)
  {
    iseed=((iseed ^ MASK) << 1) | IB1;
    return 1;
  }
  else
  {
    iseed <<= 1;
    return 0;
  }
}


//设置随机数初始种子seed的值函数,这里用机器当前时间做种子，而在FPGA硬件实现时没有机器时间，
//能否用硬件设计一个时钟计数器(16位或32位)，然后随机取计算器的制作为初始种子。或者考虑用其它更好的办法。

void GARandomSeed()
{
unsigned long int tmp;

    while(seed == 0)
	{
      tmp = time(NULL);
      for(unsigned int i=0; i < 8*sizeof(unsigned int); i++)
	     seed += (tmp & (1 << i));
    }

    idum = (long)seed;
    if (idum == 0) idum=1;
    if (idum < 0) idum = -idum;

    iseed = seed;
}


//由指定概率确定0,1翻转的函数
bool GAFlipCoin(double p)
{
  return((p == 1.0) ? true : (p == 0.0) ? false : ((GARandomFloat() <= p) ? true : false));
}


//目标函数,这里假设最优的染色体已知（放在buf里）,实际应该根据滤波器的频率响应来确定评估值，目前评估方法还没有想好，你有没有好多想法?
double objFunction1(unsigned char *p, unsigned int n)
{
double score=0.0;
unsigned int buf[24] = {1, 1, 0, 0, 1, 1, 0, 0, 1, 0, 1, 1, 0, 0, 1, 1, 1, 1, 0, 0, 1, 1, 0, 0};
unsigned int i;

    for (i=0; i < n; i++)
	{
        if (p[i] == buf[i])
	      score += 1.0;
	}

	return score;
}

//计算种群每个染色体的目标函数值
void evaluate(SGenome &s_pop)
{
double score=0.0, sum=0.0;
unsigned int i;

    for(i=0; i < s_pop.size; i++)
	{
		score = objFunction1(s_pop.data[i], s_pop.width);
		s_pop.score[i] = score;
		sum += score;
		if (i == 0)
		{
			s_pop.max_score = score;
			s_pop.min_score = score;
		}
		else
		{
            if (s_pop.max_score < score)
				s_pop.max_score = score;
			if (s_pop.min_score > score)
				s_pop.min_score = score;
		}
	}

	s_pop.ave_score = sum / s_pop.size;
}

//打印结果
void printScore(SGenome &s_pop)
{
unsigned int i, j;

    for (i=0; i < s_pop.size; i++)
	{
	    for (j=0; j < s_pop.width; j++)
	        cout << (int)s_pop.data[i][j] << "/ ";
	    cout << endl;
	}

    for (i=0; i < s_pop.size; i++)
	    cout << s_pop.score[i] << " ";
	cout << endl;
}


//采用冒泡排序对染色体根据评估值进行排序,这里要改用和谐理论排序算法来实现
void sortPop(SGenome &s_pop)
{
unsigned int i,j;
unsigned char *tmpdata = new unsigned char[s_pop.width];
double tmpfit;

    for (i= 1; i < s_pop.size; i++)
		for (j=0; j < s_pop.size-i; j++)
		{
			if (s_pop.score[j] < s_pop.score[j+1])
			{
				tmpfit = s_pop.score[j];
                memcpy(tmpdata, s_pop.data[j], s_pop.width*sizeof(unsigned char));
				s_pop.score[j] = s_pop.score[j+1];
                memcpy(s_pop.data[j], s_pop.data[j+1], s_pop.width*sizeof(unsigned char));
				s_pop.score[j+1] = tmpfit;
                memcpy(s_pop.data[j+1], tmpdata, s_pop.width*sizeof(unsigned char));
			}
		}

    delete [] tmpdata;
}


//预选择模块
void preSelect(SGenome &s_pop)
{
unsigned int i;


	if (s_pop.max_score == s_pop.min_score)
    {
	    for (i=0; i < s_pop.size; i++)
	  	    s_pop.psum[i] = (i+1)/(double)s_pop.size;
	}
    else
	{
        sortPop(s_pop);
	    s_pop.psum[0] = s_pop.score[0];
	    for (i=1; i < s_pop.size; i++)
	  	    s_pop.psum[i] = s_pop.score[i]+s_pop.psum[i-1];
		for (i=0; i < s_pop.size; i++)
		    s_pop.psum[i] /= s_pop.psum[s_pop.size-1];
	}
}

//选择模块,采用轮盘赌原则,返回的是染色体的序号
int select(double *psum, int n)
{
double cutoff;
int i, lower, upper;

      cutoff = GARandomFloat();
      lower = 0;
	  upper = n-1;
      while(upper >= lower)
	  {
          i = lower + (upper-lower)/2;
          if(psum[i] > cutoff)
              upper = i-1;
          else
              lower = i+1;
	  }
      lower = n-1 < lower ? n-1 : lower;
      lower = 0 > lower ? 0 : lower;

	  return lower;
}

//交叉模块,p1,p2为待交叉的两个染色体,c1,c2为交叉后生成的染色体,n为染色体含基因的个数
void crossOver(const unsigned char *p1, const unsigned char *p2, unsigned char* c1, unsigned char* c2, int n)
{
unsigned int momsite, momlen;
unsigned int dadsite, dadlen;

      momsite = dadsite = GARandomInt(0, n);
      momlen = dadlen = n - momsite;

      memcpy(c1, p1, momsite*sizeof(unsigned char));
      memcpy(c1+momsite, p2+dadsite, dadlen*sizeof(unsigned char));

      memcpy(c2, p2, dadsite*sizeof(unsigned char));
      memcpy(c2+dadsite, p1+momsite, momlen*sizeof(unsigned char));
}

//变异模块,p指向染色体,n为染色体含基因的个数,pmut为变异概率
void mutate(unsigned char *p, int n, double pmut)
{
int i, j;
double nmut;

    nmut = pmut * n;
    if (nmut < 1.0)
	{		// we have to do a flip test on each bit
        for (i=n-1; i >= 0; i--)
		{
            if(GAFlipCoin(pmut))
			{
                if (p[i] == 0)
	               p[i] = 1;
		        else
			       p[i] = 0;
			}
		}
	}
    else
	{				// only flip the number of bits we need to flip
        for(j=0; j<nmut; j++)
		{
           i = GARandomInt(0, n-1); // the index of the bit to flip
           if (p[i] == 0)
	          p[i] = 1;
		   else
		      p[i] = 0;
		}
	}
}

//返回最佳染色体的下标值,这里可以用和谐理论排序算法,一个时钟即可
double bestGenome(SGenome &s_pop)
{
double best = s_pop.score[0];
unsigned int index = 0;

	for (unsigned int i=1; i < s_pop.size; i++)
		if (best < s_pop.score[i])
		{
			best = s_pop.score[i];
			index = i;
		}

    return index;
}

//返回最差染色体的下标值,这里可以用和谐理论排序算法,一个时钟即可
unsigned int worstGenome(SGenome &s_pop)
{
double worst = s_pop.score[0];
unsigned int index = 0;

	for (unsigned int i=1; i < s_pop.size; i++)
		if (worst > s_pop.score[i])
		{
			worst = s_pop.score[i];
			index = i;
		}

    return index;
}

//对种群初始化
void initialize(SGenome &s_pop, int psize, int gwidth)
{
	s_pop.size = psize;
	s_pop.width = gwidth;
	s_pop.score = new double [psize];
	s_pop.psum = new double[psize];

    s_pop.data = new unsigned char * [psize];
    for (int i=0; i < psize; i++)
	{
        s_pop.data[i] = new unsigned char[gwidth];
        for(int j=0; j < gwidth; j++)
	        s_pop.data[i][j] = GARandomBit();
	}
}

//释放内存
void release(SGenome &s_pop)
{
unsigned int i;

    delete [] s_pop.score;
	delete [] s_pop.psum;

	for (i=0; i < s_pop.size; i++)
		delete [] s_pop.data[i];

	delete [] s_pop.data;
}

void main()
{
SGenome  pop, new_pop;
unsigned int i;
unsigned int gen_no, sel_no, mom, dad;
unsigned int width = 24;      //基因个数
unsigned int popsize = 100;    //种群所含的染色体数,选为10,便于用和谐理论排序算法
unsigned int ngen = 100;      //最大叠代次数
double pmut = 0.001;          //变异概率
double pcross = 0.9;          //交叉概率
unsigned int old_best, new_best;

    cout << "This program tries to len a binary pattern\n\n";

    GARandomSeed();

//initialize
	initialize(pop, popsize, width);
    initialize(new_pop, popsize, width);

//evaluate
	evaluate(pop);
	printScore(pop);


//evolve
  gen_no = 0;
  while (gen_no++ < ngen)
  {
	  preSelect(pop);   //排序
	  for (sel_no = 0; sel_no < pop.size-1; sel_no += 2)
	  {
//select,每次选择两个染色体
		  mom = select(pop.psum, popsize);
		  dad = select(pop.psum, popsize);

//crossover,进行交叉
		  if(GAFlipCoin(pcross))
		  {
              crossOver(pop.data[mom], pop.data[dad], new_pop.data[sel_no], new_pop.data[sel_no+1], pop.width);
		  }
          else
		  {
              memcpy(new_pop.data[sel_no], pop.data[mom], pop.width*sizeof(unsigned char));   //直接将父辈拷给子辈
              memcpy(new_pop.data[sel_no+1], pop.data[dad], pop.width*sizeof(unsigned char));
		  }
//mutate,对子辈种群进行变异
		  mutate(new_pop.data[sel_no], new_pop.width, pmut);
		  mutate(new_pop.data[sel_no+1], new_pop.width, pmut);
	  }

	  new_best = bestGenome(new_pop);
	  old_best = bestGenome(pop);   	  
//replace new worst with old best,若子辈的最佳不如父辈的最佳,则用父辈的最佳替换子辈的最差,这一机制很重要,可以用和谐理论排序算法实现
	  if (pop.score[old_best] > new_pop.score[new_best])
          memcpy(new_pop.data[worstGenome(new_pop)], pop.data[old_best], pop.width*sizeof(unsigned char));

//current new_pop is the pop for next generation,子辈又是下一代的父辈
	  for (i=0; i < pop.size; i++)
          memcpy(pop.data[i], new_pop.data[i], pop.width*sizeof(unsigned char));
//evaluate,计算评估函数
      evaluate(pop);
  }


  printScore(pop);

  for (i=0; i < pop.width; i++)
	  cout << (int)pop.data[new_best][i];
  cout << endl;

  release(pop);
  release(new_pop);
}