资源描述
精选优质文档-倾情为你奉上无约束优化设计实验报告力学系 型号：联想y470CPU：i5-2450M内存：2GB系统：win7-64位如下是三个目标函数（包括自定义函数）以及初值和精度选取：1.minf(x)=x12+2*x22-2*x1*x2-4*x1初值x0=1,1；精度为：0.12.minf(x, y)=x4-2*x*x*y-2*x*y+x2+2*y*y+4.5*x-4*y+4初值为（-2.5,4.25）；精度为0.013.minf(x)=x12+x22+x32初值为（3,2,1）；精度为0.01如下是运算结果：目标函数无约束方法一维搜索所需时间迭代次数极值点极值1最速下降法黄金分割法5.818（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.042（4.0,2.0）-8.0牛顿法1.242（4.0,2.0）-8.0不精确法1.353（4.0,2.0）-8.0单形替换法无0.029（4.1,2.2）-8.0目标函数无约束方法一维搜索所需时间迭代次数极值点极值2最速下降法黄金分割法20.7621（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.613（-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,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最速下降法黄金分割法2.702（0,0,0）0牛顿法0.761（0,0,0）0不精确法0.781（0,0,0）0阻尼牛顿法黄金分割法1.611（0,0,0）0牛顿法0.661（0,0,0）0不精确法0.651（0,0,0）0共轭梯度法黄金分割法3.153（0,0,0）0牛顿法0.400（0,0,0）0不精确法0.440（0,0,0）0鲍维尔法黄金分割法7.171（0,0,0）0牛顿法1.991（0,0,0）0不精确法4.080（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总结及比较：根据上面三个函数的表格可以看出：首先，从迭代时间来看，三种一维搜索方法中黄金分割法所用时间最久，牛顿法和不精确法所用时间较少，这两种方法相比较而言牛顿法所用时间更少一些。而六种无约束方法中，由于单形替代法不需要使用一维搜索方法，故迭代时间最少，紧接着在使用一维搜索的五种方法中以阻尼牛顿法迭代时间相对较少，共轭梯度法迭代时间最久；然后，从迭代次数来看，共轭梯度法往往需要较多的迭代次数，从而所需时间也最久；接着，从计算结果的精度来看，阻尼牛顿法的结果精度最高，而单形替换法的精度最低；最后，从编程来看，在编好一维搜索方法的情况下，最速下降法和阻尼牛顿法编程简单容易，而共轭梯度法、变尺度法和单形替代法需要两重循环实现，鲍威尔法和单形替换法则需要编程者对矩阵的操作能力有较高的要求，故编程较难。同时，从上面的结果也可以发现，鲍威尔法在使用不精确的一维搜索方法时，对函数1无法收敛，对函数2、3收敛到错误的结果，所以鲍威尔法是依赖于精确的一维搜索过程的，而其他几个则不依赖于精确一维搜索过程。精确的一维搜索方法通常需要花费很大的工作量，特别是当迭代点远离问题的解时，精确的求解一个一维子问题通常不是十分有效率的。因此，只要保证目标函数值在每一步都有满意的下降，使用不是非常精确的一维搜索，就可以大大节省工作量。在分析函数2的计算结果时，可以发现存在两个收敛结果，当然这两个结果都是极值，因为函数2是二元四次函数，存在多个极值。故为了验证正确性，自己曾将初始点（-2.5,4.25）调成（2.5,4.25），分别代入程序中计算，计算结果都收敛于极值为0.99的这个点上。所以，在存在多个极值点的情况下结果是和初始点的选取有关。对于单形替换法，这种方法不需要一维搜索最佳步长，故没有一维搜索方法反复地计算最佳步长的计算时间，程序运行效率快。但它的收敛条件不好选择，通过查找文献资料总结出以下三个收敛条件：1.利用最坏点函数值与最好点函数值之差判别；2. 利用相邻两次函数值差值的绝对值判别；3. 利用各点的函数值与最好点函数值之差的均方根判别。为了保证程序的执行可靠性，这三种常用的方法中自己选择了判别3，即：综上所述，阻尼牛顿法是无约束方法中最有效的方法。不仅编程简单，而且迭代次数较少，运行时间较短，结果的精度也较高。在程序的运行方面，分别设置了可变的函数选择、无约束方法选择、一维搜索方法选择、起始点、精度这五个输入，故可以在命令窗口运行主程序main，再根据提示要求分别输入这五个参数的所需值，就可以得到运行结果。程序如下：1、 主函数clear;global k;k=0;disp(1.f(x)=x12+2*x22-2*x1*x2-4*x1);disp(2.f(x,y)=x4-2*x*x*y-2*x*y+x2+2*y*y+4.5*x-4*y+4);disp(3.f(x)=x12+x22+x32);while 1 n0=input(请输入上面所想选择函数的编号（1、2、3）：); if n0=1|n0=2|n0=3 break; end disp(此次输入无效.);end disp( );disp(1.最速下降法);disp(2.阻尼牛顿法);disp(3.共轭梯度法);disp(4.鲍威尔法);disp(5.变尺度法);disp(6.单纯形法);while 1 m0=input(请输入上面所想选择无约束方法的编号（1、2、3、4、5、6）：); if m0=1|m0=2|m0=3|m0=4|m0=5|m0=6 break; end disp(此次输入无效.);end disp( );disp(1.黄金分割法);disp(2.牛顿法);disp(3.不精确一维搜索法);while 1 m1=input(请输入上面所想选择一维搜索方法的编号（1、2、3）：); if m1=1|m1=2|m1=3 break; end disp(此次输入无效.);end disp( );s=input(请输入用空格隔开的初始值坐标向量（如：1.1 2.0）：,s);x=str2num(s);x=x;disp( );while 1e=input(请输入精度(建议0.1或0.01)：); if e0 break; end disp(此次输入无效.);end disp( );disp();xx,yy=fmins(m0,m1,n0,x,e);fprintf(迭代次数为： %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=1 tic;xx,yy=zuisu(m1,n0,x,e);toc;elseif m0=2 tic;xx,yy=zuni(m1,n0,x,e);toc;elseif m0=3 tic;xx,yy=gonge(m1,n0,x,e);toc;elseif m0=4 tic;xx,yy=powell(m1,n0,x,e);toc;elseif m0=5 tic;xx,yy=bianchi(m1,n0,x,e);toc;elseif m0=6 tic;xx,yy=danxing(n0,x,e);toc;endend3、 最速法function xx,yy=zuisu(m1,n0,x,e)%UNTITLED2 此处显示有关此函数的摘要% 此处显示详细说明global k; g,ss=gra(n0);while 1 d=-double(subs(g,ss,x); a=fmin(m1,n0,x,d,e); x1=x+a*d; if norm(x1-x)e break; end x=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 1 d=-double(subs(h,ss,x)(-1)*subs(g,ss,x); a=fmin(m1,n0,x,d,e); x1=x+a*d; if norm(x1-x)e break; end x=x1; k=k+1;endxx=x1;yy=f0(n0,xx);end5、 共轭梯度法function xx,yy=gonge(m1,n0,x,e)%UNTITLED2 此处显示有关此函数的摘要% 此处显示详细说明global k; if n0=1 n=2;elseif n0=2 n=4;elseif n0=3 n=2;endg,ss=gra(n0);while 1 kk=0; d=-double(subs(g,ss,x); while 1 a=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)e break; elseif kk=n break; end beta=norm(gx1)2/norm(gx)2; d=-gx+beta*d; x=x1; k=k+1; kk=kk+1; end if norm(gx1)e break; end x=x1; k=k+1;endxx=x1;yy=f0(n0,xx);end6、 鲍威尔法function xx,yy=powell(m1,n0,x,e)%UNTITLED3 此处显示有关此函数的摘要% 此处显示详细说明global k; if n0=1 n=2;elseif n0=2 n=2;elseif n0=3 n=3;endd=zeros(n,n+1);xk=zeros(n,n+1);deta=zeros(n,1);for j=1:n d(j,j)=1;endwhile 1 xt=x; for i=1:n a=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); end xt=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 (ff3ff0)&(ff0-2*ff2+ff3)*(ff0-ff2-md)20.5*md*(ff0-ff3)2) 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:n d(:,i)=d(:,i+1); end else if ff2ff3 x=xk(:,n); else x=xk(:,n+1); end end if norm(x-xt)e break; end k=k+1;endxx=x;yy=f0(n0,xx);end7、 变尺度法function xx,yy=bianchi(m1,n0,x,e)%UNTITLED 此处显示有关此函数的摘要% 此处显示详细说明global k; if n0=1 n=2;elseif n0=2 n=2;elseif n0=3 n=3;endg,ss=gra(n0);while 1 gx=double(subs(g,ss,x); h=eye(n); kk=0; while 1 d=-h*gx; a=fmin(m1,n0,x,d,e); xk=x+a*d; if norm(xk-x)e break; end if kk=n break; end gxk=double(subs(g,ss,xk); yk=gxk-gx; sk=xk-x; h=h+(sk*sk)/(sk*yk)-(h*yk)*yk*h)/(yk*h*yk); x=xk; gx=gxk; kk=kk+1; k=k+1; end if norm(xk-x)e break; end x=xk; k=k+1;endxx=xk;yy=f0(n0,xx);end8、 单形替换法function xx,yy=danxing(n0,x,e)%UNTITLED 此处显示有关此函数的摘要% 此处显示详细说明global k; if n0=1 n=2;elseif n0=2 n=2;elseif n0=3 n=3;endf=zeros(n+5,1);xk=zeros(n,n+5);h=2*eye(n);xk(:,1)=x;for i=1:n xk(:,i+1)=x+h(:,i);endwhile 1 for i=1:n+1 f(i)=f0(n0,xk(:,i); end while 1 f(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(1); fff=f; fff(xh)=; fg=max(fff); fz=0; for i=1:n+1 fz=fz+(f(i)-f(xl)2; end fz=sqrt(fz/n); if fze break; end xkk=xk(:,1); for i=1:n xkk=xkk+xk(:,i+1); end xk(:,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)fl xk(:,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); else xk(:,xh)=xk(:,n+3); f(xh)=f(n+3); end else if f(n+3)=fh xk(:,n+3)=xk(:,xh); end xk(:,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)fh xk(:,xh)=xk(:,n+5); f(xh)=f(n+5); else for i=1:n+1 xk(:,i)=(xk(:,i)+xk(:,xl)/2; end break; end end end k=k+1; end if fzy1 h=-h; a3=a1;y3=y1; a1=a2; a2=a3;y2=y3;enda3=a2+h;y3=f_1(n0,x,d,a3);while y3=e if y1=y2 amin=a1; a1=a2; y1=y2; a2=amin+0.618*(amax-amin); y2=f_1(no,x,d,a2); else amax=a2; a2=a1; y2=y1; a1=amax-0.618*(amax-amin); y1=f_1(no,x,d,a1); endendxx=(amax+amin)/2;end12、牛顿一维法function xx=newton(n0,x,d,amax,e)%UNTITLED9 此处显示有关此函数的摘要% 此处显示详细说明syms s;z=x+s*d;if n0=1 a=z(1); b=z(2); f=a2+2*b2-2*a*b-4*a;elseif n0=2 a=z(1); b=z(2); f=a4-2*a*a*b-2*a*b+a*a+2*b*b+4.5*a-4*b+4;elseif n0=3 a=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; end x0 = x0-double(subs(diff(f,s),x0)/subs(diff(diff(f,s),s),x0); if abs(double(subs(diff(f,s),x0)e break; endendxx=x0;end13、不精确一维搜索法function alf=wolfe(n0,x,d)if n0=1 syms a b; f=a2+2*b2-2*a*b-4*a; g=gradient(f); xx=a;b;elseif n0=2 syms a b; f=a4-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=3 syms 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=double(subs(g,xx,x);while 1 xk=x+alf*d; h=-u*alf*gx*d; while fx-double(subs(f,xx,xk) h bb=alf; alf=(alf+aa)/2; h=-u*alf*gx*d; xk=x+alf*d; end gk=double(subs(g,xx,xk); if gk*d q*gx*d aa=alf; alf=min(2*alf,(alf+bb)/2); else break; endend14、求导函数function g,ss=gra(n0)%UNTITLED 此处显示有关此函数的摘要% 此处显示详细说明if n0=1 syms a b; f=a2+2*b2-2*a*b-4*a; g=gradient(f); ss=a;b;elseif n0=2 syms a b; f=a4-2*a*a*b-2*a*b+a*a+2*b*b+4.5*a-4*b+4; g=gradient(f); ss=a;b;elseif n0=3 syms a b c; f=a*a+b*b+c*c; g=gradient(f); ss=a;b;c;endend15、f0函数function yy=f0(n0,xx)%UNTITLED2 此处显示有关此函数的摘要% 此处显示详细说明if n0=1 a=xx(1); b=xx(2); yy=a2+2*b2-2*a*b-4*a;elseif n0=2 a=xx(1); b=xx(2); yy=a4-2*a*a*b-2*a*b+a*a+2*b*b+4.5*a-4*b+4;elseif n0=3 a=xx(1); b=xx(2); c=xx(3); yy=a*a+b*b+c*c;endend16、f_1函数function yy=f_1(n0,x,d,xx)syms s;z=x+s*d;if n0=1 a=z(1); b=z(2); f=a2+2*b2-2*a*b-4*a;elseif n0=2 a=z(1); b=z(2); f=a4-2*a*a*b-2*a*b+a*a+2*b*b+4.5*a-4*b+4;elseif n0=3 a=z(1); b=z(2); c=z(3); f=a*a+b*b+c*c;endyy=double(subs(f,s,xx);end专心-专注-专业
展开阅读全文