1.等式約束二次規劃問題

(1)

其中??且為對稱矩陣，,?，. 求解等式約束二次規劃問題一般有直接消元法，正交分解法和拉格朗日法。

拉格朗日法

構建Eq.(1)的拉格朗日乘子函數，得

由KKT條件可得：

寫成矩陣的形式得

系數矩陣稱為Lagrange矩陣，它是對稱不定的。如果Lagrange矩陣的逆存在

可得最優解為.如果存在，可得

如果Lagrange矩陣的逆不存在，可以求偽逆。

2.等式和不等式約束二次規劃問題

由于帶有等式約束的二次規劃問題已經有成熟的解法，所以要求解帶有不等式約束的二次規劃問題，一種很自然的想法就是將不等式約束轉換成等式約束。有效集法（Active Set Method）是求解帶不等式約束的二次規劃問題的一種經典算法。首先看下問題的幾何表述

• 可行域：目標問題的可行域是空間中的一個凸多面體（包括邊界和內部）。例如上圖中的綠色多邊形區域。
• 等值面：如果矩陣G正定，那么目標函數的等值面是空間中的一族同心超橢球面，所有橢球面相似且同向（各軸同向重合）。例如上圖中的同心橢圓曲線族。越內層的橢球面，目標函數在其上的取值越優（低）。
• 優化問題：在幾何表述下，我們的目標就是在可行域內找到達到同心橢圓曲線族最內層的點。例如上圖中，P是無約束最優解，Q是約束最優解。

顯然加入約束之后，最優解必定出現在某個或某幾個不等式約束的邊界上。

• 有效集：任給一個可行解，其有效集就是滿足等號成立的約束條件的集合。顯然，有效集必然包含所有等式約束，同時包含不等式約束的一個子集。注意這里是等號成立，也就是說，當前解讓一個不等式約束的等號成立時，這個不等式約束才會被納入有效集。
• 最優有效集：最優解的有效集稱為最優有效集。只需把最優有效集中的約束全部改寫為等式約束，并扔掉其它約束，然后用Lagrange乘子法直接求解即可。

由于約束條件數目是有限的，而最優有效集是其子集，進而最優有效集只有有限多種可能。于是很自然地就想到了下面大名鼎鼎的暴力破解法——窮舉法。

每次，我們從全體約束中選出一個子集（包含全部等式約束和部分不等式約束），稱之為本次嘗試的工作集，并求解該工作集對應的子優化問題。一個工作集對應的子優化問題是這樣定義的：在原問題中，將工作集中的約束全部改為等式約束，同時扔掉工作集之外的約束。由于這是一個等式約束問題，所以可以用Lagrange法輕松求解。之后檢查該解是否是原問題的可行解（即它是否也滿足工作集之外的約束）。如果是，則記錄下來，否則就扔掉。遍歷全部可能的子集之后，在記錄下來的所有解中，選一個使目標函數值最優的即可。

窮舉法雖然已經可以確保在有限步內獲得目標問題的解，但其運算量常常超乎想象（尤其是控制標量和不等式約束條件個數較多時），所以我們在以下兩方面尋求改進：

• 改進最優解的識別方式：利用凸二次規劃的特點，建立一組識別規則（實際上就是KKT條件），可直接判斷一個可行解是否是最優解，因此算法一旦試出最優解就立即停止，不必遍歷剩余可能性
• 改進嘗試各種可能性的順序：不再隨機選取下一次嘗試的工作集，而是通過求解子優化問題找到能使目標函數值有效下降的新工作集和迭代點。

經過以上兩項改進，就形成了效率大幅提高的有效集法。

有效集法

初始可行解：有效集法需要一個初始可行解??作為迭代的起點。由于可行解與目標函數無關，所以獲得初始可行解的方法與線性規劃完全一致，典型方法有兩階段法和大M法，，因此我們假設已經獲得一個初始可行解 x0 ，并得到有效集為W0

迭代求解。迭代的過程主要包括兩個部分，一個尋找當前有效集中的最優解，另一個是更新有效集以得到最優集。假設迭代k（k=1,2,3,...）步后，當前解為xk，工作集為wk,根據xk落在可行域上的位置，可以劃分出四種情形。

