机械优化设计惩罚函数内点法.doc

#include #include #define m 12double f(double x,double r);void jintuifa(double abmm,int n,double x0,double h,int ij,double a,double b,double r0);void hongjinfa(int n,double a,double b,double flag,double x,double r0);void baoweifa(int n,double x0,double h,double flag,double a,double b,double x,double r0);double fahansu(double tt) double ty; if(tt0) ty=-tt;else ty=0; return ty;double yuanhansu(double x) double s; /s=x0*x0+x1*x1; s=x0*x0+x1*x1+x2*x2+x3*x3; return s;double f(double x,double r) double s,t,t2; /t=1-x0; t=1-x0;t2=2-x1; / s=yuanhansu(x)-r*log(fahansu(t); s=yuanhansu(x)-r*log(fahansu(t)-r*log(fahansu(t2); return s;void jintuifa(double abmm,int n,double x0,double h,int ij,double a,double b,double r0) int i,j,z; double x1m,x2m,x3m,f1,f2,f3; double sm; for(i = 0; i n; i +) si=abiij; for(i=0;i=f1) h=(-1)*h; for(i=0;in;i+) x3i=x1i; f3=f1; for(i=0;in;i+) x1i=x2i; f1=f2; for(i=0;in;i+) x2i=x3i; f2=f3; for(i=0;in;i+) x3i=x2i+(double)h*si; f3=f(x3,r0); while(f3f2) h=2*h; for(i=0;in;i+) x1i=x2i; f1=f2; for(i=0;in;i+) x2i=x3i; f2=f3; for(i=0;in;i+) x3i=x2i+(double)h*si; f3=f(x3,r0); for(i=0;in;i+) if(x1ix3i) ai=x1i; bi=x3i; else ai=x3i; bi=x1i; void hongjinfa(int n,double a,double b,double flag,double x,double r0) int i; double x1m,x2m,f1,f2; while(1) for(i=0;in;i+) x1i=bi; f1=f(x1,r0); for(i=0;in;i+) x2i=ai; f2=f(x2,r0); if(fabs(f2-f1)/f2)=flag) for(i=0;in;i+) xi=(double)(bi+ai)/2.0; break; else for(i=0;in;i+) x1i=bi-(double)0.618*(bi-ai); for(i=0;if(x2,r0) for(i=0;in;i+) ai=x1i; else for(i=0;in;i+) bi=x2i; void baoweifa(int n,double x0,double h,double flag,double a,double b,double x,double r0) int i,j,k,r,sum,p,e; double x1m,x2m,abcm,krm,f0,f1,f2,fnm,cm,q,wo; k=0; double abmm; for(i = 0; i n; i +) for(j = 0; j n; j +) if(i=j) abij = 1; else abij =0; while(1) wo=f(x0,r0); for(p=0;pn;p+) q=f(x0,r0); jintuifa(ab,n,x0,h,p,a,b,r0); hongjinfa(n,a,b,flag,x1,r0); for(e=0;en;e+) x0e=x1e; fnp=(double)(q-f(x1,r0); k=k+1; for(i=0;in;i+) x2i=(double)2*x1i-x0i; for(i=0;in;i+) ci=(double)x1i-x0i; r=0; for(i=0;in;i+) if(fn0=f0|(double)f0-(double)2*f1+(double)f2)*(double)f0-(double)f1-(double)fn0)*(double)f0-(double)f1-(double)fn0)=0.5*fn0*(double)(f0-f2)*(double)(f0-f2) if(f1f2) for(i=0;in;i+) x0i=x1i; else for(i=0;in;i+) x0i=x2i; else int ij; ij=n-1; for(i=0;in;i+) kri=x1i; for(i = 0; i n; i +) for(j = 0; j =r) abij = abij+1; abin-1=ci; jintuifa(ab,n,kr,h,ij,a,b,r0); hongjinfa(n,a,b,flag,x1,r0); for(i=0;in;i+) x0i=x1i; if(fabs(f(x0,r0)-wo)/f(x0,r0)=flag) break; for(i=0;in;i+) xi=x0i;int main() int v,i,n; double h,flag,x0m,am,bm,xm; double fom,fxo,c,r0; c=0.5; r0=100.0; fom=100; printf(请输入维数:n); scanf(%d,&n); printf(请输入初始点:); for(i=0;iflag) fom=fxo; r0=c*r0; for(v=0;vn;v+) x0v=xv; else p