SDP半正定规划的低复杂度求解：基于块坐标下降（Block Coordinate Descent）

前言

之前的几篇博客 经典的SDR算法： 用半正定松弛法 ( Semidefinite Relaxation) 求解二次优化问题 和 经典的SDR算法（下）：SDR的具体使用细节与相关代码 中介绍了一种行之有效的 QCQP问题的求解方法。 这其中， SDP 半正定规划 是 无可避免的必由之路。 然而，传统的CVX求解方法， 如内点法等， 其复杂度为 O ( n 3.5 log ⁡ ( 1 / ϵ ) ) O\left(n^{3.5} \log (1 / \epsilon)\right) O(n3.5log(1/ϵ))， 其中 n n n 为变量维度， ϵ \epsilon ϵ 为目标精度。 可以看出， 这在现有算法中，绝不能算是低复杂度的算法。 而 SDR 本身的性能又是次优的， 这就令其实际应用显得非常尴尬。 这篇博客， 笔者希望介绍一种 基于 块坐标下降方法的 SDP 求解算法， 旨在降低计算复杂度。

主要参考文章为：《PHASE RECOVERY, MAXCUT AND COMPLEX SEMIDEFINITE PROGRAMMING》

块坐标下降

对于 块坐标下降法， 其实可以简单地理解为大型的交替优化——固定一些变量而优化另一些变量。 对于 SDP问题， 我们先给出典型的问题形式如下：

minimize ⁡ Tr ⁡ ( U M )  subject to  U ⪰ 0 \begin{array}{ll} \operatorname{minimize} & \operatorname{Tr}(U M) \\ \text { subject to } & U \succeq 0 \end{array} minimize subject to ​Tr(UM)U⪰0​ 这个限制条件非常棘手。 因此， 我们将其转化为如下的无约束问题：  minimize  Tr ⁡ ( U M ) − μ log ⁡ det ⁡ ( U ) (1) \text { minimize } \operatorname{Tr}(U M)-\mu \log \operatorname{det}(U)\tag{1}  minimize Tr(UM)−μlogdet(U)(1) 其中 μ > 0 \mu>0 μ>0 为 障碍系数。 这是基于内点法的思想：通过在目标函数中加入障碍函数，从而自然地限制解在可行集中。 注意到， 原可行集 U ⪰ 0 U\succeq 0 U⪰0 等效为： U U U的所有特征值 u 1 , ⋯   , u n u_1,\cdots, u_n u1​,⋯,un​非负。 又根据行列式的性质，有： d e t ( U ) = Π i = 1 n u i \mathrm{det}(U)=\Pi_{i=1}^{n}{u_i} det(U)=Πi=1n​ui​， 即有： log ⁡ det ⁡ ( U ) = ∑ i = 1 n l o g ( u i ) . \log \operatorname{det}(U)= \sum_{i=1}^n\mathrm{log}(u_i). logdet(U)=i=1∑n​log(ui​). 而根据 l o g \mathrm{log} log的定义域， 可知，这个障碍函数隐含了所有特征值非负的条件。 接下来，对问题进行拆解来降低求解的复杂度。 注意到存在如下等式： det ⁡ ( U ) = det ⁡ ( B ) det ⁡ ( y − x T B − 1 x ) (2) \operatorname{det}(U)=\operatorname{det}(B) \operatorname{det}\left(y-x^{T} B^{-1} x\right)\tag{2} det(U)=det(B)det(y−xTB−1x)(2) 其中 U = ( B x x T y ) (3) U=\left(\begin{array}{cc} B & x \\ x^{T} & y \end{array}\right)\tag{3} U=(BxT​xy​)(3) （注意， U U U是共轭对称矩阵， 否则特征值可能为虚数， 因此， U U U必能写成上式的形式）。 (2) 的证明放在后面的附录中， 我们先继续往下： 这样转化的目的非常明确， 我们每次优化 U U U 中的部分变量即 x x x。那么（2）已经把障碍函数中和 x x x相关的项提取出来了，即 det ⁡ ( y − x T B − 1 x ) \operatorname{det}\left(y-x^{T} B^{-1} x\right) det(y−xTB−1x)。 经此， 事实上(1)中与 x x x有关的部分，可以体现为如下的优化问题： min ⁡ x c T x − μ log ⁡ ( 1 − x T B − 1 x ) \min _{x} c^{T} x-\mu \log \left(1-x^{T} B^{-1} x\right) xmin​cTx−μlog(1−xTB−1x) 而这个问题是一个凸函数的最小化问题—— l o g ( x ) \mathrm{log}(x) log(x)为非减的凹函数， 1 − x T B − 1 x 1-x^TB^{-1}x 1−xTB−1x 为凹函数， 那么根据复合函数的结论 https://zhuyulab.blog.csdn.net/article/details/121102961， 可知 log ⁡ ( 1 − x T B − 1 x ) \log \left(1-x^{T} B^{-1} x\right) log(1−xTB−1x) 为凹。 那么 x x x 的最优解可以直接由一阶条件给出， 即目标函数对 x x x 求导为0， 即： c + 2 μ 1 − x T B − 1 x B − 1 x = 0 c + \frac{2\mu}{1-x^{T} B^{-1} x}B^{-1}x=0 c+1−xTB−1x2μ​B−1x=0 那么 x x x 的形式为 x = α B c x =\alpha Bc x=αBc （否则右边不为0）, 其中 α \alpha α 为待定常数。 代入得： c + 2 μ α 1 − α 2 c T B c c = 0 ⇒ 2 μ α 1 − α 2 c T B c = − 1 ⇒ α 2 c T B c − 2 μ α − 1 = 0 c + \frac{2\mu\alpha}{1-\alpha^2c^TBc}c=0\Rightarrow \frac{2\mu\alpha}{1-\alpha^2c^TBc}=-1\Rightarrow \alpha^2c^TBc-2\mu\alpha-1=0 c+1−α2cTBc2μα​c=0⇒1−α2cTBc2μα​=−1⇒α2cTBc−2μα−1=0 令 γ = c T B c \gamma = c^TBc γ=cTBc, 由一元二次方程根可知， α = μ − μ + γ γ \alpha=\frac{\mu-\sqrt{\mu +\gamma}}{\gamma} α=γμ−μ+γ ​​ 这里需要指出， 为什么不取另一个解 α = μ + μ + γ γ \alpha=\frac{\mu+\sqrt{\mu +\gamma}}{\gamma} α=γμ+μ+γ ​​ 呢？ 这是因为， 我们必须满足 1 − x T B − 1 x > 0 1-x^{T} B^{-1} x>0 1−xTB−1x>0。 至此， x x x 的 闭式解为： x = μ − μ + γ γ B c (4) x = \frac{\mu-\sqrt{\mu +\gamma}}{\gamma}Bc\tag{4} x=γμ−μ+γ ​​Bc(4)。

