中科大优化设计课程大作业之无约束优化试验报告

1、无约束优化设计实验报告力学系型号：联想y470CPU: i5-2450M 内存：2GB 系统:win7-64位如下是三个目标函数(包括自定义函数)以及初值和精度选取:1.minf(x)=x1A2+2*x2A2-2*x1*x2-4*x1初值x0=1,1'精度为：0.12.minf(x,y)=xA4-2*x*x*y-2*x*y+xA2+2*y*y+4.5*x-4*y+4初值为(-2.5,4.25);精度为0.013.minf(x)=x1A2+x2A2+x3A2初值为(3,2,1);精度为0.01如下是运算结果:目标函数无约束方法一维搜索所需时间迭代次数极值点极值1最速卜降法黄金分割法5.8

2、18(3.9,1.9)-8.0牛顿法2.618(3.9,1.9)-8.0不精确法2.968(3.9,1.9)-8.0阻尼牛顿法黄金分割法1.661(4.0,2.0)-8.0牛顿法1.021(4.0,2.0)-8.0不精确法0.761(4.0,2.0)-8.0共腕梯度法黄金分割法14.6524(3.9,2.0)-8.0牛顿法7.9124(4.0,2.0)-8.0不精确法11.2827(3.9,2.0)-8.0鲍维尔法黄金分割法9.532(4.0,2.0)-8.0牛顿法2.862(4.0,2.0)-8.0不精确法-变尺度法黄金分割法8.04 n2(4.0,2.0)-8.0牛顿法1.242(4.0,

3、2.0)-8.0不精确法1.353(4.0,2.0)-8.0单形替换法无0.029(4.1,2.2)-8.0目标函数无约束方法一维搜索所需时间迭代次数极值点极值2最速卜降法黄金分割法20.76P 21(1.90,3.73)0.99牛顿法3.236(1.93,3.82)0.99不精确法10.4927(1.94,3.84)0.99阻尼牛顿法黄金分割法3.503(-1.05,1.03)-0.51牛顿法1.61P 3(-1.05,1.03)-0.51不精确法1.584(-1.05,1.03)-0.51共腕梯度法黄金分割法25.0630(-1.05,1.03)-0.51牛顿法24.4035(-1.05,

4、1.03)-0.51不精确法19.8033(1.94,3.86)0.99鲍维尔法黄金分割法11.493(-1.05,1.03)-0.51牛顿法6.173(-1.05,1.03)-0.51不精确法4.751(-1.50,4.25)17.31变尺度法黄金分割法3.943(1.94,3.85)0.99牛顿法2.393(1.94,3.85)0.99不精确法2.466(1.94,3.84)0.99单形替换法无0.0116(1.92,3.81)0.99目标函数无约束方法一维搜索所需时间迭代次数极值点极值3最速卜降法r黄金分割法2.702r (0,0,0)0厂牛顿法0.761(0,0,0)0不精确法0.78

5、1(0,0,0)0阻尼牛顿法黄金分割法1.611(0,0,0)0r牛顿法0.661（0。0）0不精确法0.651（0。0）0共腕梯度法黄金分割法3.153（0。0）0牛顿法0.400（0。0）0不精确法0.440（0。0）0鲍维尔法黄金分割法7.171（0。0）0牛顿法1.991（0。0）0r不精确法4.080F (3.00,2.00,1.00)14.00变尺度法黄金分割法2.692(0,0,0)0牛顿法0.671(0,0,0)0不精确法0.771(0,0,0)0单形替换法无0.0124(-0.01,0.03,-0.08)0.01总结及比较:根据上面三个函数的表格可以看出：首先，从迭代时间来看

6、，三种一维搜索方法中黄金分割法所 用时间最久，牛顿法和不精确法所用时间较少， 这两种方法相比较而言牛顿法所用时间更少一些。 而 六种无约束方法中，由于单形替代法不需要使用一维搜索方法，故迭代时间最少，紧接着在使用一维搜索的五种方法中以阻尼牛顿法迭代时间相对较少，共腕梯度法迭代时间最久；然后，从迭代次数来看，共腕梯度法往往需要较多的迭代次数，从而所需时间也最久；接着，从计算结果的精度来看，阻 尼牛顿法的结果精度最高，而单形替换法的精度最低；最后，从编程来看，在编好一维搜索方法的情 况下，最速下降法和阻尼牛顿法编程简单容易，而共腕梯度法、变尺度法和单形替代法需要两重循环实现，鲍威尔法和单形替换法则

