數(shù)學(xué)實驗——無約束優(yōu)化.docx

數(shù)學(xué)實驗報告實驗6 無約束優(yōu)化實驗6 無約束優(yōu)化分1 黃浩 2011011743一、 實驗?zāi)康?. 掌握用MATLAB 優(yōu)化工具箱的基本用法，對不同算法進行初步分析、比較。2. 練習(xí)用無約束優(yōu)化方法建立和求解實際問題模型（包括非線性最小二乘擬合）。二、 實驗內(nèi)容1. 數(shù)學(xué)實驗第二版（問題2.1）問題敘述：取不同的初值計算非線性規(guī)劃：minz=x1x22(1-x1)21-x1-x2(1-x1)52盡可能求出所有局部極小點，進而找出全局極小點，并對不同算法（搜索方向、步長搜索、數(shù)值梯度與分析梯度等）的結(jié)果進行分析、比較。實驗過程：首先繪制這個函數(shù)的三維圖形以及等高線（程序見四.1），結(jié)果如下： 通過觀察這兩幅圖，可以得到，x2確定時，x1越負，函數(shù)值越大，x1確定時，x2絕對值越大，函數(shù)值越大。但對于x1正向偏離0的情況，并沒有很好的反映，于是擴大繪圖范圍，做出下圖（程序見四.2）：由上面兩幅圖可見，方程像是一個四角被捏起的花布，而且z的最小值為0。因此只要求解該方程的零點，即得到了方程的局部極小點，且若將原方程變形為：z=x1x22(1-x1)41-x2(1-x1)42我們?nèi)菀装l(fā)現(xiàn)，該方程的零點為：x2=0或x1=0或x1=1或x2(1-x1)4=1在求解零點之前，先針對一個零點，不妨用x1=1，x2=1，分析不同算法的優(yōu)劣。在matlab的無約束優(yōu)化中，可以使用fminumc和fminsearch兩種函數(shù)，搜索方向的算法有BFGS公式、DFP公式和最速下降法三種（書中還提到的Gill-Murray公式在matlab中已經(jīng)不再使用），步長的一維搜索有混合二次三次多項式插值和三次多項式插值兩種方法，另外，在求解函數(shù)梯度是也有數(shù)值方法和分析方法兩種。在對上述四類算法因素進行分析時，我們采用控制變量法，每次只保持一種或兩種算法因素改變，分析它的精度及效率。（一） 分析fminumc與fminsearch兩種方法的精度及效率選擇初值為x1=0.8,x2=0.8，使用fminunc和fminsearch的默認算法及控制參數(shù)，輸出結(jié)果如下（程序見四.3、四.4）：最優(yōu)解x1最優(yōu)解x2最優(yōu)值z迭代次數(shù)z調(diào)用次數(shù)fminunc9.9400e-0017.9238e-0018.0435e-0101239fminsearch1.0000e+0005.8842e-0019.5246e-05767123 因為精確解為x1=1，z=0，我們便可以比較出不同算法的精度。由上表可見，fminunc的迭代次數(shù)和目標函數(shù)的調(diào)用次數(shù)較少，效率較高，但最優(yōu)值的精度較低，與之相反的是，fminsearch的迭代次數(shù)和目標函數(shù)的調(diào)用次數(shù)較多，效率比較低，然而精度很高。（二） 比較不同搜索方向和不同步長的算法精度及效率使用fminunc函數(shù)，選擇初值為x1=0.6,x2=0.8，分別用BFGS公式、DFP公式和最速下降法三種搜索方向，混合二次三次多項式插值和三次多項式插值兩種步長算法，共2*3=6種方式進行運算，結(jié)果如下（程序見四.5）：搜索方向步長算法最優(yōu)解x1最優(yōu)解x2最優(yōu)值z迭代次數(shù)z調(diào)用次數(shù)BFGS混合二次三次多項式插值1.0072e+0007.6715e-0011.5634e-009727DFP1.0074e+0007.5404e-0011.7336e-009521最速下降法9.9290e-0017.2793e-0011.3266e-00939255BFGS三次多項式插值1.0072e+0007.6715e-0011.5634e-009727DFP1.0074e+0007.5404e-0011.7336e-009521最速下降法9.9290e-0017.2793e-0011.3266e-00939255顯然，精確解應(yīng)為x1=1，z=0.由上表可知，計算效率：DFPBFGS最速下降法，而且DFP和BFGS的迭代次數(shù)近似，最速下降法的迭代次數(shù)遠遠大于前兩者。這是因為最速下降法的收斂階數(shù)為1，開始時收斂非?？欤斀咏顑?yōu)解時，收斂就變慢了。而前兩種算法均利用到了二階導(dǎo)數(shù)項，收斂階數(shù)較高，因此總體上收斂速率較快，計算效率更高。對于解的精度，三者近似相同，沒有很大的區(qū)別。 此外，通過對比可以看出，對于混合二次三次多項式插值和三次多項式插值兩種步長算法，他們的最優(yōu)解以及迭代次數(shù)都是一樣的。但事實上不是二者沒有區(qū)別，我通過翻閱matlab的warning消息，得知LineSearchType的控制參數(shù)已經(jīng)不存在了，因而matlab只使用了某一種步長算法。（三） 比較數(shù)值方法和分析方法的精度及效率 首先編寫含梯度的函數(shù)M文件（程序見四.6），然后選擇初值x1=0.6，x2=0.8，分別使用數(shù)值方法和分析方法進行計算（程序見四.7），結(jié)果如下：搜索方向算法最優(yōu)解x1最優(yōu)解x2最優(yōu)值迭代次數(shù)z調(diào)用次數(shù)BFGS數(shù)值方法1.00720E+007.67150E-011.56340E-09727DFP1.00740E+007.54040E-011.73360E-09521最速下降法9.92900E-017.27930E-011.32660E-0939255BFGS分析方法1.00720E+007.67150E-011.56340E-0979DFP1.00740E+007.54040E-011.73360E-0957最速下降法9.92900E-017.27930E-011.32660E-093985由上表可知，分析方法的目標函數(shù)調(diào)用次數(shù)少于數(shù)值方法，在進行大規(guī)模計算的時候，分析方法不失為數(shù)值方法的一種替代。但從本例來看，分析方法對于提高精度的效果不明顯，而且對迭代次數(shù)也沒有減少。在進行完上述算法的比較之后，我們開始逐個檢驗之前得到的極小值點，選定fminunc函數(shù)和BFGS搜索方向進行分析, 取不同的初值（程序見四.8），結(jié)果如下：初值x01初值x02期望解X1期望解X2最優(yōu)解x1最優(yōu)解x2最優(yōu)值0.30.3002.2336e-0012.8437e-0061.4678e-0130.2-20-7.5875e-009-1.3095e+0005.2657e-0161.241-1.0026e+0003.9964e+0007.2181e-0100.3180.5161.0088e+0001.7980e+0011.9765e-0060.4515.60.5164.9682e-0011.5599e+0012.6331e-014在第四步實驗（即上表第四行）中，作者本來希望初值（0.3,18）能收斂到期望值（0.5,16），然而卻最終得到的最優(yōu)解為（1,18），雖然也為極值點，但搜索方向卻沒有按照預(yù)期進行。于是進行了第五步的實驗，將初值與期望值的偏差進一步縮小，結(jié)果獲得了期望的最優(yōu)解。由上表的最優(yōu)值，我們可以判定，當x2=0或x1=0或x1=1或x2(1-x1)4=1時，目標函數(shù)值為0，又因為目標函數(shù)是非負連續(xù)的，因而零點即為局部極小值，也為全局極小值。但這只是局部極小值的充分條件，僅憑此并不能找到所有的局部極小值。然而，通過對圖像的觀察，可以直觀地判斷出零點和極小值點是等價的，至于更分析上的探討，由于本人數(shù)學(xué)能力有限，無法做出全面的證明。得出結(jié)論：該函數(shù)的極小值點為：x2=0或x1=0或x1=1或x2(1-x1)4=1，極小值點即為零點，因而極小值為0.對于不同算法的優(yōu)劣，有如下四點結(jié)論：a) fminunc的迭代次數(shù)和目標函數(shù)的調(diào)用次數(shù)較少，效率較高，但最優(yōu)值的精度較低，與之相反的是，fminsearch的迭代次數(shù)和目標函數(shù)的調(diào)用次數(shù)較多，效率比較低，然而精度很高。b) 對于DFP、BFGS和最速下降法這三種搜索方向，計算效率：DFPBFGS最速下降法，而三者的計算精度在數(shù)量級上是近似的。c) 對于混合二次三次多項式插值和三次多項式插值兩種步長算法，由于matlab中已經(jīng)沒有了步長的設(shè)置參數(shù)，因此無法進行比較。d) 對于數(shù)值方法和分析方法計算梯度，結(jié)論是分析方法可以大幅減少目標函數(shù)的調(diào)用次數(shù)，但對于提高精度的貢獻很小，迭代次數(shù)也沒有區(qū)別。2. 數(shù)學(xué)實驗第二版（問題6）問題描述：有一組數(shù)據(jù)(ti,yi)(i=1,2,33)，其中ti=10(i-1)，yi由表中給出?，F(xiàn)在用這組數(shù)據(jù)擬合函數(shù)：f(x,t)=x1+x2e-x4t+x3e-x5t其中參數(shù)x，初值可選為(0.5,1.5,-1,0.01,0.02)，用GN和LM兩種方法求解。對yi作一擾動，即yi+ei，ei為(-0.05,0.05)內(nèi)的隨機數(shù)，觀察并分析迭代收斂是否會變慢。iyiiyiiyi10.844 120.718 230.478 20.908 130.685 240.467 30.932 140.658 250.457 40.936 150.628 260.448 50.925 160.603 270.438 60.908 170.580 280.431 70.881 180.558 290.424 80.850 190.538 300.420 90.818 200.522 310.414 100.784 210.506 320.411 110.751 220.490 330.406 實驗過程：在matlab中，實現(xiàn)最小二乘擬合有l(wèi)sqnonlin和lsqcurvefit，鑒于后者的代碼更清晰簡潔，因而我使用了lsqcurvefit來進行擬合。又因為在此版本的matlab中，已經(jīng)不能設(shè)置LevenbergMarquardt的控制參數(shù)，即不能手動設(shè)置非線性最小二乘擬合的下降方向，因此下面我只給出了LM法的擬合結(jié)果（程序見四.9、四.10）：x1x2x3x4x5迭代數(shù)z調(diào)用次數(shù)norm0.375411.9358-1.46470.012870.022128545.4649e-005因此，對于函數(shù)f(x,t)=x1+x2e-x4t+x3e-x5t，五個參數(shù)的取值分別為0.37541、1.9358、-1.4647、0.01287、0.02212時，與上述數(shù)據(jù)點的差別最小，擬合度最高。誤差平方和為5.4649e-005然后，對yi作一擾動，即yi+ei，ei為(-0.05,0.05)內(nèi)的隨機數(shù)，并重復(fù)四次（程序見四.11）所得結(jié)果如下：序號x1x2x3x4x5迭代數(shù)z調(diào)用次數(shù)norm10.4012.399-1.9960.0150.0239602.2390e-00220.3421.307-0.8170.0100.0258542.2705e-00230.40725.920-25.5000.0180.018503063.0588e-00240.3671.540-1.1040.0120.02611722.1883e-002由上表可以看出，擾動后輸出的誤差平方和norm數(shù)量級在0.01左右，遠遠大于未擾動之前的誤差。因此，當一個隨機小擾動加到y(tǒng)上之后，新的數(shù)據(jù)點一般會產(chǎn)生較大的偏離，致使數(shù)據(jù)點與函數(shù)的吻合性變差了。不僅如此，當加入擾動后，迭代數(shù)與目標函數(shù)調(diào)用次數(shù)也有或多或少的增加，這說明收斂速度變慢了，這也是由于數(shù)據(jù)點擬合度不高的必然結(jié)果。得出結(jié)論：對于函數(shù)f(x,t)=x1+x2e-x4t+x3e-x5t，五個參數(shù)的取值分別為0.37541、1.9358、-1.4647、0.01287、0.02212時，擬合度最高。而且在加入擾動后，擬合度和精度大大下降，收斂速度也有所變慢。三、 實驗總結(jié)本次實驗是利用matlab優(yōu)化工具箱進行無約束優(yōu)化，并進行多元函數(shù)極小值的分析及非線性含參量函數(shù)的最小二乘擬合。這一部分無論在基本原理上還是代碼上都比之前的實驗更有難度，如果說之前的實驗是灰箱的話，那么我認為本次實驗完全是“半黑”箱，各種搜索方向和步長算法完全是照搬代碼，而對其內(nèi)部的原理不甚清楚，這就導(dǎo)致在得出實驗結(jié)果之后，只能進行表面上的分析而不能從內(nèi)部原理來解釋各種算法的優(yōu)劣。當然，對于今后可能的工程應(yīng)用而言，數(shù)學(xué)的工具性高于它的理論性，通過這次實驗，我也基本掌握了函數(shù)極值的基本求解思路及非線性最小二乘的方法，對算法和精度控制的有關(guān)參數(shù)也能夠熟練應(yīng)用，達到了基本要求。四、 程序清單1. 第一題繪制第一幅三維圖和等高線x1,x2=meshgrid(-2:0.001:1,-1:0.001:1);z=(x1.*x2).2.*(1-x1).2.*(1-x1-x2.*(1-x1).5).2;mesh(x1,x2,z)pause;contour(x1,x2,z,120)2. 第一題繪制第二幅三維圖和等高線x1,x2=meshgrid(-8:0.01:10,-5:0.01:5);z=(x1.*x2).2.*(1-x1).2.*(1-x1-x2.*(1-x1).5).2;mesh(x1,x2,z)pause;contour(x1,x2,z,220)3. 第一題設(shè)置無約束優(yōu)化的函數(shù)function z=profun1(x)z=(x(1)*x(2)2*(1-x(1)2*(1-x(1)-x(2)*(1-x(1)5)2;end4. 第一題比較fminunc和fminsearch的優(yōu)劣x0=0.8,0.8;format shortex1,z1,ef1,out1=fminunc(profun1,x0)x2,z2,ef2,out2=fminsearch(profun1,x0)5. 第一題比較不同搜索方向和不同步長的算法精度及效率x0=0.6,0.8;format shorteopt=optimset(LargeScale,off,MaxFunEvals,1000);opt1=optimset(opt,HessUpdate,bfgs);opt2=optimset(opt,HessUpdate,dfp);opt3=optimset(opt,HessUpdate,steepdesc);opt4=optimset(opt,HessUpdate,bfgs,LineSearchType,cubicpoly);opt5=optimset(opt,HessUpdate,dfp,LineSearchType,cubicpoly);opt6=optimset(opt,HessUpdate,steepdesc,LineSearchType,cubicpoly);x,fv,exit,out=fminunc(profun1,x0,opt1)pausex,fv,exit,out=fminunc(profun1,x0,opt2)pausex,fv,exit,out=fminunc(profun1,x0,opt3)pausex,fv,exit,out=fminunc(profun1,x0,opt4)pausex,fv,exit,out=fminunc(profun1,x0,opt5)pausex,fv,exit,out=fminunc(profun1,x0,opt6)pause6. 第一題設(shè)置分析方法的函數(shù)function z,g=profun11(x);z=(x(1)*x(2)2*(1-x(1)2*(1-x(1)-x(2)*(1-x(1)5)2;if nargout1 g(1)=2*x(1)*x(2)2*(1-x(1)-x(2)*(1-x(1)5)*(2*x(1)2-3*x(1)+1)*(1-x(1)-x(2)*(1-x(1)5)+(x(1)3-2*x(1)2+x(1)*(5*x(2)*(1-x(1)4-1); g(2)=2*x(1)2*(1-x(1)2*x(2)*(1-x(1)-x(2)*(1-x(1)5)2-2*x(1)2*(1-x(1)7*x(2)2*(1-x(1)-x(2)*(1-x(1)5);end7. 第一題比較數(shù)值方法和分析方法計算梯度x0=0.6,0.8;format shorteopt=optimset(LargeScale,off,MaxFunEvals,1000);opt1=optimset(opt,HessUpdate,bfgs);opt2=optimset(opt,HessUpdate,dfp);opt3=optimset(opt,HessUpdate,steepdesc);opt=optimset(LargeScale,off,MaxFunEvals,1000,GradObj,on);opt4=optimset(opt,HessUpdate,bfgs);opt5=optimset(opt,HessUpdate,dfp);opt6=optimset(opt,HessUpdate,steepdesc);x1,fv1,exit,out1=fminunc(profun1,x0,opt1)x2,fv2,exit,out2=fminunc(profun1,x0,opt2)x3,fv3,exit,out3=fminunc(profun1,x0,opt3)x4,fv4,exit,out4=fminunc(profun11,x0,opt4)x5,fv5,exit,out5=fminunc(profun11,x0,opt5)x6,fv6,exit,out6=fminunc(profun11,x0,opt6)x=x1;x2;x3;x4;x5;x6fv=fv1;fv2;fv3;fv4;fv5;fv6iterations=out1.iterations;out2.iterations;out3.iterations;out4.iterations;out5.iterations;out6.iterationsfunccount=out1.funcCount;out2.funcCount;out3.funcCount;out4.funcCount;out5.funcCount;out6.funcCount8. 第一題選定fminunc函數(shù)和BFGS搜索