机械优化设计复合型法.doc

机械优化设计实验报告

(一) 题目：用复合形法求约束优化问题
；；的最优解。
基本思路：在可行域中构造一个具有K个顶点的初始复合形。对该复合形各顶点的目标函数值进行比较，找到目标函数值最大的顶点（即最坏点），然后按一定的法则求出目标函数值有所下降的可行的新点，并用此点代替最坏点，构成新的复合形，复合形的形状每改变一次，就向最优点移动一步，直至逼近最优点。
(二) 复合形法的计算步骤
1） 选择复合形的顶点数k，一般取，在可行域内构成具有k个顶点的初始复合形。
2） 计算复合形个顶点的目标函数值，比较其大小，找出最好点xL、最坏点xH、及此坏点xG..
3） 计算除去最坏点xH以外的（k-1）个顶点的中心xC。判别xC是否可行，若xC为可行点，则转步骤4）；若xC为非可行点，则重新确定设计变量的下限和上限值，即令，然后转步骤1），重新构造初始复合形。
4） 按式计算反射点xＲ，必要时改变反射系数α的值，直至反射成功，即满足式。然后xＲ以取代xＨ，构成新的复合形。
5） 若收敛条件得到满足，计算终止。约束最优解为：。
(三) 复合形法程序框图见下图：
是
否
否
否
否
是
是
是
是
否

求反射点
计算各顶点的目标函数值
f（xj）（j=1,2，……，k）
一次坏点xG代替最坏点xH

结 束

xR可行？

xC可行？

？

计算出去xH后的各顶点中心

将各定点的目标函数值和坐
标按目标函数值的大小排序
形成初始复合形的k个顶点
xj（j=1,2，……，k）
输入 n，k，ε

开 始

(四) 源程序如下：
#include<math.h>
double objfx(double x[])
{
double ff;
ff=(x[0]-5)*(x[0]-5)+4*(x[1]-6)*(x[1]-6);
return ff;
}
#include<stdio.h>
void main()
{
void comple(int n,int k,int kg,double ep,double x[],double bl[],double bu[], double xcom[][100],double *f) ;
double a[]={0,0},b[]={10,20},f,x[2],xcom[2][100];
comple(2,3,3,0.00001,x,a,b,xcom,&f);
printf(" \n\n\n 输出最优解及目标函数值：\n");
printf("\n x1=%.5f x2=%.5f f(x1,x2)=%.5f\n\n ",x[0],x[1],f);
}
/*wsd说明 n 优化模型维数 k复合形顶点数 kg 约束函数个数 ep 收敛精度 x a初始点下限 b初始点上限 xcom进程中的优化结果 *f */
#include <stdlib.h>
#include <math.h>
#include <stdio.h>
double objfx(double x[]);
void constraint(double x[],double g[]);
void constraint(double x[],double g[])
{
g[0]=64-x[0]*x[0]-x[1]*x[1];
g[1]=x[1]-x[0]-10;
g[2]=x[0]-10;
}