7、需要编程者对矩阵的操作能力有较高的要求，故编程较难。同时，从上面的结果也可以发现，鲍威尔法在使用不精确的一维搜索方法时， 对函数1无法收敛， 对函数2、3收敛到错误的结果，所以鲍威尔法是依赖于精确的一维搜索过程的，而其他几个则不依 赖于精确一维搜索过程。精确的一维搜索方法通常需要花费很大的工作量，特别是当迭代点远离问题 的解时，精确的求解一个一维子问题通常不是十分有效率的。因此，只要保证目标函数值在每一步都有满意的下降，使用不是非常精确的一维搜索，就可以大大节省工作量。在分析函数2的计算结果时，可以发现存在两个收敛结果，当然这两个结果都是极值，因为函数 2是二元四次函数，存在多个极值。故为了验

8、证正确性，自己曾将初始点(-2.5,4.25)调成(2.5,4.25), 分别代入程序中计算，计算结果都收敛于极值为 0.99的这个点上。所以，在存在多个极值点的情况 下结果是和初始点的选取有关。对于单形替换法，这种方法不需要一维搜索最佳步长，故没有一维搜索方法反复地计算最佳步长 的计算时间，程序运行效率快。但它的收敛条件不好选择，通过查找文献资料总结出以下三个收敛条 件：1.利用最坏点函数值与最好点函数值之差判别；2.利用相邻两次函数值差值的绝对值判别； 3.利用各点的函数值与最好点函数值之差的均方根判别。为了保证程序的执行可靠性，这三种常用的方法 中自己选择了判别3,即：综上所述，阻尼牛顿

