hga.cpp

hga是一种混合遗传算法c程序源代码

//计算函数的最大值

//库文件
#include<iostream.h>
#include"stdio.h"
#include"stdlib.h"
#include"time.h"
#include"math.h"

//随机数定义
#define rdint(i) (rand()%(int)(i))
#define rdft() (float)((double)rdint(16384)/(16383.0))
#define rnd(a,b) (float)(rdint((int)(b)-(int)(a)+1)+(int)(a))

//宏定义
#define POPSIZE 40
#define NVARS 4
#define Error 0.9
#define LargestNumber 50
#define MaxNumber 10000
#define Max_Largest_Number 10000
#define Simplex_Number 10
#define alfa 2


//结构体定义
struct indi
{
	double variable[NVARS];
	double fitness;
}individual[POPSIZE],parent[2*POPSIZE],bestindividual,temp[NVARS+1],midpoint,addpoint1,addpoint2,addpoint3,temp_high;

//全局变量定义
double PC=0.85;
double PM=0.1;
double Minpoint[NVARS]={1.2,0.02,3,1.05};
double Maxpoint[NVARS]={3,0.5,12,3};

//函数定义
void initiate(void);
void calculate(int flag,int number);
void crossover(void);
void selection(void);
void simplex_method(void);
void order(void);
void area_range(void);   
void determine(void);
double cal(double x1,double x2,double x3,double x4);

//区域自适应收缩
void area_range()
{
	int i,j;
	double s,h,minpoint_error[NVARS],maxpoint_error[NVARS];
	//计算分量平均值
	for(j=0;j<NVARS;j++)
	{
		s=0.0;
		for(i=0;i<2*POPSIZE;i++)
			s=s+parent[i].variable[j];
		midpoint.variable[j]=s/POPSIZE;
	}
	//计算加权平均
	s=0.0;
	for(i=0;i<2*POPSIZE;i++)
		s=s+parent[i].fitness;
	for(j=0;j<NVARS;j++)
	{
		h=0.0;
		for(i=0;i<2*POPSIZE;i++)
			h=h+(parent[i].variable[j]-midpoint.variable[j])*parent[i].fitness/s;
		midpoint.variable[j]=midpoint.variable[j]+h;
	}
	//估算方差
	for(j=0;j<NVARS;j++)
	{
		h=0.0;
		for(i=0;i<2*POPSIZE;i++)
			h=h+pow(parent[i].variable[j]-midpoint.variable[j],2.0);
		h=sqrt(h/(POPSIZE-1));
		minpoint_error[j]=midpoint.variable[j]-h;
		maxpoint_error[j]=midpoint.variable[j]+h;
	}
	//可行域的调整
	for(j=0;j<NVARS;j++)
	{
		if(minpoint_error[j]<Minpoint[j])
			minpoint_error[j]=Minpoint[j];
		if(maxpoint_error[j]>Maxpoint[j])
			maxpoint_error[j]=Maxpoint[j];
	}
	//包括最优个体
	for(j=0;j<NVARS;j++)
	{
		if(bestindividual.variable[j]<minpoint_error[j])
			minpoint_error[j]=bestindividual.variable[j];
		if(bestindividual.variable[j]>maxpoint_error[j])
			maxpoint_error[j]=bestindividual.variable[j];
	}
	//新种群的生成
	for(i=0;i<2*POPSIZE;i++)
	if(rdft()<PM)
	{
		for(j=0;j<NVARS;j++)
		{
			if(rdft()<0.9)
				parent[i].variable[j]=minpoint_error[j]+rdft()*(maxpoint_error[j]-minpoint_error[j]);
			else
				parent[i].variable[j]=Minpoint[j]+(Maxpoint[j]-Minpoint[j])*rdft();
		}
		calculate(1,i);
	}
}

