ga.cpp

遗传算法类的声明及定义

//遗传算法类cpp文件
//Edit by Wang Shimin
//2004.7.5～2004.7.30
//本文件为最新文件，修改于2006.3.12

// GA.cpp: implementation of the GA class.
//
//////////////////////////////////////////////////////////////////////

#include "stdafx.h"
#include "GA.h"
#include "math.h"
#include "fstream"

#ifdef _DEBUG
#undef THIS_FILE
static char THIS_FILE[]=__FILE__;
#define new DEBUG_NEW
#endif

//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////

GA::GA()
{
   /* 初始化计数变量和一些数值*/
	m_convergencenumber = 0;
	m_precision = 0.000001;
	m_errorflag = TRUE;
	m_oldbestfit.fitness = 0.0;
    m_bestfit.fitness = 0.0;
    m_bestfit.generation = 0;
	m_pc2 = 0.3;
	m_pc1 = 0.9;
	m_pm2 = 0.01;
	m_pm1 = 0.1;
	m_perlchrom = 16; //染色体中单个参数的长度默认为16位
    m_pcross = 0.0;     //交叉率
	m_pmutation = 0.0;  //变异率
	m_numbndload = 0;
	m_popsize = 10;
	m_maxgen = 20;
	m_currentgen = 0;
	m_pernumbytes = 1;

	/* 初始化随机种子 */
	srand((unsigned)time(NULL));
	m_randseed = double((rand()%10000)/10000.0);

	/*将指针置零*/
	p_Bndload = NULL;
	p_BAKindflag = NULL;
	p_iparameter = NULL;
	p_jparameter = NULL;
	p_calculation = NULL;
	p_observation = NULL;
	p_highboundary = NULL;             
	p_lowboundary = NULL;
	p_deltavariable = NULL;
	p_oldpop = NULL;
	p_newpop = NULL;
	m_bestfit.chrom = NULL;
	m_bestfit.variable = NULL;

	m_blErrorExit = false;
}

GA::~GA()
{
	/*释放内存空间*/
	for(int i = 0; i < m_popsize; i++)
    {
		delete[]p_oldpop[i].chrom;
		delete[]p_oldpop[i].variable;
		delete[]p_newpop[i].chrom;
		delete[]p_newpop[i].variable;
    }
	delete[]p_oldpop;
	delete[]p_newpop;
	delete[]m_bestfit.chrom;
	delete[]m_bestfit.variable;
	
	delete[]p_BAKindflag;
	delete[]p_iparameter;
	delete[]p_jparameter;
	delete[]p_highboundary;
	delete[]p_lowboundary;
	delete[]p_deltavariable;
	delete[]p_calculation;
	delete[]p_observation;
	delete[]p_Bndload;
}

/* 遗传算法初始化 */
void GA::initialize()  
{
    /* 输入遗传算法参数 */
	ReadINP();
	if(m_blErrorExit)return;

	CString csMEAFileName = mGcsFileName+_T(".MEA");
	ifstream cMEAFile;
	cMEAFile.open((LPCSTR)csMEAFileName); 
	if (cMEAFile.fail())
	{
		AfxMessageBox("Cannot open input MEA file!!!",MB_OK);
		m_blErrorExit = true;
		return;
	}
	
	cMEAFile >> m_numpoint;  /* 从MEA文件中读取观测点的个数 */
	cMEAFile.close();

    /*分配给全局数据结构空间 */
    initmalloc();

	/* 读入待反演参数 */
	ReadBAA();
	if(m_blErrorExit)return;

	/* 读入测点的观测值 */
	ReadMEA();
	if(m_blErrorExit)return;
	
 	/* 初始化种群，并统计计算结果 */
    initpop();
}

/* 读取遗传算法参数 */
void GA::ReadINP()
{
	CString csINPFileName = mGcsFileName+_T(".INP");
	ifstream cINPFile;
	cINPFile.open((LPCSTR)csINPFileName); 
	if (cINPFile.fail())
	{
		AfxMessageBox("Cannot open input INP file!!!",MB_OK);
		m_blErrorExit = true;
		return;
	}

	char tempch[256];

	cINPFile.getline(tempch,256);
	cINPFile >> m_popsize >> m_numvar >> m_pm2 >> m_maxgen
		>> m_pc2 >> m_perlchrom >> m_precisionnum;

	m_pernumbytes = m_perlchrom/8;

	cINPFile.close();

	if((m_popsize%2)==1)
	{
		AfxMessageBox(GIDS_PROMPT_NOT_SINGLE,MB_OK);
		m_blErrorExit = true;
		return;
	}
	/* 确定染色体的总长度 */
	m_sumlchrom = m_perlchrom*m_numvar ;
}

