优化设计-惩罚函数法程序c语言_精品文档.doc
《优化设计-惩罚函数法程序c语言_精品文档.doc》由会员分享，可在线阅读，更多相关《优化设计-惩罚函数法程序c语言_精品文档.doc（6页珍藏版）》请在冰豆网上搜索。
优化设计-惩罚函数法程序（c语言）
#include
#include
#definem10
floatf（floatx[],floatp）;
voidmjtf（intn,floatx0[],floath,floats[],floata[],floatb[],floatp）;
voidmhjfgf（intn,floata[],floatb[],floatflag,floatx[],floatp）;
voidmbwef（intn,floatx0[],floath,floatflag,floata[],floatb[],floatx[],floatp）;
floatmax（floatl）
{
floatresu;
if（l<=0）
resu=0;
else
resu=l*l;
returnresu;
}
floatf（floatx[],floatp）
{
floatresult,t1;
floatt2,t3,t4;
t1=（8-x[0]-x[1]）;
result=x[0]*x[0]+x[1]*x[1]-x[0]*x[1]-10*x[0]-4*x[1]+60+p*max（t1）;
returnresult;
}
floatg（floatx[]）
{
floatre;
re=x[0]*x[0]+x[1]*x[1]-x[0]*x[1]-10*x[0]-4*x[1]+60;
returnre;
}
voidmjtf（intn,floatx0[],floath,floats[],floata[],floatb[],floatp）
{
inti;
floatx1[m],x2[m],x3[m],f1,f2,f3;
for（i=0;i{
x1[i]=x0[i];
x2[i]=x0[i]+h*s[i];
}
f1=f（x1,p）;
f2=f（x2,p）;
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（x3,p）;
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（x3,p）;
}/*已完成*/
for（i=0;i{
if（x1[i]{
a[i]=x1[i];
b[i]=x3[i];
}
else
{
a[i]=x3[i];
b[i]=x1[i];
}
}
}
voidmhjfgf（intn,floata[],floatb[],floatflag,floatx[],floatp）
{
inti;
floatx1[m],x2[m],f1,f2,sum;
for（i=0;ix1[i]=b[i]-（float）0.618*（b[i]-a[i]）;
f1=f（x1,p）;
for（i=0;ix2[i]=a[i]+（float）0.618*（b[i]-a[i]）;
f2=f（x2,p）;
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（x1,p）;
}
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（x2,p）;
}
sum=0;
for（i=0;isum+=（b[i]-a[i]）*（b[i]-a[i]）;
printf（"%f\n",sum）;
}while（sqrt（sum）>flag*0.1）;/*判断是否未达到精度要求，若未达到则返回继续缩小区间*/
for（i=0;ix[i]=（float）0.5*（b[i]+a[i]）;/*已达到，输出极小点*/
}
voidmbwef（intn,floatx0[],floath,floatflag,floata[],floatb[],floatx[],floatp）
{
inti,j,k,r,c;
floatx1[m],x2[m],f0,f1,f2,fn[m],s[m][m],sum;
for（i=0;ifor（k=0;kif（i==k）
s[i][k]=1;
else
s[i][k]=0;
k=1;
c=5;
p=0.3;
while
（1）
{
for（i=0;ix1[i]=x0[i];
for（i=0;i{
mjtf（n,x1,h,s[i],a,b,p）;/*调用多维进退法子程序*/
mhjfgf（n,a,b,flag,x1,p）;/*调用多维黄金分割法子程序*/
fn[i]=f（x0,p）-f（x1,p）;/*计算函数下降值*/
}
for（i=0;ix2[i]=2*x1[i]-x0[i];
for（i=1;iif（fn[0]{
fn[0]=fn[i];
r=i;
}
else
r=0;
f0=f（x0,p）;/*计算始点、终点和映射点的函数值*/
f1=f（x1,p）;
f2=f（x2,p）;
if（f2>=f0||（f0-2*f1+f2）*（f0-f1-fn[0]）*（f0-f1-fn[0]）>=0.5*fn[0]*（f0-f2）*（f0-f2））
{/*判断是否需要换方向*//*不需要换*/
sum=0;/*计算一轮中终点与始点的距离*/
for（i=0;isum+=（x1[i]-x0[i]）*（x1[i]-x0[i]）;
if（f1<=f2）/*判断将终点还是将映射点赋给下一轮始点*/
for（i=1;ix0[i]=x1[i];
else
for（i=1;ix0[i]=x2[i];
}
else
{/*需要换*/
for（i=r;ifor（j=0;js[i][j]=s[i+1][j];
for（i=0;is[n][i]=x1[i]-x0[i];
mjtf（n,x1,h,s[n],a,b,p）;/*按新方向搜索一次*/
mhjfgf（n,a,b,flag,x1,p）;
sum=0;/*计算一轮中终点与始点的距离*/
for（i=0;isum+=（x1[i]-x0[i]）*（x1[i]-x0[i]）;
for（i=0;ix0[i]=x1[i];
}
if（sqrt（sum）<=flag）/*判断是否满足精度要求*/
break;
else
{k+=1;
p=p*c;
}
}
for（i=0;ix[i]=x1[i];
}
/*主程序*/
voidmain（）
{
inti,n;
floath,flag,x0[m],a[m],b[m],x[m],p,l;
printf（"\n<鲍威尔法>\n"）;
printf（"请输入维数:
\n"）;
scanf（"%d",&n）;
printf（"请输入初始
升级会员 