9、法是无约束方法中最有效的方法。不仅编程简单，而且迭代次数较少，运行 时间较短，结果的精度也较高。在程序的运行方面，分别设置了可变的函数选择、无约束方法选择、一维搜索方法选择、起始点、 精度这五个输入，故可以在命令窗口运行主程序 main,再根据提示要求分别输入这五个参数的所需 值，就可以得到运行结果。程序如下：1、主函数clear;global k;k=0;disp('1.f(x)=x1A2+2*x2A2-2*x1*x2-4*x1');disp('2.f(x,y)=xA4-2*x*x*y-2*x*y+xA2+2*y*y+4.5*x-4*y+4');disp(&#

10、39;3.f(x)=x1A2+x2A2+x3A2');while 1n0=input('请输入上面所想选择函数的编号(1、2、3):'); if n0=1|n0=2|n0=3break;enddisp（'此次输入无效.,）;enddisp（''）;disp（'1.最速下降法'）;disp（'2.阻尼牛顿法'）;disp（'3.共腕梯度法'）;disp（'4.鲍威尔法'）;disp（'5.变尺度法'）;disp（'6.单纯形法'）;while 1m0=i

11、nput（'请输入上面所想选择无约束方法的编号（1、2、3、4、5、6）:'）;if m0=1|m0=2|m0=3|m0=4|m0=5|m0=6break;enddisp（'此次输入无效.'）;enddisp（' '）;disp（'1.黄金分割法'）;disp（'2.牛顿法'）;disp（'3.不精确一维搜索法'）;while 1m1=input（'请输入上面所想选择一维搜索方法的编号（1、2、3）:'）;if m1=1|m1=2|m1=3break;enddisp（'此次输

12、入无效.'）;enddisp（' '）;s=input（'请输入用空格隔开的初始值坐标向量（如：1.1 2.0） ： ','s'）;x=str2num（s）;x=x'disp（' '）;while 1e=input（'请输入精度（建议0.1或0.01）：'）;if e>0break;enddisp（'此次输入无效.'）;enddisp（' '）;disp（''）;xx,yy=fmins（m0,m1,n0,x,e）;fprintf（'迭代次

13、数为：8.0fn', k）;disp（'所求极值点的坐标向量为：'）;fprintf（'%16.5fn', xx）;fprintf（'所求函数的极值为：%16.5fn', yy）;2、外部多维的调用函数function xx,yy=fmins(m0,m1,n0,x,e) %UNTITLED此处显示有关此函数的摘要 %此处显示详细说明if m0=1tic;xx,yy=zuisu(m1,n0,x,e);toc;elseif m0=2tic;xx,yy=zuni(m1,n0,x,e);toc;elseif m0=3tic;xx,yy=gonge

14、(m1,n0,x,e);toc;elseif m0=4tic;xx,yy=powell(m1,n0,x,e);toc;elseif m0=5tic;xx,yy=bianchi(m1,n0,x,e);toc;elseif m0=6tic;xx,yy=danxing(n0,x,e);toc; endend3、最速法function xx,yy=zuisu(m1,n0,x,e)%UNTITLED2此处显示有关此函数的摘要%此处显示详细说明global k;g,ss=gra(n0);while 1d=-double(subs(g,ss,x);a=fmin(m1,n0,x,d,e);x1=x+a*d;i

15、f norm(x1-x)<ebreak;endx=x1;k=k+1;endxx=x1;yy=f0(n0,xx); end4、阻尼法function xx,yy=zuni(m1,n0,x,e)%UNTITLED此处显示有关此函数的摘要%此处显示详细说明global k;g,ss=gra(n0);h=jacobian(g',ss);while 1d=-double(subs(h,ss,x)A(-1)*subs(g,ss,x);a=fmin(m1,n0,x,d,e);x1=x+a*d;if norm(x1-x)<ebreak;endx=x1;k=k+1;endxx=x1;yy=f

16、0(n0,xx);end5、共腕梯度法function xx,yy=gonge(m1,n0,x,e)%UNTITLED2此处显示有关此函数的摘要%此处显示详细说明global k;if n0=1n=2;elseif n0=2n=4;elseif n0=3n=2;endg,ss=gra(n0);while 1kk=0;d=-double(subs(g,ss,x);while 1a=fmin(m1,n0,x,d,e);x1=x+a*d;gx=double(subs(g,ss,x);gx1=double(subs(g,ss,x1);if norm(gx1)<ebreak;elseifkk=nb

17、reak;endbeta=norm(gx1)A2/norm(gx)A2;d=-gx+beta*d;x=x1;k=k+1;kk=kk+1;endif norm(gx1)<ebreak;endx=x1;k=k+1;endxx=x1;yy=f0(n0,xx);end6、鲍威尔法function xx,yy=powell(m1,n0,x,e)%UNTITLED3此处显示有关此函数的摘要%此处显示详细说明global k;if n0=1n=2;elseif n0=2n=2;elseif n0=3n=3;endd=zeros(n,n+1);xk=zeros(n,n+1);deta=zeros(n,1

18、);for j=1:nd(j,j)=1;endwhile 1xt=x;for i=1:na=fmin(m1,n0,xt,d(:,i),e);xk(:,i)=xt+a*d(:,i);deta(i)=f0(n0,xt)-f0(n0,xk(:,i);xt=xk(:,i);endxt=x;xk(:,n+1)=2*xk(:,n)-x;ff0=f0(n0,x);ff2=f0(n0,xk(:,n);ff3=f0(n0,xk(:,n+1);md=max(deta);m=find(deta=md);if (ff3<ff0)&&(ff0-2*ff2+ff3)*(ff0-ff2-md)A2&l

19、t;0.5*md*(ff0-ff3)A2)d(:,n+1)=xk(:,n)-xt;a=fmin(m1,n0,xk(:,n),d(:,n+1),e); x=xk(:,n)+a*d(:,n+1);for i=m:nd(:,i)=d(:,i+1);endelseif ff2<ff3x=xk(:,n);elsex=xk(:,n+1);endendif norm(x-xt)<ebreak;endk=k+1;endxx=x;yy=f0(n0,xx);end7、变尺度法function xx,yy=bianchi(m1,n0,x,e)%UNTITLED此处显示有关此函数的摘要%此处显示详细说明g

20、lobal k;if n0=1n=2;elseif n0=2n=2;elseif n0=3n=3;endg,ss=gra(n0);while 1gx=double(subs(g,ss,x);h=eye(n);kk=0;while 1d=-h*gx;a=fmin(m1,n0,x,d,e);xk=x+a*d;if norm(xk-x)<ebreak;endifkk=nbreak;endgxk=double(subs(g,ss,xk);yk=gxk-gx;sk=xk-x;h=h+(sk*sk')/(sk'*yk)-(h*yk)*yk'*h)/(yk'*h*yk)

21、;x=xk;gx=gxk;kk=kk+1;k=k+1;endif norm(xk-x)<ebreak;endx=xk;k=k+1;endxx=xk;yy=f0(n0,xx);end8、单形替换法function xx,yy=danxing(n0,x,e)%UNTITLED此处显示有关此函数的摘要%此处显示详细说明global k;if n0=1n=2;elseif n0=2n=2;elseif n0=3n=3;endf=zeros(n+5,1);xk=zeros(n,n+5);h=2*eye(n);xk(:,1)=x;for i=1:nxk(:,i+1)=x+h(:,i);endwhil

22、e 1for i=1:n+1f(i)=f0(n0,xk(:,i);endwhile 1f(n+2)=nan;f(n+3)=nan;f(n+4)=nan;f(n+5)=nan;fl=min(f);xll=find(f=fl);xl=xll(1);fh=max(f);xhh=find(f=fh);xh二xhh ； fff=f;fff(xh 尸口； fg=max(fff);fz=O;for i=1 :n+1fz=fz+(f(i)-f(xl)A2;endfz=sqrt(fz/n);iffz<ebreak;endxkk=xk(:,1);for i=1:nxkk=xkk+xk(:,i+1);endx

23、k(:,n+2)=(xkk-xk(:,xh)/n;xk(:,n+3)=2*xk(:,n+2)-xk(:,xh);f(n+3 尸 f0(n0,xk(:,n+3);if f(n+3)<flxk(:,n+4)=xk(:,n+2)+2*(xk(:,n+3)-xk(:,n+2);f(n+4 尸 f0(n0,xk(:,n+4);if f(n+4)<f(n+3)xk(:,xh)=xk(:,n+4);f(xh)=f(n+4);elsexk(:,xh)=xk(:,n+3);f(xh)=f(n+3);end else if f(n+3)<fg xk(:,xh)=xk(:,n+3); f(xh)=

24、f(n+3);elseif f(n+3)>=fhxk(:,n+3)=xk(:,xh);endxk(:,n+5)=xk(:,n+2)+0.5*(xk(:,n+3)-xk(:,n+2); f(n+5)=f0(n0,xk(:,n+5);if f(n+5)<fhxk(:,xh)=xk(:,n+5);f(xh)=f(n+5);elsefor i=1:n+1xk(:,i)=(xk(:,i)+xk(:,xl)/2;endbreak;endendendk=k+1;end iffz<e break;endk=k+1;endxx=xk(:,xl);yy=f0(n0,xx);end9、内部循环一维

25、调用函数function xx=fmin(m1,n0,x,d,e)x0=0;%初始步长默认为0amin,amax=range1(n0,x,d,x0);if m1=1xx=gold(n0,x,d,amin,amax,e);elseif m1=2xx=newton(n0,x,d,amax,e);elseif m1=3xx=wolfe(n0,x,d);end end10、一维搜索确定区间函数function amin,amax = range1(n0,x,d,x0) %UNTITLED5此处显示有关此函数的摘要 %此处显示详细说明h=1;a1=x0;y1=f_1(n0,x,d,a1);a2=a1+h

26、;y2=f_1(n0,x,d,a2);if y2>y1h=-h; a3=a1;y3=y1;a1=a2;a2=a3;y2=y3;end a3=a2+h;y3=f_1(n0,x,d,a3);while y3<y2h=h*2;a1=a2;a2=a3;y2=y3;a3=a2+h;y3=f_1(n0,x,d,a3);endamin=min(a1,a3);amax=max(a1,a3);end11、黄金一维法function xx=gold(no,x,d,amin,amax,e) %UNTITLED6此处显示有关此函数的摘要 %此处显示详细说明a1=amax-0.618*(amax-amin)

27、;y1=f_1(no,x,d,a1);a2=amin+0.618*(amax-amin);y2=f_1(no,x,d,a2);while abs(amax-amin)>=eif y1>=y2amin=a1;a1=a2;y1=y2;a2=amin+0.618*(amax-amin);y2=f_1(no,x,d,a2);elseamax=a2;a2=a1;y2=y1;a1=amax-0.618*(amax-amin);y1=f_1(no,x,d,a1);end endxx=(amax+amin)/2;end12、牛顿一维法function xx=newton(n0,x,d,amax,e

28、)%UNTITLED9此处显示有关此函数的摘要%此处显示详细说明syms s;z=x+s*d;if n0=1a=z(1);b=z(2);f=aA2+2*bA2-2*a*b-4*a;elseif n0=2a=z(1);b=z(2);f=aA4-2*a*a*b-2*a*b+a*a+2*b*b+4.5*a-4*b+4;elseif n0=3a=z(1);b=z(2);c=z(3);f=a*a+b*b+c*c;endx0=amax;while(1)if subs(diff(diff(f,s),s),x0)=0 break;endx0 = x0-double(subs(diff(f,s),x0)/sub

29、s(diff(diff(f,s),s),x0);if abs(double(subs(diff(f,s),x0)<e break;end end xx=x0;end13、不精确一维搜索法functionalf=wolfe(n0,x,d) if n0=1 syms a b;f=aA2+2*bA2-2*a*b-4*a;g=gradient(f);xx=a;b;elseif n0=2 syms a b;f=aA4-2*a*a*b-2*a*b+a*a+2*b*b+4.5*a-4*b+4; g=gradient(f);xx=a;b;elseif n0=3syms a b c;f=a*a+b*b+c*c;g=gradient(f);xx=a;b;c;endu=0.1; q=0.4; aa=0; bb=inf; alf=1;fx=double(subs(f,xx,x);gx=d