//复制算子
void selection()
{
	int i,j,k;
	double s;
	//排序
	for(i=0;i<(2*POPSIZE-1);i++)
		for(j=(i+1);j<2*POPSIZE;j++)
			if(parent[i].fitness>parent[j].fitness)
			{
				s=parent[i].fitness;
				parent[i].fitness=parent[j].fitness;
				parent[j].fitness=s;
				for(k=0;k<NVARS;k++)
				{
					s=parent[i].variable[k];
					parent[i].variable[k]=parent[j].variable[k];
					parent[j].variable[k]=s;
				}
			}
	//复制
	for(i=0;i<POPSIZE;i++)
	{
		individual[i].fitness=parent[i].fitness;
		for(j=0;j<NVARS;j++)
			individual[i].variable[j]=parent[i].variable[j];
	}
	//最优选择
	if(bestindividual.fitness>individual[0].fitness)
	{
		bestindividual.fitness=individual[0].fitness;
		for(j=0;j<NVARS;j++)
			bestindividual.variable[j]=individual[0].variable[j];
	}
}

//决定当前最优解
void determine()
{
	int i,j;
	//复制单纯型法的解，并转化为求解最大值
	for(i=0;i<(NVARS+1);i++)
	{
		individual[i].fitness=Max_Largest_Number-temp[i].fitness;
		for(j=0;j<NVARS;j++)
			individual[i].variable[j]=temp[i].variable[j];
	}
	//统计当前最优解
	for(i=0;i<(NVARS+1);i++)
		if(bestindividual.fitness<individual[i].fitness)
		{
			bestindividual.fitness=individual[i].fitness;
			for(j=0;j<NVARS;j++)
				bestindividual.variable[j]=individual[i].variable[j];
		}
}

//杂交算子
void crossover()
{
	int i,j,k;
	double s,radio[POPSIZE],temp1,temp2;
	//排序
	order();
	//计算概率分布
	s=0.0;
	for(i=0;i<POPSIZE;i++)
		s=s+individual[i].fitness;
	for(i=0;i<POPSIZE;i++)
		radio[i]=individual[i].fitness/s;
	for(i=1;i<POPSIZE;i++)
		radio[i]=radio[i-1]+radio[i];
	//杂交母体的选择
	for(i=0;i<2*POPSIZE;i++)
	{
		s=rdft();
		for(j=0;j<POPSIZE;j++)
			if(radio[j]>s)
			{
				parent[i].fitness=individual[j].fitness;
				for(k=0;k<NVARS;k++)
					parent[i].variable[k]=individual[j].variable[k];
				break;
			}
	}
	//单点杂交
	for(i=0;i<POPSIZE;i++)
	  if(rdft()<PC)
	  {
		for(j=0;j<NVARS;j++)
		{
			s=rdft();
			temp1=parent[2*i].variable[j]*s+(1-s)*parent[2*i+1].variable[j];
			temp2=parent[2*i].variable[j]*(1-s)+s*parent[2*i+1].variable[j];
			parent[2*i].variable[j]=temp1;
			parent[2*i+1].variable[j]=temp2;
		}
		calculate(2,2*i);
		calculate(2,2*i+1);
	  }
}

