机械优化设计C语言程序.docx

机械优化设计C语言程序.docx

《机械优化设计C语言程序.docx》由会员分享，可在线阅读，更多相关《机械优化设计C语言程序.docx（18页珍藏版）》请在冰豆网上搜索。

机械优化设计C语言程序

一、黄金分割法

二、二次插值法

三、最速下降法（阶梯法）

四、变尺度法

五、鲍威尔法

一、黄金分割法

#include

#include#definer0.618#definef（x）x*x+2*xgolden（float,float,float）;

main（）

{

floata,b,e;

printf（"\n请输入区间和收敛精度:

a,b,e\n"）;scanf（"%f,%f,%f",&a,&b,&e）;

golden（a,b,e）;

}

golden（floata,floatb,floate）

{

floaty1,y2,a1,a2,A,Y;

intn=0;

a1=b-r*（b-a）;a2=a+r*（b-a）;

y1=f（a1）;

y2=f（a2）;

printf（"黄金分割法的搜索过程:

"）;

do

%d

{printf（"\na=%f,b=%f,a1=%f,a2=%f,y1=%f,y2=%f",n,a,b,a1,a2,y1,y2）;

if（y1>=y2）

{

a=a1;a1=a2;y1=y2;a2=a+r*（b-a）;y2=f（a2）;

}else

{

b=a2;a2=a1;y2=y1;a1=b-r*（b-a）;y1=f（a1）;

n++;

}while（fabs（（b-a）/b）>=e||fabs（（y2-y1）/y2）>=e）;

A=（a+b）/2;

Y=f（A）;

printf（"\n%da=%f,b=%f,a1=%f,a2=%f",n,a,b,a1,a2）;

printf（"\n结果:

\n极值点及其函数值：

A=%f,Y=%f\n",A,Y）;}

二、二次插值法

#include"stdio.h"

#include"math.h"#include"conio.h"voidmain（）{

float*area（floata1,floatp,floata[3]）;

floatf（floatx）;

floatar,fr;

floata1=10,p=0.01,e=0.000001;

floatpa[3];

area（a1,p,pa）;

a1=pa[0];

floata2=pa[1];

floata3=pa[2];

floatf1=f（a1）;

floatf2=f（a2）;

floatf3=f（a3）;

do

{ar=（（a3*a3-a2*a2）*f1+（a1*a1-a3*a3）*f2+（a2*a2-a1*a1）*f3）;ar=ar/2/（（a3-a2）*f1+（a1-a3）*f2+（a2-a1）*f3）;

fr=f（ar）;

if（ar>a2）

{if（fr>f2）

{a3=ar;f3=fr;}

elseif（fr

{a1=a2;f1=f2;

a2=ar;f2=fr;}

else

{a3=ar;a1=a2;a2=（a1+a3）/2;

f1=f2;f3=fr;f2=f（a2）;}

}

elseif（ar

{if（fr>f2）

{a1=ar;f1=fr;}elseif（fr

a2=ar;f2=fr;}

else

{a1=ar;a3=a2;a2=（a1+a3）/2;f1=fr;f3=f2;f2=f（a2）;}

}

if（fabs（a1-a3）<=e）break;

}while

（1）;if（f2

{ar=a2;fr=f2;}printf（"\nx=%f\nf（x*）=%f",ar,fr）;

}

float*area（floata1,floatp,floata[3]）{floatf（floatx）;

floata2,f2,a3,f3,temp;

floatacc=0.00001;

floatf1=f（a1）;

floatstorep=p;while

（1）

{a2=a1+p;f2=f（a2）;

if（f2>=f1）

{if（fabs（f2-f1）

elsep=-p;

}

elsebreak;

}while

（1）

{a3=a2+p;f3=f（a3）;if（f2<=f3）break;p=2*p;a1=a2;f1=f2;a2=a3;f2=f3;

}if（a1>a3）

{temp=a1;a1=a3;a3=temp;}a[0]=a1;a[1]=a2;a[2]=a3;

returna;

}

floatf（floatx）

{

floaty=pow（fabs（x-1）,1.5）+pow（fabs（x-1）,2.7）;

returny;

}

