razor.cpp

复合形优化算法的VC2005环境下的实现实例

#include "StdAfx.h"
#include "razor.h"

razor::razor(int dimens,int num)
{
	dimension = dimens;
	vertex_num = num;
	g_function = NULL;
	a = (double)1.3;
	x = new double[dimens*num];
	l = new double[dimens];
	h = new double[dimens];
	xc = new double[dimens];
	xr = new double[dimens];
	xe = new double[dimens];
}

razor::~razor(void)  //释放申请的各种空间
{
	delete []x;delete []l;delete []h;delete []xc;delete []xe;
	g_node *temp,*del;
	temp = g_function;
	del = temp;
	while(temp!=NULL){
		temp = temp->next;
		delete del;
		del = temp;}
}

void razor::Add_g(double (*p)(double*))
{
	g_node *node,*temp1,*temp2;
	temp1 = g_function;
	temp2 = temp1;
	while(temp1!=NULL){
		temp2 = temp1;
		temp1 = temp1->next;}
	node = new g_node;
	node->g = p;
	node->next = NULL;
	if(temp1 == g_function)
		g_function = node;
	else
		temp2->next = node;
}

void razor::Init()
{
	int i,j=-1,k=0,m=0;
	cout<<"请选择初始点坐标产生方法: 1.随机选取  2.输入选取"<<endl<<"   ";
	while((j!=1)&&(j!=2))
		cin>>j;
	if(j==1)
		cout<<"请输入各自变量取值范围:";
	else
		cout<<"请输入各自变量取值范围和初始点坐标:";
	while(m==0){
		if((k!=0)&&(j==2))
			cout<<endl<<"坐标不在定义域内,请重新输入"<<endl;
		if(k>=RANDOM_H){
			cout<<"随机生成次数已超过"<<RANDOM_H<<"次,请重新确定各自变量取值范围:"<<endl;
			k = 0;}
		for(i=0;i<dimension;i++){
			if(k==0){
				cout<<endl<<" x"<<i<<"下限 = ";cin>>l[i];
				cout<<" x"<<i<<"上限 = ";cin>>h[i];	}
			if(j==2)
			{cout<<" x"<<i<<"初始值 = ";cin>>x[i];}}
		if(j==1)random_vertex(x);
		m = Is_fillin(x);
		k++;}
///*
	for(i=0;i<dimension;i++)
		cout<<x[i]<<" ";
	cout<<endl;
//*/
}
void razor::optimize()
{
	int step,num;
	Init();                  //step 1
S2:
	num = 0;
	construct_shape();       //step 2 3
S4:
	if(get_queue()<DELTA)    //step 4
		goto S10;
	num++;
	if(num>NUMLIMUP)         //如果计算次数过多,强制结束,给出当前结果
		goto S10;
S5:
	if(get_XC()==2)          //step 5
		goto S2;
S6:
	if(get_XR()==9)          //step 6
		goto S9;
	step = compare();        //step 7
	if(step==4)
		goto S4;
S8:
	step = get_XE();         //step 8
	if(step==8)
		goto S8;
	if(step==4)
		goto S4;
S9:
	step = zoom();           //step 9
	if(step==5)
		goto S5;
	else
		goto S6;
S10:
	step = over();           //结束
	if(step == 2)
		goto S2;
	return;
}

int razor::Is_fillin(double *x)
{
	int i;
	for(i=0;i<dimension;i++)
	{
		if(x[i]<l[i])return 0;
		if(x[i]>h[i])return 0;
	}
	g_node *temp;
	temp = g_function;
	while(temp!=NULL){
		if((*temp->g)(x)>0)
			return 0;
		temp = temp->next;}
	return 1;
}

void razor::random_vertex(double *x)
{
	int i;double r;
	for(i=0;i<dimension;i++){
		r = (double)(rand()%10000)/10000;
		x[i] = l[i]+(h[i]-l[i])*r;}	
}

void razor::construct_shape()
{
	a = 1.3;
	int *g,i,j,bad_num;
	double *x_c;
	g = new int[vertex_num];      //用于记录顶点是否合格,合格为1,否则为0
	x_c = new double[dimension];   //用于存储中点坐标
STEP3:
	for(i=0;i<vertex_num;i++)
		g[i] = 1;
	bad_num = 0;                  //初始顶点必然合格
	for(i=1;i<vertex_num;i++){    //随机生成余下顶点
		random_vertex(&x[dimension*i]);
		g[i] = Is_fillin(&x[dimension*i]);
		if(g[i]==0){
			bad_num++;
/*			
			cout<<i<<"ww->";
			for(j=0;j<dimension;j++)
				cout<<"\t"<<x[dimension*i+j];
			cout<<endl;
//*/
		}}
GET_C:
	if(bad_num==0)
		goto OVER;
	for(i=0;i<dimension;i++)
		x_c[i] = 0;
	for(i=0;i<vertex_num;i++){
		if(g[i]==1)
			for(j=0;j<dimension;j++)
				x_c[j] += x[i*dimension+j];}
	for(i=0;i<dimension;i++)
		x_c[i] /= vertex_num-bad_num;
/*
	cout<<"x_c->";
	for(i=0;i<dimension;i++)
		cout<<x_c[i]<<"\t";
	cout<<endl;
//*/
	if(Is_fillin(x_c)==0)
		goto STEP3;
	for(i=1;i<vertex_num;i++)
		if(g[i]==0){
			while(g[i]==0){
				for(j=0;j<dimension;j++)
					x[dimension*i+j] = x_c[j] + (double)0.5*(x[dimension*i+j] - x_c[j]);
				g[i] = Is_fillin(&x[dimension*i]);}
			bad_num--;
			goto GET_C;}
OVER:
	delete []g;delete []x_c;
///*
	for(i=0;i<vertex_num;i++){
		cout<<i<<"->";
		for(j=0;j<dimension;j++)
			cout<<"\t"<<x[dimension*i+j];
		cout<<endl;}
//*/
	return;
}

