哪位能告诉我正确的复合形算法啊?谢谢!
我倒. 手头上三本关于优化的书, 介绍的复合形法都有点不同,试着编了三个程序验证一下哪个比较好, 优化出的结果都离正确值偏差很多.
搞不懂怎么会这样!
请大家帮帮我!谢谢! 复合形法怎么会不一样的？不能理解
/* 以下为复合形法法子程序 */
#include "stdlib.h"
#include "math.h"
#include "stdio.h"
float objfx(float x[]);
void constraint(float x[],float g[]);
int gau(float x[],float g[],int kg)
{
int i;
constraint(x,g);
for(i=0;i<kg;i++)
{
if(g<0)
goto s333;
}
return 1;
s333:return 0;
}
void xcent(int n,int k,int ll,int lh,float x0[],float xcom[])
{
int i,l;
float xs;
for(i=0;i<n;i++)
{
xs=0;
for(l=0;l<ll;l++)
{
if(l!=lh)
xs=xs+xcom;
}
if(lh>-1)
x0=xs/(ll-1);
else
x0=xs/ll;
}
}
void fxse(int n,int k,float x[],float xcom[],float fxk[])
{
int l,lp,lp1,i;
float w;
for(l=0;l<k-1;l++)
for(lp=0;lp<k-l;lp++)
{
lp1=lp+1;
if(fxk<=fxk)
{
w=fxk;
fxk=fxk;
fxk=w;
for(i=0;i<n;i++)
{
x=xcom;
xcom=xcom;
xcom=x;
}
}
}
}
void complex(int n,int k,int kg,float ep,float x[],float bl[],float bu[],
float xcom[],float *f,int *nf,int *ng)
{
int i,iw,l,ll,lh,it;
float fx,fx0,sdx,fxh,fxr,alp;
float *x0=(float*)calloc(n,sizeof(float));
float *xh=(float*)calloc(n,sizeof(float));
float *xr=(float*)calloc(n,sizeof(float));
float *fxk=(float*)calloc(k,sizeof(float));
float *g=(float*)calloc(kg,sizeof(float));
s5: for(i=0;i<n;i++)
x=bl+rand()/40000.0*(bu-bl);
iw=gau(x,g,kg);
*ng=*ng+1;
if(iw==0)goto s5;
for(i=0;i<n;i++)
xcom=x;
for(l=1;l<k;l++)
for(i=0;i<n;i++)
xcom=bl+rand()/50000.0*(bu-bl);
lh=-1;
for(ll=1;ll<k;ll++)
{
xcent(n,k,ll,lh,x0,xcom);
iw=gau(x0,g,kg);
*ng=*ng+1;
if(iw==0)goto s5;
for(i=0;i<n;i++)
x=xcom;
s24: iw=gau(x,g,kg);
*ng=*ng+1;
if(iw==0)
{
for(i=0;i<n;i++)
x=x0+0.5*(x-x0);
goto s24;
}
else
{
for(i=0;i<n;i++)
xcom=x;
}
}
for(l=0;l<k;l++)
{
for(i=0;i<n;i++)
x=xcom;
fx=objfx(x);
*nf=*nf+1;
fxk=fx;
}
it=0;
s14: it=it+1;
printf("\n\n +++ ITERATION COMPUTE +++");
printf("\n ITER = %d",it);
lh=-1;
xcent(n,k,k,lh,x0,xcom);
fx0=objfx(x0);
*nf=*nf+1;
iw=gau(x0,g,kg);
*ng=*ng+1;
printf("\n Fmid = %f",fx0);
for(i=0;i<n;i++)
printf("\n X(%d)mid=%f",i,x0);
for(i=0;i<kg;i++)
printf("\n G(%d)mid=%f",i,g);
sdx=0;
for(l=0;l<k;l++)
sdx=sdx+(fx0-fxk)*(fx0-fxk);
sdx=sqrt(sdx/(float)k);
if(sdx<ep)goto s38;
fxse(n,k,x,xcom,fxk);
lh=0;
s22: fxh=fxk;
for(i=0;i<n;i++)
xh=xcom;
xcent(n,k,k,lh,x0,xcom);
iw=gau(x0,g,kg);
*ng=*ng+1;
if(iw==0)goto s36;
alp=1.3;
s12: for(i=0;i<n;i++)
xr=x0+alp*(x0-xh);
iw=gau(xr,g,kg);
*ng=*ng+1;
if(iw==0)
{
alp=alp*0.5;
goto s12;
}
fxr=objfx(xr);
*nf=*nf+1;
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=xr;
fxk=fxr;
goto s14;
s36: for(i=0;i<n;i++)
{
bl=xcom;
bu=x0;
}
goto s5;
s38: for(i=0;i<n;i++)
x=x0;
*f=objfx(x);
*nf=*nf+1;
free(x0);
free(xh);
free(xr);
free(g);
free(fxk);
}
来自:http://hi.baidu.com/mec_auto/blog/item/7783d3fe6a999f305d600850.html
页:
[1]