MATLAB黃金分割法課程論文

1、中南林業(yè)科技大學本科課程論文學 院： 理學院 專業(yè)年級： 14級信息與計算科學2班 學生姓名： 邱文林 學 號： 20144349 課 程： MATLAB程序設計教程 設計題目：基于MATLAB的黃金分割法 與拋物線插值法 指導教師： 龔志偉 2016年4月中文摘要 為了求解最優(yōu)化模型的最優(yōu)解，可使用基于MATLAB算法編程的黃金分割法與拋物線插值法，來實現求解的過程。黃金分割法是通過所選試點的函數值而逐步縮短單谷區(qū)間來搜索最優(yōu)點，利用迭代進而得出結論。拋物線插值法亦稱二次插值法，是一種多項式插值法，逐次以擬合的二次曲線的極小點，逼近原尋求函數極小點的一種方法。通過將MATLAB與最優(yōu)化問題相

2、結合，不僅可以加深對黃金分割法、拋物線插值法的基本理解和算法框圖及其步驟的全面理解，也有利于幫助我們掌握MATLAB的使用方法。關鍵詞：MATLAB，黃金分割法，拋物線插值法，最優(yōu)解，迭代英文摘要In order to solve the optimization model of the optimal solution, using MATLAB algorithm based on the golden section method and the parabola interpolation method, to realize the process of solving. The

3、golden section method is used to search the most advantage through the function value of the selected pilot, which can be used to search for the most advantage. Parabolic interpolation method, also known as the two interpolation method, is a polynomial interpolation method, successive to fit the two

4、 curve of the minimum point, the original search function to find a very small point of the method. By combining MATLAB and optimization problems can not only deepen the comprehensive understanding of the golden section method, the parabola interpolation basic understanding and block diagram of the

5、algorithm and steps, but also conducive to help us to grasp the method of using MATLAB.Key words: MATLAB, golden section method, parabolic interpolation method, optimal solution, iteration目 錄1. 黃金分割法 21.1 算法原理21.2 算法步驟21.3黃金分割法算法框圖 32. 拋物線插值法42.1 算法原理42.2 算法步驟42.3 拋物線插值法算法框圖53. 算法的MATLAB實現63.1黃金分割法程

6、序代碼 63.2 實例驗證63.3 誤差分析93.4 拋物線插值法程序代碼103.5 實例驗證113.6 誤差分析123.7 黃金分割法與拋物線插值法的對比13參考文獻15引言為了對最優(yōu)化問題進行求解，為了更好的掌握MATLAB的運用，本文主要介紹一維最優(yōu)化問題的一維搜索法黃金分割法與拋物線插值法，利用MATLAB算法將其實現。黃金分割法也叫0.618法，它是一種基于區(qū)間收縮的極小點搜索算法，當用進退法確定搜索區(qū)間后，我們只知道極小點包含于搜索區(qū)間內，但是具體是哪個點，無法得知。黃金分割法的思想很直接，既然極小點包含于搜索區(qū)間內，那么可以不斷的縮小搜索區(qū)間，就可以使搜索區(qū)間的端點逼近到極小點。

7、拋物線插值法(parabolic interpolation met-hod)亦稱二次插值法一種多項式插值法.逐次以擬合的二次曲線的極小點，逼近原尋求函數極小點的一種方法,能求得精度較高的解，但速度會比較慢。時至今日，經過Math Works公司的不斷完善，MATLAB已經發(fā)展成為適合多學科、多種工作平臺的功能強勁的大型軟件。在國外，MATLAB已經經受了多年考驗。在歐美等高校，MATLAB已經成為線性代數、自動控制理論、數理統計、數字信號處理、時間序列分析、動態(tài)系統仿真等高級課程的基本教學工具；成為攻讀學位的大學生、碩士生、博士生必須掌握的基本技能。在設計研究單位和工業(yè)部門，MATLAB被廣