double razor::get_queue()
{
	double *y;
	y = new double[vertex_num];
	int i;
	for(i=0;i<vertex_num;i++)
		y[i] = (*targetFunc)(&x[dimension*i]);
	g = 0;b = 0;s = 0;
	for(i=1;i<vertex_num;i++){
		if(y[i]<y[g])
			g = i;}
	for(i=1;i<vertex_num;i++){
		if(y[i]>y[b])
			b = i;
		if(b==0)
			s = 1;}
	for(i=1;i<vertex_num;i++){
		if(i!=b)
			if(y[i]>y[s])
				s = i;}
	cout<<y[0]<<" "<<y[1]<<" "<<y[2]<<endl;
	cout<<g<<" =  "<<y[g]<<" "<<b<<" =  "<<y[b]<<" "<<s<<" = "<<y[s]<<endl;
	cout<<"****** S4 over ******"<<endl;
	return y[b]-y[g];
}

int razor::get_XC()
{
	int i,j;
	cout<<"xc: ";
	for(i=0;i<dimension;i++){
		xc[i] = 0;
		for(j=0;j<vertex_num;j++)
			xc[i] += x[dimension*j+i];
		xc[i] -= x[dimension*b+i];
		xc[i] /= vertex_num-1;
		cout<<xc[i]<<" ";}
	cout<<endl;
	if(Is_fillin(xc)==1){
		cout<<"****** S5 -> S6 ******"<<endl;
		return 1;
	}
	else{
		for(i=0;i<dimension;i++){
			x[i] = x[dimension*g+i];
			l[i] = MINIMUM(x[i],xc[i]);
			h[i] = MAXIMUM(x[i],xc[i]);
			cout<<x[i]<<" ";}
		cout<<endl;
		cout<<"****** S5 -> S2 ******"<<endl;
		return 2;}
}

int razor::get_XR()
{
	int i;
	cout<<"xr: ";
	for(i=0;i<dimension;i++){
		xr[i] = xc[i] + a * (xc[i] - x[dimension*b+i]);
		cout<<xr[i]<<" ";}
	cout<<endl;
	if(Is_fillin(xr)==1)
	{
		cout<<"****** S6 -> S7 ******"<<endl;
		return 1;
	}
	else
	{
		cout<<"****** S6 -> S9 ******"<<endl;
		return 9;
	}
}

int razor::compare()
{
	int i;
	cout<<"f(xr):"<<(*targetFunc)(xr)<<" f(s):"<<(*targetFunc)(&x[dimension*s])<<" f(g):"<<(*targetFunc)(&x[dimension*g])<<" f(b):"<<(*targetFunc)(&x[dimension*b])<<endl;
	if((*targetFunc)(xr)<(*targetFunc)(&x[dimension*s])){
		if((*targetFunc)(xr)<(*targetFunc)(&x[dimension*g]))
		{
			cout<<"****** S7 -> S8 ******"<<endl;
			return 8;}
		else{
			for(i=0;i<dimension;i++)
				x[dimension*b+i] = xr[i];
			cout<<"****** S7 -> S4 ******"<<endl;
			return 4;}}
	else{
		if((*targetFunc)(xr)<(*targetFunc)(&x[dimension*b])){
			for(i=0;i<dimension;i++)
				x[dimension*b+i] = xr[i];
			cout<<"****** S7 -> S4 ******"<<endl;
			return 4;}
		else
		{
			cout<<"****** S7 -> S9 ******"<<endl;
			return 9;
		}
	}
	return 1;
}

int razor::get_XE()
{
	int i;
	cout<<"xe: ";
	for(i=0;i<dimension;i++){
		xe[i] = xc[i] + 2 * (xr[i] - xc[i]);
		cout<<" "<<xe[i];}
	cout<<endl;
	if(Is_fillin(xe)==1){
		if((*targetFunc)(xe)<(*targetFunc)(xr)){
			for(i=0;i<dimension;i++)
				xr[i] = xe[i];
			cout<<"****** S8 -> S8 ******"<<endl;
			return 8;}
		else{
			for(i=0;i<dimension;i++)
				x[dimension*b+i] = xr[i];
			cout<<"****** S8 -> S4 ******"<<endl;
			return 4;}}
	else{
		for(i=0;i<dimension;i++)
			x[dimension*b+i] = xr[i];
		cout<<"****** S8 -> S4 ******"<<endl;
		return 4;}
}

int razor::zoom()
{
	if(a<DELTA){
		for(int i=0;i<dimension;i++)
			x[dimension*b+i] = x[dimension*s+i];
		return 5;}
	else{
		a = a/2;
		return 6;}
}

int razor::over()
{
	FILE *file = fopen("optimize.txt","w");
	cout<<"优化结果函数值: "<<(*targetFunc)(&x[dimension*g])<<endl;
	fprintf(file,"targetFunc=%f\n",(*targetFunc)(&x[dimension*g]));
	cout<<"各自变量值:"<<endl;
	int i;
	for(i=0;i<dimension;i++){
		cout<<"\tx"<<i<<" = "<<x[dimension*g+i];
		fprintf(file,"x%d=%f\n",i,x[dimension*g+i]);}
	cout<<endl;
	fclose(file);
	cout<<"是否以当前结果为初始继续优化? 1.是  other.退出"<<endl;
	int j;
	cin>>j;
	if(j==1)
	{
		for(i=0;i<dimension;i++)
			x[i] = x[dimension*g+i];
		return 2; 
	}
	return 0;
}