C++ 语言片段如何用其他语言也实现？

cyandlhb 2009-01-15 11:36:12

/*基于复合形法的圆柱弹簧优化设计*/

#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[i]<0)
goto s333;
}
return 1;
s333:return 0;}

void xcent (int n, int k, int ll ,int lh, float x0[ ], float xcom[ ][6])/*求平均值*/
{ 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[i][l];
}
if(lh>-1)
x0[i]=xs/(ll-1);
else
x0[i]=xs/ll;
}
}
void fxse(int n,int k,float x[],float xcom[][6],float fxk[])/*对函数值fxk进行从小到大的排序*/
{int l,lp,lp1,i;
float w;
/*printf("\n fxk[0]=%f,fxk[1]=%f,fxk[2]=%f,fxk[3]=%f,fxk[4]=%f,fxk[5]=%f",fxk[0],fxk[1],fxk[2],fxk[3],fxk[4],fxk[5]);*/

for(l=0;l<k-1;l++)
for(lp=0;lp<k-1;lp++)
{
lp1=lp+1;
if(fxk[lp]<=fxk[lp1])
{w=fxk[lp];
fxk[lp]=fxk[lp1];
fxk[lp1]=w;
for(i=0;i<n;i++)
{x[i]=xcom[i][lp];
xcom[i][lp]=xcom[i][lp1];
xcom[i][lp1]=x[i];
}
}
}
}

/*add function*/
void complex(int n, int k,int kg,float ep,float x[],float bl[],float bu[],float xcom[][6],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(n,sizeof(float));
float *g=(float*)calloc(n,sizeof(float));
s5:for(i=0;i<n;i++)
x[i]=bl[i]+rand()/40000.0*(bu[i]-bl[i]);
iw=gau(x,g,kg);
*ng=*ng+1;
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=0;ll<k-1;ll++)
{xcent(n,k,ll+1,lh,x0,xcom);

iw=gau(x0,g,kg);
*ng=*ng+1;
if(iw==0)goto s5;
for(i=0;i<n;i++)
x[i]=xcom[i][ll+1];

s24: iw=gau(x,g,kg);
*ng=*ng+1;
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);*nf=*nf+1;fxk[l]=fx;
printf("\n x[0]=%f,x[1]=%f,x[2]=%f,x[3]=%f",x[0],x[1],x[2],x[3]);
printf("\n fxk[l]=%f",fxk[l]);
}
it=0;
s14:it=it+1;
printf("\n\n +++compute +++\n");;
printf("\n iter = %d ",it);
lh=-1;
/*add function */
xcent(n,k,k,lh,x0,xcom);
/*printf("\n x0[0]=%f,x0[1]=%f,x0[2]=%f,x0[3]=%f",x0[0],x0[1],x0[2],x0[3]);*/
fx0=objfx(x0);
printf("\n fxk=%f",fx0);
*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[i]);
for(i=0;i<n;i++)
printf("\n G(%d)mid= %f",i,g[i]);
sdx=0;
for(l=0;l<k;l++)
sdx=sdx+(fx0-fxk[1])*(fx0-fxk[1]);
sdx=sqrt(sdx/(float)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);

*ng=*ng+1;
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);
*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[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);*nf=*nf+1;
/*free(x0);free(xh);free(xr);free(g);free(fxk);*/
}

float objfx(float x[])/*目标函数*/

{return 0.192457e-4*(x[1]+2)*x[0]*x[0]*x[2];
}

void constraint(float x[],float g[])/*约束条件*/
{g[0]=350.0-(163.0*float(pow(x[0],-2.86))*float(pow(x[2],0.86)));
g[1]=0.4e-2*pow(x[0],-4)*x[1]*pow(x[2],3)-10.0;
g[2]=3.7*x[2]-0.44e-2*pow(x[0],-4)*pow(x[2],3)-(x[1]+1.50)*x[0];
g[3]=0.356e6*x[0]*pow(x[1],-1)*pow(x[2],-2)-375;
g[4]=pow(x[0],-1)*x[2]-4.0;
g[5]=x[0]-1.0;
g[6]=4.0-x[0];
g[7]=x[1]-4.50;
g[8]=50.0-x[1];
g[9]=x[2]-10.0;
g[10]=30.0-x[2];
}
void main()
{

float x[3],bl[3],bu[3],xcom[3][6],f;
int n=3,nf=0,ng=0;
bl[0]=1.0;bl[1]=4.5;bl[2]=10.0;
bu[0]=4.0;bu[1]=50.0;bu[2]=30.0;
/*clrscr();*/
complex(n,6,11,1.0e-1,x,bl,bu,xcom,&f,&nf,&ng);
printf("\n x[0]=%f,x[1]=%f,x[2]=%f",x[0],x[1],x[2]);
printf("\n minF(x)=%f",f);
printf("\n nf=%d,ng=%d",nf,ng);

getchar();
}

...全文