rtcgps.cpp

基于GPS的多目标优化及动态多目标优化源代码 里面 包含MOP, DMOP的程序

	{  
		uniform_crossover (&m_Pop[i],&m_Pop[i + 1]);
		sbx (&m_Pop[i],&m_Pop[i + 1]);
	}
	
	/* do mutation */
	for(i = 0; i < PopSize; i++)
	{
		mutation(&m_Pop[i]);
	}
	
	/* do evaluation */
	for(i = 0; i < PopSize; i++)
	{
		SetFitness(m_Pop[i]);
	}
	
	return (0);
}


double drand(double range)
{
	double j;
	j=(range * (double) rand() / (RAND_MAX + 1.0));
	return (j);
}

int RTCGPS::mutation(CSolution *ind)
{
	int i;
	for (i = 0; i < CSolution::GetGeneLength(); i++)
	{
		double eta = 20;
		double u = drand(1);//CRandom::Random(0.0,0.99999999);
		double delta = 0;
		double x = ind->Gene[i];
		double lb = 0;    /* lower bound of variable i */
		double ub = 1;    /* upper bound of variable i */
		double diff = ub - lb;  /* range of variable i */
		double maxmut0 = x - lb;
		double maxmut = ub - x;
		double delta_max = maxmut0 / diff;
		if (maxmut0 > maxmut)
		{
			delta_max = maxmut / diff;
		}
		
		if (u < 0.5)
		{
			double b =  2*u + (1-2*u)*(pow(1-delta_max,(eta+1)));
			delta = pow(b,(1.0/(eta+1))) - 1.0;
		}
		else
		{
			double b = 2*(1-u) + 2*(u-0.5)*(pow(1-delta_max,(eta+1)));
			delta = 1.0 - pow(b,(1.0/(eta+1)));
		}
		if (delta > delta_max)  /* machine accuracy problem */
			delta = delta_max;
		else if (delta < -delta_max)
			delta = -delta_max;
		
		ind->Gene[i] = x + delta * diff;
	}
	
	return (0);
}

int RTCGPS::mutation(CSolution *ind, int GeneID)
{
	int i = GeneID;
	double eta = 20;
	double u = drand(1);//CRandom::Random(0.0,0.99999999);
	double delta = 0;
	double x = ind->Gene[i];
	double lb = 0;    /* lower bound of variable i */
	double ub = 1;    /* upper bound of variable i */
	double diff = ub - lb;  /* range of variable i */
	double maxmut0 = x - lb;
	double maxmut = ub - x;
	double delta_max = maxmut0 / diff;
	if (maxmut0 > maxmut)
	{
		delta_max = maxmut / diff;
	}
	
	if (u < 0.5)
	{
		double b =  2*u + (1-2*u)*(pow(1-delta_max,(eta+1)));
		delta = pow(b,(1.0/(eta+1))) - 1.0;
	}
	else
	{
		double b = 2*(1-u) + 2*(u-0.5)*(pow(1-delta_max,(eta+1)));
		delta = 1.0 - pow(b,(1.0/(eta+1)));
	}
	if (delta > delta_max)  /* machine accuracy problem */
		delta = delta_max;
	else if (delta < -delta_max)
		delta = -delta_max;
	
	ind->Gene[i] = x + delta * diff;
	
	return (0);
}

int RTCGPS::uniform_crossover(CSolution *ind1, CSolution *ind2)
{
	int i;
	
	for (i = 0; i <CSolution::GetGeneLength(); i++)
	{
		if (CRandom::Random(0.0,1.0) <= 0.5) /* switch variable */
		{
			double x = ind2->Gene[i];
			ind2->Gene[i] = ind1->Gene[i];
			ind1->Gene[i] = x;
		} 
	}  
	
	return (0);
}



