惩罚函数法

#include #include #define m 10 //f函数 float f(float x[],float p); //进退法 void mjtf(int n,float x0[],float h,float s[],float a[],float b[],float p); void mhjfgf(int n,float a[],float b[],float flag,float x[],float p); //鲍威尔法 void mbwef(int n,float x0[],float h,float flag,float a[],float b[],float x[],float p); float max(float l) { float resu; if(l<=0) resu=0; else resu=l*l;

#include 
#include 
#define m 10
//f函数
float f（float x[]float p）;
//进退法
void mjtf（int nfloat x0[]float hfloat s[]float a[]float b[]float p）; void mhjfgf（int nfloat a[]float b[]float flagfloat x[]float p）;
//鲍威尔法
void mbwef（int nfloat x0[]float hfloat flagfloat a[]float b[]float x[]float p）;

float max（float l）
{
float resu;
if（l<=0）
resu=0;
else
resu=l*l;
return resu;
}

float f（float x[]float p）
{
float resultt1;
float t2t3t4;
t1=（8-x[0]-x[1]）;
result=（x[0]-2）*（x[0]-2）+（x[1]-3）*（x[1]-3）+（x[2]-4）*（x[2]-4）+p*max（t1）;
return result;
}

float g（float x[]）
{
float re;
re=（x[0]-2）*（x[0]-2）+（x[1]-3）*（x[1]-3）+（x[2]-4）*（x[2]-4）;
return re;
}

void mjtf（int nfloat x0[]float hfloat s[]float a[]float b[]float p） {
int i;
float x1[m]x2[m]x3[m]f1f2f3;
for（i=0;i{
 x1[i]=x0[i];
x2[i]=x0[i]+h*s[i];
}
f1=f（x1p）;
f2=f（x2p）;
if（f2>=f1）             /*判断搜索方向*/
{                      /*搜索方向为反向，转身*/
h=（-1）*h;
for（i=0;ix3[i]=x1[i];
f3=f1;
for（i=0;ix1[i]=x2[i];
f1=f2;
for（i=0;ix2[i]=x3[i];
f2=f3;
}
for（i=0;ix3[i]=x2[i]+h*s[i];
f3=f（x3p）;
while（f3{
h=2*h;
for（i=0;ix1[i]=x2[i];
f1=f2;
for（i=0;ix2[i]=x3[i];
f2=f3;
for（i=0;ix3[i]=x2[i]+h*s[i];
f3=f（x3p）;
}
for（i=0;i{
if（x1[i]{
a[i]=x1[i];
b[i]=x3[i];
}
else
{                                       /*搜索方向为正向*/ /*计算第三试点*/ /*判断是否未完成搜索*/ /*未完成，继续搜索*/ /*已完成*/ /*输出初始搜索区间*/
a[i]=x3[i];
b[i]=x1[i];
}
}
}


void mhjfgf（int nfloat a[]float b[]float flagfloat x[]float p） {
int i;
float x1[m]x2[m]f1f2sum;
for（i=0;ix1[i]=b[i]-（float）0.618*（b[i]-a[i]）;
f1=f（x1p）;
for（i=0;ix2[i]=a[i]+（float）0.618*（b[i]-a[i]）;
f2=f（x2p）;
do
{
if（f1<=f2）           /*判断消去区间*/ { /*消去右*/
for（i=0;ib[i]=x2[i];
for（i=0;ix2[i]=x1[i];
f2=f1;
for（i=0;ix1[i]=b[i]-（float）0.618*（b[i]-a[i]）;
f1=f（x1p）;
}
else
{                                    /*消去左*/
for（i=0;ia[i]=x1[i];
for（i=0;ix1[i]=x2[i];
f1=f2;
for（i=0;ix2[i]=a[i]+（float）0.618*（b[i]-a[i]）;
f2=f（x2p）;
}
sum=0;
for（i=0;isum+=（b[i]-a[i]）*（b[i]-a[i]）;
}while（sqrt（sum）>fl