/*为数组数据变量分配空间 */
void GA::initmalloc()
{
	//当重新输入参数的时候调用此函数，释放原来的动态内存
	if(p_oldpop)
	{
		delete[]p_oldpop ;
		p_oldpop = NULL ;
	}
	
	if(p_newpop)
	{
		delete[]p_newpop ;
		p_newpop = NULL;
	}

	/* 分配给当前代和新一代种群内存空间 */
	p_oldpop=new individual[m_popsize];
	p_newpop=new individual[m_popsize];

	/* 分配给观测数据和计算数据内存空间 */
	if(p_calculation)
		delete[]p_calculation;
	
	if(p_observation)
		delete[]p_observation;

	p_calculation = new double[m_numpoint];
	p_observation = new double[m_numpoint];

	/* 为存储识别待反演参数性质的数组分配内存空间 */

	if(p_BAKindflag)
		delete[]p_BAKindflag;

	if(p_iparameter)
		delete[]p_iparameter;

	if(p_jparameter)
		delete[]p_jparameter;

	p_BAKindflag = new int[m_numvar];
	p_iparameter = new int[m_numvar];
	p_jparameter = new int[m_numvar];

	/* 为存储各个参数上下限和精度的数组分配内存空间 */
	if(p_highboundary)
		delete[]p_highboundary;
	
	if(p_lowboundary)
		delete[]p_lowboundary;
	
	if(p_deltavariable)
		delete[]p_deltavariable;
	
	p_highboundary = new double[m_numvar];
	p_lowboundary = new double[m_numvar];
	p_deltavariable = new double[m_numvar];
	
	/* 分配给染色体内存空间 */
	m_sumnumbytes = m_sumlchrom/8;
	
	for(int i = 0; i < m_popsize; i++)
	{
		p_oldpop[i].chrom = new BYTE[m_sumnumbytes];
		p_newpop[i].chrom = new BYTE[m_sumnumbytes];
		m_bestfit.chrom = new BYTE[m_sumnumbytes];

		p_oldpop[i].variable = new double[m_numvar];
		p_newpop[i].variable = new double[m_numvar];
		p_oldpop[i].objvalue = new double[m_numpoint];
		p_newpop[i].objvalue = new double[m_numpoint];
		m_bestfit.variable = new double[m_numvar];
		m_bestfit.objvalue = new double[m_numpoint];
	}
}

/* 从BAA文件中读入待反演参数的性质标志以及上下限 */
void GA::ReadBAA()
{
	char tempchar[256];
	int  index;

	CString csBaaFileName = mGcsFileName+_T(".BAA");
	ifstream cBaaFile;
	cBaaFile.open((LPCSTR)csBaaFileName); 
	if (cBaaFile.fail())
	{
		AfxMessageBox("Cannot open input BAA file!!!",MB_OK);
		m_blErrorExit = true;
		return;
	}
	

	cBaaFile.getline(tempchar,256);
	
	//读入各个参数的上下限，并计算出各个参数的精度
	for(int i=0; i<m_numvar; i++)
	{
		cBaaFile >> index;
		cBaaFile >> p_BAKindflag[i];
		cBaaFile >> p_iparameter[i];
		cBaaFile >> p_jparameter[i];
		cBaaFile >> p_lowboundary[i];
		cBaaFile >> p_highboundary[i];
		cBaaFile.getline(tempchar,256);
		p_deltavariable[i] = (p_highboundary[i]-p_lowboundary[i])/(pow(2.0,(double)m_perlchrom)-1.0);
	}

	cBaaFile >> m_numbndload;

	if(m_numbndload>0)
	{
		if(p_Bndload)
		{
			delete[]p_Bndload ;
			p_Bndload = NULL ;
		}

		p_Bndload = new bndload[m_numbndload];

		for(int j=0; j<m_numbndload; j++)
		{
			cBaaFile >> p_Bndload[j].index >> p_Bndload[j].wAddStage
				>> p_Bndload[j].wAddStep >> p_Bndload[j].iNode1
				>> p_Bndload[j].iNode2 >> p_Bndload[j].iCoordFlag
				>> p_Bndload[j].BeginLoad1 >> p_Bndload[j].BeginLoad2
				>> p_Bndload[j].BeginLoad3 >> p_Bndload[j].EndLoad1
				>> p_Bndload[j].EndLoad2 >> p_Bndload[j].EndLoad3;

			cBaaFile.getline(tempchar,256);
		}
	}

	cBaaFile.close();
}