至此 对于 问题（1） 中的优化， 我们可以对其每一列（事实上是每一列去掉一个元素， 以契合(3)中 x x x的形式）元素进行依次优化，来求解 SDP问题。

最后给出伪代码： 不得不说有点抽象， 我简单解释一下： 这里 X X X 就是 (1) 中的变量 U U U。 i c i_c ic​即列索引中去掉 i i i后的集合。 如 X i c , i c k X_{i^{c}, i^{c}}^{k} Xic,ick​ 就是指把 X X X 矩阵的第 i i i行第 i i i列去除后的结果。 X i c , i X_{i^c,i} Xic,i​ 其实就对应 每次的优化变量 x x x。 step 4 的结果可以看做是 (4) 的 简化形式。

证明：

det ⁡ ( U ) = det ⁡ ( B ) det ⁡ ( y − x T B − 1 x ) \operatorname{det}(U)=\operatorname{det}(B) \operatorname{det}\left(y-x^{T} B^{-1} x\right) det(U)=det(B)det(y−xTB−1x)

根据分块矩阵乘法，有： det ⁡ [ A B C D ] = det ⁡ ( [ A O C D − C A − 1 B ] [ I A − 1 B O I ] ) = det ⁡ ( A ) ⋅ det ⁡ ( D − C A − 1 B ) \begin{aligned} \operatorname{det}\left[\begin{array}{ll} \boldsymbol{A} & \boldsymbol{B} \\ \boldsymbol{C} & \boldsymbol{D} \end{array}\right] &=\operatorname{det}\left(\left[\begin{array}{cc} \boldsymbol{A} & \boldsymbol{O} \\ \boldsymbol{C} & \boldsymbol{D}-\boldsymbol{C} \boldsymbol{A}^{-1} \boldsymbol{B} \end{array}\right]\left[\begin{array}{cc} \boldsymbol{I} & \boldsymbol{A}^{-1} \boldsymbol{B} \\ \boldsymbol{O} & \boldsymbol{I} \end{array}\right]\right) \\ &=\operatorname{det}(\boldsymbol{A}) \cdot \operatorname{det}\left(\boldsymbol{D}-\boldsymbol{C A}^{-1} \boldsymbol{B}\right) \end{aligned} det[AC​BD​]​=det([AC​OD−CA−1B​][IO​A−1BI​])=det(A)⋅det(D−CA−1B)​ 得证。