//单纯型法
void simplex_method()
{
	int i,j,k,flag,high_flag,second_high_flag,low_flag;
	double s,reflect=1.0,extend=2.0,reduce=0.5,contract=0.5;
	//得到单纯型法的初始种群，变换函数为求最小值
	for(i=0;i<(NVARS+1);i++)
	{
		temp[i].fitness=Max_Largest_Number-individual[i].fitness;
		for(j=0;j<NVARS;j++)
			temp[i].variable[j]=individual[i].variable[j];
	}
	//在迭代次数内进行单纯型法寻优
	for(i=0;i<Simplex_Number;i++)
	{
		//寻找最高点
		high_flag=0;
		for(k=1;k<(NVARS+1);k++)
			if(temp[high_flag].fitness>temp[k].fitness)
				high_flag=k;
		//寻找次高点
		if(high_flag=0)
			second_high_flag=1;
		else
			second_high_flag=0;
		for(k=0;k<(NVARS+1);k++)
			if(k!=high_flag)
				if(temp[second_high_flag].fitness>temp[k].fitness)
					second_high_flag=k;
		//寻找最差点
		low_flag=0;
		for(k=1;k<(NVARS+1);k++)
			if(temp[low_flag].fitness<temp[k].fitness)
				low_flag=k;
		//确定中心位置
		for(j=0;j<NVARS;j++)
		{
			s=0.0;
		    for(k=0;k<(NVARS+1);k++)
				if(k!=high_flag)
					s=s+temp[k].variable[j];
			s=s-temp[high_flag].variable[j];
			midpoint.variable[j]=s/NVARS;
		}
		calculate(3,1);
		//反射操作
		for(j=0;j<NVARS;j++)
		{
			addpoint1.variable[j]=midpoint.variable[j]+reflect*(midpoint.variable[j]-temp[high_flag].variable[j]);
			if(addpoint1.variable[j]>Maxpoint[j])
				addpoint1.variable[j]=2*Maxpoint[j]-addpoint1.variable[j];
			if(addpoint1.variable[j]<Minpoint[j])
				addpoint1.variable[j]=2*Minpoint[j]-addpoint1.variable[j];
		}
		calculate(4,1);
		//扩展操作
		if(addpoint1.fitness<temp[low_flag].fitness)
		{
			for(j=0;j<NVARS;j++)
			{
				addpoint2.variable[j]=midpoint.variable[j]+extend*(addpoint1.variable[j]-midpoint.variable[j]);
				if(addpoint2.variable[j]>Maxpoint[j])
				    addpoint2.variable[j]=2*Maxpoint[j]-addpoint2.variable[j];
			    if(addpoint2.variable[j]<Minpoint[j])
				    addpoint2.variable[j]=2*Minpoint[j]-addpoint2.variable[j];
			}
			calculate(5,1);
			flag=1;
		}
		else if(addpoint1.fitness<=temp[second_high_flag].fitness)
		{
			for(j=0;j<NVARS;j++)
				temp[high_flag].variable[j]=addpoint1.variable[j];
			temp[high_flag].fitness=addpoint1.fitness;
			flag=3;
		}
		else
		{
			//压缩操作
			if(temp[high_flag].fitness>addpoint1.fitness)
			{
				temp_high.fitness=addpoint1.fitness;
				for(j=0;j<NVARS;j++)
					temp_high.variable[j]=addpoint1.variable[j];
			}
			else
			{
				temp_high.fitness=temp[high_flag].fitness;
				for(j=0;j<NVARS;j++)
					temp_high.variable[j]=temp[high_flag].variable[j];
			}
			for(j=0;j<NVARS;j++)
			{
				addpoint3.variable[j]=midpoint.variable[j]+reduce*(temp_high.variable[j]-midpoint.variable[j]);
				if(addpoint3.variable[j]>Maxpoint[j])
				   addpoint3.variable[j]=2*Maxpoint[j]-addpoint3.variable[j];
			    if(addpoint3.variable[j]<Minpoint[j])
				   addpoint3.variable[j]=2*Minpoint[j]-addpoint3.variable[j];
			}
			calculate(6,1);
			flag=2;
		}
		if(flag==1)
		{
			if(addpoint2.fitness<addpoint1.fitness)
			{
				temp[high_flag].fitness=addpoint2.fitness;
				for(j=0;j<NVARS;j++)
					temp[high_flag].variable[j]=addpoint2.variable[j];
			}
			else
			{
				temp[high_flag].fitness=addpoint1.fitness;
				for(j=0;j<NVARS;j++)
					temp[high_flag].variable[j]=addpoint1.variable[j];
			}
		}
		if(flag==2)
		{
			if(addpoint3.fitness<=temp_high.fitness)
			{
				temp[high_flag].fitness=addpoint3.fitness;
				for(j=0;j<NVARS;j++)
					temp[high_flag].variable[j]=addpoint3.variable[j];
			}
			else
			{
				for(k=0;k<(NVARS+1);k++)
				{
					for(j=0;j<NVARS;j++)
					{
						temp[k].variable[j]=temp[k].variable[j]+contract*(temp[low_flag].variable[j]-temp[k].variable[j]);
						if(temp[k].variable[j]>Maxpoint[j])
							temp[k].variable[j]=2*Maxpoint[j]-temp[k].variable[j];
						if(temp[k].variable[j]<Minpoint[j])
				            temp[k].variable[j]=2*Minpoint[j]-temp[k].variable[j];
					}
					calculate(7,k);
				}
			}
		}
	}
}