三、最速下降法

#include

#include

floatfun1（floatx1,floatx2）

{floata;

a=2*x1;

returna;

}

floatfun2（floatx1,floatx2）

{floatb;

b=50*x2;

returnb;

}

floatfun3（floatx1,floatx2）

{floaty;

y=x1*x1+25*x2*x2;

returny;

}

main（）

{floatt,e=0.01,x1=0,x2=1,a,b,y,m;

a=fun1（x1,x2）;

b=fun2（x1,x2）;

m=sqrt（a*a+b*b）;

while（m>e）

{t=（a*a+b*b）/（2*a*a+50*b*b）;x1=x1-a*t;x2=x2-b*t;a=fun1（x1,x2）;b=fun2（x1,x2）;m=sqrt（a*a+b*b）;

}

y=fun3（x1,x2）;

printf（"Theminis%f",y）;

}

四、DFP变尺度法

#include

#include

#include

#include#definett0.01

#defineff1.0e-6

#defineac1.0e-6

#definead1.0e-6

#definen2

doubleia;

doublefny（double*x）

{

doublex1=x[0],x2=x[1];

doublef;f=x1*x1+2*x2*x2-4*x1-2*x1*x2;

returnf;

}

double*iterate（double*x,doublea,double*s）{

double*x1;

inti;

x1=（double*）malloc（n*sizeof（double））;for（i=0;i

x1[i]=x[i]+a*s[i];

returnx1;

}

doublefunc（double*x,doublea,double*s）

{double*x1;doublef;x1=iterate（x,a,s）;f=fny（x1）;

returnf;

}

voidfinding（doublea[3],doublef[3],double*xk,double*s）{

doublet=tt;

inti;

doublea1,f1;

a[0]=0;f[0]=func（xk,a[0],s）;

for（i=0;;i++）{

a[1]=a[0]+t;f[1]=func（xk,a[1],s）;

if（f[1]=ad）

{

t=-t;

a[0]=a[1];f[0]=f[1];

}

else

{

if（ia==1）return;//break

t=t/2;ia=1;

}

}

for（i=0;;i++）

{a[2]=a[1]+t;

f[2]=func（xk,a[2],s）;

if（f[2]>f[1]）break;

t=2*t;

a[0]=a[1];f[0]=f[1];

a[1]=a[2];f[1]=f[2];

}

if（a[0]>a[2]）

{

a1=a[0];

f1=f[0];

a[0]=a[2];

f[0]=f[2];

a[2]=a1;

f[2]=f1;

return;

}

doublelagrange（double*xk,double*ft,double*s）

{

inti;

doublea[3],f[3];

doubleb,c,d,aa;

finding（a,f,xk,s）;

for（i=0;;i++）

{

if（ia==1）{aa=a[1];*ft=f[1];break;}

d=（pow（a[0],2）-pow（a[2],2））*（a[0]-a[1]）-（pow（a[0],2）-pow（a[1],2））*（a[0]-a[2]）;if（fabs（d）==0）break;

c=（（f[0]-f[2]）*（a[0]-a[1]）-（f[0]-f[1]）*（a[0]-a[2]））/d;

if（fabs（c）==0）break;

b=（（f[0]-f[1]）-c*（pow（a[0],2）-pow（a[1],2）））/（a[0]-a[1]）;

aa=-b/（2*c）;

*ft=func（xk,aa,s）;

if（fabs（aa-a[1]）<=ad）{if（*ft>f[1]）aa=a[1];break;}

if（aa>a[1]）

{

if（*ft>f[1]）{a[2]=aa;f[2]=*ft;}

elseif（*ft

elseif（*ft==f[1]）

{a[2]=aa;a[0]=a[1];

f[2]=*ft;f[0]=f[1];

a[1]=（a[0]+a[2]）/2;

f[1]=func（xk,a[1],s）;

}

}

else

{

if（*ft>f[1]）{a[0]=aa;f[0]=*ft;}

elseif（*ft

elseif（*ft==f[1]）{a[0]=aa;a[2]=a[1];

f[0]=*ft;f[2]=f[1];

a[1]=（a[0]+a[2]）/2;

f[1]=func（xk,a[1],s）;

}

}

}

if（*ft>f[1]）{*ft=f[1];aa=a[1];}

returnaa;

}

double*gradient（double*xk）

{

double*g,f1,f2,q;

inti;

g=（double*）malloc（n*sizeof（double））;

f1=fny（xk）;

for（i=0;i

{q=ff;

xk[i]=xk[i]+q;f2=fny（xk）;g[i]=（f2-f1）/q;xk[i]=xk[i]-q;

}

returng;

}

double*bfgs（double*xk）

{

doubleu[n],v[n],h[n][n],dx[n],dg[n],s[n];

doubleaa,ib;

double*ft,*xk1,*g1,*g2,*xx,*x0=xk;

doublefi;

inti,j,k;

ft=（double*）malloc（sizeof（double））;

xk1=（double*）malloc（n*sizeof（double））;

for（i=0;i

{

s[i]=0;

for（j=0;j

{

h[i][j]=0;

if（j==i）h[i][j]=1;

}

}

g1=gradient（xk）;

fi=fny（xk）;x0=xk;

for（k=0;k

if（ia==1）{xx=xk;break;}

ib=0;for（i=0;i

s[i]+=-h[i][j]*g1[j];aa=lagrange（xk,ft,s）;xk1=iterate（xk,aa,s）;g2=gradient（xk1）;

for（i=0;i=ac）&&（fabs（g2[i]-g1[i]）>=ac））{ib=ib+1;}

if（ib==0）{xx=xk1;break;}fi=*ft;

if（k==n-1）

{intj;

xk=xk1;for（i=0;i

if（j==i）h[i][j]=1;

}g1=g2;k=-1;

}

else

{

intj;

doublea1=0,a2=0;for（i=0;i

for（i=0;i

{intj;

u[i]=0;v[i]=0;for（j=0;j

{u[i]=u[i]+dg[j]*h[j][i];v[i]=v[i]+dg[j]*h[i][j];

}

}

for（j=0;j

{a1+=dx[j]*dg[j];a2+=v[j]*dg[j];

}if（fabs（a1）!

=0）

{a2=1+a2/a1;for（i=0;i

for（j=0;j

}xk=xk1;g1=g2;

}

}if（*ft>fi）{*ft=fi;xx=xk;}xk=x0;