• 情形1：?若xk??滿足KKT條件，則它就是最優解，算法停止。
• 情形2：若 xk??在當前工作集（等式約束集）下已達最優，但工作集中存在至少一個約束，它原來是不等式約束，而且當 xk??朝該不等式約束的可行一側移動時，可以使目標函數值進一步下降，那么我們令迭代點不動（即 xk+1=xk并從工作集中去掉那個約束（即將它還原放松成不等式約束）。
• 情形3：若 xk在當前工作集（等式約束集）下還可以改進，但步長不能走滿（否則會走出原始問題的可行域），那么就取 xk?延改進方向恰走到原始可行域邊界處的位置為下一個迭代點 xk+1?，同時把碰到的邊界處對應的不等式約束加入工作集。情形3的例子請參見下圖中x3?
• 情形4：若 xk在當前工作集（等式約束集）下還可以改進，并且步長能夠走滿（即改進到最優仍未出原始可行域），那么就取這個最優點作為下一個迭代點 xk+1，并保持工作集不變。情形4的例子請參見下圖中 x1和 x4。

下面是使用有效集求解帶不等式約束優化問題的經典例子。

G = eye(2,2)*2; d = [-2 -5]'; A = [-1 2; 1 2; 1 -2; -1 0; 0 -1]; b = [2;6;2;0;0]; Aeq = [-1 1]; beq = [0.03]; MaxIterations = 100; Tolerance = 1.0e-6; x0 = [0;1]; workset0= [3;5]; x = active_set_method(x0,G,d, A, b, Aeq, beq, MaxIterations, Tolerance) % options = optimoptions(@quadprog,'Algorithm','active-set','MaxIterations',500); y = quadprog(G, d, A, b, Aeq, beq, [], [], x0, options)function x = active_set_method(x0,G,d, Ai, bi, Aeq, beq, MaxIterations, Tolerance) %min 1/2*x'Gx+d'x %st. A*x <= b % Aeq*x = bne = size(Aeq,1);ni = size(Ai,1);workset = zeros(ni,1);x = x0;%給定初始解for i = 1:ni%將滿足初始解的不等式都加入有效集if abs(Ai(i,:)*x0 -bi(i)) <= Toleranceworkset(i) = 1;endendfor k = 1:MaxIterations%迭代求解 a = [Aeq; Ai(workset == 1,:)];b_ = [beq; bi(workset == 1,:)];[x1, lambda]=kkt(G, d, a, b_);s = x1 - x;if norm(s,2) < Tolerance%得到的解是本次有效集中的最優解x = x1; min_lambda = 0;if length(lambda) > ne[min_lambda,index] = min(lambda(ne+1:end));endif min_lambda >= 0 %滿足KKT條件，得到最優解break;else%尋找新的有效集for i=1:niif workset(i)&&sum(workset(1:i))==indexworkset(i)=0;break;endend endelse%得到的解不是本次有效集中的最優解，繼續優化[alpha, index] = Alpha(workset, x, s, Ai, bi);x = x + alpha*s;if alpha <1workset(index) = 1;endendendfunction [s, lambda]=kkt(G, g, Ae, be)ws = size(Ae,1);kkt_A = [G, Ae'; Ae, zeros(ws,ws)];kkt_B = [-g;be];slambda = pinv(kkt_A)*kkt_B; %kkt_A不一定可逆，需要求偽逆dim = size(G,2);s = slambda(1:dim);lambda = slambda(dim+1:end);endfunction [alpha, index] = Alpha(workset, x, s, Ai, bi)outset = find(workset == 0);Atp = Ai(outset,:)*s;Atp_index = Atp>0;out_plus_ind = outset(Atp_index);alphatset = (bi(out_plus_ind) - Ai(out_plus_ind,:)*x)./Atp(Atp_index);[alpha, index] = min([alphatset;1]);if alpha < 1index = out_plus_ind(index);endend end

參考：

https://zhuanlan.zhihu.com/p/29525367

https://zhuanlan.zhihu.com/p/26514613

https://blog.csdn.net/dymodi/article/details/50358278

https://zhuanlan.zhihu.com/p/50906541

總結

以上是生活随笔為你收集整理的二次规划--积极集法（active set method）的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。