實驗一工業(yè)生產(chǎn)的產(chǎn)量問題

1、實驗七 生產(chǎn)計劃中的產(chǎn)量問題實驗?zāi)康摹?介紹無約束最優(yōu)化方法的一些基本概念。2了解幾種常見的無約束優(yōu)化問題的求解方法，如迭代算法、最速下降法(梯度法)、牛頓法(Newt on)、擬牛頓法。3學(xué)習(xí)掌握用MATLAB優(yōu)化工具箱中的命令來求解無約束優(yōu)化問題。實驗內(nèi)容】某公司生產(chǎn)一種產(chǎn)品有甲、 乙兩個品牌， 試討論產(chǎn)銷平衡下的最大利潤。 所謂產(chǎn)銷平衡 指公司的產(chǎn)量等于市場上的銷量。 利潤既取決于銷量和單件價格， 也依賴于產(chǎn)量和單件成本。按照市場規(guī)律，甲種品牌的價格 pi固然會隨其銷量Xi的增長而降低；同時乙品牌銷量X2的 增長也會使甲的價格有稍微下降，根據(jù)實際情況，可以確定價格與銷量成線性關(guān)系，即p1

2、 = 300 2.35 x1 0.09 x2乙的價格卩2遵循同樣的規(guī)律，有p2 = 480 0.14 x1 2.98 x2甲品牌的成本會隨著其產(chǎn)量的增長而降低，按實際情況可假設(shè)為負(fù)指數(shù)關(guān)系，即有q1 = 38 e 0.023X1 + 1 16乙品牌的成本q2遵循同樣的規(guī)律，有0.018 x2q2 = 94 e2 + 145試確定甲、乙兩種品牌的產(chǎn)量，使公司獲得的總利潤最大?！緦嶒灉?zhǔn)備】 無約束最優(yōu)化方法是指在沒有約束條件限制下， 求多變量實值函數(shù)極值的方法。 無約束 最優(yōu)化問題的數(shù)學(xué)表達(dá)式為min f (x) , x =( x1, x2，Xn) Rn(1)一般f為非線性函數(shù)，x是n維實變量，實

3、際上這是一個多元函數(shù)無條件極值 問題。由于一個求極大值問題， 可以添加負(fù)號的方式轉(zhuǎn)化為求極小值問題， 因此通常只討論求極小 值問題。應(yīng)該注意的是，極值問題的解，即極值點，都是 局部最優(yōu)解 ，全局最優(yōu)解 只能從局 部最優(yōu)解的比較中得到。如何求解無約束最優(yōu)化問題( 1)的最優(yōu)解呢？一般是采用迭代方法，即先選擇一個初 始點， 再尋找該點處的下降方向，我們稱為搜索方向。在該方向上求極小點，得到一個新的 點。這個新的點要優(yōu)于原來的點， 即新點處的目標(biāo)函數(shù)值小于原來點處的目標(biāo)函數(shù)值。然后在新點處再尋找下降方向和該方向上求極小點，如此下去，最終得到最優(yōu)點。因此， 求解無約束最優(yōu)化問題需解決兩個問題： 一是在

4、某些方向上的一維極小點，我們也稱為一維搜索； 另一個是尋找某些點處的下降方向， 這是無約束最優(yōu)化方法中最重要的一 個問題。我們先了解第一個問題最常用的搜索方法。1. 求一維極小的二點二次插值方法設(shè)d(k)是點x(k)處的一個搜索方向，要在該方向上尋優(yōu)問題，轉(zhuǎn)化為求一維函數(shù)()=f ( x(k) + d(k)(2)的求極值問題。最常用的一維搜索方法是插值法，即用某些點的函數(shù)值(或?qū)?shù)值)構(gòu)造插值函數(shù)，用插值函數(shù)的極小點來近似函數(shù)()的極小點。這里介紹一種有效的插值方法，稱為二點二次插值方法，即用二點處的函數(shù)值和一個點處的導(dǎo)數(shù)值構(gòu)造二次函數(shù)，反復(fù)用二次函數(shù)的極小點來逼近函數(shù)()的極小點。算法1二點

