1、學 院：數(shù)學與統(tǒng)計學院班 級：碩2041班姓 名：王彭學 號：指導教師：阮小娥同 組 人：錢東東最優(yōu)化方法課程實驗報告求解二次規(guī)劃問題的拉格朗日及有效集方法求解二次規(guī)劃問題的拉格朗日及有效集方法摘要二次規(guī)劃師非線性優(yōu)化中的一種特殊情形，它的目標函數(shù)是二次實函數(shù)，約束函數(shù)都是線性函數(shù)。由于二次規(guī)劃比較簡單，便于求解（僅次于線性規(guī)劃），并且一些非線性優(yōu)化問題可以轉(zhuǎn)化為求解一些列的二次規(guī)劃問題，因此二次規(guī)劃的求解方法較早引起人們的重視，稱為求解非線性優(yōu)化的一個重要途徑。二次規(guī)劃的算法較多，本文僅介紹求解等式約束凸二尺規(guī)劃的拉格朗日方法以及求解一般約束凸二次規(guī)劃的有效集方法。關(guān)鍵字：二次規(guī)劃，拉格朗日

2、方法，有效集方法?！灸夸洝空? 1 -1 等式約束凸二次規(guī)劃的解法- 3 -1.1 問題描述- 3 -1.2 拉格朗日方法求解等式約束二次規(guī)劃問題- 3 -1.2.1 拉格朗日方法的推導- 3 -1.2.2 拉格朗日方法的應用- 4 -2 一般凸二次規(guī)劃問題的解法- 5 -2.1 問題描述- 5 -2.2 有效集法求解一般凸二次規(guī)劃問題- 6 -2.2.1 有效集方法的理論推導- 6 -2.2.2 有效集方法的算法步驟- 9 -2.2.3 有效集方法的應用- 10 -3 總結(jié)與體會- 11 -4 附錄- 11 -4.1 拉格朗日方法的matlab程序- 11 -4.2 有效集方法的Matla

3、b程序- 11 -1 等式約束凸二次規(guī)劃的解法1.1 問題描述我們考慮如下的二次規(guī)劃問題(1.1)其中對稱正定，行滿秩，。1.2 拉格朗日方法求解等式約束二次規(guī)劃問題1.2.1 拉格朗日方法的推導首先寫出拉格朗日函數(shù)：，(1.2)令，得到方程組將上述方程組寫成分塊矩陣形式：(1.3)我們稱傷處方程組的系數(shù)矩陣為拉格朗日矩陣。下面的定理給出了線性方程組(1.1)有唯一解的充分條件。定理1 設對稱正定，行滿秩。若在問題(1.1)的解處滿足二階充分條件，即則線性方程組(1.4)的系數(shù)矩陣非奇異，即方程組（1.4）有唯一解。其中，方程組(1.4)為（1.1）對應的齊次方程組：（1.4）.下面，我們來推

4、導方程(1.3)的求解公式。根據(jù)定理1，拉格朗日矩陣必然是非奇異的，故可設其逆為.由恒等式可得于是由上述四個等式得到矩陣的表達式因此，由(1.3)可得解得表達式其中，分別由(1.5),(1.6),(1.7)給出。下面給出和的另一種等價表達式。設是問題(1.1)的任一可行點，即滿足。而在此點處目標函數(shù)的梯度為，利用和，可將(1.8)改寫為1.2.2 拉格朗日方法的應用（1）拉格朗日方法的Matlab程序見附錄。（2）利用拉格朗日方法求解下列問題：解 容易寫出利用Matlab程序求解該問題可以結(jié)果如下：2 一般凸二次規(guī)劃問題的解法2.1 問題描述考慮一般二次規(guī)劃其中是階對稱陣。記，下面的定理給出了

5、問題(2.1)的一個最優(yōu)性充要條件。定理2 是二次規(guī)劃問題(2.1)的局部極小點當且僅當（1）存在，使得(2)記則對于任意的，均有.容易發(fā)現(xiàn)，問題(2.1)是凸二次規(guī)劃的充要條件是半正定。此時，定理2的第二部分自然滿足。注意到凸優(yōu)化問題的局部極小點也是全局極小點的性質(zhì)，我們有下面的定理：定理3 是凸二次規(guī)劃的全局極小點的充要條件是滿足條件，即存在，使得2.2 有效集法求解一般凸二次規(guī)劃問題2.2.1 有效集方法的理論推導首先引入下面的定理，它是有效集方法理論基礎。定理4 設是一般凸二次規(guī)劃問題的全局極小點，且在處的有效集為，則也是下列等式約束凸二次規(guī)劃的全局極小點。從上述定理可以發(fā)現(xiàn)，有效集方

6、法的最大難點是事先一般不知道有效集，因此只有想辦法構(gòu)造一個集合序列去逼近它，即從初始點出發(fā)，計算有效集，解對應的等式約束子問題。重復這一做法，得到有效集序列使之，以獲得原問題的最優(yōu)解?；谏鲜龆ɡ?，我們分4步來介紹有效集方法的算法原理和實施步驟。第一步，形成子問題并求出搜索方向.設是問題(2.1)的一個可行點，據(jù)此確定相應的有效集，其中求解相應的子問題上述問題等價于其中設求出問題(2.4)的全局極小點為是對應的拉格朗日乘子。第二步，進行線性搜索確定步長因子.假設，分兩種情形討論。(1) 若是問題(2.1)的可行點，即則令(2) 若不是問題(2.1)的可行點，則通過線性搜索求出下降最好的可行點。

