中科大《優(yōu)化設(shè)計(jì)》課程大作業(yè)之無約束優(yōu)化實(shí)驗(yàn)報(bào)告(共17頁)

精選優(yōu)質(zhì)文檔-傾情為你奉上無約束優(yōu)化設(shè)計(jì)實(shí)驗(yàn)報(bào)告力學(xué)系 型號：聯(lián)想y470CPU：i5-2450M內(nèi)存：2GB系統(tǒng)：win7-64位如下是三個(gè)目標(biāo)函數(shù)（包括自定義函數(shù)）以及初值和精度選取：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如下是運(yùn)算結(jié)果：目標(biāo)函數(shù)無約束方法一維搜索所需時(shí)間迭代次數(shù)極值點(diǎn)極值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目標(biāo)函數(shù)無約束方法一維搜索所需時(shí)間迭代次數(shù)極值點(diǎn)極值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目標(biāo)函數(shù)無約束方法一維搜索所需時(shí)間迭代次數(shù)極值點(diǎn)極值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總結(jié)及比較：根據(jù)上面三個(gè)函數(shù)的表格可以看出：首先，從迭代時(shí)間來看，三種一維搜索方法中黃金分割法所用時(shí)間最久，牛頓法和不精確法所用時(shí)間較少，這兩種方法相比較而言牛頓法所用時(shí)間更少一些。而六種無約束方法中，由于單形替代法不需要使用一維搜索方法，故迭代時(shí)間最少，緊接著在使用一維搜索的五種方法中以阻尼牛頓法迭代時(shí)間相對較少，共軛梯度法迭代時(shí)間最久；然后，從迭代次數(shù)來看，共軛梯度法往往需要較多的迭代次數(shù)，從而所需時(shí)間也最久；接著，從計(jì)算結(jié)果的精度來看，阻尼牛頓法的結(jié)果精度最高，而單形替換法的精度最低；最后，從編程來看，在編好一維搜索方法的情況下，最速下降法和阻尼牛頓法編程簡單容易，而共軛梯度法、變尺度法和單形替代法需要兩重循環(huán)實(shí)現(xiàn)，鮑威爾法和單形替換法則需要編程者對矩陣的操作能力有較高的要求，故編程較難。同時(shí)，從上面的結(jié)果也可以發(fā)現(xiàn)，鮑威爾法在使用不精確的一維搜索方法時(shí)，對函數(shù)1無法收斂，對函數(shù)2、3收斂到錯(cuò)誤的結(jié)果，所以鮑威爾法是依賴于精確的一維搜索過程的，而其他幾個(gè)則不依賴于精確一維搜索過程。精確的一維搜索方法通常需要花費(fèi)很大的工作量，特別是當(dāng)?shù)c(diǎn)遠(yuǎn)離問題的解時(shí)，精確的求解一個(gè)一維子問題通常不是十分有效率的。因此，只要保證目標(biāo)函數(shù)值在每一步都有滿意的下降，使用不是非常精確的一維搜索，就可以大大節(jié)省工作量。在分析函數(shù)2的計(jì)算結(jié)果時(shí)，可以發(fā)現(xiàn)存在兩個(gè)收斂結(jié)果，當(dāng)然這兩個(gè)結(jié)果都是極值，因?yàn)楹瘮?shù)2是二元四次函數(shù)，存在多個(gè)極值。故為了驗(yàn)證正確性，自己曾將初始點(diǎn)（-2.5,4.25）調(diào)成（2.5,4.25），分別代入程序中計(jì)算，計(jì)算結(jié)果都收斂于極值為0.99的這個(gè)點(diǎn)上。所以，在存在多個(gè)極值點(diǎn)的情況下結(jié)果是和初始點(diǎn)的選取有關(guān)。對于單形替換法，這種方法不需要一維搜索最佳步長，故沒有一維搜索方法反復(fù)地計(jì)算最佳步長的計(jì)算時(shí)間，程序運(yùn)行效率快。但它的收斂條件不好選擇，通過查找文獻(xiàn)資料總結(jié)出以下三個(gè)收斂條件：1.利用最壞點(diǎn)函數(shù)值與最好點(diǎn)函數(shù)值之差判別；2. 利用相鄰兩次函數(shù)值差值的絕對值判別；3. 利用各點(diǎn)的函數(shù)值與最好點(diǎn)函數(shù)值之差的均方根判別。為了保證程序的執(zhí)行可靠性，這三種常用的方法中自己選擇了判別3，即：綜上所述，阻尼牛頓法是無約束方法中最有效的方法。不僅編程簡單，而且迭代次數(shù)較少，運(yùn)行時(shí)間較短，結(jié)果的精度也較高。在程序的運(yùn)行方面，分別設(shè)置了可變的函數(shù)選擇、無約束方法選擇、一維搜索方法選擇、起始點(diǎn)、精度這五個(gè)輸入，故可以在命令窗口運(yùn)行主程序main，再根據(jù)提示要求分別輸入這五個(gè)參數(shù)的所需值，就可以得到運(yùn)行結(jié)果。程序如下：1、 主函數(shù)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('請輸入上面所想選擇函數(shù)的編號（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('請輸入用空格隔開的初始值坐標(biāo)向量（如：1.1 2.0）：','s');x=str2num(s);x=x'disp(' ');while 1e=input('請輸入精度(建議0.1或0.01)：'); if e>0 break; end disp('此次輸入無效.');end disp(' ');disp('');xx,yy=fmins(m0,m1,n0,x,e);fprintf('迭代次數(shù)為： %8.0fn', k);disp('所求極值點(diǎn)的坐標(biāo)向量為：');fprintf(' %16.5fn', xx);fprintf('所求函數(shù)的極值為： %16.5fn', yy);2、 外部多維的調(diào)用函數(shù)function xx,yy=fmins(m0,m1,n0,x,e)%UNTITLED 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明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 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明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 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明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 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明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 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明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 (ff3<ff0)&&(ff0-2*ff2+ff3)*(ff0-ff2-md)2<0.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 ff2<ff3 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 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明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 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明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 fz<e 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)<fg xk(:,xh)=xk(:,n+3); f(xh)=f(n+3); 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 fz<e break; end k=k+1;endxx=xk(:,xl);yy=f0(n0,xx);end9、 內(nèi)部循環(huán)一維調(diào)用函數(shù)function xx=fmin(m1,n0,x,d,e)x0=0;%初始步長默認(rèn)為0amin,amax=range1(n0,x,d,x0);if m1=1 xx=gold(n0,x,d,amin,amax,e);elseif m1=2 xx=newton(n0,x,d,amax,e);elseif m1=3 xx=wolfe(n0,x,d);endend10、一維搜索確定區(qū)間函數(shù)function amin,amax = range1(n0,x,d,x0)%UNTITLED5 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明h=1;a1=x0;y1=f_1(n0,x,d,a1);a2=a1+h;y2=f_1(n0,x,d,a2);if y2>y1 h=-h; a3=a1;y3=y1; a1=a2; a2=a3;y2=y3;enda3=a2+h;y3=f_1(n0,x,d,a3);while y3<y2 h=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 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明a1=amax-0.618*(amax-amin);y1=f_1(no,x,d,a1);a2=amin+0.618*(amax-amin);y2=f_1(no,x,d,a2);while abs(amax-amin)>=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 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明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、求導(dǎo)函數(shù)function g,ss=gra(n0)%UNTITLED 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明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函數(shù)function yy=f0(n0,xx)%UNTITLED2 此處顯示有關(guān)此函數(shù)的摘要% 此處顯示詳細(xì)說明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函數(shù)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專心-專注-專業(yè)