//个体排序
void order()
{
	int i,j,k;
	double s;
	for(i=0;i<(POPSIZE-1);i++)
		for(j=(i+1);j<POPSIZE;j++)
			if(individual[i].fitness>individual[j].fitness)
			{
				s=individual[i].fitness;
				individual[i].fitness=individual[j].fitness;
				for(k=0;k<NVARS;k++)
				{
					s=individual[i].variable[k];
					individual[i].variable[k]=individual[j].variable[k];
					individual[j].variable[k]=s;
				}
			}
}

//适应值函数计算
void calculate(int flag,int number)
{
	double s=0;
	double x1=0,x2=0,x3=0,x4=0;
	int i=0;	
	if(flag==1)
	{
		x1=individual[number].variable[0];
		x2=individual[number].variable[1];
		x3=individual[number].variable[2];
		x4=individual[number].variable[3];
		s=cal(x1,x2,x3,x4);
		individual[number].fitness=Max_Largest_Number-(Max_Largest_Number-s);
	}
	if(flag==2)
	{
		x1=parent[number].variable[0];
		x2=parent[number].variable[1];
		x3=parent[number].variable[2];
        x4=parent[number].variable[3];
		s=cal(x1,x2,x3,x4);
        parent[number].fitness=Max_Largest_Number-(Max_Largest_Number-s);
	}
	if(flag==3)
	{
		x1=midpoint.variable[0];
		x2=midpoint.variable[1];
		x3=midpoint.variable[2];
		x4=midpoint.variable[3];
		s=cal(x1,x2,x3,x4);
		midpoint.fitness=Max_Largest_Number-(Max_Largest_Number-s);
	}
	if(flag==4)
	{
		x1=addpoint1.variable[0];
		x2=addpoint1.variable[1];
		x3=addpoint1.variable[2];
		x4=addpoint1.variable[3];
		s=cal(x1,x2,x3,x4);
		addpoint1.fitness=Max_Largest_Number-(Max_Largest_Number-s);
	}
	if(flag==5)
	{
		x1=addpoint2.variable[0];
		x2=addpoint2.variable[1];
		x3=addpoint2.variable[2];
		x4=addpoint2.variable[3];
		s=cal(x1,x2,x3,x4);
		addpoint2.fitness=Max_Largest_Number-(Max_Largest_Number-s);
	}
	if(flag==6)
	{
		x1=addpoint3.variable[0];
		x2=addpoint3.variable[1];
		x3=addpoint3.variable[2];
		x4=addpoint3.variable[3];
		s=cal(x1,x2,x3,x4);
		addpoint3.fitness=Max_Largest_Number-(Max_Largest_Number-s);
	}
	if(flag==7)
	{
		x1=temp[number].variable[0];
		x2=temp[number].variable[1];
		x3=temp[number].variable[2];
        x4=temp[number].variable[3];
		s=cal(x1,x2,x3,x4);
		temp[number].fitness=Max_Largest_Number-(Max_Largest_Number-s);
    }
	
}
double cal(double x1,double x2,double x3,double x4)
{
        double s=0,g1=0,g2=0,g3=0,g4=0,g5=0,g51=0,g6=0;
			
		double alfat=1;	
  	
		s=(2*7.9-0.4*8.9)*x3+4*7.9*x1+0.5*7.9*x4*x4-7.9/(2*x4*x4)+8.9*0.4*x3*x1*x1;    
        s=s*3.1416*6.4*6.4*6.4/(x2*x2*x2);
        
        g1=x4-1;
		g1=g1*1000000;
        if(g1>0)
          g1=0;
          
        g2=x1-x4;
		g2=g2*1000000;
        if(g2>0)
          g2=0;
          
        g3=0;
        
		g4=3.1416*x4*x4/(x2*x2)*(1+x4*x4/(4*x1*x1));
        g4=g4+32*(x4-1/x4)/(1-1/(x4*x4)+4*x2)/(1-1/(x4*x4)+4*x2);
        g4=g4-0.9;
        g4=g4*1.0764*x3*x3*x1*(x1-1)/(x2*x2*x2*(x1+1)*1000000);
		g4=g4*1000000;
        if(g4>0)
           g4=0;
           
		g51=(x1-1)*x1/((x1+1)*x2);
		g51=sqrt(g51);
		g51=1.4112*x3*g51/1000;
		g5=x1*x1-1;
		g5=sqrt(g5);
		g5=x1+g5;
		g5=log(g5);
		g5=3.1416*x3/g5;
		g5=g5+8*x4*(x4*x4-1)/(x4*x4+4*x2*x4*x4-1)+(3*x4*x4*x4+x4*x4-1)/(x4*x4)+0.9*x2+3.1416/x2*(x4*x4+x4*x4*x4*x4/(4*x1*x1));
		g5=1.4-g5*g51;
		g5=g5*1000000;
		if(g5>0)
			g5=0;
			
        g6=85-2.3125*(x1*x1-1)/(x1*x2);
		g6=g6*1000000;
        if(g6>0)
          g6=0;		
        
        s=s-alfat*(g1+g2+g3+g4+g6);
		s=s/1000000;
		if (s<0)
		s=100;
		return s;
}