int RTCGPS::sbx(CSolution *ind1, CSolution *ind2)
{
	int i;
	
	for (i = 0; i <  CSolution::GetGeneLength(); i++)
	{
		if ((ind1->Gene[i]>1) ||(ind1->Gene[i]<0)) ind1->Gene[i] = CRandom::Random(0.0,1.0);
		if ((ind2->Gene[i]>1) ||(ind2->Gene[i]<0)) ind2->Gene[i] = CRandom::Random(0.0,1.0);
		
		if (drand(1) <= 1.0)  
		{
			double di = 15; //eta_recombination; /* distribution index */
			int bounded = 1;
			double lb = 0;    /* lower bound of variable i */
			double ub = 1;    /* upper bound of variable i */	     
			double u = drand(1);
			double b0=0, b1=0;   /* spread factors */
			double x0 = ind1->Gene[i];
			double x1 = ind2->Gene[i];
			double bl=0, bu=0, p_bl=0, p_bu=0, bll=0, buu=0, blll=0, buuu=0;
			double dx = 0;
			double u0=0, u1=0;
			
            /* calculate spread factor(s) b0, b1 */ 
            if (bounded == 1)
            {
                dx = fabs(x1-x0);   /* difference of x values */
                if (dx > 0)
                {
                    bl = (x0 + x1 - 2*lb) / dx;
                    bu = (2*ub - x0 - x1) / dx;
                    bll = (x0 + x1 - 2*(x0-lb)) / dx;
                    buu = (2*(ub-x1)-x0-x1) / dx;
                    if (x0 < x1)
                    {
                        blll = 1 + 2 * (x0 - lb) / dx;
                        buuu = 1 + 2 * (ub - x1) / dx;
                    }
                    else
                    {
                        blll = 1 + 2 * (x1 - lb) / dx;
                        buuu = 1 + 2 * (ub-x0) / dx;
					}
					
					bl = blll; /* take Deb's version (numerically stable) */
                    bu = buuu;
					
                    /* switch off symmetric recombination to avoid
					* getting stuck on a line where one parameter
					* equals an extreme value. */
                    if (1)//(use_symmetric_recombination)
                    {
						if (bl < bu)  /* symmetric distribution, like Deb */
							bu = bl;
						else
							bl = bu;
						assert (b0 == b1);
                    }
                    assert(bl > 0 && bu > 0);
                    p_bl = 1 - 1/(2*pow(bl,di+1));
                    p_bu = 1 - 1/(2*pow(bu,di+1));
                }
                else
                {
                    p_bl = 1;
                    p_bu = 1;
                }
                u0 = u*p_bl;
                u1 = u*p_bu;
                if (u0<=0.5)
                    b0 = pow(2*u0,1/(di+1));
                else
                    b0 = pow(0.5/(1-u0),1/(di+1));
                if (u1<=0.5)
                    b1 = pow(2*u1,1/(di+1));
                else
                    b1 = pow(0.5/(1-u1),1/(di+1));
                assert(dx==0 || (b0<=bl && b1<=bu)); /* machine accuracy */
            }
            else
            {
                if (u<=0.5)
                    b0 = pow(2*u,1/(di+1));
                else
                    b0 = pow(0.5/(1-u),1/(di+1));
                b1 = b0;
            }
			
            if (x0 < x1)
            {
                ind1->Gene[i] = 0.5*(x0+x1 + b0*(x0-x1));
                ind2->Gene[i] = 0.5*(x0+x1 + b1*(x1-x0));
            }
            else
            {
                ind1->Gene[i] = 0.5*(x0+x1 + b1*(x0-x1));
                ind2->Gene[i] = 0.5*(x0+x1 + b0*(x1-x0));
            }
		}
	}  
	
	return (0);
	
}


int RTCGPS::sbx(CSolution *ind1, CSolution *ind2, int GeneID)
{
	int i = GeneID;
	
	if (drand(1) <= 1.0)  
	{
		double di = 15; //eta_recombination; /* distribution index */
		int bounded = 1;
		double lb = 0;    /* lower bound of variable i */
		double ub = 1;    /* upper bound of variable i */	     
		double u = drand(1);
		double b0=0, b1=0;   /* spread factors */
		double x0 = ind1->Gene[i];
		double x1 = ind2->Gene[i];
		double bl=0, bu=0, p_bl=0, p_bu=0, bll=0, buu=0, blll=0, buuu=0;
		double dx = 0;
		double u0=0, u1=0;
		
		/* calculate spread factor(s) b0, b1 */ 
		if (bounded == 1)
		{
			dx = fabs(x1-x0);   /* difference of x values */
			if (dx > 0)
			{
				bl = (x0 + x1 - 2*lb) / dx;
				bu = (2*ub - x0 - x1) / dx;
				bll = (x0 + x1 - 2*(x0-lb)) / dx;
				buu = (2*(ub-x1)-x0-x1) / dx;
				if (x0 < x1)
				{
					blll = 1 + 2 * (x0 - lb) / dx;
					buuu = 1 + 2 * (ub - x1) / dx;
				}
				else
				{
					blll = 1 + 2 * (x1 - lb) / dx;
					buuu = 1 + 2 * (ub-x0) / dx;
				}
				
				bl = blll; /* take Deb's version (numerically stable) */
				bu = buuu;
				
				/* switch off symmetric recombination to avoid
				* getting stuck on a line where one parameter
				* equals an extreme value. */
				if (1)//(use_symmetric_recombination)
				{
					if (bl < bu)  /* symmetric distribution, like Deb */
						bu = bl;
					else
						bl = bu;
					assert (b0 == b1);
				}
				assert(bl > 0 && bu > 0);
				p_bl = 1 - 1/(2*pow(bl,di+1));
				p_bu = 1 - 1/(2*pow(bu,di+1));
			}
			else
			{
				p_bl = 1;
				p_bu = 1;
			}
			u0 = u*p_bl;
			u1 = u*p_bu;
			if (u0<=0.5)
				b0 = pow(2*u0,1/(di+1));
			else
				b0 = pow(0.5/(1-u0),1/(di+1));
			if (u1<=0.5)
				b1 = pow(2*u1,1/(di+1));
			else
				b1 = pow(0.5/(1-u1),1/(di+1));
			assert(dx==0 || (b0<=bl && b1<=bu)); /* machine accuracy */
		}
		else
		{
			if (u<=0.5)
				b0 = pow(2*u,1/(di+1));
			else
				b0 = pow(0.5/(1-u),1/(di+1));
			b1 = b0;
		}
		
		if (x0 < x1)
		{
			ind1->Gene[i] = 0.5*(x0+x1 + b0*(x0-x1));
			ind2->Gene[i] = 0.5*(x0+x1 + b1*(x1-x0));
		}
		else
		{
			ind1->Gene[i] = 0.5*(x0+x1 + b1*(x0-x1));
			ind2->Gene[i] = 0.5*(x0+x1 + b0*(x1-x0));
		}
	}
	
	return (0);
	
}