7、注意到目標函數(shù)是凸二次函數(shù)，那么這一點應該在可行域的邊界上達到。因此只要求出滿足可行條件的最大步長即可。當時，對于任意的，都有和，此時，不受限制。當時，即第個約束是嚴格的不等式約束，此時要求滿足，即注意到上式右端非正，故當時，上式恒成立。而當時，由上式可解得故有合并第(1)(2)可得.第三步，修正.當時，有效集不變，即.而當時，故，因此在處增加了一個有效約束，即.第四步，考慮的情形。此時是問題(2.3)的全局極小點。若這是對應的不等式約束的拉格朗日乘子均為非負，則也是問題(2.1)的全局極小點，迭代終止。否則，如果對應的不等式約束的拉格朗日乘子有負的分量，那么需要重新尋找一個下降可行方向。設，

8、現(xiàn)在要求一個下降可行方向，滿足且，為簡便計算，按下述方式選取：即另一方面，注意到是子問題(2.3)的全局極小點，故有即其中從而，由(2.6)知于是有上式表明，由(2.6)確定的是一個下降可行方向。因此，令，則修正后的子問題的全局極小點必然是原問題的一個下降可行方向。2.2.2 有效集方法的算法步驟經(jīng)過上面的分析和推導，我們現(xiàn)在可寫出有效集方法的算法步驟：（1）選取初值。給定初始可行點，令.（2）解子問題。確定相應的有效集.求解子問題得極小點和拉格朗日乘子向量.若轉(zhuǎn)步驟（4）；否則，轉(zhuǎn)步驟（3）。（3）檢驗終止準則。計算拉格朗日乘子其中令若，則是全局極小點，停算。否則，若，則令，轉(zhuǎn)步驟（2）。（

9、4）確定步長.令，其中令（5）若，則令；否則，若，則令，其中滿足（6）令，轉(zhuǎn)步驟(1)。2.2.3 有效集方法的應用（1）有效集法的Matlab程序見附錄。（2）用有效集方法求解下列二次規(guī)劃問題：解 首先確定有關(guān)數(shù)據(jù)：利用Matlab計算可得結(jié)果如下：3 總結(jié)與體會通過本次實驗，筆者對求解等式約束凸二次規(guī)劃問題的拉格朗日方法及一般約束條件下凸二次規(guī)劃問題的有效集方法有了較深的認識和了解。感謝阮老師的悉心教誨和指導，使得筆者對最優(yōu)化課程中的理論推導、算法步驟及編程都比較熟悉，對以后的科研工作有很好的指導和借鑒意義。4 附錄4.1 拉格朗日方法的matlab程序(1) 拉格朗日算法函數(shù)%本程序用拉

10、格朗日方法求解燈飾約束條件的二次規(guī)劃問題。function x,lam,fval=qlag(H,A,b,c)%功能：用拉格朗日方法求解等式約束二次規(guī)劃：% min f(x)=0.5*xHx+cx， s.t. Ax=b%輸入：H,c分別是目標函數(shù)的矩陣和向量，A,b分別是約束條件中的矩陣和向量%輸出：(x,lam)是KT點，fval是最優(yōu)值IH=inv(H);AHA=A*IH*A;IAHA=inv(AHA);AIH=A*IH;G=IH-AIH*IAHA*AIH;B=IAHA*AIH;C=-IAHA;x=B*b-G*c;lam=B*c-C*b;fval=0.5*x*H*x+c*x; (2) 拉格朗

11、日算法求解等式約束凸二次規(guī)劃問題主程序：H=2 -2 0;-2 4 0;0 0 2;c=0 0 1;A=1 1 1;2 -1 1;b=4 2;x,lam,fval=qlag(H,A,b,c) 4.2 有效集方法的Matlab程序（1）有效集方法函數(shù)%本程序主要適用于求解一般約束條件下的凸二次規(guī)劃問題function x,lamk,exitflag,output=qpact(H,c,Ae,be,Ai,bi,x0)%功能：用有效集方法解一般約束二次規(guī)劃問題%min f(x)=0.5*x*H*x+c*x,% s.t. a_i*x-b_i=0,(i=1,l),% a_i*x-b_i=0,(i=l+1,

12、m)%輸入：x0是初始點，H,c分別是目標函數(shù)二次矩陣和向量；%Ae=（a_1,.,a_l）,be=(b_1,.,b_l);%Ai=(a_l+1,.,a_m),bi=(b_l+1,.,b_m).%輸出：x是最優(yōu)解，lambda是對應的乘子向量；output是結(jié)構(gòu)變量% 輸出極小值f(x),迭代次數(shù)k等信息，exitflag是算法終止類型%初始化epsilon=1.0e-9;err=1.0e-6;k=0;x=x0;n=length(x);kmax=1.0e3;ne=length(be);ni=length(bi);lamk=zeros(ne+ni,1);index=ones(ni,1);for

13、i=1:ni if Ai(i,:)*xbi(i)+epsilon index(i)=0; endend%算法主程序while k0 Aee=Ae; end for j=1:ni if index(j)0 Aee=Aee;Ai(j,:); end end gk=H*x+c; m1,n1=size(Aee); dk,lamk=qsubp(H,gk,Aee,zeros(m1,1); if norm(dk)ne y,jk=min(lamk(ne+1:length(lamk); end if y=0 exitflag=0; else exitflag=1; for i=1:ni if index(i)&(ne+sum(index(1:i)=jk index(i)=0; break; end end end k=k+1; else exitflag=1; %求步長 alpha=1.0;tm=1.0; for i=1:ni if (index(i)=0)&(Ai(i,:)*dk0) tm1=(bi(i)-Ai(i,:)*x)/(Ai(i,:)*dk); if tm1tm tm=tm1; ti=i; end end end alpha=min(alpha,tm); x=x+alpha*dk; %修正有效集 if tm0 rb=Ae*ginvH*c+be; lambda=