(5.3.10)-2.5 多維非線性規(guī)劃

腳本——多維非線性規(guī)劃(ppt1，ppt2)同學(xué)們，你們好，這節(jié)課我們來學(xué)習(xí)數(shù)學(xué)建模中的多維非線性規(guī)劃問題。(ppt3)今天我們講的主要是無約束多維極值問題的最優(yōu)化方法，首先，我們一起了解這類問題的背景分析。(ppt4)（動畫1）隨著社會的發(fā)展，最優(yōu)化方法受到普遍的關(guān)注，廣泛的應(yīng)用于物流運(yùn)輸網(wǎng)絡(luò)、公司企業(yè)經(jīng)營管理等領(lǐng)域，利用最優(yōu)化方法解決問題，以追求效益和利潤最大化目標(biāo)，以最大程度滿足自身期望。（動畫2）左圖是全國肯德基門店的數(shù)目分布示意圖，其中主要調(diào)研了每個地區(qū)消費(fèi)者購買情況和人口分布，進(jìn)而合理規(guī)劃門店的數(shù)量。（動畫3）右圖展現(xiàn)的是電商貨物倉庫中心，也涉及了優(yōu)化問題模型的求解，目的是為了設(shè)計(jì)最優(yōu)的解決方案。由于上述問題考慮了諸多因素，這些因素很復(fù)雜，因此蘊(yùn)含了多個變量之間的非線性關(guān)系。這就需要用到我們今天所講的知識，求解多維非線性規(guī)劃問題模型的優(yōu)化方法。(ppt5)首先，我們將介紹相關(guān)的數(shù)學(xué)基礎(chǔ)知識。(ppt6)（動畫1）凸集的定義為：設(shè)S為n維歐氏空間R^n中的一個集合。若對任意兩點(diǎn)x^((1))，x^((2))∈S及每個實(shí)數(shù)lammda∈[0,1]，有l(wèi)ammda*x^((1))+(1-lammda)*x^((2))∈S則稱S為凸集。lammda*x^((1))+(1-lammda)*x^((2))稱為x^((1))和x^((2))的凸組合。（動畫2）如左圖是凸集示意圖，形象地說就是在區(qū)域中找兩點(diǎn)，連接這兩點(diǎn)的線段還是屬于這個區(qū)域。（動畫3）右圖是非凸集的示意圖。在區(qū)域中找兩點(diǎn)，連接這兩點(diǎn)的線段，或者線段的部分不屬于這個區(qū)域。(ppt7)下面我們介紹向量函數(shù)的一階偏導(dǎo)的定義。（動畫1）先來看如何建立模型）先來看定義為：設(shè)f：R^n→R^m是連續(xù)函數(shù)，對于x_0∈R^n和R^n空間的標(biāo)準(zhǔn)基{e_1,e_2,…,e_n}，則函數(shù)f(x_0+t*e_j)在t=0的導(dǎo)數(shù)(若存在的話)稱為f在點(diǎn)x_0關(guān)于x_j的一階偏導(dǎo)數(shù)(j=1,2,…,n)。（動畫2）因此，如果f(x)=[(f_1(x)，…，f_m(x))]^T，則有?f/(?x_j)(x_0)=[((?f_1)/(?x_j)(x_0)，…，(?f_m)/(?x_j)(x_0))]^T。(ppt8)（動畫1，2）雅可比矩陣的定義為：向量函數(shù)f=[f_1,f_2,…,f_n]^T：R^n→R^m，在x_0的導(dǎo)數(shù)矩陣記為L∈R^(m×n)，則稱該矩陣為向量函數(shù)f在x_0點(diǎn)的Jacobi矩陣，記為Df(x_0)，具體形式如下：看背板(ppt9)梯度的定義為：（動畫1）如果f：R^n→R是可微的，注意這里f是標(biāo)量函數(shù)，那么函數(shù)（動畫2）?f(x)=[(?f/(?x_1)(x)，…，?f/(?x_n)(x))]^T=Df(x)^T，（動畫3）稱為f的梯度。梯度是從一個從R^n到R^n的函數(shù)，如果繪制梯度向量?f(x)，其起點(diǎn)為點(diǎn)x，箭頭表示方向，如此可將梯度表示為向量場。(ppt10)黑塞矩陣的定義為：（動畫1）當(dāng)f(x)在x點(diǎn)存在二階偏導(dǎo)時，函數(shù)f在點(diǎn)x的黑塞（Hessen）矩陣定義為梯度算子作用在f(x)的梯度上。具體如下所示。當(dāng)二階偏導(dǎo)連續(xù)時，它是一個對稱矩陣。(ppt11)進(jìn)一步，我們介紹兩個基本定理。（動畫1）定理1設(shè)f(x)具有連續(xù)的一階偏導(dǎo)，x^*為f(x)的局部極值，則（動畫2）?f(x^*)=Df(x^*)^T?((?f(x^?))/(?x_1),(?f(x^?))/(?x_2),(?f(x^))/(?x_n))^T=0，即函數(shù)f(x)的梯度等于零。(動畫3)定理2設(shè)f(x)在x∈R^n中有連續(xù)的二階偏導(dǎo)數(shù)，且梯度?f(x^?)等于零，黑塞（Hessen）矩陣?^2f(x^?)正定，則x^*為f(x)的嚴(yán)格局部極小值點(diǎn)。(ppt12)（動畫1）下面我引入多維極值問題模型(ppt13)我們考慮的問題類為：（動畫1）無約束多維極值問題的一般表達(dá)式為：minf(x),x∈R^n，其中，x為向量，f(x)為標(biāo)量函數(shù)。（動畫2）下面在介紹最優(yōu)化方法求解問題時，將要用到問題模型的導(dǎo)數(shù)(即多元函數(shù)的梯度)，有些算法甚至需要多元函數(shù)的雅可比矩陣。因此，這一類問題在搜索過程中需要用到的梯度，故稱為梯度方法。(ppt14)首先，我們介紹的是最速下降法。（動畫1）一個實(shí)值可微函數(shù)在某點(diǎn)處函數(shù)值增加最快的方向，正交于經(jīng)過該點(diǎn)的函數(shù)水平集，即在梯度方向上，自變量的細(xì)微變動所導(dǎo)致的目標(biāo)函數(shù)值的增加幅度要超過其他任意方向。（動畫2）因此，我們可知梯度方向就是函數(shù)增加最快的方向。反之，梯度負(fù)方向就是減小最快的方向。如圖所示?f(x_0)方向?yàn)閤_0梯度方向。(ppt15)最速下降法的主要思想為：（動畫1）第一，最速下降法的搜索方向是目標(biāo)函數(shù)的負(fù)梯度方向，目的是搜索到目標(biāo)函數(shù)的最低點(diǎn)。第二，每次迭代，從迭代點(diǎn)x^((k))出發(fā)，沿著梯度負(fù)方向-?f(x^((k)))開展一維搜索，直到找到步長的最優(yōu)值，確定新的迭代點(diǎn)x^((k+1))。第三,最速下降法的相鄰搜索方向是正交的。（動畫2）如圖所示，最速下降法搜索的方向，在點(diǎn)x_0處第一步搜索的方向是①，走的步長太小，效率低，太長，走過頭，可能函數(shù)值就不一定是減少了，從示意圖可以看出，走到x1這個位置剛剛好。（動畫3）第二步搜索的方向如②所示（動畫4）第三步搜索的方向如③所示，依次類推，我們就可找到極小值點(diǎn)。(ppt16)下面，具體學(xué)習(xí)最速下降法的求解步驟：（看背板）（動畫1）①給定初始點(diǎn)x^((0))和允許誤差epsilon>0，令k=0。（動畫2）②計(jì)算?f(x^((k)))，若‖?f(x^((k)))‖≤epsilon，即梯度接近于零時，則停止迭代，極小點(diǎn)就是x^?=x^((k))；否則，轉(zhuǎn)入③。（動畫3）③求min_(lammda≥0)f(x^((k))-lammda?f(x^((k))))，得到極小值lammda_k。（動畫4）④令x^((k+1))=x^((k))-lammda_k?f(x^((k)))，用k+1代替k轉(zhuǎn)入②，繼續(xù)迭代。(ppt17)我們介紹的第二個方法是牛頓法，其主要思想為（動畫1）利用二次型函數(shù)來局部近似目標(biāo)函數(shù)f，并求解近似函數(shù)的極小值點(diǎn)作為下一個迭代點(diǎn)，并且將-[?^2f(x^((k)))]^(-1)?f(x^((k)))作為搜索方向，迭代公式為（看背板）x^((k+1))=x^((k))-[?^2f(x^((k)))]^(-1)?f(x^((k))).（動畫2）算法的基本步驟如下：（看背板）（動畫3）①給定初始點(diǎn)x^((0))和精度epsilon>0，令k=0。（動畫4）②若‖?f(x^((k)))‖<epsilon，則停止迭代，最小點(diǎn)就是x^*≈x^((k)否則轉(zhuǎn)③。（動畫5）③計(jì)算[?^2f(x^((k)))]^(-1)，令p^((k))=-[?^2f(x^((k)))]^(-1)?f(x^((k)))。（動畫6）④令x^((k+1))=x^((k))+p^((k))，令k=k+1，轉(zhuǎn)②。(ppt18)下面我們總結(jié)牛頓法優(yōu)缺點(diǎn)（動畫1，動畫2）優(yōu)點(diǎn)為：第一，針對利用牛頓法求解無約束多維極值問題，如果目標(biāo)函數(shù)是正定的二次函數(shù)，則用牛頓法經(jīng)過一次迭代就可以達(dá)到最優(yōu)點(diǎn)。但我們引入的最速下降法做不到。（動畫3）第二，如果初值點(diǎn)離極值點(diǎn)近時，有好的收斂性。最速下降法只有一階收斂率，而牛頓法具有二階收斂率，收斂速度更快。（動畫4，動畫5）缺點(diǎn)為：第一，對優(yōu)化問題的條件加強(qiáng)了，并且黑塞矩陣非正定時，牛頓法確定的搜索方向不一定是目標(biāo)函數(shù)值下降的方向。（動畫6）第二，某些情況下，即使黑塞矩陣正定，牛頓法也不具下降特征，比如初始值遠(yuǎn)離目標(biāo)函數(shù)的極值的時候。(ppt19)針對牛頓法的優(yōu)缺點(diǎn)，我們下面提出了以下兩種改進(jìn)的方法：（動畫1,2,3）改進(jìn)方法一的具體形式為：（看背板）類似我們引入搜索步長。x^((k+1))=x^((k))-alpha_k[?^2f(x^((k)))]^(-1)?f(x^((k)))，其中alpha_k=argmin_{(alpha_k≥0)}f(x^((k))-alpha_k[?^2f(x^((k)))]^(-1)?f(x^((k))))。（動畫4,5）改進(jìn)的方法二的具體形式為x^((k+1))=x^((k))-[?^2f(x^((k)))+mu_kI]^(-1)?f(x^((k)))（動畫6）mu_k≥0，當(dāng)mu_k足夠大時，保證修正的黑塞矩陣是正定的，一定能夠保證搜索是下降方法；mu_k→0,該方法趨近于牛頓法；mu_k→+∞,該方法趨近于較小步長牛頓法是的特性。(ppt20)（動畫1）最速下降法：計(jì)算簡單，但收斂慢。牛頓法：收斂快，但計(jì)算和存儲大。（動畫2,3,4）（圖）因此，我們想著找出一種算法，同時具備這兩種方法的優(yōu)點(diǎn)，那么就是我們學(xué)習(xí)的第三種優(yōu)化方法——共軛梯度法。（動畫5）因此，我們提出了共軛梯度法，它的收斂速度介于最速下降法與牛頓法之間，這種算法能夠具有牛頓法收斂速度快的優(yōu)點(diǎn)，又有最速下降法計(jì)算簡單的優(yōu)點(diǎn)。(ppt21)什么是共軛方向呢？其中，共軛方向的定義如下，（動畫1）定義：Q為n×n的對稱實(shí)矩陣，對于方向d^((1)),d^((2)),…,d^((m)),對所有的i≠j,有d^(i)TQd^((j))=0，則它們是關(guān)于Q共軛的。（動畫2）共軛梯度法的主要思想是：（動畫3）第一，共軛梯度法不需要預(yù)先給定Q共軛方向，而是隨著迭代的進(jìn)行不斷的產(chǎn)生Q共軛方向，搜索步長由一維極值算法確定。（動畫4）第二，每一次的迭代中，利用上一個搜索方向和目標(biāo)函數(shù)在當(dāng)前的迭代點(diǎn)的梯度向量之間的線性組合構(gòu)造一個新方向，與前面已經(jīng)產(chǎn)生的搜索方向組成Q的共軛方向。(ppt22)以求解二次型問題為例總結(jié)共軛梯度法基本的算法步驟為：（看背板）（動畫1）①給定初始點(diǎn)x^((0))和精度ε>0；計(jì)算?f(x^((0))),若‖?f(x^((0)))‖≤epsilon，則停止迭代，最小點(diǎn)就是x^?≈x^((0))；否則，轉(zhuǎn)入②。（動畫2）②取p^((0))=-?f(x^((0)))，令k=0。（動畫3）③用一維搜索法求t_k，使得min_{t≥0}f(x^((k))+tp^((k)))，令x^((k+1))=x^((k))+t_k*p^((k)),轉(zhuǎn)④。（動畫4）④若‖?f(x^((k+1)))‖≤epsilon,停止迭代，最小點(diǎn)就是x^?≈x^((k+1))，否則轉(zhuǎn)⑤。（動畫5）⑤若k+1=n，令x^((0))=x^((n)),轉(zhuǎn)②，否則轉(zhuǎn)⑥。（動畫6）⑥令p^((k+1))=-?f(x^((k+1)))+λ_k*p^((k)),λ_k=‖?f(x^((k+1)))‖^2/‖?f(x^((k)))‖^2，令k=k+1，轉(zhuǎn)③。(ppt23)（動畫1）最速下降法和共軛梯度法的區(qū)別：（動畫2）1.最速下降法和共軛梯度法都是由線性方程組的解法發(fā)展而來，但共軛梯度法具有二次終止性，程序編制也簡單，實(shí)際應(yīng)用中很廣泛。（動畫3）2.一般情況下，共軛梯度法的性能優(yōu)于最速下降法，但比牛頓法差。（動畫4）3.共軛梯度法不需要存儲n×n的矩陣，也不需要對矩陣求逆。(ppt24)下面我們將這些方法應(yīng)到具體的案例中(ppt25)（動畫1）目標(biāo)函數(shù)為二次型函數(shù)：（看背板）f(x_1,x_2,x_3)=3/2*x_1^2+2*x_2^2+3/2*x_3^2+x_1*x_3+2*x_2*x_3-3*x_1-x_3，利用共軛梯度法求解極小值點(diǎn)，初始點(diǎn)為x^((0))=[0,0,0]^T。（動畫2）函數(shù)f可以寫為f(x)=1/2x^TQx-x^Tb,（動畫3）其中，矩陣Q和b矩陣具體形式是(ppt26)下面將展示詳細(xì)的求解過程（看背板）（動畫1）函數(shù)f(x)的梯度為：g(x)=?f(x)=Qx-b=[x_1+x_3-3,4*x_2+2*x_3,x_1+2*x_2+3*x_3-1)]（動畫2）因此，g^((0))=[-3,0,-1]^T,d^((0))=-g^((0)),α_0=-(g^((0))^T*d^((0)))/(d^(0)T*Q*d^((0)))=10/36=0.2778,（動畫3）新的迭代點(diǎn)為：x^((1))=x^((0))+α_0*d^((0))=[0.8333,0,0.2778]^T，可得g^((1))=?f(x^((1)))=[-0.2222,0.5556,0.6667]^T，β_0=(g^(1)T*Q*d^((0)))/(d^(0)T*Q*d^((0)))=0.08025,（動畫4）對應(yīng)的搜索方向?yàn)椋篸^((1))=-g^((1))+β_0*d^((0))=[0.4630,-0.5556,-0.5864]^T,步長為：α_1=-(g^(1)T*d^((1)))/(d^(1)T*Q*d^((1))