//算法初始化
void initiate()
{
	int i,j;
	//个体初始化
	for(i=0;i<POPSIZE;i++)
	{
		individual[i].fitness=-1.0;
		for(j=0;j<NVARS;j++)
			individual[i].variable[j]=0.0;
	}
	//杂交母体初始化
	for(i=0;i<2*POPSIZE;i++)
	{
		parent[i].fitness=-1.0;
		for(j=0;j<NVARS;j++)
			parent[i].variable[j]=0.0;
	}
	//最优个体初始化
	bestindividual.fitness=-1.0;
	for(j=0;j<NVARS;j++)
		bestindividual.variable[j]=0.0;
	//个体初始化
	for(i=0;i<POPSIZE;i++)
	{
		for(j=0;j<NVARS;j++)
			individual[i].variable[j]=Minpoint[j]+(Maxpoint[j]-Minpoint[j])*rdft();
		calculate(1,i);
	}
	//最优个体选择
	bestindividual.fitness=individual[0].fitness;
	for(j=0;j<NVARS;j++)
		bestindividual.variable[j]=individual[0].variable[j];
	for(i=1;i<POPSIZE;i++)
		if(bestindividual.fitness<individual[i].fitness)
		{
			bestindividual.fitness=individual[i].fitness;
			for(j=0;j<NVARS;j++)
				bestindividual.variable[j]=individual[i].variable[j];
		}
}

//主程序
void main()
{
	int i,j,sum=0;
	double average=0.0,total=0.0;
	srand(time(NULL));
	//程序运行次数
	for(j=0;j<LargestNumber;j++)
	{
		//程序初始化
		initiate();
		//每次运行的最大迭代次数
		for(i=0;i<MaxNumber;i++)
		{
			//个体排序
		    order();
			//单纯型法
		    simplex_method();
			//确定当前最优解
			determine();
			//判断是否收敛
		    if((Max_Largest_Number-bestindividual.fitness)<Error)
			{
				for(int k=0;k<4;k++)
					cout<<bestindividual.variable[k]<<"  ";
				printf("%d,%e\t",i,Max_Largest_Number-bestindividual.fitness);
				cout<<endl;
				sum++;
				total=total+i;
				break;
			}
		    else
			{
				//杂交算子
		        crossover();
			    //变异算子
		        area_range();
			    //选择算子
		        selection();
			}
		}
		for(int k=0;k<4;k++)
					cout<<bestindividual.variable[k]<<"  ";
				printf("%d,%e\t",i,Max_Largest_Number-bestindividual.fitness);
				cout<<endl;
	}
	printf("\n");
	printf("total sum is %d\n",sum);
	printf("average convergence is %f\n",total/sum);
	system("pause");
}