机械优化设计上机报告
班 级
姓 名 学 号成 绩 指导教师
日 期
机械〔茅以升〕101
: 1004010510
: 张 迎 辉 : 2021.11.20
:
:
:
1 «一维搜索方法»上机实践报告
1、写出所选择的一维搜索算法的差不多过程、原理〔可附流程图说明〕。 〔一〕进退法 1. 算法原理
进退法是用来确定搜索区间〔包含极小值点的区间〕的算法,其理论依据是:f(x)为单谷函数〔只有一个极值点〕,且[a,b]为其极小值点的一个搜索区间,关于任意x1,x2[a,b],假如fx1fx2,那么[a,x2]为极小值的搜索区间,假如fx1fx2,那么[x1,b]为极小值的搜索区间。
因此,在给定初始点x0,及初始搜索步长h的情形下,第一以初始步长向前搜索一步,运算fx0h。
(1) 假如fx0fx0h
那么可知搜索区间为[x,x0h],其中x待求,为确定x,后退一步运算f(x0h),为缩小系数,且01,直截了当找到合适的*,使得f(x0*h)fx0,从而确定搜索区间[x0*h,x0h]。
(2) 假如fx0fx0h
那么可知搜索区间为[x0,x],其中x待求,为确定x,前进一步运算f(x0h),为放大系数,且1,明白找到合适的*,使得fx0hf(x0*h),从而确定搜索区间
[x0,x0*h]。
2. 算法步骤
用进退法求一维无约束问题minf(x),xR的搜索区间〔包含极小值点的区间〕的差不多算法步骤如下:
(1) 给定初始点x(0),初始步长h0,令hh0,x(1)x(0),k0; (2) 令x(4)x(1)h,置kk1;
(3) 假设fx(4)fx(1),那么转步骤〔4〕,否那么转步骤〔5〕;
(4) 令x(2)x(1),x(1)x(4),fx(2)fx(1),fx(1)fx(4),令h2h,转步骤〔2〕; (5) 假设k1,那么转步骤〔6〕否那么转步骤〔7〕; (6) 令hh,x(2)x(4),fx(2)fx(4),转步骤〔2〕;
(7) 令x(3)x(2),x(2)x(1),x(1)x(4),停止运算,极小值点包含于区间
[x(1),x(3)]或[x(3),x(1)]
〔二〕黄金分割法
1、黄金分割法差不多思路:
黄金分割法适用于[a,b]区间上的任何单股函数求极小值问题,对函数除要求〝单谷〞外不做其他要求,甚至能够不连续。因此,这种方法的适应面专门广。黄金分割法也是建立在区间消去法原理基础上的试探方法,即在搜索区间[a,b]内适当插入两点a1,a2,并运算其函数值。a1,a2将区间分成三段,应用函数的单谷性质,通过函数值大小的比较,删去其中一段,是搜索区间得以缩小。然后再在保留下来的区间上作同样的处理,如此迭代下去,是搜索区间无限缩小,从而得到极小点的数值近似解。
2 黄金分割法的差不多原理
一维搜索是解函数极小值的方法之一,其解法思想为沿某一方向求目标函数的极小值点。一维搜索的解法专门多,那个地点要紧采纳黄金分割法〔0.618法〕。该方法用不变的区间缩短率0.618代替斐波那契法每次不同的缩短率,从而能够看成是斐波那契法的近似,实现起来比较容易,也易于人们所同意。
图1
黄金分割法是用于一元函数f(x)在给定初始区间[a,b]内搜索极小点α*的一种方法。它是优化运算中的经典算法,以算法简单、收敛速度平均、成效较好而著称,是许多优化算法的基础,但它只适用于一维区间上的凸函数[6],即只在单峰区间内才能进行一维寻优,其收敛效率较低。其差不多原理是:依照〝去劣存优〞原那么、对称原那么、以及等比收缩原那么来逐步缩小搜索区间[7]。具体步骤是:在区间[a,b]内取点:a1 ,a2 把[a,b]分为三段。假如f(a1)>f(a2),令a=a1,a1=a2,a2=a+r*(b-a);假如f(a1) 3 程序流程如下: 4 实验所编程序框图 给定a=-3,b=5,收敛精度ε=0.001 开始 |(b-a)/b|<ε和 |〔y2-y1〕/y2|<ε? a2=a+r*(b-a) y2=f(a2) a1=b-r*(b-a) y1=f(a1) a=a1 a1=a2 y1=y2 b=a2 a2=a1 y2=y1 是 否 y1>=y2 r=0.618 a1=b-r*(b-a) y1=f(a1) a2=a+r*(b-a) y2=f(a2) 算例1:minf(x)= x*x+2*x (1)C++程序如下: #include double calc(double *a,double *b,double e,int *n) { double x1,x2,s; if(fabs(*b-*a)<=e) s=f((*b+*a)/2); else { x1=*b-0.618*(*b-*a); x2=*a+0.618*(*b-*a); if(f(x1)>f(x2)) *a=x1; else *b=x2; *n=*n+1; s=calc(a,b,e,n); } return s; } main() { double s,a,b,e; int n=0; scanf(\"%lf %lf %lf\s=calc(&a,&b,e,&n); printf(\"a=%lf,b=%lf,s=%lf,n=%d\\n\} 2、程序运行结果: 是 否 算例2:minf=x^2-10*x+36 理论最优解:x*=5.0,f〔x*〕=11.0 〔1〕MATLAB程序清单: function f=myfun_yi(x) f=x^2-10*x+36 >> fminbnd(@myfun_yi,1,12) 〔2〕运行结果: >> fminbnd(@myfun_yi,1,12) f = 11.0407 f = 18.8309 f = 12.9691 f = 11 f = 11.0000 f = 11.0000 ans = 5 〔3〕结果分析:由迭代程序f=11.0,ans=5,与理论结果相等 算例3:minf=x^4-5*x^3+4*x^2-6*x+60 理论最优解:x*=3.2796,f〔x*〕=22.6590 〔1〕MATLAB程序清单: function f=myfun_yi(x) f=x^4-5*x^3+4*x^2-6*x+60 >> fminbnd(@myfun_yi,1,12) 〔2〕运行结果: >> fminbnd(@myfun_yi,1,12) f = 165.3948 f = 1.5836e+03 f = 24.8730 f = 35.9194 f = 23.9089 f = 22.7621 f = 31.7507 f = 22.6673 f = 22.6594 f = 22.6590 f = 22.6590 f = 22.6590 f = 22.6590 ans = 3.2796 〔3〕结果分析:由迭代程序得f =22.659,ans =3.2796,与理论最优解相等 2 «无约束优化搜索方法»上机实践报告 1、写出所选择的无约束优化搜索算法的差不多过程、原理〔可附流程图说明〕。 鲍威尔改进方法 鲍威尔〔Powell〕法是直截了当利用函数值来构造共轭方向的一种方法 在鲍威尔差不多算法中,每一轮迭代都用连结始点和终点所产生出的搜索方向去替换原向量组中的第一个向量,而不管它的〝好坏〞,这是产生向量组线性相关的缘故所在。 在改进的算法中第一判定原向量组是否需要替换。假如需要替换,还要进一步判定原向量组中哪个向量最坏,然后再用新产生的向量替换那个最坏的向量,以保证逐次生成共轭方向。 2、程序运算结果分析:中间各步骤的结果分析及与理论运算结果分析对比。 算例1:minf=4*(x(1)-5)^2+(x(2)-6)^2 初始点:x0=[8;9],f(x0)=45 最优解:x*=[5;6],f(x*)=0 (1)MATLAB程序清单: function f=myfun_wuyueshu(x) f=4*(x(1)-5)^2+(x(2)-6)^2 >> [x,fval]=fminunc(@myfun_wuyueshu,x0) (2)运行结果: f = 45 Warning: Gradient must be provided for trust-region algorithm; using line-search algorithm instead. > In fminunc at 367 f = 45.0000 f = 45.0000 f = 23.5625 f = 23.5625 f = 23.5625 f = 2.6958 f = 2.6958 f = 2.6958 f = 1.3788 f = 1.3788 f = 1.3788 f = 0.0054 f = 0.0054 f = 0.0054 f = 6.4975e-05 f = 6.4973e-05 f = 6.4975e-05 f = 6.1579e-09 f = 6.1522e-09 f = 6.1443e-09 f = 1.7876e-12 f = 1.8627e-12 f = 1.5586e-12 Local minimum found. Optimization completed because the size of the gradient is less than the default value of the function tolerance. 5.0000 6.0000 fval = 1.7876e-12 (3)结果分析:由迭代程序得x =[ 5.0000; 6.0000],fval =1.7876e-12,与理论最优解相等。 算例2:minf=(x(1)^2+x(2)-11)^2+(x(1)+x(2)^2-7)^2 初始点:x0=[1;1],f(x0)=106 最优解:x*=[3;2],f(x*)=0 (1)MATLAB程序清单: function f=myfun_wuyueshu(x) f=(x(1)^2+x(2)-11)^2+(x(1)+x(2)^2-7)^2 >> [x,fval]=fminunc(@myfun_wuyueshu,x0) (2)运行结果: >> x0=[1;1] x0 = 1 1 >> [x,fval]=fminunc(@myfun_wuyueshu,x0) f = 106 Warning: Gradient must be provided for trust-region algorithm; using line-search algorithm instead. > In fminunc at 367 f = 106.0000 f = 106.0000 f = 29.5430 f = 29.5430 f = 29.5430 f = 1.7450e+04 f = 1.7450e+04 f = 1.7450e+04 f = 90.3661 f = 90.3661 f = 90.3661 f = 0.3575 f = 0.3575 f = 0.3575 f = 0.0179 f = 0.0179 f = 0.0179 f = 0.0064 f = 0.0064 f = 0.0064 f = 1.0048e-06 f = 1.0044e-06 f = 1.0049e-06 f = 4.8639e-09 f = 4.8567e-09 f = 4.8781e-09 f = 5.2125e-12 f = 5.8703e-12 f = 5.7870e-12 Local minimum found. Optimization completed because the size of the gradient is less than the default value of the function tolerance. 3.0000 2.0000 fval = 5.2125e-12 (3)结果分析:由迭代程序得x=[3;2],fval = 5.2125e-12,与理论最优解相等 算例3:ff=x[0]*x[0]+2*x[1]*x[1]-4*x[0]-2*x[0]*x[1]; 〔1〕鲍威尔改进算法C++程序清单: #include \"stdio.h\" #include \"stdlib.h\" #include \"math.h\" double objf(double x[]) {double ff; ff=x[0]*x[0]+2*x[1]*x[1]-4*x[0]-2*x[0]*x[1]; return(ff); } void jtf(double x0[ ],double h0,double s[ ],int n,double a[ ],double b[ ]) {int i; 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 for (i=0;i f1=f2; f2=f3; } for(;;) {h=2. *h; for(i=0;i f1=f2; f2=f3; } } if(h<0. ) for(i=0;i double gold(double a[],double b[],double eps,int n,double xx[]) { int i; double f1,f2,*x[2],ff,q,w; for(i=0;i<2;i++) x[i]=(double*)malloc (n*sizeof(double)); for(i=0;i f1=objf(x[0]); f2=objf(x[1]); do {if(f1>f2) {for(i=0;i f1=f2; for(i=0;i {for(i=0;i f2=f1; for(i=0;i for(i=0;i }while(w>eps); for(i=0;i for(i=0;i<2;i++) free(x[i]); return(ff); } double oneoptim(double x0[],double s[],double h0,double epsg,int n,double x[]) {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); } double powell(double p[],double h0,double eps,double epsg,int n,double x[]) {int i,j,m; double *xx[4],*ss,*s; double f,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 xx[i]=(double *)malloc(n*sizeof(double)); for (i=0;i } f0=f1=objf(x); dlt=-1; for (j=0;j 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 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 for (i=0;i {for (i=0;i f=oneoptim(xx[0],s,h0,epsg,n,x); for(i=0;i void main() {double p[]={1,1}; double ff,x[2],x1,x2,f; ff=powell(p,0.3,0.001,0.0001,2,x); printf(\"shuchuzuiyoujie:\\n\"); x1=x[1];x2=x[2];f=ff; printf(\"x1=%f,x2=%f,f=%f\\n\getchar(); } 〔2〕运行结果为: 3«约束优化搜索方法»上机实践报告 1、写出所选择的约束优化搜索算法的差不多过程、原理〔可附流程图说明〕。 2、程序运算结果分析:中间各步骤的结果分析及与理论运算结果分析对比。 算例1:minf=〔x(1)-2)^2+(x(2)-1)^2; s.t g1(x)=x(1)^2-x(2)<=0 g2(x)=x(1)+x(2)-2 <=0 初始点:x0=[3;3],f(x0)=5 最优解:x*=[1;1] f(x*)=1 (1)MATLAB程序清单: function f=myfun_constrain(x) f=〔x(1)-2)^2+(x(2)-1)^2; function [c,ceq]=mycon(x) c=[x(1)^2-x(2);x(1)+x(2)-2 ] ceq=[] >> [x,fval]=fmincon(@myfun_constrain,x0,A,b,[],[],[],[],@mycon) (2)运行结果: >> A=[-1,0;0,-1] b=[0;0] x0=[3;3] A = -1 0 0 -1 b = 0 0 x0 = 3 3 >> [x,fval]=fmincon(@myfun_constrain,x0,A,b,[],[],[],[],@mycon) Warning: The default trust-region-reflective algorithm does not solve problems with the constraints you have specified. FMINCON will use the active-set algorithm instead. For information on applicable algorithms, see Choosing the Algorithm in the documentation. > In fmincon at 486 f = 5 c = 6 4 ceq = [] f = 5.0000 c = 6.0000 4.0000 ceq = [] f = 5.0000 c = 6.0000 4.0000 ceq = [] f = 2.0000 c = 1.0000 -1.0000 ceq = [] f = 2.0000 c = 1.0000 -1.0000 ceq = [] f = 2.0000 c = 1.0000 -1.0000 ceq = [] f = 1.0000 c = 1.0e-15 * 0.9992 0.4441 ceq = [] f = 1.0000 c = 1.0e-07 * 0.2980 0.1490 ceq = [] f = 1.0000 c = 1.0e-07 * -0.1490 0.1490 ceq = [] Local minimum found that satisfies the constraints. Optimization completed because the objective function is non-decreasing in feasible directions, to within the default value of the function tolerance, and constraints are satisfied to within the default value of the constraint tolerance. Active inequalities (to within options.TolCon = 1e-06): lower upper ineqlin ineqnonlin 1 2 x = 1.0000 1.0000 fval = 1.0000 (3)结果分析:由迭代程序得x =[1.0000; 1.0000],fval =1.0000,与理论最优解相等 算例2.minf=1000-x(1)^2-2*x(2)^2-x(3)^2-x(1)*x(2)-x(1)*x(3); S.t g1(x)=x(1)^2+x(2)^2+x(3)^2-25<=0 g2(x)8*x(1)+14*x(2)+7*x(3)-56<=0 g3(x)=-x(1)<=0 g4(x)=-x(2)<=0 g5(x)=-x(3)<=0 初始点:x0=[2;2;2],f(x0)=976 最优解:x*=[3.512;0.217;3.552],f(x*)=961.715 〔1〕MATLAB程序清单: function f=myfun_constrain(x) f=1000-x(1)^2-2*x(2)^2-x(3)^2-x(1)*x(2)-x(1)*x(3); function [c,ceq]=mycon(x) c=[x(1)^2+x(2)^2+x(3)^2-25;8*x(1)+14*x(2)+7*x(3)-56] ceq=[] >> [x,fval]=fmincon(@myfun_constrain,x0,A,b,[],[],[],[],@mycon) 〔2〕运行结果 >> A=[-1,0,0;0,-1,0;0,0,-1] b=[0;0;0] x0=[2;2;2] A = -1 0 0 0 -1 0 0 0 -1 b = 0 0 0 x0 = 2 2 2 >> [x,fval]=fmincon(@myfun_constrain,x0,A,b,[],[],[],[],@mycon) Warning: The default trust-region-reflective algorithm does not solve problems with the constraints you have specified. FMINCON will use the active-set algorithm instead. For information on applicable algorithms, see Choosing the Algorithm in the documentation. > In fmincon at 486 c = -13 2 ceq = [] c = -13.0000 2.0000 ceq = [] c = -13.0000 2.0000 ceq = [] c = -13.0000 2.0000 ceq = [] c = -5.9320 0 ceq = [] c = -5.9320 0.0000 ceq = [] c = -5.9320 0.0000 ceq = [] c = -5.9320 0.0000 ceq = [] c = 1.1562 0.0000 ceq = [] c = 1.1562 0.0000 ceq = [] c = 1.1562 0.0000 ceq = [] c = 1.1562 0.0000 ceq = [] c = 31.5713 0.0000 ceq = [] c = 8.1851 0.0000 ceq = [] c = 8.1851 0.0000 ceq = [] c = 8.1851 0.0000 ceq = [] c = 8.1851 0.0000 ceq = [] c = 3.2553 0.0000 ceq = [] c = 3.2553 0.0000 ceq = [] c = 3.2553 0.0000 ceq = [] c = 3.2553 0.0000 ceq = [] c = 1.0789 -0.0000 ceq = [] c = 1.0789 0.0000 ceq = [] c = 1.0789 0.0000 ceq = [] c = 1.0789 0.0000 ceq = [] c = 0.4793 0.0000 ceq = [] c = 0.4793 0.0000 ceq = [] c = 0.4793 0.0000 ceq = [] c = 0.4793 0.0000 ceq = [] c = 0.0039 -0.0000 ceq = [] c = 0.0039 0.0000 ceq = [] c = 0.0039 0.0000 ceq = [] c = 0.0039 0.0000 ceq = [] c = 1.0e-06 * 0.4884 0.0000 ceq = [] c = 1.0e-06 * 0.8559 0.4187 ceq = [] c = 1.0e-06 * 0.4948 0.2086 ceq = [] c = 1.0e-06 * 0.8644 0.3705 ceq = [] Local minimum possible. Constraints satisfied. fmincon stopped because the predicted change in the objective function is less than the default value of the function tolerance and constraints are satisfied to within the default value of the constraint tolerance. Active inequalities (to within options.TolCon = 1e-06): lower upper ineqlin ineqnonlin 1 2 x = 3.5120 0.2170 3.5523 fval = 961.7152 >> 〔3〕结果分析:由迭代程序得x=[ 3.5120; 0.2170; 3.5523],fval =961.7152,与理论结果相等。 4«平面连杆机构中再现运动规律的优化设计»上机实践报告 1、写出所针对此实际的优化问题,所确定的优化目标函数和约束条件。 例8.5 设计一曲柄摇杆机构〔如图9所示〕,要求曲柄l1从0转到090时,摇杆l3的转角最正确再现的运动规律: E0承诺在45135范畴内变化。 202, 且l11,l45,0为极位角,其传动角3 2、MATLAB程序清单: function f=objfun(x) n=30; L1=1; L2=5; fx=0; fa0=acos(((L1+x(1))^2-x(2)^2+L2^2)/(2*(L1+x(1))*L2)); %Çú±ú³õʼ½Ç pu0=acos(((L1+x(1))^2-x(2)^2-L2^2)/(2*(x(2))*L2)); %Ò¡¸Ë³õʼ½Ç for i=1:n fai=fa0+0.5*pi*i/n; pui=pu0+2*(fai-fa0)^2/(3*pi); ri=sqrt(L1^2+L2^2-2*L1*L2*cos(fai)); alfai=acos((ri^2+x(2)^2-x(1)^2)/(2*ri*x(2))); betai=acos((ri^2+L2^2-L1^2)/(2*ri*L2)); if fai>0&&fai<=pi psi=pi-alfai-betai; elseif fai>pi&&fai<=2*pi psi=pi-alfai+betai; end fx=fx+(pui-psi)^2; end f=fx; function [c,ceq]=nonlinconstr(x) L1=1;L2=5; m=45*pi/180; n=135*pi/180; c(1)=x(1)^2+x(2)^2-(L1-L2)^2-2*x(1)*x(2)*cos(m);%×îС´«¶¯½ÇÔ¼Êø c(2)=-x(1)^2-x(2)^2+(L1+L2)^2+2*x(1)*x(2)*cos(m);%×î´ó´«¶¯½ÇÔ¼Êø ceq=[]; options=optimset('LargeScale','on','TolFun',1e-6); >> [x,fval]=fmincon(@objfun,x0,a,b,[],[],lb,ub,@nonlinconstr,options) 〔2〕运行结果: x0=[4.0;4.0]; lb=[1;1]; ub=[10;10]; a=[-1 -1 ;1 -1 ;-1 1 ]; b=[-6;4;4]; options=optimset('LargeScale','on','TolFun',1e-6); >> [x,fval]=fmincon(@objfun,x0,a,b,[],[],lb,ub,@nonlinconstr,options) Warning: The default trust-region-reflective algorithm does not solve problems with the constraints you have specified. FMINCON will use the active-set algorithm instead. For information on applicable algorithms, see Choosing the Algorithm in the documentation. > In fmincon at 486 No feasible solution found. fmincon stopped because the size of the current search direction is less than twice the default value of the step size tolerance but constraints are not satisfied to within the default value of the constraint tolerance. 6.6622 6.6622 fval = 0.5720 因篇幅问题不能全部显示,请点此查看更多更全内容