机械优化设计C语言程序

四、DFP变尺度法

#include #include #include #include

#define tt 0.01 #define ff 1.0e-6 #define ac 1.0e-6 #define ad 1.0e-6 #define n 2

double ia;

double fny(double *x) {

double x1=x[0],x2=x[1];

double f;

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

return f; }

double * iterate(double *x,double a,double *s) {

double *x1; int i;

x1=(double *)malloc(n*sizeof(double)); for(i=0;i

double func(double *x,double a,double *s) {

double *x1; double f;

x1=iterate(x,a,s); f=fny(x1);

return f; }

void finding(double a[3],double f[3],double *xk,double *s) {

double t=tt; int i;

double a1,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; }

double lagrange(double *xk,double *ft,double *s) { int i;

double a[3],f[3]; double b,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;}

else if(*ft

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;}

else if(*ft

} } }

if(*ft>f[1]) {*ft=f[1];aa=a[1];} return aa; }

double *gradient(double *xk) {

double *g,f1,f2,q; int i;

g=(double*)malloc(n*sizeof(double)); f1=fny(xk); for(i=0;i

xk[i]=xk[i]+q; f2=fny(xk); g[i]=(f2-f1)/q; xk[i]=xk[i]-q; }

return g; }

double * bfgs(double *xk) {

double u[n],v[n],h[n][n],dx[n],dg[n],s[n]; double aa,ib;

double *ft,*xk1,*g1,*g2,*xx,*x0=xk; double fi; int i,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

ib=0;

if(ia==1) { xx=xk; break; } ib=0;

for(i=0;i

for(i=0;i

if((fabs(g2[i])>=ac)&&(fabs(g2[i]-g1[i])>=ac)) {ib=ib+1;}

if(ib==0) { xx=xk1; break; } fi=*ft; if(k==n-1) { int j; xk=xk1;

for(i=0;i

h[i][j]=0;

if(j==i) h[i][j]=1; }

g1=g2; k=-1; }

else { int j;

double a1=0,a2=0; for(i=0;i

dg[i]=g2[i]-g1[i]; dx[i]=xk1[i]-xk[i]; }

for(i=0;i