, 刘红卫1, 谢迪21. 西安电子科技大学数学与统计学院, 西安 710126;
2. 兰州理工大学理学院, 兰州 730050
2015年03月18日 收稿; 2015年09月15日 收修改稿
基金项目: 国家自然科学基金(61072144,61179040)和中央高校基本科研业务费专项基金(K50513100007)资助
通信作者: E-mail:982402513@qq.com
摘要: 为降低半定规划(SDP)问题的迭代复杂度,并且有更好的数值实验结果,提出一种新的宽邻域上的齐次不可行内点算法.半定规划的KKT条件是单调互补问题(MCP),通过构造齐次模型(HMCP)以及提出新的宽邻域来解这个齐次模型,得到半定规划问题的最优解.这种算法容易判定原问题是否可行.在NT方向,证明迭代点在新的宽邻域内是收敛的,且迭代复杂度为O(√nlogL),其中n是SDP问题的维数,L=Tr(X0S0)/ε,其中ε是需要的精度,(X0,S0)是迭代起始点.这个复杂度比一般的半定规划不可行算法的迭代复杂度低.提供了数值实验,证明此算法比其他不可行算法具有更好的数值实验结果.
关键词: 齐次不可行内点算法单调互补问题半定规划
A homogeneous infeasible-interior-point algorithm for semidefinite programming
WU Yue1
, LIU Hongwei1, XIE Di2 1. School of Mathematics and Statistics, Xidian University, Xi'an 710126, China;
2. School of Sciences, Lanzhou University of Technology, Lanzhou 730050, China
Abstract: We propose a homogeneous infeasible-interior-point algorithm for semidefinite programming(SDP) in a new wide neighborhood in order to achieve low iteration complexity and get better experiement numerical results. Complementarity problem (MCP) is the KKT condition of SDP. We sovle MCP by constructing a homogeneous MCP model (HMCP) and proposing a new wide neighborhood. Then we derive the optimal solution of SDP. This algorithm can be easily used to determine whether SDP is feasible or not. At the direction of NT, we prove that the iteration poin is convergent in new wide neighborhood and the iteration complexity is O(√nlogL), where n is the dimension of SDP and L=Tr(X0S0)/ε with ε being the required precision and (X0,S0) the initial point. This algorithm has lower complexity degree than other algorithms for SDP. The numerical experiment is provided. We have proved that this algorithm is better than other infeasible-interior-point algorithms in numerical experimentel results.
Key words: homogeneous infeasible-starting interior-point algorithmmonotone complementarity problemsemidefinite programming
半定规划(SDP)的算法包括内点算法和非内点算法,从初始点是否可行的角度,内点算法分为可行算法和不可行算法.它们之间存在这样一个矛盾: 可行内点算法理论复杂度小,但是实验效果不好; 不可行内点算法实验效果好,但是理论复杂度很大.目前,求解半定规划问题的可行内点算法已经很成熟,所以,不可行内点算法的研究越来越广泛,目的是在较好的实验效果基础上降低其理论复杂度.目前SDP不可行内点算法最好的迭代复杂度是O(n1.5logL)[1],其中n是SDP问题的维数,L=Tr(X0S0)/ε，ε是需要的精度,(X0,S0)是迭代起始点,最好的实验结果已在文献[2]中给出.直接的算法已经很难取得更好的复杂度和实验结果,本文采用一个HMCP模型作为一个间接桥梁,提出一个新的不可行内点算法,并且采用新的宽邻域,比较起其他算法,此算法可以将复杂度降低到O(nlogL)， 并且实验效果比文献[2]中的实验结果更好,另外一个优势是可以判别原问题是否是可行问题.所以,此算法在内点算法中是有效的.
Andersan和Ye[3]研究线性规划的齐次算法,本文考虑半定规划的齐次算法.
考虑半定规划原问题和对偶问题的标准形式:
minTr(CX)
(P) s.t.Tr(AiX)=bi，i=1,…m，
X≥0，
maxbTy
| $\left( D \right) s.t.\sum\limits_{i=1}^{m}{{{y}_{i}}}{{A}_{i}}+S=C$ |
其最优性条件是一个单调互补问题(MCP):
Tr(AiX)=bi，X≥0，i=1,…m，
(MCP)
XS=0.
以上问题等价于互补问题:
minTr(XS)
(MCP) s.t.Tr(AiX)=bi，X≥0，i=1,…m，
令A=(vec(A1),…,vec(Am))T，
h(y,vec(X))=Avec(X)－b，
g(y,vec(X))=vec(C)－ATy，
| $f\left( y,vec\left( X \right) \right)=\left[ \begin{matrix} h\left( y,vec\left( X \right) \right) \\ g\left( y,vec\left( X \right) \right) \\\end{matrix} \right].$ |
注: 1) vec(M)表示矩阵M的各列依次排列成的列向量; Rm表示m维欧式空间; Sn+表示n阶实对称半正定矩阵集; Sn表示n阶实对称矩阵集.
2) 若?f是f的雅克比矩阵,则?f是半正定矩阵[4].
那么，问题MCP等价于矩阵形式:
| $\begin{align} & min vec{{\left( X \right)}^{T}}vec\left( S \right) \\ & s.t.\left[ \begin{matrix} 0 \\ vec\left( S \right) \\\end{matrix} \right]=\left[ \begin{matrix} Avec\left( X \right)-b \\ vec\left( C \right)-{{A}^{T}}y \\\end{matrix} \right]=f\left( y,vec\left( X \right) \right), \\ & S\ge 0,X\ge 0. \\ \end{align}$ |
1) 使迭代复杂度降低到O(nlogL).
2) 采用新的宽邻域,与其他SDP不可行算法相比,此算法具有较好的实验效果.
3) 判断原问题是否可行: 如果问题MCP有解,那么算法将会求得原问题最优解; 如果问题MCP无解,那么可以通过算法和齐次模型来证明.
本文中符号“?”表示2个矩阵的克罗内克积,即若A是一个m×n阶的矩阵,B是一个p×q阶的矩阵,则克罗内克积A?B是一个mp×nq的分块矩阵,表示为:
| $A\otimes B=\left[ \begin{matrix} {{a}_{11}}B & \cdots & {{a}_{1n}}B \\ \vdots & \ddots & \vdots \\ {{a}_{m1}}B & \cdots & {{a}_{mn}}B \\\end{matrix} \right].$ |
minTr(XS)+τκ
(HMCP) s.t.∑mi=1yiAi+S=τC，S≥0，τ≥0，
Tr(AiX)=τbi，X≥0，i=1,…m，
－Tr(CX)+bTy－κ=0，κ≥0.
类似于问题MCP等价形式的分析,有
| $\begin{align} & \tau f(y/\tau ,vec\left( X \right)/\tau )=\left[ \begin{matrix} \tau h(y/\tau ,vec\left( X \right)/\tau ) \\ \tau g(y/\tau ,vec\left( X \right)/\tau ) \\\end{matrix} \right] \\ & =\left[ \begin{matrix} Avec\left( X \right)-\tau b \\ \tau vec\left( C \right)-{{A}^{T}}y \\\end{matrix} \right]; \\ & -(y,vec\left( X \right))f(y/\tau ,~vec\left( X \right)/\tau ) \\ & =-{{y}^{T}}h(y/\tau ,vec\left( X \right)/\tau )\text{ }- \\ & vec{{\left( X \right)}^{T}}g(y/\tau ,vec\left( X \right)/\tau ) \\ & =-vec{{\left( C \right)}^{T}}vec\left( X \right)+{{b}^{T}}y, \\ \end{align}$ |
| $\begin{align} & min~vec{{\left( X \right)}^{T}}vec\left( S \right)+\tau \kappa \\ & s.t.\left[ \begin{matrix} 0 \\ vec\left( S \right) \\ \kappa \\\end{matrix} \right]=\left[ \begin{matrix} \tau f(y/\tau ,vec\left( X \right)/\tau ) \\ -(y,vec\left( X \right))f(y/\tau ,vec\left( X \right)/\tau ) \\\end{matrix} \right] \\ & S\ge 0,\tau \ge 0,X\ge 0,\kappa \ge 0. \\ \end{align}$ |
| $\psi (y,vec\left( X \right),\tau ):=\left[ \begin{align} & \tau f(y/\tau ,vec\left( X \right)/\tau ) \\ & -(y,vec\left( X \right))f(y/\tau ,vec\left( X \right)/\tau ) \\ \end{align} \right]$ | (1) |
| $\begin{align} & \Delta \psi = \\ & \left[ \begin{matrix} \Delta f & f-\frac{1}{\tau } \Delta f \\ {{f}^{T}}-\frac{1}{\tau } (y,vec\left( X \right)) \Delta f & -\frac{1}{\tau }(y/\tau ,vec\left( X \right)/\tau ) \Delta f \\\end{matrix} \right] \\ \end{align}$ | (2) |
引理1.1 定义ψ如(1)所示,则如下结论成立:
(a) (y,vec(X),τ)ψ(y,vec(X),τ)=0.
(b) Δψ(y,vec(X),τ)(y,vec(X),τ)T=ψ(y,vec(X),τ).
(c) (y,vec(X),τ)Δψ(y,vec(X),τ)=－ψ(y,vec(X),τ)T.
(d)Δf是定义在Rm×Rn2上的半正定矩阵,Δψ是定义在Rm×Rn2×R上的半正定矩阵.
证明 由ψ的定义和(2)易得(a),(b),(c).
因为f是从Rm×Sn+→Rm×Rn2的连续单调映射,所以?f是半正定矩阵[4]. 将下式展开可得:
| $\begin{align} & (dy,vec(dX),d\tau ) \nabla \psi (dy,vec(dX),d\tau ) \\ & {{(dy,vec(dX),d\tau )}^{T}} \\ & =[(dy,vec(dX))-d\tau (y,vec\left( X \right))/\tau ]\cdot \\ & \nabla f{{[(dy,vec(dX))-d\tau (y,vec\left( X \right))/\tau ]}^{T}} \\ & \ge 0 \\ \end{align}$ |
定理1.1
(a) 若f是从Rm+n2到Rm+n2的连续单调映射,则ψ是从Rm+n2+1到Rm+n2+1的连续单调映射.
(b) 问题HMCP是渐进可行的,并且每一个渐进可行点是一个渐进互补解.
(c) 设问题HMCP的极大互补解是(vec(X*),y*,vec(S*),τ*,κ*)，问题MCP有最优解当且仅当τ*>0成立,这时,(vec(X*)/τ*,y*/τ*,vec(S*)/τ*)是问题MCP的一个互补解.
(d) 当κ*>0时,问题MCP是强不可行的.
证明
(a) 任取(y1,vec(X1),τ1)，(y2,vec(X2),τ2)∈Rm+n2+1，因为 (y1/τ1,vec(X1)/τ1)，(y2/τ2,vec(X2)/τ2)∈Rm+n2， 由引理1.1得
| $\begin{align} & ~({{y}_{1}}-{{y}_{2}},vec({{X}_{1}})-vec({{X}_{2}}),{{\tau }_{1}}-{{\tau }_{2}})\cdot \\ & (\psi ({{y}_{1}},vec({{X}_{1}}),{{\tau }_{1}})-\psi ({{y}_{2}},vec({{X}_{2}}),{{\tau }_{2}})) \\ & ={{\tau }_{1}}{{\tau }_{2}}{{(({{y}_{1}}/{{\tau }_{1}},vec({{X}_{1}})/{{\tau }_{1}})-({{y}_{2}}/{{\tau }_{2}},vec({{X}_{2}})/{{\tau }_{2}}))}^{T}}\cdot \\ & (f({{y}_{1}}/{{\tau }_{1}},vec({{X}_{1}})/{{\tau }_{1}})- \\ & f({{y}_{2}}/{{\tau }_{2}},vec({{X}_{2}})/{{\tau }_{2}})) \\ & \ge 0. \\ \end{align}$ |
(b) 取vec(X)t=(1/2)te1(X≥0)，yt=(1/2)te2，τt=(1/2)t≥0，κt=(1/2)t≥0，vec(S)t=(1/2)te1(S≥0).
其中e1为每个元素全为1的n2维列向量,e2为每个元素全为1的m维列向量,当t→∞时,有
| $\begin{align} & Avec{{\left( X \right)}^{t}}-{{\tau }^{t}}b={{\left( 1/2 \right)}^{t}}(A{{e}_{1}}-b)\to 0, \\ & vec{{\left( S \right)}^{t}}-{{\tau }^{t}}vec\left( C \right)+{{A}^{T}}{{y}^{t}}={{\left( 1/2 \right)}^{t}}\cdot \\ & [{{e}_{1}}-vec\left( C \right)+{{A}^{T}}{{e}_{2}}]\to 0, \\ & {{\kappa }^{t}}+vec{{\left( C \right)}^{T}}vec{{\left( X \right)}^{T}}-{{b}^{T}}{{y}^{t}}={{\left( 1/2 \right)}^{t}}\cdot \\ & [1+vec{{\left( C \right)}^{T}}{{e}_{1}}-{{b}^{T}}{{e}_{2}}]\to 0, \\ \end{align}$ |
| $\begin{align} & vec{{\left( {\tilde{S}} \right)}^{T}}vec\tilde{X}+\tilde{\tau }\tilde{\kappa }=(\tilde{y},vec\left( {\tilde{X}} \right),\tilde{\tau }) \\ & {{(0,vec\left( {\tilde{S}} \right),\tilde{\kappa })}^{T}}~ \\ & =(\tilde{y},vec\left( {\tilde{X}} \right),\tilde{\tau })\psi (\tilde{y},vec\left( {\tilde{X}} \right),\tilde{\tau }) \\ & =0, \\ \end{align}$ |
(c) “?”若τ*>0，因为(vec(X*),y*,vec(S*),τ*,κ*)是问题HMCP的极大互补解,则
| $\begin{align} & vec({{S}^{*}})={{\tau }^{*}}g({{y}^{*}}/{{\tau }^{*}},vec({{X}^{*}})/{{\tau }^{*}}), \\ & vec{{({{X}^{*}})}^{T}}vec({{S}^{*}})=0, \\ & h({{y}^{*}}/{{\tau }^{*}},vec({{X}^{*}})/{{\tau }^{*}})=0, \\ \end{align}$ | (3) |
| $所以vec({{S}^{*}})/{{\tau }^{*}}=g({{y}^{*}}/{{\tau }^{*}},vec({{X}^{*}})/{{\tau }^{*}}),$ | (4) |
| $vec{{({{X}^{*}})}^{T}}vec({{S}^{*}})/{{({{\tau }^{*}})}^{2}}=0,$ | (5) |
“?”假设问题MCP的一个解为
又τ*是问题HMCP的一个极大互补解,且τ*≥0， 则问题HMCP的每一个极大互补解中有τ*>0 (否则,τ≡0， 矛盾).
(d) 因为(vec(X*);y*;vec(S*);τ*;κ*)是问题HMCP的极大互补解,则
| $Tr({{X}^{*}}{{S}^{*}})+{{\tau }^{*}}{{\kappa }^{*}}=0,$ | (6) |
| $Avec({{X}^{*}})-{{\tau }^{*}}b=0,$ | (7) |
| ${{\tau }^{*}}vec\left( C \right)-{{A}^{T}}{{y}^{*}}=vec({{S}^{*}}),$ | (8) |
| $-Tr(C{{X}^{*}})+{{b}^{T}}{{y}^{*}}={{\kappa }^{*}},$ | (9) |
| $Avec({{X}^{*}})=0,$ | (10) |
| $-{{A}^{T}}{{y}^{*}}=vec({{S}^{*}}),$ | (11) |
| $Tr(C{{X}^{*}})-{{b}^{T}}{{y}^{*}}<0,$ | (12) |
当Tr(CX*)<0成立时,由式(10)得,(P)有改良射线,则(D)强不可行,MCP强不可行.
当bTy*>0成立时,由式(11)得,(D)有改良射线,则(P)强不可行,MCP强不可行.证毕.
2 齐次模型的中心路径由定理1.1知,通过找问题HMCP的一个极大互补解来解问题MCP. 为了方便,取起始点X0=I,S0=I,τ0=1,κ0=1,y0=0，则X0S0=I,τ0κ0=1.
定义集合和差向量:
| $\begin{align} & {{H}_{++}}={{S}^{n}}_{++}\times {{S}^{n}}_{++}\times {{R}_{++}}\times {{R}_{++}}\times {{R}^{m}}, \\ & {{R}_{D}}:=-\tau C+\underset{i=1}{\overset{m}{\mathop{\sum }}}\,{{A}_{i}}^{T}{{y}_{i}}+S, \\ & {{R}_{P}}:=\tau b-Avec\left( X \right), \\ & {{R}_{G}}:=vec{{\left( C \right)}^{T}}vec\left( X \right)-{{b}^{T}}y+\kappa . \\ \end{align}$ |
| $\begin{align} & C\left( \mu \right):=\left\{ \left( X,S,\tau ,\kappa ,y \right)\in {{H}_{++}}:\text{ }\left[ \begin{matrix} XS & 0 \\ 0 & \tau \kappa \\\end{matrix} \right] \right.= \\ & \left. \theta I,\left[ \begin{matrix} {{R}_{P}} \\ {{R}_{D}} \\ {{R}_{G}} \\\end{matrix} \right]=\theta \left[ \begin{matrix} R_{P}^{0} \\ R_{D}^{0} \\ R_{G}^{0} \\\end{matrix} \right],\mu >0,0<\theta \le 1 \right\} \\ \end{align}$ | (13) |
定理2.1
(a) 对任意的
(b) 若起始点为(X0=I,S0=I,τ0=1,κ0=1,y0=0),则对任意的0<θ≤1， 中心路径(13)的解是存在唯一的.
(c) 中心路径是一条有界轨迹.
(d) θ→0时,存在唯一极限点(X(0),S(0),τ(0),κ(0),y(0))是问题HCMP的一个极大互补解.
证明
| $\begin{align} & \left( a \right) R_{D}^{0}=-{{\tau }^{0}}C+\underset{i=1}{\overset{m}{\mathop{\sum }}}\,{{A}_{i}}^{T}{{y}^{0}}_{i}+{{S}^{0}}, \\ & R_{P}^{0}=\tau b-Avec({{X}^{0}}), \\ & R_{G}^{0}=vec{{\left( C \right)}^{T}}vec({{X}^{0}})-{{b}^{T}}{{y}^{0}}+{{\kappa }^{0}}. \\ \end{align}$ |
| $\begin{align} & \left[ \begin{matrix} 0 \\ vec\left( S \right) \\ \kappa \\\end{matrix} \right]=\left[ \begin{matrix} 0 & A & -b \\ -{{A}^{T}} & 0 & vec\left( C \right) \\ {{b}^{T}} & -vec{{\left( C \right)}^{T}} & 0 \\\end{matrix} \right]\cdot \\ & \left[ \begin{matrix} y \\ vec\left( X \right) \\ \tau \\\end{matrix} \right]+\left[ \begin{matrix} \theta R_{P}^{0} \\ \theta vec(R_{D}^{0}) \\ \theta R_{G}^{0} \\\end{matrix} \right]. \\ \end{align}$ |
(c) 为了方便,我们令a=(y;vec(X);τ)，b=(0;vec(S);κ)，r=(Rp;vec(RD);RG). 则有r=b－ψ(a)=θr0，r0=b0－ψ(a0).
因为ψ单调,则(a0-a)T(ψ(a0)－ψ(a))≥0，又aTψ(a)=0，(a0)Tψ(a0)=0，则－(a0)Tψ(a)-aTψ(a0)≥0. 所以有如下不等式:
| $\begin{align} & {{a}^{T}}{{r}^{0}}={{a}^{T}}{{b}^{0}}-{{a}^{T}}\psi ({{a}^{0}}) \\ & ={{a}^{T}}{{b}^{0}}+{{b}^{T}}{{a}^{0}}-{{b}^{T}}{{a}^{0}}-{{a}^{T}}\psi ({{a}^{0}}) \\ & ={{a}^{T}}{{b}^{0}}+{{b}^{T}}{{a}^{0}}-{{({{a}^{0}})}^{T}}\cdot \\ & (\theta {{r}^{0}}+\psi \left( a \right))-{{a}^{T}}\psi ({{a}^{0}}) \\ & ={{a}^{T}}{{b}^{0}}+{{b}^{T}}{{a}^{0}}-\theta {{({{a}^{0}})}^{T}}{{r}^{0}}- \\ & {{({{a}^{0}})}^{T}}\psi \left( a \right)-{{a}^{T}}\psi ({{a}^{0}}) \\ & \ge {{a}^{T}}{{b}^{0}}+{{b}^{T}}{{a}^{0}}-\theta {{({{a}^{0}})}^{T}}{{r}^{0}} \\ & ={{a}^{T}}{{b}^{0}}+{{b}^{T}}{{a}^{0}}-\theta {{({{a}^{0}})}^{T}}({{b}^{0}}-\psi ({{a}^{0}})) \\ & ={{a}^{T}}{{b}^{0}}+{{b}^{T}}{{a}^{0}}-\theta {{({{a}^{0}})}^{T}}{{b}^{0}} \\ \end{align}$ |
| $\begin{align} & \theta {{a}^{T}}{{r}^{0}}={{a}^{T}}\left( b-\psi \left( a \right) \right)={{a}^{T}}b \\ & =\theta \left( n+1 \right)=\theta {{({{a}^{0}})}^{T}}{{b}^{0}}. \\ \end{align}$ |
| ${{a}^{T}}{{b}^{0}}+{{b}^{T}}{{a}^{0}}\le \left( 1+\theta \right){{({{a}^{0}})}^{T}}{{b}^{0}},$ |
| $\begin{align} & {{(y;vec\left( X \right);\tau )}^{T}}(0;vec({{S}^{0}});{{\kappa }^{0}})+(0;vec\left( S \right);\text{ } \\ & \kappa {{)}^{T}}({{y}^{0}};vec({{X}^{0}});{{\tau }^{0}})\le \left( 1+\theta \right)({{y}^{0}};vec({{X}^{0}}); \\ & {{\tau }^{0}}{{)}^{T}}~(0;vec({{S}^{0}});{{\kappa }^{0}}) \\ \end{align}$ |
(d) 为方便,记
3 具有O(nlogL)复杂性的齐次不可行算法下面用齐次不可行内点算法计算在中心路径附近产生的迭代点.
3.1 计算搜索方向对单调互补问题MCP采用牛顿法，可得以下线性方程组:
| $\begin{align} & \left[ \begin{matrix} 0 & vec\left( C \right) & -{{A}^{T}} \\ -vec{{\left( C \right)}^{T}} & 0 & {{b}^{T}} \\ A & -b & 0 \\\end{matrix} \right]~\left[ \begin{matrix} vec(dX) \\ d\tau \\ dy \\\end{matrix} \right]- \\ & \left[ \begin{matrix} vec(dS) \\ d\kappa \\ 0 \\\end{matrix} \right]=\eta \left[ \begin{matrix} vec(R_{D}^{k}) \\ R_{G}^{k} \\ R_{P}^{k} \\\end{matrix} \right] \\ \end{align}$ | (14a) |
| ${{X}^{k}}dS+dX{{S}^{k}}=\gamma {{\mu }^{k}}I-{{X}^{k}}{{S}^{k}}$ | (14b) |
| ${{\kappa }^{k}}d\tau +{{\tau }^{k}}d\kappa =\gamma {{\mu }^{k}}-{{\tau }^{k}}{{\kappa }^{k}}$ | (14c) |
从式(14)中能得到对称的dS，但是不能保证得到对称的dX，为了得到对称的dX， 使用Zhang[9]提出的方法,将(14b)替换为
| ${{H}_{P}}({{X}^{k}}dS+dX{{S}^{k}})={{H}_{P}}(\gamma {{\mu }^{k}}I-{{X}^{k}}{{S}^{k}}).$ | (14b1) |
| ${{H}_{P}}\left( M \right)=12[PM{{P}^{-1}}+{{P}^{-T}}{{M}^{T}}{{P}^{T}}],M\in {{R}^{n\times n}}.$ |
在文献[9]中Zhang证明了若P是一个非奇异矩阵,则
| ${{H}_{P}}\left( M \right)=\gamma {{\mu }^{k}}I\Leftrightarrow M=\gamma {{\mu }^{k}}I,$ |
| ${{H}_{P}}({{X}^{k}}dS+dX{{S}^{k}})=\gamma {{\mu }^{k}}I-{{H}_{P}}({{X}^{k}}{{S}^{k}}).$ | (14b2) |
定义(Xk0,yk,Sk0,τk>0,κk>0)是当前第k个迭代点,并且第k+1个迭代点为
(Xk+1,yk+1,Sk+1,τk+1,κk+1)=(Xk,yk,Sk,τk,κk)+α(dX,dy,dS,dτ,dκ).
其中,α是搜索步长,(dX,dτ,dy,dS,dκ)为搜索方向.
引用文献[15]给出的矩阵以及向量正部和负部的概念,本文采用如下系统计算搜索方向
| $\begin{align} & (dX,d\tau ,dy,dS,d\kappa ): \\ & \left[ \begin{matrix} 0 & vec\left( C \right) & -{{A}^{T}} \\ -vec{{\left( C \right)}^{T}} & 0 & {{b}^{T}} \\ A & -b & 0 \\\end{matrix} \right]~\left[ \begin{matrix} vec(dX) \\ d\tau \\ dy \\\end{matrix} \right]- \\ & \left[ \begin{matrix} vec(dS) \\ d\kappa \\ 0 \\\end{matrix} \right]=\eta \left[ \begin{matrix} vec(R_{D}^{k}) \\ R_{G}^{k} \\ R_{P}^{k} \\\end{matrix} \right], \\ \end{align}$ | (15a) |
| $\begin{align} & {{H}_{P}}({{X}^{k}}dS+dX{{S}^{k}})={{[\gamma {{\mu }^{k}}I-{{H}_{P}}({{X}^{k}}{{S}^{k}})]}^{-}}+ \\ & \sqrt{n+1}{{[\gamma {{\mu }^{k}}I-{{H}_{P}}({{X}^{k}}{{S}^{k}})]}^{+}}, \\ \end{align}$ | (15b) |
| $\begin{align} & {{\kappa }^{k}}d\tau +{{\tau }^{k}}d\kappa ={{[\gamma {{\mu }^{k}}-{{\tau }^{k}}{{\kappa }^{k}}]}^{-}}+\sqrt{n+1}\cdot \\ & {{[\gamma {{\mu }^{k}}-{{\tau }^{k}}{{\kappa }^{k}}]}^{+}}. \\ \end{align}$ | (15c) |
下面用克罗内克积的形式等价表示(15b):
| $\begin{align} & \left( 15b \right)\Leftrightarrow [P({{X}^{k}}dS+dX{{S}^{k}}){{P}^{-1}}+ \\ & {{P}^{-T}}({{S}^{k}}dX+dS{{X}^{k}}){{P}^{T}}]/2 \\ & ={{[\gamma {{\mu }^{k}}I-{{H}_{P}}({{X}^{k}}{{S}^{k}})]}^{-}}+ \\ & \sqrt{n+1}{{[\gamma {{\mu }^{k}}I-{{H}_{P}}({{X}^{k}}{{S}^{k}})]}^{+}} \\ & \Leftrightarrow \text{ }{{E}^{k}}vec(dX)+{{F}^{k}}vec(dS)=vec(R_{C}^{k}), \\ \end{align}$ | (15b1) |
| $\begin{align} & E:=({{P}^{-T}}S\otimes P+P\otimes {{P}^{-T}}S)/2, \\ & F:=(PX\otimes {{P}^{-T}}+{{P}^{-T}}\otimes PX)/2, \\ & R_{C}^{k}:={{[\gamma {{\mu }^{k}}I-{{H}_{P}}({{X}^{k}}{{S}^{k}})]}^{-}}+ \\ & \sqrt{n+1}{{[\gamma {{\mu }^{k}}I-{{H}_{P}}({{X}^{k}}{{S}^{k}})]}^{+}}. \\ \end{align}$ |
下面求解由系统(15a),(15b1),(15c)得到的迭代搜索方向(dX,dy,dS,dτ,dκ).
由(15b1)得,
| $vec(dS)={{({{F}^{k}})}^{-1}}[vec({{R}^{k}}_{C})-{{E}^{k}}vec(dX)],$ | (16) |
| $d\kappa ={{({{\tau }^{k}})}^{-1}}({{r}_{C}}^{k}-{{\kappa }^{k}}d\tau ),$ | (17) |
将式(16)、式(17)代入式(15a)减少牛顿方程,得牛顿搜索方向(dX,dy,dS,dτ,dκ)为如下系统的解:
| $\begin{align} & \left[ \begin{matrix} {{({{F}^{k}})}^{-1}}{{E}^{k}} & vec\left( C \right) & -{{A}^{T}} \\ -vec{{\left( C \right)}^{T}} & {{({{\tau }^{k}})}^{-1}}{{\kappa }^{k}} & {{b}^{T}} \\ A & -b & 0 \\\end{matrix} \right]\left[ \begin{matrix} vec(dX) \\ d\tau \\ dy \\\end{matrix} \right]= \\ & \left[ \begin{matrix} \eta vec({{R}^{k}}_{D})+{{({{F}^{k}})}^{-1}}vec(R_{C}^{k}) \\ \eta R_{G}^{k}+{{({{\tau }^{k}})}^{-1}}{{r}_{C}}^{k} \\ \eta R_{P}^{k} \\\end{matrix} \right]. \\ \end{align}$ | (18) |
引理3.1 若K有如上定义,则
| ${{K}^{-1}}=\left[ \begin{matrix} {{Q}^{-1}}-{{Q}^{-1}}{{A}^{T}}{{(A{{Q}^{-1}}{{A}^{T}})}^{-1}}A{{Q}^{-1}} & {{Q}^{-1}}{{A}^{T}}{{(A{{Q}^{-1}}{{A}^{T}})}^{-1}} \\ {{(A{{Q}^{-1}}{{A}^{T}})}^{-1}}A{{Q}^{-1}} & -{{(A{{Q}^{-1}}{{A}^{T}})}^{-1}} \\\end{matrix} \right],$ |
证明 易验证K－1K=KK－1=Im+n2.□
由上述可知,牛顿搜索方向(18)等价于
| $\begin{align} & \left[ \begin{matrix} K & (vec\left( C \right);-b) \\ {{(-vec\left( C \right);b)}^{T}} & {{({{\tau }^{k}})}^{-1}}{{\kappa }^{k}} \\\end{matrix} \right]\cdot \\ & \left[ \begin{matrix} (vec(dX);dy) \\ d\tau \\\end{matrix} \right]=\left[ \begin{matrix} {\tilde{r}} \\ \eta (R_{G}^{k})+{{({{\tau }^{k}})}^{-1}}{{r}_{C}}^{k} \\\end{matrix} \right], \\ \end{align}$ | (19) |
定义p和q是如下线性系统的解:
| $Kp=(vec\left( C \right);-b);\text{ }Kq=\tilde{r}$ |
(第二行)
| $\begin{align} & {{({{\tau }^{k}})}^{-1}}{{\kappa }^{k}}d\tau +{{(-vec\left( C \right);b)}^{T}}(vec(dX);dy) \\ & =\eta (R_{G}^{k})+{{({{\tau }^{k}})}^{-1}}{{r}_{C}}^{k}, \\ \end{align}$ | (20) |
| $K(vec(dX);dy)+(vec\left( C \right);-b)d\tau =\tilde{r}$ |
又由引理3.1知,K是非奇异矩阵,则
| $(vec(dX);dy)=q-pd\tau ,$ | (21) |
(τk)－1κkdτ+(－vec(C);b)T(q-pdτ)
=η(RkG)+(τk)－1rCk,
即
| $\begin{align} & [{{({{\tau }^{k}})}^{-1}}{{\kappa }^{k}}-(-vec\left( C \right),b)p]d\tau \\ & =\eta R_{G}^{k}+{{({{\tau }^{k}})}^{-1}}{{r}_{C}}^{k}-(-vec\left( C \right),b)q, \\ \end{align}$ | (22) |
因为搜索方向系统的解是存在唯一的,所以dτ可以从式(22)中获得,进而可以从式(21)中得到vec(dX),dy，从式(16),式(17)中得到vec(dS),dκ.
通过以上分析,计算搜索方向的主要工作在于计算向量p和q.
定义
| $K\left( p;q \right)=\bar{K}\left( p;-q \right)=(vec\left( C \right);-b;\tilde{r}).$ |
3.2 不可行性和对偶间隙的下降关系本文将限制关系γ=τ1，τ1≤1/4，β≤1/2成立.由文献[15]引理3.5知,若(X,y,S,τ,κ)∈N(τ1,β)，τ1≤1/4，β≤1/2， 有如下不等式成立:
| $\begin{align} & \left( \tau -1 \right){{\mu }^{k}}\left( n+1 \right)\le Tr(R_{C}^{k})+r_{C}^{k}\le \\ & \left( \tau +\beta \tau /\sqrt{n+1}-1 \right){{\mu }^{k}}\left( n+1 \right)<0. \\ \end{align}$ |
| $\eta =\frac{-Tr({{R}^{k}}_{C})+r_{C}^{k}}{\left( n+1 \right){{\mu }^{k}}}.$ | (23) |
引理3.2 算法更新后,迭代点的不可行性有如下结论：
(a) Rk+1P=(1－αη)RkP.
(b) Rk+1D=(1－αη)RkD.
(c) Rk+1G=(1－αη)RkG.
证明 (a) 根据RP,RD,RG的定义及(15a)可得,
Rk+1P=(τk+αdτ)b-A(vec(Xk)+αvec(dX))
=τkb-Avec(Xk)+α[bdτ-Avec(dX)]
=RkP-αηRkP=(1-αη)RkP.
同理,可得(b),(c).
引理3.3 若τ1≤1/4，β≤1/2， 则由系统(15)定义的方向满足
vec(dX)Tvec(dS)+dτdκ=0.
证明 由系统(15)得
| $\begin{align} & vec{{(dX)}^{T}}vec(dS)+d\tau d\kappa \\ & =vec{{(dX)}^{T}}[-{{A}^{T}}dy+vec\left( C \right)d\tau -\eta vec(R_{D}^{k})]+ \\ & d\tau [-vec{{\left( C \right)}^{T}}vec(dX)+{{b}^{T}}dy-\eta R_{G}^{k}] \\ & =[-vec{{(dX)}^{T}}{{A}^{T}}+d\tau {{b}^{T}}]dy- \\ & vec{{(dX)}^{T}}vec(R_{D}^{k})\eta -d\tau R_{G}^{k}\eta ~ \\ & =-[{{(R_{P}^{k})}^{T}}dy+vec{{(dX)}^{T}}vec(R_{D}^{k})+d\tau R_{G}^{k}]\eta . \\ \end{align}$ | (24) |
| $\begin{align} & {{({{R}_{P}})}^{T}}y+vec{{\left( X \right)}^{T}}vec({{R}_{D}})+\tau {{R}_{G}} \\ & =vec{{\left( X \right)}^{T}}vec\left( S \right)+\tau \kappa \\ \end{align}$ |
| $\begin{align} & {{({{y}^{k}})}^{T}}{{R}^{k}}_{P}+vec{{({{X}^{k}})}^{T}}vec(R_{D}^{k})+{{\tau }^{k}}R_{G}^{k} \\ & =vec{{({{X}^{k}})}^{T}}vec({{S}^{k}})+{{\tau }^{k}}{{\kappa }^{k}} \\ \end{align}$ | (25) |
| $\begin{align} & {{({{y}^{k}}+dy)}^{T}}R_{P}^{k}+{{(vec({{X}^{k}})+vec(dX))}^{T}}\cdot \\ & vec(R_{D}^{k})+({{\tau }^{k}}+d\tau )R_{G}^{k} \\ & ={{(vec({{X}^{k}})+vec(dX))}^{T}}(vec({{S}^{k}})+ \\ & vec(dS))+({{\tau }^{k}}+d\tau )({{\kappa }^{k}}+d\kappa ). \\ \end{align}$ | (26) |
| $\begin{align} & [vec{{({{X}^{k}})}^{T}}vec({{S}^{k}})+{{\tau }^{k}}{{\kappa }^{k}}]\left( 1-\eta \right)+ \\ & [{{(R_{P}^{k})}^{T}}dy+vec{{(dX)}^{T}}vec(R_{D}^{k})+d\tau R_{G}^{k}]\cdot \\ & \left( 1-\eta \right) \\ & =[vec{{({{X}^{k}})}^{T}}vec({{S}^{k}})+{{\tau }^{k}}{{\kappa }^{k}}]+ \\ & [vec{{(dX)}^{T}}vec(dS)+d\tau d\kappa ]+ \\ & vec{{({{X}^{k}})}^{T}}vec(dS)+vec{{(dX)}^{T}}vec({{S}^{k}})+ \\ & {{\tau }^{k}}d\kappa +d\tau {{\kappa }^{k}}. \\ \end{align}$ |
| $\begin{align} & {{(R_{P}^{k})}^{T}}dy+vec{{(dX)}^{T}}vec(R_{D}^{k})+d\tau R_{G}^{k} \\ & =[vec{{({{X}^{k}})}^{T}}vec({{S}^{k}})+{{\tau }^{k}}{{\kappa }^{k}}]\eta + \\ & [vec{{({{X}^{k}})}^{T}}vec(dS)+vec{{(dX)}^{T}}vec({{S}^{k}})+ \\ & {{\tau }^{k}}d\kappa +d\tau {{\kappa }^{k}}]. \\ \end{align}$ | (27) |
| $\begin{align} & vec{{({{X}^{k}})}^{T}}vec(dS)+vec{{(dX)}^{T}}vec({{S}^{k}})+{{\tau }^{k}}d\kappa + \\ & d\tau {{\kappa }^{k}}=Tr(R_{C}^{k})+r_{C}^{k}. \\ \end{align}$ | (28) |
| $\matrix{ \matrix{ vec{\left( {dX} \right)^T}vec\left( {dS} \right) + d\tau d\kappa \hfill \cr = - [\left( {vec{{\left( {{X^k}} \right)}^T}vec\left( {{S^k}} \right) + {\tau ^k}{\kappa ^k}} \right)\eta + \hfill \cr} \hfill \cr \matrix{ Tr(R_C^k) + r_C^k\eta ] \hfill \cr = - [\left( {n + 1} \right){\mu ^k}\eta + Tr(R_C^k) + r_C^k] \hfill \cr \eta = 0 \hfill \cr} \hfill \cr } $ |
引理3.4 若τ1≤1/4，β≤1/2， 则算法更新后对偶间隙满足如下等式
μk+1=(1－αη)μk.
证明
| $\begin{align} & {{\mu }^{k+1}}=[vec{{({{X}^{k+1}})}^{T}}vec({{S}^{k+1}})+{{\tau }^{k+1}}{{\kappa }^{k+1}}]/\left( n+1 \right) \\ & =[{{(vec({{X}^{k}})+\alpha vec(dX))}^{T}}(vec({{S}^{k}})+ \\ & \alpha vec(dS))+({{\tau }^{k}}+\alpha d\tau )({{\kappa }^{k}}+\alpha d\kappa )]/\left( n+1 \right) \\ & =\{[vec{{({{X}^{k}})}^{T}}vec({{S}^{k}})+{{\tau }^{k}}{{\kappa }^{k}}]+ \\ & \alpha [vec{{({{X}^{k}})}^{T}}vec(dS)+ \\ & vec{{(dX)}^{T}}vec({{S}^{k}})+{{\tau }^{k}}d\kappa +d\tau {{\kappa }^{k}}]+ \\ & {{\alpha }^{2}}[vec{{(dX)}^{T}}vec(dS)+d\tau d\kappa ]\}/\left( n+1 \right) \\ & ={{\mu }^{k}}+\alpha Tr({{R}^{k}}_{C}+{{r}^{k}}_{C})/\left( n+1 \right) \\ & =\left( 1-\alpha \eta \right){{\mu }^{k}} \\ \end{align}$ |
引理3.2和引理3.4说明当
3.3 齐次不可行内点算法上一节已限制尺度矩阵:
P(X,S):={P∈Sn++:PXSP－1∈Sn}.
本文取P=W1/2， 得到的方向称为NT方向. 其中矩阵
定义F++=Sn++×Rm×Sn++， SDP中常用的一个宽邻域是－∞邻域,定义为:
N－∞(1－τ):={(X,y,S)∈F++:λmin(XS)≥τμ}.
其中,0<τ<1，μk=Tr(XkSk)/(n+1).
后来,Ai和Zhang[11]提出求解线性互补问题的新的宽邻域,Li和Terlaky[12]将其推广到半定规划,Liu等[1]又将其推广到对称锥规划,刘新泽[13]将这种新的宽邻域推广到关于1范数的宽邻域. 本文采用Li和Terlaky[12]提出的新的宽邻域:
N(τ1,τ2)={(X,y,S)∈F++:
‖[τ1μI－X12SX12]+‖F≤(τ1－τ2)μ}.
其中,0<τ2<τ1<1.
文献[16]已证明,下面关系成立:
| $\begin{align} & N_{\infty }^{-}(1-{{\tau }_{1}})N({{\tau }_{1}},{{\tau }_{2}})N_{\infty }^{-}(1-{{\tau }_{2}}), \\ & \forall 0<{{\tau }_{2}}<{{\tau }_{1}}<1. \\ \end{align}$ |
| $\begin{align} & N\left( {{\tau }_{1}},\beta \right)=\{\left( X,y,S \right)\in {{F}_{++}}: \\ & \|{{\left[ {{\tau }_{1}}\mu I-{{X}^{\frac{1}{2}}}S{{X}^{\frac{1}{2}}} \right]}^{+}}{{\|}_{F}}\le \beta {{\tau }_{1}}\mu \} \\ \end{align}$ |
算法3.1
输入精度ε>0，以及参数τ1≤1/4，β≤1/2.
初始点满足(X0,y0,S0,τ0,κ0)∈N(τ1,β)，并假设k:=0，μ0=(Tr(X0S0)+τ0κ0)/(n+1).
步骤1 如果μk≤εμ0， 则停止.
步骤2 计算NT尺度矩阵P=W1/2，由(23)计算参数η.
步骤3 由系统(15)计算方向(dXk,dyk,dSk,dτk,dκk).
步骤4 计算使(X(αk),y(αk),S(αk),τ(αk),κ(αk))∈N(τ1,β)
成立的步长
步骤5 令
3.4 算法的复杂性分析本节在给出一些引理之后,将证明算法在宽邻域内是收敛的,且迭代复杂度为O(nlogL).
将问题元素尺度化如下:
| $\begin{align} & \left( \hat{X},\hat{y},\hat{S} \right)=\left( PXP,y,{{P}^{\text{-}1}}S{{P}^{\text{-}1}} \right), \\ & (\hat{C},{{{\hat{A}}}_{i}},{{{\hat{b}}}_{i}})=({{P}^{\text{-}1}}C{{P}^{\text{-}1}},{{P}^{\text{-}1}}Ai{{P}^{\text{-}1}},{{b}_{i}})\text{,} \\ & (d\hat{X},d\hat{y},d{{{\hat{S}}}^{)=(}}PdXP,dy,{{P}^{\text{-}1}}dS{{P}^{\text{-}1}}). \\ \end{align}$ |
| $\begin{align} & P(X,S)=\left\{ P\in S_{++}^{n}:{{(PXS{{P}^{\text{-}1}})}^{T}}=PXS{{P}^{\text{-}1}} \right\} \\ & =\left\{ P\in S_{++}^{n}:\hat{X}\hat{S}=\hat{S}\hat{X} \right\} \\ \end{align}$ |
| $(15b)\text{和}(15c)\Leftrightarrow {{{\hat{X}}}^{k}}d\hat{S}+d\hat{X}{{{\hat{S}}}^{k}}=R_{C}^{k}\text{,}$ | (29) |
| $\Leftrightarrow {{H}_{P}}({{X}^{k}}dS+dX{{S}^{k}})=R_{C}^{k},$ | (30) |
| $\begin{align} & {{{\hat{R}}}_{C}}={{\left[ \gamma \mu I-\hat{X}\hat{S} \right]}^{-}}+\sqrt{n+1}{{\left[ \gamma \mu I-\hat{X}\hat{S} \right]}^{+}}, \\ & \mu =Tr\left( XS \right)/\left( n+1 \right). \\ \end{align}$ |
| ${{{\hat{E}}}^{k}}vec(d\hat{X})+{{{\hat{F}}}^{k}}vec(d\hat{S})=vec(\hat{R}_{C}^{k}).$ | (31) |
| ${{{\hat{E}}}^{=}}(\hat{S}\otimes I+I\otimes \hat{S})/2\text{,}{{{\hat{F}}}^{=}}(\hat{X}\otimes I+I\otimes {{{\hat{X}}}^{)}}/2.$ | (32) |
引理3.5[16] 设U,V∈Sn， 则下面不等式成立:
‖(U+V)+‖F≤‖U++V+‖F≤
‖U+‖F+‖V+‖F.
引理3.6[16] 设(X,y,S)∈N(τ1,β)，矩阵
| $\begin{align} & \left\| {{{\hat{E}}}^{-1}} \right.vec[\gamma \mu I\text{-}\hat{X}\hat{S}){{\left. {{]}^{-}} \right\|}^{2}}\le (n+1)\mu \text{,} \\ & \left\| {{{\hat{E}}}^{-1}} \right.vec[\gamma \mu I\text{-}\hat{X}\hat{S}){{\left. {{]}^{+}} \right\|}^{2}}\le \beta {{\tau }_{1}}\mu . \\ \end{align}$ |
| $\begin{align} & \|{{(B{{A}^{T}})}^{-1/2}}Au{{\|}^{2}}+\|{{(B{{A}^{T}})}^{-1/2}}Bv{{\|}^{2}}+ \\ & 2{{u}^{T}}v=\|{{(B{{A}^{T}})}^{-1/2}}r{{\|}^{2}}. \\ \end{align}$ |
| $\frac{\alpha }{\sqrt{n+1}}{{\left\| {{\left[ d\hat{X}d\hat{S} \right]}^{-}} \right\|}_{F}}\le \beta {{\tau }_{1}}\mu \left( \alpha \right).$ |
| $\begin{align} & {{\left\| {{\left[ d\hat{X}d\hat{S} \right]}^{-}} \right\|}_{F}}\le {{\left\| d\hat{X}d\hat{S} \right\|}_{F}}\le {{\left\| d\hat{X} \right\|}_{F}}\le {{\left\| d\hat{S} \right\|}_{F}} \\ & =\left\| vec(d\hat{X}) \right\|\left\| vec(d\hat{S}) \right\| \\ & \le \frac{1}{2}\left( \left\| vec(d\hat{X}) \right\|+\left\| vec(d\hat{S}) \right\| \right), \\ \end{align}$ |
| $\begin{align} & {{\left\| vec(d\hat{X}) \right\|}^{2}}+{{\left\| vec(d\hat{S}) \right\|}^{2}} \\ & \le {{\left\| {{{\hat{E}}}^{-1}}vec{{\left[ \gamma \mu I-\hat{X}\hat{S} \right]}^{-}} \right\|}^{2}}+ \\ & n+1{{\left\| {{{\hat{E}}}^{-1}}vec{{\left[ \gamma \mu I-\hat{X}\hat{S} \right]}^{+}} \right\|}^{2}} \\ & \le \left( 1+\beta {{\tau }_{1}} \right)\left( n+1 \right)\mu , \\ \end{align}$ |
| $\begin{align} & \text{又} Tr(\hat{X}(\alpha )\hat{S}(\alpha )) \\ & =Tr(\hat{X}\hat{S})+\alpha [Tr({{\tau }_{1}}\mu I\text{-}\hat{X}\hat{S})+ \\ & (\sqrt{n+1}\text{-}1)Tr{{({{\tau }_{1}}\mu I\text{-}\hat{X}\hat{S})}^{+}}]. \\ & \mu (\alpha )=\mu +\alpha [({{\tau }_{1}}-1)\mu + \\ & \frac{\sqrt{n+1}-1}{n+1}Tr{{({{\tau }_{1}}\mu I\text{-}\hat{X}\hat{S})}^{+}}] \\ & \ge \mu +\alpha ({{\tau }_{1}}-1)\mu \\ & \ge \mu +{{\tau }_{1}}\mu \text{-}\mu ={{\tau }_{1}}\mu . \\ \end{align}$ |
引理3.9 假设(X,y,S)∈N(τ1,β)，矩阵
(X(α),y(α),S(α))∈N(τ1,β).
证明 因为
| $\begin{align} & {{\left\| {{\left[ {{\tau }_{1}}\mu (\alpha )I\text{-}(\hat{X}\hat{S}+\alpha {{{\hat{R}}}_{C}}) \right]}^{+}} \right\|}_{F}}\le \\ & (1\text{-}\alpha \sqrt{n+1})\beta {{\tau }_{1}}\mu (\alpha ). \\ \end{align}$ |
| $\begin{align} & {{\left\| {{\left[ {{\tau }_{1}}\mu (\alpha )I\text{-X}{{(\alpha )}^{1/2}}\text{S}(\alpha )\text{X}{{(\alpha )}^{1/2}} \right]}^{+}} \right\|}_{F}} \\ & ={{\left\| {{\left[ {{\tau }_{1}}\mu (\alpha )I\text{-}\hat{X}(\alpha )\hat{S}(\alpha ) \right]}^{+}} \right\|}_{F}} \\ & ={{\left\| {{\left[ {{\tau }_{1}}\mu (\alpha )I\text{-}(\hat{X}+\alpha d\hat{X})\hat{S}(\hat{S}+\alpha d\hat{S}) \right]}^{+}} \right\|}_{F}} \\ & ={{\left\| {{\left[ {{\tau }_{1}}\mu (\alpha )I\text{-}(\hat{X}\hat{S}+\alpha {{{\hat{R}}}_{C}}) \right]}^{+}} \right\|}_{F}}+ \\ & {{\alpha }^{2}}{{\left\| {{\left[ d\hat{X}d\hat{S} \right]}^{-}} \right\|}_{F}} \\ & \le \left( 1-\alpha \sqrt{n+1} \right)\beta {{\tau }_{1}}\mu \left( \alpha \right)+\alpha \sqrt{n+1}\beta {{\tau }_{1}}\mu \left( \alpha \right) \\ & =\beta {{\tau }_{1}}\mu \left( \alpha \right). \\ \end{align}$ |
引理3.10[15] 若τ1≤1/4，β≤1/2， 则
μ(α)≤(1－ξα)μ，
其中ξ=1－τ1－βτ1.
定理3.1 当使用NT搜索方向时,算法4.1的迭代复杂度是
证明 由引理3.10得μk≤(1－ξα)μk－1≤(1－ξα)k－1μ0，即Tr(XkSk)≤(1－ξα)kTr(X0S0).
所以,要使Tr(XkSk)≤ε成立,即算法终止,只需保证
(1－ξα)kTr(X0S0)≤ε
即可,即
| $k\log \left( 1-\xi \alpha \right)\le -\log \frac{Tr{{X}^{0}}{{S}^{0}}}{\varepsilon }.$ |
所以,
| $\begin{align} & k\ge \frac{1}{\xi \alpha }\frac{logTr({{X}^{0}}{{S}^{0}})}{\varepsilon }=\frac{(1+\beta {{\tau }_{1}})\sqrt{n+1}}{2\beta {{\tau }_{1}}^{2}\xi }\cdot \\ & \frac{logTr({{X}^{0}}{{S}^{0}})}{\varepsilon } \\ \end{align}$ |
4 数值实验本节将比较算法3.1和Rangarajan[2]中算法3.2. 所有测试是在Windows XP系统下,使用MATLAB R2010(b),测试一些随机产生的半定规划问题,包括随机半定规划(表 1),最大割问题(表 2),教育测试问题(表 3),模极小化问题(表 4),关于4个问题的具体理论详见文献[17].
Table 1
 表 1 随机半定规划Table 1 Random SDP