/*产生0到1 之间均匀分布的随机数*/
double GA::randomperc(double *seed) 
{
  
 //r双精度实型变量指针。该指针指向的单元存放随机数的种子；返回一个新值。但是本函数返回的不是r
	int m;
    double s,u,v;
	double result;

    s = 65536.0;
	u = 2053.0;
	v = 13849.0;
    m = (int)(*seed/s);
	*seed = *seed-m*s;
    *seed = u*(*seed)+v;
	m = (int)(*seed/s);
    *seed = *seed-m*s;
	result = *seed/s;
    return(result);
}

/* 随机初始化种群 */
void GA::initpop()  
{

    int i, j, k,l;
    BYTE mask = 1;

    for(i = 0; i < m_popsize; i++)
    {	
		p_oldpop[i].fitness = 0.0;

        for(k = 0; k < m_sumnumbytes; k++)
        {
            p_oldpop[i].chrom[k] = 0;
            
			for(j = 1; j <= 8; j++)
            {
				p_oldpop[i].chrom[k] = p_oldpop[i].chrom[k]<<1;
				if(flip(0.5))
					p_oldpop[i].chrom[k] = p_oldpop[i].chrom[k]|mask;
            }
        }
		
    	/* 解码得到编码所对应的实数值 */
        decode(&(p_oldpop[i]));
		/* 计算适应度函数值 */
		
		WriteToFem(&(p_oldpop[i]));
	    if(m_blErrorExit)return;
		
		objectfunc();
	    if(m_blErrorExit)return;
		
		ReadFromFem(p_calculation);
	    if(m_blErrorExit)return;

		for(l=0; l<m_numpoint; l++)
		{
			p_oldpop[i].objvalue[l] = p_observation[l]-p_calculation[l];
			p_oldpop[i].fitness += pow(p_oldpop[i].objvalue[l],2.0);
		}
    }
	
	Keepbest(p_oldpop);
}

/* 以一定概率产生0或1 */
bool GA::flip(double prob)
{
	if(randomperc(&m_randseed) <= prob)
        return(1);
    else
        return(0);
}

/* 解码 */
void GA::decode(individual *critter)  
{
    BYTE mask=1;
    int bitpos;   //染色体中的每位在参数解码时的位数
    BYTE tp;
    double  bitpow ;
    int i, j, k, l;
	for(i=0; i<m_numvar; i++)
	{
		critter->variable[i] = 0.0;
		
		for(j=0; j<m_pernumbytes; j++)
		{
			k = i*m_pernumbytes+j;
			tp = critter->chrom[k];
			for(l=0; l<8; l++)
			{
				bitpos = j*8+l;

				if((tp&mask) == 1)
				{
					bitpow = pow(2.0,(double) bitpos);
					critter->variable[i] = critter->variable[i] + bitpow;
				}
				tp = tp>>1;
			}
		}
		
		critter->variable[i] = p_lowboundary[i]+critter->variable[i]*p_deltavariable[i];
	}
}

/* 计算有限元值 */
void GA::objectfunc() 
{
	if(EXCUTE_NUMERICAL_ANALYSIS()==false)
		m_blErrorExit=true;
}

/* 遗传演化 */
void GA::generation()
{
  int i,mate1, mate2,j = 0;
  double fitness;


  Keepbest(p_oldpop);
  
  //保留上一步最佳适应值
  m_oldbestfit.fitness = m_bestfit.fitness; 

  /* 选择, 交叉, 变异 */
  do
  {
      /* 挑选交叉配对 */
      mate1 = select();
      mate2 = select();

	  /* 调整交叉概率 */
	  fitness=min(p_oldpop[mate1].fitness,p_oldpop[mate2].fitness);
	  m_pcross=adjust_cp(fitness);
	  
      /* 交叉 */
      crossover(p_oldpop[mate1].chrom, p_oldpop[mate2].chrom, p_newpop[j].chrom, p_newpop[j+1].chrom);

	  /********************** 第j个个体******************/

	  /* 解码并计算个体适应值 */
      decode(&(p_newpop[j]));
	  WriteToFem(&(p_newpop[j]));
	  if(m_blErrorExit)return;

	  objectfunc();
	  if(m_blErrorExit)return;
	  
	  ReadFromFem(p_calculation);
	  if(m_blErrorExit)return;