returnxx;

}

voidmain（）

{

intk;

double*xx,f;

doublexk[n]={1,1};xx=bfgs（xk）;

f=fny（xx）;printf（"\n\nTheOptimalDesignResultIs:

\n"）;for（k=0;k

{printf（"\n\tx[%d]*=%f",k+1,xx[k]）;}printf（"\n\tf*=%f",f）;

getch（）;

五、鲍威尔法

#include"stdio.h"

#include"stdlib.h"

#include"math.h"

doubleobjf（doublex[]）

{doubleff;ff=x[0]*x[0]+2*x[1]*x[1]-4*x[0]-2*x[0]*x[1];return（ff）;

}

n,doublea[],doubleb[]）

voidjtf（doublex0[],doubleh0,doubles[],int{inti;

double*x[3],h,f1,f2,f3;

for（i=0;i<3;i++）

x[i]=（double*）malloc（n*sizeof（double））;

h=h0;for（i=0;i

f1=objf（x[0]）;for（i=0;i

*（x[1]+i）=*（x[0]+i）+h*s[i];f2=objf（x[1]）;

if（f2>=f1）

{h=-h0;

for（i=0;i

*（x[2]+i）=*（x[0]+i）;

f3=f1;for（i=0;i

{*（x[0]+i）=*（x[1]+i）;

*（x[1]+i）=*（x[2]+i）;

}

f1=f2;

f2=f3;

}

for（;;）

{h=2.*h;

for（i=0;i

*（x[2]+i）=*（x[1]+i）+h*s[i];

f3=objf（x[2]）;if（f2

break;

else

{for（i=0;i

{*（x[0]+i）=*（x[1]+i）;

*（x[1]+i）=*（x[2]+i）;

}

f1=f2;

f2=f3;

}

}

if（h<0.）

for（i=0;i

{a[i]=*（x[2]+i）;

b[i]=*（x[0]+i）;

}

else

for（i=0;i

{a[i]=*（x[0]+i）;

b[i]=*（x[2]+i）;

}for（i=0;i<3;i++）free（x[i]）;

}

doublegold（doublea[],doubleb[],doubleeps,intn,doublexx[]）{

inti;

doublef1,f2,*x[2],ff,q,w;for（i=0;i<2;i++）

x[i]=（double*）malloc（n*sizeof（double））;for（i=0;i

{*（x[0]+i）=a[i]+0.618*（b[i]-a[i]）;*（x[1]+i）=a[i]+0.382*（b[i]-a[i]）;

}f1=objf（x[0]）;f2=objf（x[1]）;

do{if（f1>f2）

{for（i=0;i

*（x[0]+i）=*（x[1]+i）;

}f1=f2;

for（i=0;i

*（x[1]+i）=a[i]+0.382*（b[i]-a[i]）;

f2=objf（x[1]）;

}

else

{for（i=0;i

{a[i]=*（x[1]+i）;

*（x[1]+i）=*（x[0]+i）;

}

f2=f1;

for（i=0;i

*（x[0]+i）=a[i]+0.618*（b[i]-a[i]）;f1=objf（x[0]）;

}

q=0;

for（i=0;i

q=q+（b[i]-a[i]）*（b[i]-a[i]）;

w=sqrt（q）;

}while（w>eps）;

for（i=0;i

xx[i]=0.5*（a[i]+b[i]）;

ff=objf（xx）;

for（i=0;i<2;i++）

free（x[i]）;

return（ff）;

}