int gau(double x[],double g[],int kg)
{
int i;
constraint(x,g);
for(i=0;i<kg;i++)
{
if(g[i]>0)
goto s33;
}
return 1;
s33:return 0;
}
void xcent(int n,int k,int ll,int lh,double x0[],double xcom[][100])
{
int i,l;
double xs;
for(i=0;i<n;i++)
{
xs=0;
for(l=0;l<ll;l++)
{
if(l!=lh)
xs=xs+xcom[i][l];
}
if(lh>-1)
x0[i]=xs/(ll-1);
else
x0[i]=xs/ll;
}
}
void fxse(int n,int k,double x[],double xcom[][100],double fxk[])
{
int l,lp,lp1,i;
double temp;
for(l=0;l<k-1;l++)
for(lp=0;lp<k-l;lp++)
{
lp1=lp+1;
if(fxk[lp]<=fxk[lp1])
{
temp=fxk[lp];
fxk[lp]=fxk[lp1];
fxk[lp1]=temp;
for(i=0;i<n;i++)
{
x[i]=xcom[i][lp];
xcom[i][lp]=xcom[i][lp1];
xcom[i][lp1]=x[i];
}
}
}
}
void comple(int n,int k,int kg,double ep,double x[],double bl[],double bu[], double xcom[][100],double *f)
{
int i,iw,l,ll,lh,it;
double fx,fx0,sdx,fxh,fxr,alp;
double *x0=(double*)calloc(n,sizeof(double));
double *xh=(double*)calloc(n,sizeof(double));
double *xr=(double*)calloc(n,sizeof(double));
double *fxk=(double*)calloc(k,sizeof(double));
double *g=(double*)calloc(kg,sizeof(double));
s5: for(i=0;i<n;i++)
x[i]=bl[i]+rand()/40000.0*(bu[i]-bl[i]);
iw=gau(x,g,kg);
if(iw==0)
goto s5;
for(i=0;i<n;i++)
xcom[i][0]=x[i];
for(l=1;l<k;l++)
for(i=0;i<n;i++)
xcom[i][l]=bl[i]+rand()/50000.0*(bu[i]-bl[i]);
lh=-1;
for(ll=1;ll<k;ll++)
{
xcent(n,k,ll,lh,x0,xcom);
iw=gau(x0,g,kg);
if(iw==0)
goto s5;
for(i=0;i<n;i++)
x[i]=xcom[i][ll+1];
s24: iw=gau(x,g,kg);
if(iw==0)
{
for(i=0;i<n;i++)
x[i]=x0[i]+0.5*(x[i]-x0[i]);
goto s24;
}
else
{
for(i=0;i<n;i++)
xcom[i][ll+1]=x[i];
}
}
for(l=0;l<k;l++)
{
for(i=0;i<n;i++)
x[i]=xcom[i][l];
fx=objfx(x);
fxk[l]=fx;
}
it=0;
s14: it=it+1;
printf("\n\n\n N= %2d\n",it);
lh=-1;
xcent(n,k,k,lh,x0,xcom);
fx0=objfx(x0);
iw=gau(x0,g,kg);
for(i=0;i<n;i++)
printf(" X(%d)min=%.5f ",i+1,x0[i]);
printf(" Fmin=%.5f\n",fx0);
for(i=0;i<kg;i++)
printf( " G(%d)min=%.5f",i,g[i]);
sdx=0;
for(l=0;l<k;l++)
sdx=sdx+(fx0-fxk[l])*(fx0-fxk[l]);
sdx=sqrt(sdx/(double)k);
if(sdx<ep) goto s38;
fxse(n,k,x,xcom,fxk);
lh=0;
s22: fxh=fxk[lh];
for(i=0;i<n;i++)
xh[i]=xcom[i][lh];
xcent(n,k,k,lh,x0,xcom);
iw=gau(x0,g,kg);
if(iw==0)goto s36;
alp=1.3;
s12: for(i=0;i<n;i++)
xr[i]=x0[i]+alp*(x0[i]-xh[i]);
iw=gau(xr,g,kg);
if(iw==0)
{
alp=alp*0.5;
goto s12;
}
fxr=objfx(xr);
if(fxr>=fxh)
{
if(alp>1.0e-4)
{
alp=alp*0.5;
goto s12;
}
lh=lh+1;
if(lh<3)goto s22;
}
for(i=0;i<n;i++)
xcom[i][lh]=xr[i];
fxk[lh]=fxr;
goto s14;
s36: for(i=0;i<n;i++)
{
bl[i]=xcom[i][k];
bu[i]=x0[i];
}
goto s5;
s38: for(i=0;i<n;i++)
x[i]=x0[i];
*f=objfx(x);
free(x0);
free(xh);
free(xr);
free(g);
free(fxk);
}
(五)程序运行结果

8