8、泛用于科學研究和解決各種具體問題。1. 黃金分割法1.1 算法原理黃金分割法的基本原理：黃金分割法又稱0.618法，它是通過不斷縮短搜索區(qū)間的長度來尋求一維函數的極小點。這種方法的原理是：在搜索區(qū)間a, b內按如下規(guī)則對稱地取兩點：計算它們的函數值 f1=f(a1), f2=f(a2) ,比較它們的大小，結果有兩種可能：1.2 算法步驟（1）f1>f2,如圖1所示，極小點必在a1, b內，消去區(qū)間a, a1), 令a=a1，產生新區(qū)間a, b,到此區(qū)間縮短了一次。值得注意的是新區(qū)間的a1點與原區(qū)間的a2點重合，可令a1=a2,這樣可少找一個新點和節(jié)省一次函數值計算。（2）f1<=f

9、2,極小點必在a, a2內，消去區(qū)間 (a2,b,令b=a2,產生新區(qū)間a ,b,到此區(qū)間縮短了一次。同樣新區(qū)間a2點與原區(qū)間的a1點重合，可令a2=a1,f2<=f1。（3）當縮短的新區(qū)間長度小于等于某一精度，即b-a<=時，取a* =(a1+a2)/2。1.3黃金分割法算法框圖給定 a,b,ea+0.382(b-a) a1, f(a1) f1a+0.618(b-a) a2, f(a2) f2NYf1 > f2 ?a=a1 , a1 =a2 , f1 =f2a2 =a+0.618(b-a)f2 = f(a2)b=a2 ,a2 =a1 , f2= f1a1=a+0.382(b

10、-a)f1 = f(a1)結束輸出： a*=(a+b)/2NYb a < e?2. 拋物線插值法2.1 算法原理: 拋物線也叫二次插值法，其理論依據為二次多項式可以在最優(yōu)點附近較好的逼近函數的形狀，做法是在函數f(x)的最優(yōu)點附近取三個構造點，然后用這三個點構造一條拋物線，把這條拋物線的極值點作為函數f(x)的極值點的近似。每次構造一條拋物線后，拋物線的極值點就可作為一個新的構造點，新的構造點與原來的三個構造點經過某種算法，得到下一步拋物線逼近的三個構造點，這就是拋物線法的算法過程。 2.2 算法步驟：用拋物線求無約束問題minf(x),xR的算法步驟如下： 確定初始區(qū)間a,b,選定初始

11、插值內點t0(a,b)及精度e>0,令a0=a,b0=b,k=0; 求二次插值多項式的極小點：tk+1= . (2.10) 若f(tk+1) f(tk),則當|tk+1 -tk|e時，停止迭代輸出tk+1;否則轉; 若f(tk+1) f(tk),則當|tk+1 -tk|e時，停止迭代輸出tk;否則轉; 若tk+1 tk,令ak+1=ak,bk+1=tk,tk+1=tk+1;置k=k+1轉,否則令ak+1=tk,bk+1=bk,tk+1=tk+1;置k=k+1轉。 若tk+1 tk,令 ak+1=tk+1,bk+1=tk,tk+1=tk;置k=k+1轉,否則令ak+1=ak,bk+1=tk

12、+1,tk+1=tk;置k=k+1轉。初始插值內點可以取搜索區(qū)間a, b的中點。2.3拋物線插值法算法框圖開始確定tk,ak,bk,要求f(ak) f(tk),f(bk) f(tk)按（2.10）計算tk+1t*=tk+1 ,f*=f(tk+1)|tk+1 -tk|eYN輸出t*,f*Ntk+1tk ?結束Yf(tk+1) f(tk)?f(tk+1) f(tk) ?ak=tk+1,f(ak)=f(tk+1)YYbk=tk, tk=tk+1,f(bk)=f(tk),f(tk)=f(tk+1)ak=tk, tk=tk+1,f(ak)=f(tk),f(tk)=f(tk+1)bk =tk+1 ,f(b

13、k)=f(tk+1)3. 算法的MATLAB實現3.1 黃金分割法程序代碼在MATLAB中編程實現的黃金分割法函數為: golden。功能：用黃金分割法求解一維函數的極值。調用格式：tmin=golden(f,a,b,e)其中， f=目標函數； a=極值區(qū)間左端點； b=極值區(qū)間右端點； e=精度； tmin=目標取最小值時的自變量值；黃金分割法的MATLAB程序代碼如下：主程序：syms t a b e ;a=input('搜索區(qū)間第一點 a=');b=input('搜索區(qū)間第二點 b=');e=input('搜索精度ne=');disp(&#

14、39;需求的最優(yōu)化函數 f=f(t),tmin=golden(f,a,b,e)'); m文件：function tmin=golden(f,a,b,e)format long;k=0;a1=b-0.618*(b-a);a2=a+0.618*(b-a);while b-a>e y1=subs(f,a1); y2=subs(f,a2);if y1>y2 a=a1; a1=a2; y1=y2; a2=a+0.618*(b-a);else b=a2; a2=a1; y2=y1; a1=b-0.618*(b-a);endk=k+1;endtmin=(a+b)/2;fmin=subs(

15、f,tmin)fprintf('k=n');disp(k);x=-3:0.001:5;y=x.*(x+2);plot(x,y,'k-',tmin,fmin,'bp');3.2 實例驗證minf(t)=t*(t+2), 已知初始單谷區(qū)間a,b=-3,5,按精度e=0.001計算。代入到主程序代碼：搜索區(qū)間的第一點 a=-3搜索區(qū)間的第二點 b=5搜索精度e=0.001需求的最優(yōu)化函數 f=f(t),tmin=golden(f,a,b,e)>> f=t*(t+2) f = t*(t + 2) >> tmin=golden(f,

16、a,b,e)k= 19tmin =>> fmin=tmin*(tmin+2)fmin = -0.999999985524973目標函數的曲線圖如圖1-1所示:注： 圖中所標識的點為黃金分割法所求的極小點。 圖1-1 函數f(t)=t(t+2)的曲線圖3.3 誤差分析:1. f(t)=t*(t+2) 的精確解如下：f=1,2,0;p=polyder(f);t1=roots(p)f1=polyval(f,t1)t1 = -1f1 = -1極小點誤差： e1=(tmin-t1)/t1=(-(-1)/(-1) 1e-4極小值誤差： e2=(f1-fmin)/f1=( -1 +0.99999

17、9985524973)/(-1) 1e-8結論：在精度e=0.001的情況下，通過黃金分割法經過19次迭代，求得：tmin=,與精確值的誤差e1=1e-4<e=0.001,滿足題目要求，黃金分割法是正確的，也可預見，隨著精度e的提高，tmin的值會越來越接近極小點，最終會等于精確值-1.由圖1-1也可以看出，圖中所標識的位置幾乎就在t=-1的位置上。3.4 拋物線插值法程序代碼在MATLAB中編程實現的拋物線法函數為: minPWX。功能：用拋物線插值法求解一維函數的極值。調用格式：x,minf=minPWX(f,a,b,e)其中， f:目標函數； A:初始搜索區(qū)間左端點； b:初始搜索

18、區(qū)間右端點； e:精度； x:目標取最小值時的自變量值； minf:目標函數的最小值。拋物線法的MATLAB程序代碼如下：function x,minf=minPWX(f,a,b,e)%目標函數：f;%初始搜索區(qū)間左端點：a;%初始搜索區(qū)間右端點：b;%精度：e;%目標函數取最小值時的自變量值: x;%目標函數的最小值: minf format long;if nargin = 3 e=1.0e-3;endt0=(a+b)/2;k=0;tol=1; while tol>e fa=subs(f,findsym(f),a); %區(qū)間左端點函數值 fb=subs(f,findsym(f),b)

19、; %區(qū)間右端點函數值 ft0=subs(f,findsym(f),t0); %內插點函數值 tu=fa*(b2-t02)+fb*(t02-a2)+ft0*(a2-b2); td=fa*(b-t0)+fb*(t0-a)+ft0*(a-b); t1=tu/2/td;%插值多項式的極小點 ft1=subs(f,findsym(f),t1); %插值多項式的極小值 tol=abs(t1-t0); if ft1<=ft0 if t1<=t0 b=t0; t0=t1; else a=t0; t0=t1; end k=k+1; else if t1<=t0 a=t1; else b=t1; end k=k+1; endendx = t1;minf= subs(f,findsym(f),x);format short; 3.5 實例驗證 用拋物線法求函數f(t)=t*(t+2)，已知初始單谷區(qū)間a,b=-3,5,按精度e=0.001計算。解： 在MATLAb 命令窗口中輸入：>> syms t;>> f=t*(t+2);>> tmin,fmin=minPWX(f,-3,5)tmin =-1fmin = -1目標函數的曲線圖如圖1-2所示:注： 圖中O所標識的點為拋物線插值法