doubleoneoptim（doublex0[],doubles[],doubleh0,doubleepsg,intn,doublex[]）{double*a,*b,ff;

a=（double*）malloc（n*sizeof（double））;

b=（double*）malloc（n*sizeof（double））;jtf（x0,h0,s,n,a,b）;

ff=gold（a,b,epsg,n,x）;

free（a）;

free（b）;return（ff）;

}

doublepowell（doublep[],doubleh0,doubleeps,doubleepsg,intn,doublex[]）{inti,j,m;

double*xx[4],*ss,*s;

doublef,f0,f1,f2,f3,fx,dlt,df,sdx,q,d;

ss=（double*）malloc（n*（n+1）*sizeof（double））;s=（double*）malloc（n*sizeof（double））;

for（i=0;i

{for（j=0;j<=n;j++）

*（ss+i*（n+1）+j）=0;

*（ss+i*（n+1）+i）=1;

}

for（i=0;i<4;i++）

xx[i]=（double*）malloc（n*sizeof（double））;for（i=0;i

*（xx[0]+i）=p[i];

for（;;）

{for（i=0;i

x[i]=*（xx[1]+i）;

}f0=f1=objf（x）;dlt=-1;

for（j=0;j

{for（i=0;i

}f=oneoptim（xx[0],s,h0,epsg,n,x）;df=f0-f;

if（df>dlt）

{dlt=df;

m=j;

}

}

sdx=0.;

for（i=0;i

free（s）;

for（i=0;i<4;i++）free（xx[i]）;return（f）;

}

for（i=0;i

*（xx[2]+i）=x[i];

f2=f;

for（i=0;i

fx=objf（x）;

f3=fx;q=（f1-2*f2+f3）*（f1-f2-dlt）*（f1-f2-dlt）;d=0.5*dlt*（f1-f3）*（f1-f3）;

if（（f3

{if（f2<=f3）

for（i=0;i

*（xx[0]+i）=*（xx[2]+i）;

else

for（i=0;i

}

else

{for（i=0;i

}f=oneoptim（xx[0],s,h0,epsg,n,x）;for（i=0;i

for（j=m+1;j<=n;j++）

for（i=0;i

}

}

}

voidmain（）

{doublep[]={1,1};

doubleff,x[2];ff=powell（p,0.3,0.001,0.0001,2,x）;

printf（"输出最优点及其目标函数值：

\n"）;printf（"x[0]=%f,x[1]=%f,f