5、二次插值法(1) 取初始點 1 (如(如 =0.1 )，置精度要求(2)則置(3)(4)計算如果(5)計算(6)，則停止計算(1 = 0),初始步長 (如，并計算1 =( 1),;否則置 =2 = ( 2 )2 1 ),則置=21( 2 J'1 ( 2 1)作為極小點)；和步長縮減因子(1)，轉(zhuǎn)(3)否則置(0, 1)= ，轉(zhuǎn)2.最速下降法前面介紹了一維搜索的二點二次插值方法, 來看看兩個概念。2)F面討論如何選擇搜索方向的問題,我們先定義1稱n維向量(X1X2)T為函數(shù)f(x)在x處的梯度，記為 f(x)Xn f (X)=(X1定義2設(shè)d是任意的單位向量，若極限 limn 0Xnf

6、(X ")一空存在，則稱該極限為(3)函數(shù)f (X)在X處沿方向d的一階方向?qū)?shù)，簡稱為方向?qū)?shù)，記為f(X)，即d_f(x) = lim f(X ad) f(X)(4dn 0a最速下降法的基本思想：選取一點x1作為初始點，計算該點的梯度f(x)，求該點處的最速下降方向，即令 d1 = f (X1)，再沿d1方向前進(jìn)，尋找該方向上的極小點，得 到點x2，再計算 f(x2)，令d2 = f(x2)，沿d2方向前進(jìn)，得到點x3，如此下去 具體算法如下：算法2最速下降法(1) 取初始點x1，置精度要求為，置k = 1(2) 若| f (xk) | v ，則停止計算(xk為無約束問題的最優(yōu)解

7、)；否則置dk = f (xk)(3) 一維搜索，求一維問題min (a) = f ( xk ad k)ak 1kk得 ak，置 x = x ad 。(4) 置k = k +1，轉(zhuǎn)步驟(2)在算法中，為精度要求，即當(dāng)梯度接近于0時，我們就認(rèn)為達(dá)到極小點，終止計算。這樣做的目的是避免算法產(chǎn)生死循環(huán)，算法中的一維搜索可用算法1來求解。最速下降法是一種最基本的算法，它在最優(yōu)化方法中占有重要地位。其優(yōu)點是工作量小， 存儲變量較少，初始點要求不高；缺點是收斂慢，效率不高，有時達(dá)不到最優(yōu)解。下面介紹一種簡單而直觀的方法 Newton法。3. 牛頓法22 f將f (x)的二階導(dǎo)數(shù)記作f(X)=()是一個n

8、n矩陣，稱黑塞(Hessian )XiXj矩陣(簡記H(X)。當(dāng)f (x)二階偏導(dǎo)數(shù)連續(xù)時，H(X)是對稱矩陣。Newton法的基本思想是用一個二次函數(shù)去近似目標(biāo)函數(shù)f(x)，然后精確地求出這個二次函數(shù)的極小點，并把這個極小點作為所示函數(shù)的極小點x的近似值。k 1*k設(shè)x 是f(x)極小點x的第k 1次近似，將f (X)在x點作二階Taylor展開，略去 高于二階的項，用dk = xk 1 xk得到f (xk 1) = f (xk dk)f(xk) + f (xk) dk + 1 d 2 f (xk) d k(5)為了求dk使f(xk dk)最小，只需將(5)式右端對dk求導(dǎo)，并令其為零，可得

9、 f (xk) + 2 f(xk) d k = 0(6)稱為牛頓方程。當(dāng)黑塞矩陣2 f (xk)正定時，可解得dk 一 ( 2f(xk)-1 f(xk)( 7)稱為牛頓方向。其迭代公式為k 1 k2 k -1kx = x - ( f(x ) f(x )(8)4. 擬牛頓法牛頓法的優(yōu)點是在接近最優(yōu)解時具有2階收斂性，缺點是計算復(fù)雜，而且會出現(xiàn)病態(tài)、不正定等情況。為克服這些問題，同時又保持較快的收斂，禾U用第k和第k 1步得到的xk、k 1kk 1k 12kx 、 f(x ) . f (x )根據(jù)擬牛頓條件構(gòu)造一正定矩陣G近似代替f(x )或H k 1近似代替(2 f (xk 1)-1，得到的迭代

10、算法稱為擬牛頓法。其中有兩個常用的著名的迭代公式為BFGS公式和DFP公式，可參見其他最優(yōu)化方法書籍或MATLAB幫助。5. MATLAB求解無約束優(yōu)化問題的主要命令(1) MATLAB5.2及以下版本使用的命令x = fmin( ' fun ' , x1 , x2 )求一元函數(shù) y = f( x )在x1 , x2 內(nèi)的極小值x = fmin( ' fun ' , x1 , x2 , optio ns )同上，參數(shù) optio ns 的定義由表 1給出x , opti ons = fmin( .)同上，同時返回參數(shù)optio ns 的值fmin函數(shù)采用黃金分割

11、法和拋物線插值法，fun可直接用函數(shù)表達(dá)式表示，也可以是用M文件定義的函數(shù)名。建立函數(shù)形式的M-文件格式如下，文件名fun.m (般M函數(shù)文件名取函數(shù)名)function f = fun( x )f = f( x )x = fmins( ' fun ' , x0 )求多元函數(shù)(1)以x0為迭代初值的局部極小值x = fmins( ' fun ' , x0 , options )同上，參數(shù) options 的定義由表 1 給出x , optio ns = fmin s(.)同上，同時返回參數(shù) optio ns 的值fmins函數(shù)采用Nelder - Meade單純

12、形搜索法，fun可直接用函數(shù)表達(dá)式表示，也 可以是用M文件定義的函數(shù)名。x = fminu( ' fun ' , xO )x = fminu( ' fun ' , x0 , opti ons )x , opti ons = fminu( ' fun ' , x0 ,.)的值fminu函數(shù)為無約束優(yōu)化提供了三種算法，由 提供了二種算法，由optio ns (7)控制。這里,求多元函數(shù)(1)以x0為迭代初值的局部極小值同上，參數(shù)options 的定義由表1給出同時返回參數(shù)optio ns同上,opti ons (6) fminu必須先用控制，為步長一

13、維搜索 M-文件定義函數(shù) fun。(2) MATLAB5.3以上版本使用的命令fminbnd替代fmin求解一元函數(shù)極值，使用格式、搜索算法與之相同，x , Fval =fminbnd( ' fun ' , x0 ,.)同時返回解x處的函數(shù)值，而不是參數(shù)options的返回值fminsearch替代fmins求解多元函數(shù)(1)極值，使用格式、搜索算法與之相同fminunc替代fminu求解多元函數(shù)(1)極值，使用格式、搜索算法與之相同有關(guān)上述命令的詳細(xì)信息和使用方式可在幫助文件中了解，或在命令框里輸入help“命令名”查閱。在MATLAB5.3以上的各種最新版本中fmin、fm

14、ins和fminu命令仍然有用，但在MATLAB勺未來版本將可能刪除這些命令。(3) 命令中參數(shù)optio ns的有關(guān)定義在大多數(shù)MATLAB化命令函數(shù)中有一個控制參數(shù)options，它是一個有18個分量的向量，包含了在優(yōu)化程序中要用到的參數(shù)，以便在計算最優(yōu)值時控制精度要求、輸出形式、搜索算法、迭代次數(shù)、步長等等是。在對優(yōu)化命令函數(shù)的第一次調(diào)用時，向量options會自動使用缺省值；當(dāng)然在調(diào)用前，可以對optio ns的某些分量進(jìn)行賦值，以達(dá)到控制要求。表1 參數(shù)options各分量說明序號功能定義分量說明缺省值含義1輸出形式options(1)=1 ，有中間結(jié)果輸出 options(1)=

15、1,給出警告信息0，無中間結(jié)果輸出2解x(k)的精度(k)options(2)設(shè)置x 的精度|x(k1)(k) xx(k)|< 1e 43函數(shù)f (k)的精度options(3)設(shè)置f (k)的精度1)f (k)f (k)-< 1e 44函數(shù)g(k)的精度options(4)設(shè)置g (k)的精度，由 命令 constr、minmax使用終止判 據(jù)g(k1) g(k)-W 彳a 7(k) g5主要算法options(5)=1 ，高斯-牛頓法0，LM方法6搜索方向算法options(6)=1 ,擬牛頓法DFP公式 options(6)=2，最速下降法0，擬牛頓法的BFGS公式7步長一維

16、搜索算法options=1，三次多項式插值0,混合的二次和三次多項式8函數(shù)值輸岀options©)輸出極值點處的函數(shù) 值9梯度檢查options(9)=1，最初幾步，梯度將與有限微分計算的結(jié)果比較0,不作檢查10函數(shù)計算次數(shù)option(10)輸出函數(shù)計算次數(shù)11梯度計算次數(shù)option(11)輸出梯度計算次數(shù)12約束梯度計算次數(shù)option(12)輸出約束梯度計算次 數(shù)13等式約束的個數(shù)option(13)確定等式約束數(shù)目0，等式約束個數(shù)為014最大迭代次數(shù)option(14)輸入迭代的最大次數(shù)100 n，n為自變量個數(shù)200 n (fmins )500 n (fmin )15目標(biāo)

17、優(yōu)化由 attgoal 命令使用，option(15) 輸入須達(dá)目的的目標(biāo)個數(shù)016差分步長最小值用差分計算梯度時步長的下限1e 817差分步長最大值用差分計算梯度時步長的上限0.118步長參數(shù)，一般取1或更小【實驗方法與步驟】1. MATLAB化工具箱中解無約束問題的命令和參數(shù)options的基本用法下面用例題來予以說明例 1 求解 min f(x) = ex 5 x , 1< x< 2輸入命令：>> x=fmin( 'exp(x) - 5*x' , 1 , 2 ) , f=exp( x ) - 5*xx =1.6094f =-3.0472同時，該問題

18、可以用 fminbnd求解，得到的解答一樣輸入命令：>> x , fval=fmi nbn d( 'exp( x ) - 5*x' , 1 , 2 )x =1.6094fval =-3.0472例2用擬牛頓法的DFP公式求解min f (x) = 4x" + x； x3x2，初值為(1, 2)，要 求精度為1e 7,給出函數(shù)計算次數(shù)以及函數(shù)值首先建立M-文件，文件名取函數(shù)名fun.mfunction f = fun( x )f = 4 * x( 1 )A2 + x( 2 )A2 - x( 1 )A3 * x( 2 )輸入命令：>> x0=1,2

19、;>> opt( 2 )=1e-6;%設(shè)置自變量要求的精度>> opt( 3 )=1e-6;%設(shè)置函數(shù)所要求的精度>> opt( 6 )=1;%搜索算法用擬牛頓法的DFP公式求解>> x, opt=fminu( ' fun ', x0, opt);%返回參數(shù) options 的值給向量 opt可得到結(jié)果x =1.0e-008 *-0.34930.7566參數(shù)options的值返回給向量opt，全部元素略。>> n=opt( 10 );%函數(shù)計算次數(shù)賦值給nn =35>> y=opt( 8 );%在極值點處函

20、數(shù)值賦值給yy =1.0606e-0162. 引例問題的分析與求解由題設(shè)可知，甲品牌產(chǎn)品單件獲利為p1 q1，乙品牌產(chǎn)品單件獲利為P2 q2，由產(chǎn) 銷平衡原理，所有產(chǎn)品的銷量即為產(chǎn)量，則甲、乙兩種產(chǎn)品總獲利為Z(X1 ,X2)=( Pi qi) Xi +( P2 q2) x?容易看出，原問題實際上轉(zhuǎn)化為求二元函數(shù)z(x1 ,X2)的極大值，為用MATLAB優(yōu)化工具箱中的fminunc求解，需將其轉(zhuǎn)化為求函數(shù)y(x1, x2) = z(x1, x2)的最小值。為確定初始值，先忽略成本，并令p1價格中x2項的較小系數(shù)0.09和p2中x1項較小的 系數(shù)0.14等于零(因為它們對價格的作用比較微弱，暫時可忽略不計)，則確定初值問題轉(zhuǎn)化為求Z(Xi,X2)=( 300 2.35 Xi) Xi +( 480 2.98 X?) X：的極值，很容易可以求得xi =300 = 63.83 , x2 = 480= 80.54，我們用它作為原2 2.3522.98問