void RTCGPS::RunStep6()
{
	   CSolution m1(TRUE);
	   CSolution m2(TRUE);
	   BOOL IsParetoPoint;
	   BOOL isBetter ;
	   int iCount;
	   
	   for (iCount=0;iCount< RTCGPS::PopSize-1 ;iCount+=2)
	   {
		   m1.Copy(&m_Pop[iCount]);
		   m2.Copy(&m_Pop[iCount+1]);
		   uniform_crossover(&m1,&m2); 
		   sbx(&m1,&m2);
		   mutation(&m1);
		   mutation(&m2);
		   SetFitness(m1);
		   SetFitness(m2);
		   IsParetoPoint = m_Elist.AddSolution(m1);
		   isBetter =Better(m1,m_Pop[iCount],m_Pop[iCount].type );
           if (IsParetoPoint) 
		   {
			   double f1 =(m1.Objective[CSolution::OBJectiveNum]);
			   f1+=  m1.Objective[m_Pop[iCount].type];
               double f2 = (m_Pop[iCount].Objective[CSolution::OBJectiveNum]);
			   f2+= m_Pop[iCount].Objective[m_Pop[iCount].type];
			   if ( f1 < f2 ) 
			   {
				   ReplaceSolution(iCount,m1);
			   }
		   }
           if (isBetter)  
		   {
			   ReplaceSolution(iCount,m1);
		   }
		   IsParetoPoint = m_Elist.AddSolution(m2);
		   isBetter =Better(m2,m_Pop[iCount+1],m_Pop[iCount+1].type );
           if (IsParetoPoint) 
		   {
			   double f1 =(m2.Objective[CSolution::OBJectiveNum]);
			   f1+=  m1.Objective[m_Pop[iCount+1].type];
               double f2 = (m_Pop[iCount+1].Objective[CSolution::OBJectiveNum]);
			   f2+= m_Pop[iCount+1].Objective[m_Pop[iCount+1].type];
			   if ( f1 < f2 ) 
			   {
				   ReplaceSolution(iCount+1,m2);
			   }
		   }
           if (isBetter)  
		   {
			   ReplaceSolution(iCount+1,m1);
		   }
	   }
}

void RTCGPS::SetFitness(CSolution &Ind)
{
	ComputeFitness(Ind);
	double f1=0;
	//	 double weigh;
	double sum=0;
	for (int i=0; i< CSolution::OBJectiveNum;i++)
	{
		//        if (newSolution.Objective[i]==0) AssignFitness[i]=0;
		//	    else
		//		 if (i != CSolution::OBJectiveNum){
		//		 weigh = CRandom::Random(0.0, (double) 1/(CSolution::OBJectiveNum-1));
		//		 sum +=weigh;
		//		 }
		//		else
		//		{
		//			weigh = 1- sum;
		//		}
		//			f1+= weigh *atan(PI * (2*Ind.Objective[i] - 2*m_Elist.Baseline[i] -m_Elist.Scale[i])   / (2*m_Elist.Scale[i])) ;
		//			f1+= atan(PI * (2*Ind.Objective[i] - 2*m_Elist.Baseline[i] -m_Elist.Scale[i])   / (2*m_Elist.Scale[i])) ;
		//        if (_isnan(f1)) 
		//			f1 =0;
		f1+=Ind.Objective[i]/ m_Elist.Scale[i];
		
	}
	Ind.Objective[CSolution::OBJectiveNum] = f1;
}

void RTCGPS::SaveBest(CSolution *pInd)
{
	   for (int icnt2=0;icnt2 < CSolution::OBJectiveNum+1;icnt2++)
	   {
		   if (Better(*pInd,BestSolutions[icnt2], icnt2) )
		   {
			   BestSolutions[icnt2].Copy(pInd);
		   }
	   }
}

BOOL RTCGPS::WSMBetter(CSolution *Ind1, CSolution *Ind2)
{
	double f1=0, f2 =0;
	double sum=0;
	for (int iobj =0; iobj <CSolution::OBJectiveNum-1;iobj++  )
	{
		double weight = CRandom::Random(0.0, 1.0/CSolution::OBJectiveNum);
		sum+=weight;
		f1 +=weight * (Ind1->Objective[iobj]);
		f2 += weight * (Ind2->Objective[iobj]);
	}
	f1+= (1-sum) *Ind1->Objective[iobj];
	f2 += (1-sum) * Ind2->Objective[iobj];
	if (f1 < f2) return  TRUE;
	else return FALSE;
}