A block symmetric Gauss-Seidel decomposition theorem for convex composite quadratic programming and its applications

TL;DR

This paper establishes a decomposition theorem for the block symmetric Gauss-Seidel method, linking it to convex quadratic programming, and extends it to solve convex composite quadratic problems with inexact and accelerated variants.

Contribution

The paper introduces the block sGS decomposition theorem and extends the classical block sGS method to convex composite quadratic programming, including inexact and accelerated versions.

Findings

Exact solution of quadratic programming via block sGS cycles

Extension to convex composite quadratic programming with nonsmooth terms

Achieves accelerated convergence rate of O(1/k^2) with inexact computation

Abstract

For a symmetric positive semidefinite linear system of equations $\mathcal{Q} {\bf x} = {\bf b}$, where ${\bf x} = (x_1,\ldots,x_s)$ is partitioned into $s$ blocks, with $s \geq 2$, we show that each cycle of the classical block symmetric Gauss-Seidel (block sGS) method exactly solves the associated quadratic programming (QP) problem but added with an extra proximal term of the form $\frac{1}{2} \| {\bf x}-{\bf x}^k \|_{\mathcal T}^2$, where ${\mathcal T}$ is a symmetric positive semidefinite matrix related to the sGS decomposition and ${\bf x}^k$ is the previous iterate. By leveraging on such a connection to optimization, we are able to extend the result (which we name as the block sGS decomposition theorem) for solving a convex composite QP (CCQP) with an additional possibly nonsmooth term in $x_1$, i.e., $\min\{ p(x_1) + \frac{1}{2}\langle {\bf x},\, \mathcal{Q} {\bf x} \rangle -\langle {\bf b},\, {\bf x}\rangle\}$, where $p(\cdot)$ is a proper closed convex function. Based on the block sGS decomposition theorem, we are able to extend the classical block sGS method to solve a CCQP. In addition, our extended block sGS method has the flexibility of allowing for inexact computation in each step of the block sGS cycle. At the same time, we can also accelerate the inexact block sGS method to achieve an iteration complexity of $O(1/k^2)$ after performing $k$ block sGS cycles. As a {fundamental} building block, the block sGS decomposition theorem has played a key role in various recently developed algorithms such as the inexact semiproximal {ALM/ADMM} for linearly constrained multi-block convex composite conic programming (CCCP), and the accelerated block coordinate descent method for multi-block CCCP.

Equations

$$\displaystyle{\cal Q}\mbox{\boldmath{$x$}}=\mbox{\boldmath{$b$}},$$

$$\mathcal{Q}$$

$$\displaystyle\mbox{(QP)}\quad\min\Big{\{}q(\mbox{\boldmath{$x$}}):=\frac{1}{2}\langle\mbox{\boldmath{$x$}},\,{\cal Q}\mbox{\boldmath{$x$}}\rangle-\langle\mbox{\boldmath{$b$}},\,\mbox{\boldmath{$x$}}\rangle\mid\mbox{\boldmath{$x$}}\in{\cal X}\Big{\}}.$$

$$\mathcal{Q}$$

$$\displaystyle{\cal U}=\left[\begin{array}[]{cccc}\mbox{\bf 0}&Q_{1,2}&\dots&Q_{1,s}\\ &\ddots&&\vdots\\ &&\ddots&Q_{s-1,s}\\[5.0pt] &&&\mbox{\bf 0}\end{array}\right],\quad{\cal D}=\left[\begin{array}[]{cccc}Q_{1,1}&&\\[3.0pt] &Q_{2,2}&&\\ &&\ddots&\\ &&&Q_{s,s}\end{array}\right],$$

$$\displaystyle\widehat{{\cal Q}}(\mbox{\boldmath{$x$}}^{k+1}-\mbox{\boldmath{$x$}}^{k})=\mbox{\boldmath{$b$}}-{\cal Q}\mbox{\boldmath{$x$}}^{k},$$

$$\displaystyle\mbox{(CCQP)}\quad\min\Big{\{}F(\mbox{\boldmath{$x$}})\;:=\;p(x_{1})+\frac{1}{2}\langle\mbox{\boldmath{$x$}},\,{\cal Q}\mbox{\boldmath{$x$}}\rangle-\langle\mbox{\boldmath{$b$}},\,\mbox{\boldmath{$x$}}\rangle\mid\mbox{\boldmath{$x$}}=[x_{1};\ldots;x_{s}]\in{\cal X}\Big{\}},$$

$$\displaystyle\min\Big{\{}p(x_{1})+\frac{1}{2}\langle\mbox{\boldmath{$x$}},\,{\cal Q}\mbox{\boldmath{$x$}}\rangle-\langle\mbox{\boldmath{$b$}},\,\mbox{\boldmath{$x$}}\rangle+\frac{1}{2}\|\mbox{\boldmath{$x$}}-\mbox{\boldmath{$x$}}^{k}\|_{{\cal T}_{\cal Q}}^{2}\mid\mbox{\boldmath{$x$}}\in{\cal X}\Big{\}},$$

$${\mathcal{T}}_{\mathcal{Q}}=\mathcal{U}{\mathcal{D}}^{-1}{\mathcal{U}}^{\ast }.$$

$$\displaystyle\mbox{\boldmath{$x$}}^{+}:=\mbox{argmin}_{\mbox{\boldmath{$x$}}\in{\cal X}}\;\Big{\{}p(x_{1})+q(\mbox{\boldmath{$x$}})+\frac{1}{2}\|\mbox{\boldmath{$x$}}-\bar{}\mbox{\boldmath{$x$}}\|_{{\cal T}_{\cal Q}}^{2}-\langle\mbox{\boldmath{$x$}},\,\Delta(\mbox{\boldmath{$\delta^{\prime}$}},\mbox{\boldmath{$\delta$}})\rangle\Big{\}},$$

$$\displaystyle\Delta(\mbox{\boldmath{$\delta^{\prime}$}},\mbox{\boldmath{$\delta$}}):=\mbox{\boldmath{$\delta$}}+{\cal U}{\cal D}^{-1}(\mbox{\boldmath{$\delta$}}-\mbox{\boldmath{$\delta^{\prime}$}}).$$

$$x_{\ge i}=(x_i;\dots ;x_s),x_{\le i}=(x_1;\dots ;x_i),i=1,\dots ,s.$$

$$\hat{\mathcal{Q}}:=\mathcal{Q}+{\mathcal{T}}_{\mathcal{Q}}=(\mathcal{D}+\mathcal{U}){\mathcal{D}}^{-1}(\mathcal{D}+{\mathcal{U}}^{\ast })\succ 0.$$

$$x_i^{\prime }:=$$

$$=$$

