趙 寒,何洪津
(杭州電子科技大學理學院,浙江 杭州 310018)
特征值互補問題(Eigenvalue Complementarity Problem,ECP)是源于工程力學中的一類重要問題[1]。近年來,此類問題已被推廣至張量形式,并稱為張量特征值互補問題[2-3]。在特征值互補問題的內嵌矩陣為嚴格協正矩陣條件下,Júdice等[4]證明了特征值互補問題至少有一個解,且可等價轉化為一個變分不等式問題。特別地,若相關矩陣為對稱的,特征值互補問題還可轉化為單純形約束的瑞利商極大化問題。從轉化后的等價形式上分析,瑞利商極大化問題為一般的約束優化模型,可采用譜投影梯度算法(Spectral Projected Gradient, SPG)[5]進行求解。目前,針對特征值互補問題提出的眾多有效算法中,譜投影梯度算法是最受歡迎的方法之一,并被推廣至求解張量特征值互補問題[6-7]。但是,譜投影梯度算法并未充分利用單純形約束的特殊結構,每次需調用子程序來實現單純形集合的投影。另外,譜投影梯度算法需調用線搜索尋找合適步長以保證目標函數值有一定的增長。對于大規模問題,線搜索無疑會增加計算成本。近年來,針對可分離結構優化問題設計的交替方向乘子法(Alternating Direction Method of Multipliers,ADMM)在大規模圖像處理、機器和統計學習等領域取得了巨大成功[8]。為了提高特征值互補問題的求解效率,本文引入輔助變量,將單純形約束進行分離,使得瑞利商優化模型變成可分的優化問題。在此基礎上,采用交替方向乘子法對新模型進行求解,使其子問題都具有封閉解,達到有效規避單純形約束的投影計算和線搜索循環過程的目的。
特征值互補問題是指尋找1個實數μ>0和非零向量x∈Rn,使得
w=(μB-C)x,w≥0,x≥0,xTw=0
(1)


(2)
文獻[10]的命題9指出優化模型(5)的穩定點即為對稱特征值互補問題的解。此結論提供了一種求解對稱特征值互補問題的有效途徑。而且,文獻[10]的命題10揭示:若C為嚴格協正矩陣(即對任意n維向量x≥0且x≠0,都有xTCx>0)時,對稱特征值互補問題至少有1個解。本文假設對稱特征值互補問題的解集非空。
不難發現,瑞利商極大化模型(2)的單純形約束可表示為2個簡單凸集的交集形式,即,Δn={x∈Rn|eTx=1,x≥0}=X∩Y,其中

(3)
因此,通過引入輔助變量x=y,則問題(2)等價為:

(4)
上述問題是一類具有可分結構的簡單約束優化模型。因此,本文采用交替方向乘子法對問題(4)進行求解。首先引入增廣拉格朗日函數:
其中,λ∈Rn為拉格朗日乘子,β>0為等式約束的罰參數。給定第k步迭代點(yk,λk),交替方向乘子法的迭代格式為:
(5)
迭代格式(5)中的xk+1,yk+1子問題都是約束優化問題。但不難發現,xk+1子問題具有封閉解:
(6)



(7)
綜上,可得求解對稱特征值互補問題的交替方向乘子法,主要步驟如下。

初始化:給定初始迭代點(y0,λ0)∈YY×RRn,選擇β>0和γ>0。(1)若迭代點滿足某個停止準則,算法停止,否則轉步驟2;(2)按式(6)更新xk+1;(3)按式(7)更新yk+1;(4)更新λk+1=λk-β(xk+1-yk+1)。
從本文算法不難發現,每一個子問題都具有封閉迭代形式。相比以往的算法,其迭代格式簡單易實現。對于本文算法的停機準則,采用
來判斷算法是否停止,然后再根據返回結果驗證其是否為原問題的解。當然,也可根據優化問題的一階最優性條件設計停機準則,如文獻[5]中的停機準則。
由于優化模型(4)的目標函數為2個二次函數的商,根據文獻[11]可知μ(y)為半代數函數,則μ(y)滿足Kurdyka-Lojasiewicz不等式。又因單純形約束集合為有界閉集,可知交替方向乘子法產生的迭代序列是有界的。因此,由文獻[12]的推論3.1可知,{(xk,yk,λk)}收斂到模型(4)的KKT點,即原對稱特征值互補問題的解點。因證明過程完全類似文獻[12]的證明,本文省略算法的收斂性證明。
本文通過2個算例來驗證本文算法求解對稱特征值互補問題的可行性。算法執行過程中,統一設定線性化參數γ=1。根據計算經驗發現,交替方向乘子法取固定參數β時,對不同規模問題有一定影響。因此本文采用簡單的調節準則調節算法,即每迭代10步,β值擴大0.2倍。所有試驗運行環境為Windows 10操作系統下的MATLAB R2016a,配置為Intel Core i7-7700HQ CPU, 8GB內存的筆記本電腦。
算例1取B矩陣為單位矩陣,C矩陣為對稱正定矩陣,滿足C=MMT,

由文獻[10]可知,算例1和算例2都至少有1個解,這保證了問題解集的非空性。特別地,根據文獻[10]可知,算例1有唯一解,且解恰是C矩陣的最大特征值λmax(C)。


表1 2種算法求解算例1的數值結果
針對算例1,進一步討論不同初始迭代點對交替方向乘子法的影響。本文選取2組初始迭代點y0=(1,…,1),λ0=(0,…,0)和隨機的y0=rand(n,1),λ0=(0,…,0)進行分析,結果如表2所示。實驗結果表明,初始迭代點的選擇對本文算法影響較小,即使較大規模問題,本文算法依然能夠在較短的計算時間內成功求解。

表2 不同初始迭代點對改進的交替方向乘子法求解算例1的影響
為了更進一步測試算法的收斂性,針對構造的算例2對每一個維度情況進行100次實驗,并觀察100次實驗成功的次數Sn。表3列出了本文算法與SPG算法在停機精度取ε=10-6,罰參數β為103時的計算結果。實驗結果表明,本文算法能在更短的時間內求解對稱的特征值互補問題,其成功率可達100%。

表3 2種算法求解算例2的數值結果
本文針對對稱特征值互補問題提出了一種改進的交替方向乘子法。在計算較大規模問題時,改進算法用時較少,且初始迭代點的選取表現較穩定。但是,部分實際問題不一定滿足本文的對稱性要求,下一步將探討非對稱特征值互補問題的求解。