| 表 1 随机半定规划Table 1 Random SDP | |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Table 2
表 2 最大割问题Table 2 Max-Cut problem
| 表 2 最大割问题Table 2 Max-Cut problem | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Table 3
表 3 教育测试问题Table 3 Educational testing problem
| 表 3 教育测试问题Table 3 Educational testing problem | |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Table 4
表 4 模极小化问题Table 4 Norm minimization problem
| 表 4 模极小化问题Table 4 Norm minimization problem | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
算法3.1中定义:
normp=‖τb-AX‖F，
normd=‖τC-ATy-S‖F.
算法3.2中定义:
normp=‖b-AX‖F，
normd=‖C-ATy-S‖F.
对算法选取最优参数为τ1=0.05，β=0.01. 使用NT尺度化方向,对同一个m和n， 运行10次取平均值结果列在下表.从表 1—表 4可以看出算法3.1的平均迭代次数比算法3.2减少了38.35%.虽然这2个程序是粗糙的,但整体上不影响2个算法的比较.
5 结论通过以上讨论,算法迭代后得到τ和κ的值,结合定理1.1可得知可行点的存在性,即原问题是可行问题或者不可行问题.若原问题是可行问题,结合定理1.1可以得到原问题的最优解.除此结论还可得到如下结论: 采用NT尺度化矩阵和新的宽邻域得到的复杂度比一般算法的复杂度小,同时也具有更好的实验结果,从而证明本文研究的算法是有效的.
本文提出的是齐次不可行内点算法,此外,也可以研究齐次可行内点算法,具体过程与本文十分相似,进而可以求解辅助问题HMCP,计算迭代复杂度,测试实验结果,与本文所提出的算法进行比较,具体可以参考文献[16].
参考文献
| [1] | Liu H, Yang X, Liu C. A new wide neighborhood primal-dual infeasible-interior-point method for symmetric cone programming[J].Journal of Optimization Theory and Applications, 2013, 158(3):796–815.DOI:10.1007/s10957-013-0303-y |
| [2] | Rangarajan B K. Polynomial convergence of infeasible-interior-point methods over symmetric cones[J].Siam Journal on Optimization, 2006, 16(4):1211–1229.DOI:10.1137/040606557 |
| [3] | Andersen E D, Ye Y. On a homogeneous algorithm for the monotone complementarity problem[J].Mathematical Programming, 1995, 84(2):375–399. |
| [4] | Cottle R W, Pang J S, Stone R E. The linear complementarity problem[J].Computer Science and Scientific Computing, 1992, 1:132–133. |
| [5] | Güler O. Existence of interior points and interior paths in nonlinear monotone complementarity problems[J].Mathematics of Operations Research, 1993, 18(1):128–147.DOI:10.1287/moor.18.1.128 |
| [6] | Kojima M, Mizuno S, Yoshise A. A convex property of monotone complementarity problems . Department of Information Sciences, Tokyo Institute of Technology, Tokyo, Japan:1993. |
| [7] | 韩继业, 修乃华, 戚厚铎. 非线性互补理论与算法[M].上海: 上海科学技术出版社, 2006. |
| [8] | Klerk E. Aspects of semidefinite programming[M].New York, Boston, Dordrecht, London, Moscow: Kluwer Academic Publishers, 2004. |
| [9] | Zhang Y. On extending some primal-dual interior-point algorithms from linear programming to semidefinite programming[J].Siam Journal on Optimization, 1998, 8(2):365–386.DOI:10.1137/S1052623495296115 |
| [10] | Shida M, Shindoh S, Kojima M. Existence and uniqueness of search directions in interior-point algorithms for the SDP and the monotone SDLCP[J].Siam Journal on Optimization, 1998, 8(2):387–396.DOI:10.1137/S1052623496300611 |
| [11] | Ai W, Zhang S. An O( n L) iteration primal-dual path-following method, based on wide neighborhoods and large updates, for monotone LCP[J].Siam Journal on Optimization, 2005, 16(2):400–417.DOI:10.1137/040604492 |
| [12] | Li Y, Terlaky T. A new class of large neighborhood path-following interior point algorithms for semidefinite optimization with O n log Tr(X0S0) ε iteration complexity[J].Siam Journal on Optimization, 2010, 20(6):2853–2875.DOI:10.1137/080729311 |
| [13] | 刘新泽. 对称锥互补问题若干内点算法的复杂性研究 . 西安: 西安电子科技大学, 2014. |
| [14] | Jin S, Ariyawansa K A, Zhu Y. Homogeneous self-dual algorithms for stochastic semidefinite programming[J].Journal of Optimization Theory and Applications, 2012, 155(3):1073–1083.DOI:10.1007/s10957-012-0110-x |
| [15] | 杨喜美. 对称锥规划的宽邻域内点算法研究 . 西安: 西安电子科技大学, 2014. |
| [16] | 刘长河. 锥规划中若干内点算法的复杂性研究 .西安: 西安电子科技大学, 2012.http://cdmd.cnki.com.cn/article/cdmd-10701-1013114233.htm |
| [17] | Todd M J, Toh K C, T R H. On the Nesterov-Todd direction in semidefinite programming[J].Siam Journal on Optimization, 1997, 8(3):769–796. |