$$\left\{\begin{aligned} x_{1}^{+}={}&\mathop{\rm argmin}_{x_{1}\in{\cal X}_{1}}\;p(x_{1})+q(x_{1};x^{\prime}_{\geq 2})-\langle\delta_{1},\,x_{1}\rangle,\\[5.0pt] x_{i}^{+}={}&\mathop{\rm argmin}_{x_{i}\in{\cal X}_{i}}\;p(x_{1}^{+})+q(x^{+}_{\leq i-1};x_{i};x^{\prime}_{\geq i+1})-\langle\delta_{i},\,x_{i}\rangle\\[5.0pt] ={}&Q_{ii}^{-1}\big{(}b_{i}+\delta_{i}-\mbox{$\sum_{j=1}^{i-1}$}Q_{ji}^{*}x_{j}^{+}-\mbox{$\sum_{j=i+1}^{s}$}Q_{ij}x^{\prime}_{j}\big{)},\quad i=2,\ldots,s.\end{aligned}\right.$$

$$\mathcal{Q}+{\mathcal{T}}_{\mathcal{Q}}=\mathcal{D}+\mathcal{U}+{\mathcal{U}}^{\ast }+\mathcal{U}{\mathcal{D}}^{-1}{\mathcal{U}}^{\ast }=(\mathcal{D}+\mathcal{U}){\mathcal{D}}^{-1}(\mathcal{D}+{\mathcal{U}}^{\ast }).$$

$$x_1^{\prime }={\mathrm{argmin}}_{x_1\in {\mathcal{X}}_1}p(x_1)+q(x_1;x_{\ge 2}^{\prime })-⟨{\delta }_1^{\prime },x_1⟩={\mathrm{argmin}}_{x_1\in {\mathcal{X}}_1}p(x_1)+q(x_1;x_{\ge 2}^{\prime })-⟨{\delta }_1,x_1⟩=x_1^{+}.$$

$$\displaystyle Q_{11}x_{1}^{\prime}=b_{1}-\gamma_{1}+\delta_{1}^{\prime}-\mbox{$\sum_{j=2}^{s}$}Q_{ij}x^{\prime}_{j},$$

$$\displaystyle Q_{11}x_{1}^{+}=b_{1}-\gamma_{1}+\delta_{1}-\mbox{$\sum_{j=2}^{s}$}Q_{ij}x^{\prime}_{j},$$

$$({\cal D}+{\cal U})\mbox{\boldmath{$x$}}^{\prime}=\mbox{\boldmath{$b$}}-\mbox{\boldmath{$\gamma$}}+\mbox{\boldmath{$\delta$}}^{\prime}-{\cal U}^{*}\bar{}\mbox{\boldmath{$x$}},$$

$$({\cal D}+{\cal U}^{*})\mbox{\boldmath{$x$}}^{+}=\mbox{\boldmath{$b$}}-\mbox{\boldmath{$\gamma$}}+\mbox{\boldmath{$\delta$}}-{\cal U}\mbox{\boldmath{$x$}}^{\prime}.$$

$$\displaystyle({\cal D}+{\cal U}^{*})\mbox{\boldmath{$x$}}^{+}={}$$

$$=$$

$$({\cal Q}+{\cal T}_{{\cal Q}})\mbox{\boldmath{$x$}}^{+}=\mbox{\boldmath{$b$}}-\mbox{\boldmath{$\gamma$}}+{\cal T}_{\cal Q}\bar{}\mbox{\boldmath{$x$}}+\Delta(\mbox{\boldmath{$\delta$}}^{\prime},\mbox{\boldmath{$\delta$}}).$$

$$\xi\leq\|{\cal D}^{-1/2}(\mbox{\boldmath{$\delta$}}-\mbox{\boldmath{$\delta$}}^{\prime})\|+\|\widehat{{\cal Q}}^{-1/2}\mbox{\boldmath{$\delta$}}^{\prime}\|.$$

$${\hat{\mathcal{Q}}}^{-1}=(\mathcal{D}+{\mathcal{U}}^{\ast }{)}^{-1}\mathcal{D}(\mathcal{D}+\mathcal{U}{)}^{-1}=(\mathcal{D}+{\mathcal{U}}^{\ast }{)}^{-1}{\mathcal{D}}^{1/2}{\mathcal{D}}^{1/2}(\mathcal{D}+\mathcal{U}{)}^{-1},$$

$$\xi =$$

$$x_1^{+}$$

$$\widehat{\Theta}_{j}:=\left[\begin{array}[]{c}Q_{1,j}\\ \vdots\\ Q_{j-1,j}\end{array}\right]Q_{j,j}^{-1}\;\left[Q_{1,j}^{*},\;\dots,\;Q_{j-1,j}^{*}\right]$$

$${\Theta}_{j}:=\left[\begin{array}[]{c}Q_{1,2}\\ 0\\ \vdots\\ 0\end{array}\right]Q_{2,2}^{-1}\;\left[Q_{1,2}^{*},\;0,\;\dots,\;0\right]+\dots+\left[\begin{array}[]{c}Q_{1,j}\\ \vdots\\ Q_{j-1,j}\\ 0\end{array}\right]Q_{j,j}^{-1}\;\left[Q_{1,j}^{*},\;\dots,\;Q_{j-1,j}^{*},\;0\right].$$

Full text

Keywords: Convex composite quadratic programming, block symmetric Gauss-Seidel, Schur complement, augmented Lagrangian method

AMS subject classifications: 90C06, 90C20, 90C25, 65F10

1 Introduction

It is well known that the classical block symmetric Gauss-Seidel (block sGS) method [1, 7, 13, 23] can be used to solve a symmetric positive semidefinite linear system of equations {\cal Q}\mbox{\boldmath{x}}=\mbox{\boldmath{b}} where \mbox{\boldmath{x}}=(x_{1};\ldots;x_{s}) is partitioned into $s$ blocks with $s\geq 2$. We are particularly interested in the case when $s>2$. In this paper, we show that each cycle of the classical block sGS method exactly solves the corresponding convex quadratic programming (QP) problem but added with an extra proximal term depending on the previous iterate (say \mbox{\boldmath{x}}^{k}). Through such a connection to optimization, we are able to extend the result (which we name as the block sGS decomposition theorem) to a convex composite QP (CCQP) with an additional possibly nonsmooth term in $x_{1}$, and subsequently extend the classical block sGS method to solve a CCQP. We can also extend the classical block sGS method to the inexact setting, where the underlying linear system for each block of the new iterate \mbox{\boldmath{x}}^{k+1} need not be solved exactly. Moreover, by borrowing ideas in the optimization literature, we are able to accelerate the classical block sGS method and provide new convergence results. More details will be given later.

Assume that ${\cal X}_{i}=\mathbb{R}^{n_{i}}$ for $i=1,\ldots,s$, and ${\cal X}={\cal X}_{1}\times\cdots\times{\cal X}_{s}$, where $s\geq 2$ is a given integer. Consider the following symmetric positive semidefinite block linear system of equations:

[TABLE]

where \mbox{\boldmath{x}}=[x_{1};\;\ldots;\;x_{s}]\in{\cal X}, \mbox{\boldmath{b}}=[b_{1};\;\ldots;\;b_{s}]\in{\cal X}, and

[TABLE]

with $Q_{i,j}\in\mathbb{R}^{n_{i}\times n_{j}}$ for $1\leq i,j\leq s$. It is well known that (1) is the optimality condition for the following unconstrained QP:

[TABLE]

Note that even though our problem is phrased in the matrix-vector setting for convenience, one can consider the setting where each ${\cal X}_{i}$ is a real $n_{i}$-dimensional inner product space and $Q_{i,j}$ is a linear map from ${\cal X}_{i}$ to ${\cal X}_{j}$. Throughout the paper, we make the following assumption:

Assumption 1**.**

${\cal Q}$* is symmetric positive semidefinite and each diagonal block $Q_{i,i}$ is symmetric positive definite for $i=1,\ldots,s$.

From the following decomposition of ${\cal Q}$:

[TABLE]

where

[TABLE]

the classical block sGS iteration in numerical analysis is usually derived as a natural generalization of the classical pointwise sGS for solving a symmetric positive definite linear system of equations, and the latter is typically derived as a fixed-point iteration for the sGS matrix splitting based on (7); see for example [23, Sec. 4.1.1], [13, Sec. 4.5]. Specifically, the block sGS fixed-point iteration in the third normal form (in the terminology used in [13]) reads as follows:

[TABLE]

where $\widehat{{\cal Q}}=({\cal D}+{\cal U}){\cal D}^{-1}({\cal D}+{\cal U}^{*})$.

In this paper, we give a derivation of the classical block sGS method (17) from the optimization perspective. By doing so, we are able to extend the classical block sGS method to solve a structured CCQP problem of the form:

[TABLE]

where $p:{\cal X}_{1}\rightarrow(-\infty,\infty]$ is a proper closed convex function such as $p(x_{1})=\|x_{1}\|_{1}$ or $p(x_{1})=\delta_{\mathbb{R}^{n_{1}}_{+}}(x_{1})$ (the indicator function of $\mathbb{R}^{n_{1}}_{+}$ defined by $\delta_{\mathbb{R}^{n_{1}}_{+}}(x_{1})=0$ if $x_{1}\in\mathbb{R}^{n_{1}}_{+}$ and $\delta_{\mathbb{R}^{n_{1}}_{+}}(x_{1})=\infty$ otherwise). Our specific contributions are described in the next few paragraphs. We note that the main results presented here are parts of the thesis of the first author [18].

First, we establish the key result of the paper, the block sGS decomposition theorem, which states that each cycle of the block sGS method, say at the $k$th iteration, corresponds exactly to solving (18) with an additional proximal term \frac{1}{2}\|\mbox{\boldmath{x}}-\mbox{\boldmath{x}}^{k}\|^{2}_{{\cal T}_{\cal Q}} added to its objective function, i.e.,

[TABLE]

where ${\cal T}_{\cal Q}={\cal U}{\cal D}^{-1}{\cal U}^{*}$, and \|\mbox{\boldmath{x}}\|_{{\cal T}_{\cal Q}}^{2}=\langle\mbox{\boldmath{x}},\,{\cal T}_{\cal Q}\mbox{\boldmath{x}}\rangle. It is clear that when $p(\cdot)\equiv 0$, the problem (18) is exactly the QP (6) associated with the linear system (1). Therefore, we can interpret the classical block sGS method as a proximal-point minimization method for solving the QP (6), and each cycle of the classical block sGS method solves exactly the proximal subproblem (19) associated with the QP (6). As far as we are aware of, this is the first time in which the classical block sGS method (17) (and also the pointwise sGS method) is derived from an optimization perspective.

Second, we also establish a factorization view of the block sGS decomposition theorem and show its equivalence to the Schur complement based (SCB) reduction procedure proposed in [17] for solving a recursively defined variant of the proximal subproblem (19). The SCB reduction procedure in [17] is derived by inductively finding an appropriate proximal term to be added to the objective function of (18) so that the block variables $x_{s},x_{s-1},\ldots,x_{2}$ can be eliminated in a sequential manner and thus ending with a minimization problem involving only the variable $x_{1}$. In a nutshell, we show that the SCB reduction procedure sequentially eliminates the blocks (in the reversed order starting from $x_{s}$) in the variable $x$ of the proximal subproblem (19) by decomposing the proximal term \frac{1}{2}\|\mbox{\boldmath{x}}-\mbox{\boldmath{x}}^{k}\|^{2}_{{\cal T}_{\cal Q}} also in a sequential manner. In turn, each of the reduction step corresponds exactly to one step in a cycle of the block sGS method.

Third, based on the block sGS decomposition theorem, we are able to extend the classical block sGS method for solving the QP (6) to solve the CCQP (18), and each cycle of the extended block sGS method corresponds precisely to solving the proximal subproblem (19). Our extension of the block sGS method has thus overcome the limitation of the classical method by allowing us to solve the nonsmooth CCQP which often arises in practice, for example, in semidefinite programming where $p(x_{1})=\delta_{\mathbb{S}^{n_{1}}_{+}}(x_{1})$ and $\mathbb{S}^{n_{1}}_{+}$ is the cone of $n_{1}\times n_{1}$ symmetric positive semidefinite matrices. Moreover, our extension also allows the updates of the blocks to be inexact. As a consequence, we also obtain an inexact version of the classical block sGS method, where the iterate \mbox{\boldmath{x}}^{k+1} need not be computed exactly from (17). We should emphasize that the inexact block sGS method is potentially very useful when a diagonal block, say $Q_{i,i}$, in (5) is large and the computation of \mbox{\boldmath{x}}^{k+1}_{i} must be done via an iterative solver rather than a direct solver. Note that even for the linear system (17), our systematic approach (in section 4) to derive the inexact extension of the classical block sGS method appears to be new. The only inexact variant of the classical block sGS method for (17) with a convergence proof we are aware of is the pioneering work of Bank et al. in [3]. In [3], the authors showed that by modifying the diagonal blocks in ${\cal D}$, the linear system involved in each block can be solved by a given fixed number of pointwise sGS cycles.

Fourth, armed with the optimization interpretation of each cycle of the block sGS method, it becomes easy for us to adapt ideas from the optimization literature to establish the iteration complexity of O(\|\mbox{\boldmath{x}}^{0}-\mbox{\boldmath{x}}^{*}\|_{{\widehat{{\cal Q}}}}^{2}/k) for the extended block sGS method as well as to accelerate it to obtain the complexity of O(\|\mbox{\boldmath{x}}^{0}-\mbox{\boldmath{x}}^{*}\|_{{\widehat{{\cal Q}}}}^{2}/(k+1)^{2}), after running for $k$ cycles, where \mbox{\boldmath{x}}^{*} is an optimal solution for (18). Just as in the classical block sGS method, we can obtain a linear rate of convergence for our extended inexact block sGS method under the assumption that ${\cal Q}$ is positive definite. With the help of an extensive optimization literature on the linear convergences of proximal gradient methods, we are further able to relax the positive definiteness assumption on ${\cal Q}$ to a mild error bound assumption on the function $F$ in (18) and derive at least R-linear convergence results for our extended block sGS method. The error bound assumption in fact holds automatically for many interesting applications, including the important case when $p(\cdot)$ is a piecewise linear-quadratic function. We note that there is active research in studying the convergence of proximal gradient methods for a convex composite minimization problem of the form \min\{f(\mbox{\boldmath{x}})+g(\mbox{\boldmath{x}})\mid\mbox{\boldmath{x}}\in{\cal X}\}, with $f$ being a smooth convex function and $g$ a proper closed convex function whose proximal map is easy to compute; see for example [24] and the references therein. In each iteration of a typical proximal gradient method, a simple proximal term \frac{L}{2}\|\mbox{\boldmath{x}}-\bar{\mbox{\boldmath{x}}}\|^{2}, where $L$ is a Lipschitz constant for the gradient of $f$, is added to the objective function. Our extended block sGS method for (CCQP) differs from those proximal gradient methods in the literature in that the proximal term we add comes from the sophisticated positive semidefinite linear operator associated with the sGS decomposition of ${\cal Q}.$

Recent research works in [6, 16, 17, 25, 26] have shown that our block sGS decomposition theorem for the CCQP (18) can play an essential role in the design of efficient algorithms for solving various convex optimization problems such as convex composite quadratic semidefinite programming problems. Indeed, the block sGS decomposition based ADMM algorithms designed in [6, 25, 26] have found applications in various recent papers such as [2, 10, 15]. Our experiences have shown that the inexact block sGS cycle can provide the much needed computational efficiency when one is designing an algorithm based on the framework of the proximal augmented Lagrangian (ALM) or proximal alternating direction method of multipliers (ADMM) for solving important classes of large scale convex composite optimization problems. As a concrete illustration of the application of our block sGS decomposition theorem, we will briefly describe in section 5 on how to utilize the theorem in the design of the proximal augmented Lagrangian method for solving a linearly constrained convex composite quadratic programming problem.

The idea of sequentially updating the blocks of a multi-block variable, either in the Gauss-Seidel fashion or the successive over-relaxation (SOR) fashion, has been incorporated into quite a number of optimization algorithms [5] and in solving nonlinear equations [22]. Indeed the Gauss-Seidel (also known as the block coordinate descent) approach for solving optimization problems has been considered extensively; we refer the readers to [4, 12] for the literature review on the recent developments, especially for the case where $s>2$. Here we would like to emphasize that even for the case of an unconstrained smooth convex minimization problem \min\{f(\mbox{\boldmath{x}})\mid\mbox{\boldmath{x}}\in{\cal X}\}, whose objective function f(\mbox{\boldmath{x}}) (not necessarily strongly convex) has a Lipschitz continuous gradient of modulus $L$, it is only proven recently in [4] that the block coordinate (gradient) descent method is globally convergent with the iteration complexity of $O(Ls/k)$ after $k$ cycles, where $s$ is the number of blocks. When f(\mbox{\boldmath{x}}) is the quadratic function in (6), the block coordinate descent method is precisely the classical block Gauss-Seidel (GS) method. In contrast to the block sGS method, each iteration of the block GS method does not appear to have an optimization equivalence. Despite the extensive work on the Gauss-Seidel approach for solving convex optimization problems, surprisingly, little is known about the symmetric Gauss-Seidel approach for solving the same problems except for the recent paper [25] which utilized our block sGS decomposition theorem to design an inexact accelerated block coordinate descent method to solve a problem of the form \min\{p(x_{1})+f(\mbox{\boldmath{x}})\mid\mbox{\boldmath{x}}\in{\cal X}\}.

The remaining parts of the paper are organized as follows. The next section is devoted to the block sGS decomposition theorem for the CCQP (18). In section 3, we present a factorization view of the block sGS theorem and prove its equivalence to the SCB reduction procedure proposed in [17, 18]. In the following section, we derive the block sGS method from an optimization perspective and extend it to solve the CCQP (18). The convergence results for our extended block sGS method are also presented in this section. In section 5, the application of our block sGS decomposition theorem is demonstrated in the design of a proximal augmented Lagrangian method for solving a linearly constrained convex composite quadratic programming problem. The extension of the classical block symmetric SOR method for solving (18) is presented in section 6. We conclude our paper in the final section.

We end the section by giving some notation. For a symmetric matrix ${\cal Q}$, the notation ${\cal Q}\succeq 0$ (${\cal Q}\succ 0$) means that the matrix ${\cal Q}$ is symmetric positive semidefinite (definite). The spectral norm of ${\cal Q}$ is denote by $\|{\cal Q}\|_{2}.$

2 Derivation of the block sGS decomposition theorem for (18)

In this section, we present the derivation of one cycle of the block sGS method for (18) from the optimization perspective as mentioned in the introduction.

Recall the decomposition of ${\cal Q}$ in (7), ${\cal U}$,${\cal D}$ in (16) and the sGS linear operator defined by

[TABLE]

Given \bar{}\mbox{\boldmath{x}}\in{\cal X}, corresponding to problem (18), we consider solving the following subproblem

[TABLE]

where \mbox{\boldmath{\delta^{\prime}}},\mbox{\boldmath{\delta}}\in{\cal X} are two given error vectors with $\delta_{1}^{\prime}=\delta_{1}$, and

[TABLE]

We note that the vectors \mbox{\boldmath{\delta^{\prime}}},\mbox{\boldmath{\delta}} need not be known a priori. We should view \mbox{\boldmath{x}}^{+} as an approximate solution to (21) without the perturbation term \langle\mbox{\boldmath{x}},\,\Delta(\mbox{\boldmath{\delta^{\prime}}},\mbox{\boldmath{\delta}})\rangle. Once \mbox{\boldmath{x}}^{+} has been computed, the associated error vectors can then be obtained, and \mbox{\boldmath{x}}^{+} is then the exact solution to the perturbed problem (21).

The following theorem shows that \mbox{\boldmath{x}}^{+} can be computed by performing exactly one cycle of the block sGS method for (18). In particular, if $p(x_{1})\equiv 0$ and \mbox{\boldmath{\delta}}^{\prime}=0=\mbox{\boldmath{\delta}}, then the computation of \mbox{\boldmath{x}}^{+} corresponds exactly to one cycle of the classical block sGS method. For the proof, we need to define the following notation for a given \mbox{\boldmath{x}}=(x_{1};\ldots;x_{s}),

[TABLE]

We also define $x_{\geq s+1}=\emptyset$.

Theorem 1** (sGS Decomposition).**

Assume that ${\cal Q}\succeq 0$ and the self-adjoint linear operators ${\cal Q}_{ii}$ are positive definite for all $i=1,\ldots,s$. Then, it holds that

[TABLE]

For $i=s,\ldots,2,$ suppose that we have computed $x_{i}^{\prime}\in{\cal X}_{i}$ defined by

[TABLE]

Then the optimal solution \mbox{\boldmath{x}}^{+} for (21) can be computed exactly via the following steps:

[TABLE]

Proof.

Since ${\cal D}\succ 0$, we know that ${\cal D}$, ${\cal D}+{\cal U}$ and ${\cal D}+{\cal U}^{*}$ are all nonsingular. Then, (23) can easily be obtained from the following observation

[TABLE]

Next we show the equivalence between (21) and (25). By noting that $\delta_{1}=\delta_{1}^{\prime}$ and ${\cal Q}_{11}\succ 0$, we can define $x_{1}^{\prime}$ as follows:

[TABLE]

The optimality conditions corresponding to $x_{1}^{\prime}$ and $x_{1}^{+}$ in (27) can be written as

[TABLE]

where $\gamma_{1}\in\partial p(x_{1}^{\prime})\equiv\partial p(x_{1}^{+})$. Simple calculations show that (28a) together with (24) can equivalently be rewritten as

[TABLE]

where \mbox{\boldmath{\gamma}}=(\gamma_{1};0;\ldots,0)\in{\cal X}, while (25) can equivalently be recast as

[TABLE]

By substituting \mbox{\boldmath{x}}^{\prime}=({\cal D}+{\cal U})^{-1}(\mbox{\boldmath{b}}-\mbox{\boldmath{\gamma}}+\mbox{\boldmath{\delta}}^{\prime}-{\cal U}^{*}\bar{}\mbox{\boldmath{x}}) into the above equation, we obtain that

[TABLE]

which, together with (26), (22) and the definition of ${\cal T}_{{\cal Q}}$ in (20), implies that

[TABLE]

In the above, we have used the fact that $({\cal D}+{\cal U}){\cal D}^{-1}{\cal U}({\cal D}+{\cal U})^{-1}={\cal U}{\cal D}^{-1}$. By noting that (29) is in fact the optimality condition for (21) and ${\cal Q}+{\cal T}_{{\cal Q}}\succ 0$, we have thus obtained the equivalence between (21) and (25). This completes the proof of the theorem. ∎

We shall explain here the roles of the error vectors \mbox{\boldmath{\delta}}^{\prime} and $\delta$ in the above block sGS decomposition theorem. There is no need to choose these error vectors in advance. We emphasize that $x_{i}^{\prime}$ and $x_{i}^{+}$ obtained from (24) and (25) should be viewed as approximate solutions to the minimization problems without the terms involving $\delta_{i}^{\prime}$ and $\delta_{i}$. Once these approximate solutions have been computed, they would generate $\delta_{i}^{\prime}$ and $\delta_{i}$ automatically. With these known error vectors, we know that the computed approximate solutions are the exact solutions to the minimization problems in (24) and (25).

The following proposition is useful in estimating the error term \Delta(\mbox{\boldmath{\delta}}^{\prime},\mbox{\boldmath{\delta}}) in (21).

Proposition 1**.**

Denote $\widehat{{\cal Q}}:={\cal Q}+{\cal T}_{{\cal Q}}$, which is positive definite. Let \xi=\|\widehat{\cal Q}^{-1/2}\Delta(\mbox{\boldmath{\delta}}^{\prime},\mbox{\boldmath{\delta}})\|. It holds that

[TABLE]

Proof.

Recall that $\widehat{\cal Q}=({\cal D}+{\cal U}){\cal D}^{-1}({\cal D}+{\cal U}^{*})$. Thus, we have

[TABLE]

which, together with the definition of \Delta(\mbox{\boldmath{\delta}}^{\prime},\mbox{\boldmath{\delta}}) in (22), implies that

[TABLE]

The desired result then follows. ∎

Theorem 1 shows that instead of solving the QP subproblem (21) directly with an $N$-dimensional variable $x$, where $N=\sum_{i=1}^{s}n_{i}$, the computation can be decomposed into $s$ pieces of smaller dimensional problems involving only the variable $x_{i}$ for each $i=1,\ldots,s$. Such a decomposition is obviously highly useful for dealing with a large scale CCQP of the form (18) when $N$ is very large. The benefit is especially important because the computation of $x_{i}$ for $i=2,\ldots,s$ involves only solving linear systems of equations. Of course, one would still have to solve a potentially difficult subproblem involving the variable $x_{1}$ due to the presence of the possibly nonsmooth term $p(x_{1})$, i.e.,

[TABLE]

where $c_{1}$ is a known vector depending on the previously computed $x_{s}^{\prime},\ldots,x_{2}^{\prime}$. However, in many applications, $p(x_{1})$ is usually a simple nonsmooth function such as $\|x_{1}\|_{1}$, $\|x_{1}\|_{\infty}$, or $\delta_{\mathbb{R}^{n_{1}}_{+}}(x_{1})$ for which the corresponding subproblem is not difficult to solve. As a concrete example, suppose that $Q_{11}=I_{n_{1}}$. Then $x_{1}^{+}={\rm Prox}_{p}(c_{1})$ and the Moreau-Yosida proximal map ${\rm Prox}_{p}(c_{1})$ can be computed efficiently for various nonsmooth function $p(\cdot)$ including the examples just mentioned. In fact, one can always make the subproblem easier to solve by (a) adding an additional proximal term $\frac{1}{2}\|x_{1}-\bar{x}_{1}\|^{2}_{J_{1}}$ to (21), where $J_{1}=\mu_{1}I_{n_{1}}-Q_{11}$ with $\mu_{1}=\|Q_{11}\|_{2}$; and (b) modifying the sGS operator to ${\cal U}\widehat{{\cal D}}^{-1}{\cal U}^{*}$, where $\widehat{{\cal D}}={\cal D}+\mbox{diag}(J_{1},0,\ldots,0)$. With the additional proximal term involving $J_{1}$, the subproblem corresponding to $x_{1}$ then becomes

[TABLE]

In fact, more generally, one can also modify the other diagonal blocks in ${\cal D}$ to make the linear systems involved easier to solve by adding the proximal term \frac{1}{2}\|\mbox{\boldmath{x}}-\bar{\mbox{\boldmath{x}}}\|^{2}_{\mbox{diag}(J_{1},J_{2},\ldots,J_{s})} to (21), where $J_{i}\succeq 0$, $i=1,\ldots,s$ are given symmetric matrices. Correspondingly, the sGS linear operator for the proximal term to be added to the problem (21) then becomes ${\cal T}_{{\cal Q}+\mbox{diag}(J_{1},\ldots,J_{s})}={\cal U}\widehat{{\cal D}}^{-1}{\cal U}^{*},$ where $\widehat{{\cal D}}={\cal D}+\mbox{diag}(J_{1},J_{2},\ldots,J_{s})$, and $\widehat{{\cal Q}}$ in (23) becomes $\widehat{{\cal Q}}={\cal Q}+\mbox{diag}(J_{1},\ldots,J_{s})+{\cal U}\widehat{{\cal D}}^{-1}{\cal U}^{*}$. There are many suitable choices for $J_{i}$, $i=2,\ldots,s.$ A conservative choice would be $J_{i}=\|Q_{ii}\|_{2}I_{n_{i}}-Q_{ii}$, in which case the linear system to be solved has its coefficient matrix given by $Q_{ii}+J_{i}=\|Q_{ii}\|_{2}I_{n_{i}}$. Another possible choice of $J_{i}$ is the sGS linear operator associated with the matrix $Q_{ii}$, in which case the linear system involved has its coefficient matrix given by $Q_{ii}+{\cal T}_{Q_{ii}}$ and its solution can be computed by using one cycle of the sGS method. The latter choice has been considered in [3] for its variant of the classical block sGS method. Despite the advantage of simplifying the linear systems to be solved, one should note that the price to pay for adding the extra proximal term \frac{1}{2}\|\mbox{\boldmath{x}}-\bar{\mbox{\boldmath{x}}}\|^{2}_{\mbox{diag}(J_{1},\ldots,J_{s})} is worsening the convergence rate of the overall block sGS method.

3 A factorization view of the block sGS decomposition

theorem and its equivalence to the SCB reduction procedure

In this section, we present a factorization view of the block sGS decomposition theorem and show its equivalence to the Schur complement based (SCB) reduction procedure developed in [17, 18].

Let ${\Theta}_{1}$ be the zero matrix in $\mathbb{R}^{n_{1}\times n_{1}}$ and $N_{1}:=n_{1}$. For $j=2,\ldots,s$, let $N_{j}:=\sum_{i=1}^{j}n_{i}$ and define $\widehat{\Theta}_{j}\in\mathbb{R}^{N_{j-1}\times N_{j-1}}$ and ${\Theta}_{j}\in\mathbb{R}^{N_{j}\times N_{j}}$ as follows:

[TABLE]

and

[TABLE]

Then, the above definitions indicate that, for $2\leq j\leq s$,

[TABLE]

In [17, 18], the SCB reduction procedure corresponding to problem (18) is derived through the construction of the above self-adjoint linear operator ${\Theta}_{s}$ on ${\cal X}$. Now we recall the key steps in the SCB reduction procedure derived in the previous work. For $j=1,\ldots,s$, define

[TABLE]

Note that $\widehat{\Theta}_{j}=R_{j}Q_{j,j}^{-1}R_{j}^{*}$. It is easy to show that

[TABLE]

By first solving the inner minimization problem with respect to $x_{s}$, we get the solution as a function of $x_{1},\ldots,x_{s-1}$ as follows:

[TABLE]

And the minimum value is given by

[TABLE]

Thus (32) reduces to a problem involving only the variables $x_{1},\ldots,x_{s-1}$, which, up to a constant, is given by

[TABLE]

where $x_{s}^{\prime}=Q_{s,s}^{-1}(b_{s}-R_{s}^{*}\bar{x}_{\leq s-1}).$ Observe that (41) has exactly the same form as (32). By repeating the above procedure to sequentially eliminate the variables $x_{s-1},\ldots,x_{2}$, we will finally arrive at a minimization problem involving only the variable $x_{1}$. Once that minimization problem is solved, we can recover the solutions for $x_{2},\ldots,x_{s}$ in a sequential manner.

Now we will prove the equivalence between the block sGS decomposition theorem and the SCB reduction procedure in the subsequent analysis by proving that ${\Theta}_{s}={\cal T}_{\cal Q}$, where ${\cal T}_{\cal Q}$ is given in (20). For $j=2,\dots,s$, define the block matrices $\widehat{\cal V}_{j}\in\mathbb{R}^{N_{j}\times N_{j}}$ and ${\cal V}_{j}\in\mathbb{R}^{N\times N}$ by

[TABLE]

where $I_{n_{j}}$ is the $n_{j}\times n_{j}$ identity matrix. Note that $\widehat{\cal V}_{s}={\cal V}_{s}$. Given $j\geq 2$, we have, by simple calculations, that for any $k<j$,

[TABLE]

From (30) and (43), we have that

[TABLE]

Lemma 1**.**

Let ${\cal U}$ and ${\cal D}$ be given in (16). It holds that

[TABLE]

Proof.

It can be verified directly that

[TABLE]

The second equality follows readily from the first. ∎

In the proof of the next lemma, we will make use of the well known fact that for given symmetric matrices $A,C$ such that $C\succ 0$ and $M:=A-BC^{-1}B^{*}\succ 0$, we have that

[TABLE]

Theorem 2**.**

It holds that

[TABLE]

Proof.

By using (54), for $j=2,\ldots,s$, we have that

[TABLE]

where

[TABLE]

Thus, from (44), we know that for $2\leq j\leq s$,

[TABLE]

For $2\leq j\leq s$, by (31), we have that

[TABLE]

and consequently,

[TABLE]

Thus, by recalling the definitions of $\widehat{\cal V}_{j}$ and ${\cal V}_{j}$ in (42) and using (56), we obtain through simple calculations that

[TABLE]

Thus, by using the fact that ${\cal Q}_{1}+{\Theta}_{1}=Q_{1,1}$, we get

[TABLE]

By Lemma 1, it follows that

[TABLE]

where the last equation follows from (23) in Theorem 1. Since ${\cal Q}_{s}={\cal Q}$, we know that

[TABLE]

This completes the proof of the theorem. ∎

4 An extended block sGS method for solving the CCQP (18)

With the block sGS decomposition theorem (Theorem 1) and Proposition 1, we can now extend the classical block sGS method to solve the CCQP (18). The detail steps of the algorithm for solving (18) are given as follows.

Algorithm 1: An sGS based inexact proximal gradient method for (18).

Input \widetilde{}\mbox{\boldmath{x}}^{1}=\mbox{\boldmath{x}}^{0}\in{\rm dom}(p)\times\mathbb{R}^{n_{2}}\times\ldots\times\mathbb{R}^{n_{s}}, $t_{1}=1$ and a summable sequence of nonnegative numbers $\{\epsilon_{k}\}$. For $k=1,2,\ldots$, perform the following steps in each iteration.

Step 1.

Compute

\mbox{\boldmath{x}}^{k}=\mbox{argmin}_{\mbox{\boldmath{x}}\in{\cal X}}\;\Big{\{}p(x_{1})+q(\mbox{\boldmath{x}})+\frac{1}{2}\|\mbox{\boldmath{x}}-\widetilde{}\mbox{\boldmath{x}}^{k}\|_{{\cal T}_{\cal Q}}^{2}-\langle\mbox{\boldmath{x}},\,\Delta(\,\widetilde{}\mbox{\boldmath{\delta}}^{k},\mbox{\boldmath{\delta}}^{k})\rangle\Big{\}},

(57)

via the sGS decomposition procedure described in Theorem 1, where \widetilde{}\mbox{\boldmath{\delta}}^{k},\,\mbox{\boldmath{\delta}}^{k}\in{\cal X} are error vectors such that

\max\{\|\widetilde{}\mbox{\boldmath{\delta}}^{k}\|,\|\mbox{\boldmath{\delta}}^{k}\|\}\leq\frac{\epsilon_{k}}{t_{k}}.

(58)

Step 2.

Choose $t_{k+1}$ such that $t_{k+1}^{2}-t_{k+1}\leq t_{k}^{2}$ and set $\beta_{k}=\frac{t_{k}-1}{t_{k+1}}$. Compute

\widetilde{}\mbox{\boldmath{x}}^{k+1}=\mbox{\boldmath{x}}^{k}+\beta_{k}(\mbox{\boldmath{x}}^{k}-\mbox{\boldmath{x}}^{k-1}).

We have the following iteration complexity convergence results for Algorithm 1.

Proposition 2**.**

*Suppose \mbox{\boldmath{x}}^{*} is an optimal solution of problem (18). Let \{\mbox{\boldmath{x}}^{k}\} be the sequence generated by Algorithm 1. Define $M=2\|{\cal D}^{-1/2}\|_{2}+\|\widehat{{\cal Q}}^{-1/2}\|_{2}$.

(a) If $t_{k+1}=\frac{1+\sqrt{1+4t_{k}^{2}}}{2}$ for all $k\geq 1$, it holds that*

[TABLE]

*where $\widehat{{\cal Q}}={\cal Q}+{\cal T}_{\cal Q}$ and $\bar{\epsilon}_{k}=2M\sum_{i=1}^{k}\epsilon_{i}$.

(b) If $t_{k}=1$ for all $k\geq 1$, it holds that*

[TABLE]

where $\tilde{\epsilon}_{k}=4M\sum_{i=1}^{k}i\epsilon_{i}$.

Proof.

(a) The result can be proved by applying Theorem 2.1 in [14]. In order to apply the theorem, we need to verify that the error \mbox{\boldmath{e}}:=\mbox{\boldmath{\gamma}}+{\cal Q}\mbox{\boldmath{x}}^{k}-\mbox{\boldmath{b}}+{\cal T}_{\cal Q}(\mbox{\boldmath{x}}^{k}-\widetilde{\mbox{\boldmath{x}}}^{k}), where \mbox{\boldmath{\gamma}}=(\gamma_{1};0;\ldots;0) and $\gamma_{1}\in\partial p(x_{1}^{k})$, incurred for solving the subproblem (without the perturbation term \Delta(\,\widetilde{}\mbox{\boldmath{\delta}}^{k},\mbox{\boldmath{\delta}}^{k})) in Step 1 inexactly is sufficiently small. From Theorem 1, we know that

[TABLE]

The theorem is proved via Theorem 2.1 in [14] if we can show that \|\widehat{{\cal Q}}^{-1/2}\Delta(\,\widetilde{}\mbox{\boldmath{\delta}}^{k},\mbox{\boldmath{\delta}}^{k})\|\leq M\frac{\epsilon_{k}}{t_{k}}. But from (58) and Proposition 1, we have that

[TABLE]

thus the required inequality indeed holds true, and the proof is completed.

(b) There is no straightforward theorem for which we can apply to prove the result, we will provide the proof in the Appendix. ∎

Remark 1**.**

It is not difficult to show that if $p(\cdot)\equiv 0$, $t_{k}=1$, and \mbox{\boldmath{\delta}}^{k}=\widetilde{}\mbox{\boldmath{\delta}}^{k}=0 for all $k\geq 1$, then Algorithm 1 exactly coincides with the classical block sGS method (17); and if \mbox{\boldmath{\delta}}^{k}, \widetilde{}\mbox{\boldmath{\delta}}^{k} are allowed to be non-zero but satisfy the condition (58) for all $k\geq 1$, then we obtain the inexact extension of the classical block sGS method.

Remark 2**.**

Proposition 2 shows that the classical block sGS method for solving (1) can be extended to solve the convex composite QP (18). It also demonstrates the advantage of interpreting the block sGS method from the optimization perspective. For example, one can obtain the $O(1/k)$ iteration complexity result for the classical block sGS method without assuming that ${\cal Q}$ is positive definite. To the best of our knowledge, such a complexity result for the classical block sGS is new. More importantly, inexact and accelerated versions of the block sGS method can also be derived for (1).

Remark 3**.**

In solving (57) via the sGS decomposition procedure to satisfy the error condition (58), let \mbox{\boldmath{x}}^{\prime}=[x^{\prime}_{1};\ldots;x^{\prime}_{s}] be the intermediate solution computed during the backward GS sweep (in Theorem 1) and the associated error vector be \widetilde{\mbox{\boldmath{\delta}}}^{k}=[\widetilde{\delta}^{k}_{1};\ldots;\widetilde{\delta}^{k}_{s}]. In the forward GS sweep, one can often save computations by using the computed $x^{\prime}_{i}$ to estimate $x^{k+1}_{i}$ for $i\geq 2$, and the resulting error vector will be given by \delta^{k}_{i}=\widetilde{\delta}^{k}_{i}+\mbox{\sum_{j=1}^{i-1}}Q_{ji}^{*}(x_{j}^{k+1}-\widetilde{x}^{k}_{j}). If we have that

[TABLE]

where \rho=\frac{c}{\sqrt{s}}\|\widetilde{\mbox{\boldmath{\delta}}}^{k}\| and $c>0$ is some given constant, then clearly $\|\delta^{k}_{i}\|^{2}\leq 2\|\widetilde{\delta}^{k}_{i}\|^{2}+2\rho^{2}$. When all the error components $\|\delta^{k}_{i}\|^{2}$ satisfy the previous bound for $i=1,\ldots,s$, regardless of whether $x^{k+1}_{i}$ is estimated from $x^{\prime}_{i}$ or computed afresh, we get \|\mbox{\boldmath{\delta}}^{k}\|\leq\sqrt{2(1+c^{2})}\|\widetilde{\mbox{\boldmath{\delta}}}^{k}\|. Consequently the error condition (58) can be satisfied with a slightly larger error tolerance $\sqrt{2(1+c^{2})}\,\epsilon_{k}/t_{k}$. It is easy to see that one can use the condition in (59) to decide whether $x^{k+1}_{i}$ can be estimated from $x^{\prime}_{i}$ without contributing a large error to \|\mbox{\boldmath{\delta}}^{k}\| for each $i=2,\ldots,s$.

Besides the above iteration complexity results, one can also study the linear convergence rate of Algorithm 1. Indeed, just as in the case of the classical block sGS method, the convergence rate of our extended inexact block sGS method for solving (18) can also be established when ${\cal Q}\succ 0$. The precise result is given in the next theorem.

Theorem 3**.**

Suppose that the relative interior of the domain of $p$, ${\rm ri}({\rm dom}(p))$, is non-empty, ${\cal Q}\succ 0$ and $t_{k}=1$ for all $k\geq 1$. Then

[TABLE]

where ${\cal B}=I-{\widehat{{\cal Q}}}^{-1/2}{\cal Q}{\widehat{{\cal Q}}}^{-1/2}$, and $M$ is defined as in Proposition 2. Note that $0\preceq{\cal B}\prec I$.

Proof.

For notational convenience, we let \Delta^{j}=\Delta(\widetilde{\mbox{\boldmath{\delta}}}^{j},\mbox{\boldmath{\delta}}^{j}) in this proof.

Define $E_{1}:{\cal X}\to{\cal X}_{1}$ by E_{1}(\mbox{\boldmath{x}})=x_{1} and $\widehat{p}:{\cal X}\to(-\infty,\infty]$ by \widehat{p}(\mbox{\boldmath{x}})=p(E_{1}{\widehat{{\cal Q}}}^{-1/2}\mbox{\boldmath{x}}). Since ${\widehat{{\cal Q}}}\succ 0$, it is clear that ${\rm range}(E_{1}{\widehat{{\cal Q}}}^{-1/2})={\cal X}_{1}$ and hence ${\rm ri}({\rm dom}(p))\cap{\rm range}(E_{1}{\widehat{{\cal Q}}}^{-1/2})\not=\emptyset$. By [21, Theorem 23.9], we have that

[TABLE]

From the optimality condition of \mbox{\boldmath{x}}^{j}, we have that

[TABLE]

where \mbox{\boldmath{\gamma}}^{j}=(\gamma_{1}^{j};0;\ldots;0) with $\gamma_{1}^{j}\in\partial p(x_{1}^{j})$. Let \hat{\mbox{\boldmath{x}}}^{j}={\widehat{{\cal Q}}}^{1/2}\mbox{\boldmath{x}}^{j} and \hat{\mbox{\boldmath{x}}}^{j-1}={\widehat{{\cal Q}}}^{1/2}\mbox{\boldmath{x}}^{j-1}. Then we have that

[TABLE]

Similarly if \mbox{\boldmath{x}}^{*} is an optimal solution of (18), then we have that

[TABLE]

By using the nonexpansive property of ${\rm Prox}_{\widehat{p}}$, we have that

[TABLE]

By applying the above inequality sequentially for $j=k,k-1,\ldots,1$, we get the required result in (60). ∎

Remark 4**.**

In fact, one can weaken the positive definiteness assumption of ${\cal Q}$ in the above theorem and still expect a linear rate of convergence. As a simple illustration, we only discuss here the exact version of Algorithm 1, i.e., \widetilde{}\mbox{\boldmath{\delta}}^{k}=\mbox{\boldmath{\delta}}^{k}=0, under the error bound condition [19, 20] on $F$ which holds automatically if $p$ is a convex piecewise quadratic/linear function such as $p(x_{1})=\|x_{1}\|_{1}$, $p(x_{1})=\delta_{\mathbb{R}^{n_{1}}_{+}}$ or if ${\cal Q}\succ 0$. When $t_{k}=1$ for all $k\geq 1$, one can prove that \{F(\mbox{\boldmath{x}}^{k})\} converges at least Q-linearly and \{\mbox{\boldmath{x}}^{k}\} converges at least R-linearly to an optimal solution of problem (18) by using the techniques developed in [19, 20]. Interested readers may refer to [28, 30] for more details. For the accelerated case, with the additional fixed restarting scheme incorporated in Algorithm 1, both the R-linear convergences of \{F(\mbox{\boldmath{x}}^{k})\} and \{\mbox{\boldmath{x}}^{k}\} can be obtained from [29, Corollary 3.8].

5 An illustration on the application of the block sGS decomposition theorem in designing

an efficient proximal ALM

In this section, we demonstrate the usefulness of our block sGS decomposition theorem as a building block for designing an efficient proximal ALM for solving a linearly constrained convex composite QP problem given by

[TABLE]

where ${\cal P}$ is a positive semidefinite linear operator on ${\cal X}$, ${\cal A}:{\cal X}\to{\cal Y}$ is a given linear map, and \mbox{\boldmath{g}}\in{\cal X}, $d\in{\cal Y}$ are given data. Here ${\cal X}$ and ${\cal Y}$ are two finite dimensional inner product spaces. Specifically, we show how the block sGS decomposition theorem given in Theorem 1 can be applied within the proximal ALM. We must emphasize that our main purpose here is to briefly illustrate the usefulness of the block sGS decomposition theorem but not to focus on the proximal ALM itself. Indeed, simply being capable of handling the nonsmooth function $p(\cdot)$ has already distinguished our approach from other approaches of using the sGS technique in optimization algorithms, e.g, [8, 9], where the authors incorporated the pointwise sGS splitting as a preconditioner within the Douglas–Rachford splitting method for a convex-concave saddle point problem.

In depth analysis of various recently developed ADMM-type algorithms and accelerated block coordinate descent algorithms employing the block sGS decomposition theorem as a building block can be found in [6, 15, 16, 17, 25]. Thus we shall not elaborate here again on the essential role played by the block sGS decomposition theorem in the design of those algorithms.

Although the problem (62) looks deceivingly simple, in fact it is a powerful model which includes the important class of standard convex quadratic semidefinite programming (QSDP) in the dual form given by

[TABLE]

where $h\in\mathbb{R}^{p}$, $C\in\mathbb{S}^{n}$ are given data, ${\cal B}:\mathbb{S}^{n}\to\mathbb{R}^{p}$ is a given linear map that is assumed to be surjective, ${\cal H}:\mathbb{S}^{n}\to\mathbb{S}^{n}$ is a self-adjoint positive semidefinite linear operator, and ${\cal W}\subseteq\mathbb{S}^{n}$ is any subspace containing ${\rm Range}({\cal H})$, the range space of ${\cal H}$. Here $\mathbb{S}^{n}$ denotes the space of $n\times n$ symmetric matrices and $\mathbb{S}^{n}_{+}$ denotes the cone of symmetric positive semidefinite matrices in $\mathbb{S}^{n}$. One can obviously express the QSDP problem (63) in the form of (62) by defining \mbox{\boldmath{x}}=(Z;\xi;W), $p(Z)=\delta_{\mathbb{S}^{n}_{+}}(Z)$, ${\cal P}=\mbox{diag}(0,0,{\cal H})$, and ${\cal A}=({\cal I},{\cal B}^{*},{\cal H})$.

We begin with the augmented Lagrangian function associated with (62):

[TABLE]

where $\sigma>0$ is a given penalty parameter and $y\in{\cal Y}$ is the multiplier associated with the equality constraint. The template for a proximal ALM is given as follows. Given ${\cal T}\succeq 0$, \mbox{\boldmath{x}}^{0}\in{\cal X} and $y^{0}\in{\cal Y}$. Perform the following steps in each iteration.

Step 1.

Compute

[TABLE]

where \mbox{\boldmath{b}}=\mbox{\boldmath{g}}+{\cal A}^{*}(\sigma d-y^{k}).

Step 2.

Compute y^{k+1}=y^{k}+\tau\sigma({\cal A}\mbox{\boldmath{x}}^{k}-d), where $\tau\in(0,2)$ is the step-length.

It is clear that the subproblem (65) has the form given in (18). Thus, one can apply the block sGS decomposition theorem to efficiently solve the subproblem if we choose ${\cal T}={\cal T}_{{\cal P}+\sigma{\cal A}^{*}{\cal A}}$, i.e., the sGS operator associated with ${\cal Q}:={\cal P}+\sigma{\cal A}^{*}{\cal A}$. For the QSDP problem (63) with ${\cal W}:={\rm Range}({\cal H})$, we have that

[TABLE]

and that the subproblem (65) can be efficiently solved by one cycle of the extended block sGS method explicitly as follows, given the iterate $(Z^{k},\xi^{k},{\cal H}W^{k})$ and multiplier $y^{k}$.

Step 1a.

Compute ${\cal H}W^{\prime}$ as the solution of $(\sigma^{-1}{\cal I}+{\cal H}){\cal H}W^{\prime}=\sigma^{-1}b_{W}-{\cal H}Z^{k}-{\cal H}{\cal B}^{*}\xi^{k},$ where $b_{W}={\cal H}(\sigma C-Y^{k})$.

Step 1b.

Compute $\xi^{\prime}$ from ${\cal B}{\cal B}^{*}\xi^{\prime}=\sigma^{-1}b_{\xi}-{\cal B}Z^{k}-{\cal B}{\cal H}W^{\prime},$ where $b_{\xi}=h+{\cal B}(\sigma C-Y^{k})$.

Step 1c.

Compute Z^{k+1}=\mbox{argmin}\Big{\{}\delta_{\mathbb{S}^{n}_{+}}(Z)+\frac{\sigma}{2}\|Z+{\cal B}^{*}\xi^{\prime}+{\cal H}W^{\prime}-\sigma^{-1}b_{Z}\|^{2}\Big{\}}, where $b_{Z}=\sigma C-Y^{k}$.

Step 1d.

Compute $\xi^{k+1}$ from ${\cal B}{\cal B}^{*}\xi^{k+1}=\sigma^{-1}b_{\xi}-{\cal B}Z^{k+1}-{\cal B}{\cal H}W^{\prime}.$

Step 1e.

Compute ${\cal H}W^{k+1}$ from $(\sigma^{-1}{\cal I}+{\cal H}){\cal H}W^{k+1}=\sigma^{-1}b_{W}-{\cal H}Z^{k+1}-{\cal H}{\cal B}^{*}\xi^{k+1}.$

From the above implementation, one can see how simple it is for one to apply the block sGS decomposition theorem to solve the complicated subproblem (65) arising from QSDP. Note that in Step 1a and Step 1e, we only need to compute ${\cal H}W^{\prime}$ and ${\cal H}W^{k+1}$, respectively, and we do not need the values of $W^{\prime}$ and $W^{k+1}$ explicitly. Here, for simplicity, we only write down the exact version of a proximal ALM by using our exact block sGS decomposition theorem. Without any difficulty, one can also apply the inexact version of the block sGS decomposition theorem to derive a more practical inexact proximal ALM for solving (62), say when the linear systems involved are large scale and have to be solved by a Krylov subspace iterative method.

6 Extension of the classical block symmetric SOR method for solving (18)

In a way similar to what we have done in section 2, we show in this section that the classical block symmetric SOR (block sSOR) method can also be interpreted from an optimization perspective.

Given a parameter $\omega\in[1,2)$, the $k$th iteration of the classical block sSOR method in the third normal form is defined by

[TABLE]

where

[TABLE]

$\tau=1/\omega$, and $\rho=2\tau-1$. Note that for $\omega\in[1,2)$, we have that $\tau\in(1/2,1]$ and $\rho\in(0,1]$. We should mention that the classical block sSOR method is typically not derived in the form given in (67), see for example [13, p.117], but one can show with some algebraic manipulations that (67) is an equivalent reformulation.

Denote

[TABLE]

In the next proposition, we show that ${\cal W}$ can be decomposed as the sum of ${\cal Q}$ and ${\cal T}_{\rm sSOR}$. Similar to the linear operator ${\cal T}_{\cal Q}$ in section 2, ${\cal T}_{\rm sSOR}$ is the key ingredient which enables us to derive the block sSOR method from the optimization perspective, and to extend it to solve the CCQP (18).

Proposition 3**.**

Let $\omega\in[1,2)$, and denote $\tau=1/\omega\in(1/2,1]$, $\rho=2\tau-1$. It holds that

[TABLE]

Proof.

Let $\bar{\tau}:=\tau-\frac{1}{2}>0$ and $\overline{{\cal U}}={\cal U}+\frac{1}{2}{\cal D}$. Note that $\rho=2\bar{\tau}$ and

[TABLE]

Now

[TABLE]

From here, we get the required expression for $\widehat{{\cal Q}}_{\omega}$ in (68). ∎

Given two error tolerance vectors $\delta$ and \mbox{\boldmath{\delta}}^{\prime} with $\delta_{1}=\delta^{\prime}_{1}$, let

[TABLE]

Given \bar{}\mbox{\boldmath{x}}\in{\cal X}, similar to Theorem 1, one can prove without much difficulty that the optimal solution of the following minimization subproblem

[TABLE]

can be computed by performing exactly one cycle of the block sSOR method. In particular, when $p(\cdot)\equiv 0$ and \mbox{\boldmath{\delta}}=\mbox{\boldmath{\delta}}^{\prime}=\mbox{\bf 0}, the optimal solution to (69) can be computed by (67), i.e., set \bar{}\mbox{\boldmath{x}}=\mbox{\boldmath{x}}^{k}, then \mbox{\boldmath{x}}^{k+1} obtained from (67) is the optimal solution to (69). By replacing ${\cal T}_{\cal Q}$ and $\Delta(\cdot,\cdot)$ in Algorithm 1 with ${\cal T}_{\rm sSOR}$ and $\Delta_{\rm sSOR}(\cdot,\cdot)$, respectively, one can obtain a block sSOR based inexact proximal gradient method for solving (18) and the convergence results presented in Proposition 2 and Theorem 3 still remain valid with $\widehat{{\cal Q}}$ replaced by $\widehat{{\cal Q}}_{\omega}$.

Remark 5**.**

For the classical pointwise sSOR method, it was shown in [13, Theorem 4.8.14] that if there exist positive constants $\gamma$ and $\Gamma$ such that

[TABLE]

then its convergence rate is $\|I-{\cal Q}^{1/2}\widehat{{\cal Q}}_{\omega}^{-1}{\cal Q}^{1/2}\|_{2}\leq 1-\frac{2\bar{\tau}}{\bar{\tau}^{2}/\gamma+\bar{\tau}+\Gamma/4},$ where $\bar{\tau}=1/\omega-1/2.$ Interestingly, for the convergence rate of our block sSOR method in Theorem 3, we also have a similar estimate given by

[TABLE]

In order to minimize the upper bound, we can choose $\omega_{*}=2/(1+\sqrt{\gamma\Gamma})$ and then we get

[TABLE]

7 Conclusion

In this paper, we give an optimization interpretation that each cycle of the classical block sGS method is equivalent to solving the associated multi-block convex QP problem with an additional proximal term. This equivalence is fully characterized via our block sGS decomposition theorem. A factorization view of this theorem and its equivalence to the SCB reduction procedure are also established. The classical block sGS method, viewed from the optimization perspective via the block sGS decomposition theorem, is then extended to the inexact setting for solving a class of multi-block convex composite QP problems involving nonsmooth functions. Moreover, we are able to derive $O(1/k)$ and $O(1/k^{2})$ iteration complexities for our inexact block sGS method and its accelerated version, respectively. These new interpretations and convergence results, together with the incorporation of the (inexact) sGS decomposition techniques in the design of efficient algorithms for core optimization problems in [6, 15, 16, 17, 25], demonstrate the power and usefulness of our simple yet elegant block sGS decomposition theorem. We believe this decomposition theorem will be proven to be even more useful in solving other optimization problems and beyond.

Appendix: Proof of part (b) of Proposition 2

To begin the proof, we state the following lemma from [24].

Lemma 2**.**

Suppose that $\{u_{k}\}$ and $\{\lambda_{k}\}$ are two sequences of nonnegative scalars, and $\{s_{k}\}$ is a nondecreasing sequence of scalars such that $s_{0}\geq u_{0}^{2}$. Suppose that for all $k\geq 1$, the inequality $u_{k}^{2}\leq s_{k}+2\sum_{i=1}^{k}\lambda_{i}u_{i}$ holds. Then for all $k\geq 1$, $u_{k}\leq\bar{\lambda}_{k}+\sqrt{s_{k}+\bar{\lambda}_{k}^{2}},$ where $\bar{\lambda}_{k}=\sum_{i=1}^{k}\lambda_{i}$.

Proof.

In this proof, we let \Delta^{j}=\Delta(\tilde{\mbox{\boldmath{\delta}}}^{j},\mbox{\boldmath{\delta}}^{j}). Note that under the assumption that $t_{j}=1$ for all $j\geq 1$, \widetilde{\mbox{\boldmath{x}}}^{j}=\mbox{\boldmath{x}}^{j-1}. Note also that from (58), we have that $\|{\widehat{{\cal Q}}}^{-1/2}\Delta^{j}\|\leq M\epsilon_{j}$, where $M$ is given as in Proposition 2.

From the optimality of \mbox{\boldmath{x}}^{j} in (57), one can show that

[TABLE]

Let \mbox{\boldmath{e}}^{j}=\mbox{\boldmath{x}}^{j}-\mbox{\boldmath{x}}^{*}. By setting \mbox{\boldmath{x}}=\mbox{\boldmath{x}}^{j-1} and \mbox{\boldmath{x}}=\mbox{\boldmath{x}}^{*} in (70), we get

[TABLE]

By multiplying $j-1$ to (71) and combining with (72), we get

[TABLE]

where a_{j}=2j[F(\mbox{\boldmath{x}}^{j})-F(\mbox{\boldmath{x}}^{*})] and b_{j}=\|\mbox{\boldmath{e}}^{j}\|_{\widehat{{\cal Q}}}. Note that the last inequality follows from (72) with $j=1$ and some simple manipulations. To summarize, we have $b_{j}^{2}\leq b_{0}^{2}+2\sum_{i=1}^{j}2Mi\epsilon_{i}b_{i}$. By applying Lemma 2, we get

[TABLE]

where $\bar{\lambda}_{j}=\sum_{i=1}^{j}\lambda_{i}$ with $\lambda_{i}=2Mi\epsilon_{i}$. Applying the above result to (73), we get

[TABLE]

From here, the required result in Part (b) of Proposition 2 follows. ∎