Overlapping Schwarz Decomposition for Constrained Quadratic Programs

Structured optimization problems arise in many applications such as trajectory planning and model predictive control [1], multistage stochastic programming [2, 3], optimization with embedded partial differential equations (PDEs) [4], and in optimization of energy infrastructures [5, 6, 7]. Decomposition schemes that exploit the problem structure have been used to overcome scalability limits of centralized schemes [8, 9, 10, 11, 12]. A wide range of decomposition schemes have been proposed in the literature such as Lagrangian decomposition [13], the alternating direction method of multipliers (ADMM) [14], and Jacobi/Gauss-Seidel methods [15]. The basic tenet behind such algorithms is to decompose the original problem into subproblems and to coordinate subproblem solutions by using primal-dual information. These decomposition schemes are different from waterfilling methods [16] in that the constraints are not necessarily satisfied at each iteration within the algorithm. Moreover, these schemes are advantageous in that they can handle a wide variety of objective/constraint formulations and in that they have well-understood convergence properties. However, a disadvantage of these schemes is that convergence can be rather slow [17].

Overlapping Schwarz decomposition schemes have been recently used to solve structured optimization problems, and they have demonstrated to outperform popular schemes such as ADMM and Jacobi/Gauss-Seidel in certain problem classes [18]. Schwarz algorithms were originally developed for the parallel solution of linear algebra systems arising in PDEs, but such schemes can also be used to handle general linear systems and optimization problems by exploiting their underlying algebraic topology [19]. As the name suggests, overlapping Schwarz algorithms decompose the full problem (the underlying graph) into subproblems that are defined over overlapping subdomains. Solutions for the subproblems are computed in parallel and convergence is enforced by updating primal-dual information in the overlapping regions. In the context of quadratic programs (QP) that are unconstrained and convex, we have shown that the convergence rate of Schwarz algorithms improves exponentially with the size of the overlapping region [18]. Overlapping Schwarz schemes provide a bridge between fully decentralized Jacobi/Gauss-Seidel algorithms (no overlap) and centralized algorithms (the overlap is the entire domain).

This paper presents an overlapping Schwarz algorithm for constrained QPs. We analyze the convergence properties of the algorithm and derive an explicit relationship between its convergence rate and the size of the overlap. This result extends existing convergence results for unconstrained QPs [18] to equality- and inequality-constrained QPs. In particular, we show that the algorithm converges with sufficiently large overlap and that the convergence rate improves exponentially with the size of overlap. This convergence result relies on a property called exponential decay of sensitivity [20, 21, 22, 23, 18]. The property states that the sensitivity of the primal-dual solution at a given node decays exponentially with respect to the distance from the perturbation. Such a property has been established for discrete-time [20, 21] and continuous-time [22, 23] optimal control problems (the graph is a line) and for unconstrained QPs (general graph) [18]. This paper establishes the exponential sensitivity property for constrained QPs over general graphs, thus making the results broadly applicable.

The paper is organized as follows. In the remainder of this section we introduce basic notation and the problem under study. In Section II we introduce the overlapping Schwarz algorithm. In Section III we present the main theoretical results. We first analyze the sensitivity of the solution of structured QPs against parametric perturbations and then use the results to establish convergence conditions for the algorithm. Numerical results are given in Section IV.

Notation. The set of real numbers and the set of integers are denoted by $\mathbb{R}$ and $\mathbb{I}$, respectively, and we define $\mathbb{I}_{a:b}:=\mathbb{I}\cap[a,b]$, $\mathbb{I}_{>0}:=\mathbb{I}\cap(0,\infty)$, $\mathbb{R}_{>0}:=(0,\infty)$, and $\overline{\mathbb{R}}:=\mathbb{R}\cup\{\infty\}$. By default, vectors are assumed to be column vectors, and we use the syntax $(M_{1},\cdots,M_{n}):=[M_{1}^{\top}\,\cdots\,M_{n}^{\top}]^{\top}$, $\{M_{i}\}_{i\in\mathcal{I}}:=(M_{i_{1}},\cdots,M_{i_{m}})$, and $\{M_{i,j}\}_{i\in\mathcal{I},j\in\mathcal{J}}:=\{\{M_{i,j}^{\top}\}_{j\in\mathcal{J}}^{\top}\}_{i\in\mathcal{I}}$, where $\mathcal{I}=\{i_{1}<\cdots<i_{m}\}$ and $\mathcal{J}=\{j_{1}<\cdots<j_{n}\}$. Furthermore, $v[i]$ is the $i$th component of $v$, $M[i,j]$ is the $(i,j)$th component of $M$, $v[\mathcal{I}]:=\{v[i]\}_{i\in\mathcal{I}}$, and $M[\mathcal{I},\mathcal{J}]:=\{M[i,j]\}_{i\in\mathcal{I},j\in\mathcal{J}}$. Vector 2-norms and induced 2-norms are denoted by $\|\cdot\|$. For matrices $A$ and $B$, $A\succeq B$ indicates that $A-B$ is positive semi-definite while $A\geq B$ represents a componentwise inequality.

Setting. We consider a (potentially infinite) parent graph $\mathcal{G}(\mathcal{V},\mathcal{E})$, where $\mathcal{V}$ is the set of nodes and $\mathcal{E}$ is the set of edges. We also consider the finite node subset $U\subseteq\mathcal{V}$ and the following QP:

Here $x_{i}\in\mathbb{R}^{r_{i}}$ are the decision variables for node $i$; $\lambda^{\mathcal{E}}_{i}\in\mathbb{R}^{m^{\mathcal{E}}_{i}}$ and $\lambda^{\mathcal{I}}_{i}\in\mathbb{R}^{m^{\mathcal{I}}_{i}}$ are the dual variables; $Q_{i,j}\in\mathbb{R}^{r_{i}\times r_{j}}$, $A^{\mathcal{E}}_{i,j}\in\mathbb{R}^{m^{\mathcal{E}}_{i}\times r_{j}}$, $A^{\mathcal{I}}_{i,j}\in\mathbb{R}^{m^{\mathcal{I}}_{i}\times r_{j}}$, $f_{i}\in\mathbb{R}^{r_{i}}$, $g^{\mathcal{E}}_{i}\in\mathbb{R}^{m^{\mathcal{E}}_{i}}$, and $g^{\mathcal{I}}_{i}\in\mathbb{R}^{m^{\mathcal{I}}_{i}}$ are the data; and $N_{U}[X]:=N_{U}(X)\cup X$, where $N_{U}(X):=\{j\in U\setminus X:\exists i\in X\text{ such that }\{i,j\}\in\mathcal{E}\}$ and the argument is considered as a singleton if $X$ is a single node. We define $A_{i,j}:=(A^{\mathcal{E}}_{i,j},A^{\mathcal{I}}_{i,j})$, $g_{i}:=(g^{\mathcal{E}}_{i},g^{\mathcal{I}}_{i})$, $\lambda_{i}:=(\lambda^{\mathcal{E}}_{i},\lambda^{\mathcal{I}}_{i})$, $z_{i}:=(x_{i},\lambda_{i})$, $d_{i}:=(f_{i},g_{i})$, $m_{i}:=m^{\mathcal{E}}_{i}+m^{\mathcal{I}}_{i}$, and $n_{i}:=r_{i}+m_{i}$. We assume that $Q_{i,j}=Q_{j,i}^{\top}$.

An equivalent problem can be written in a compact form:

Here, $\boldsymbol{x}_{U}:=\{x_{i}\}_{i\in U}$; $\boldsymbol{\lambda}^{\mathcal{E}}_{U}:=\{\lambda^{\mathcal{E}}_{i}\}_{i\in U}$; $\boldsymbol{\lambda}^{\mathcal{I}}_{U}:=\{\lambda^{\mathcal{I}}_{i}\}_{i\in U}$; $\boldsymbol{\lambda}_{U}:=\{\lambda_{i}\}_{i\in U}$; $\boldsymbol{z}_{U}:=\{z_{i}\}_{i\in U}$; $\boldsymbol{f}_{U}:=\{f_{i}\}_{i\in U}$; $\boldsymbol{g}^{\mathcal{E}}_{U}:=\{g^{\mathcal{E}}_{i}\}_{i\in U}$; $\boldsymbol{g}^{\mathcal{I}}_{U}:=\{g^{\mathcal{I}}_{i}\}_{i\in U}$; $\boldsymbol{g}_{U}:=\{g_{i}\}_{i\in U}$; $\boldsymbol{d}_{U}:=\{d_{i}\}_{i\in U}$; $\boldsymbol{Q}_{U}=\{Q_{i,j}\}_{i,j\in U}$, $\boldsymbol{A}^{\mathcal{E}}_{U}:=\{A^{\mathcal{E}}_{i,j}\}_{i,j\in U}$; $\boldsymbol{A}^{\mathcal{I}}_{U}:=\{A^{\mathcal{I}}_{i,j}\}_{i,j\in U}$; $\boldsymbol{A}_{U}:=\{A_{i,j}\}_{i,j\in U}$; $r_{U}:=\sum_{i\in U}r_{i}$; $m_{U}=\sum_{i\in U}m_{i}$; and $n_{U}=\sum_{i\in U}n_{i}$. The problem is denoted as the parametric form $P_{U}(\boldsymbol{d}_{U})$.

This section introduces the Schwarz algorithm for the solution of $P_{V}(\boldsymbol{d}_{V})$ (referred to as the full problem) with $V\subseteq\mathcal{V}$. We consider a non-overlapping partition $\{V_{k}\subseteq V\}_{k=1}^{K}$ of $V$ and an overlapping partition $\{W_{k}\subseteq V\}_{k=1}^{K}$ of $V$ such that $V_{k}\subseteq W_{k}$ holds for $k\in\mathbb{I}_{1:K}$. We call $V_{1},\cdots,V_{K}$ original subdomains and $W_{1},\cdots,W_{K}$ expanded subdomains. The Schwarz algorithm is defined below.

Here, we use a syntax that can be applied to any $U\subseteq\mathcal{V}$: $\boldsymbol{H}_{U}:=\{H_{i,j}\}_{i,j\in U}$; $\boldsymbol{H}_{-U}:=\{H_{i,j}\}_{i\in U,j\in N_{V}(U)}$; $H_{i,j}:=[Q_{i,j}\;A_{j,i}^{\top};A_{i,j}\;0]$; and $\boldsymbol{z}_{-U}:=\{z_{i}\}_{i\in N_{V}(U)}$. Furthermore, $\boldsymbol{z}^{*}_{U}(\cdot)$ is the primal-dual solution mapping of the parametric optimization problem $P_{U}(\cdot)$; $T_{U_{1}\leftarrow U_{2}}:=\{T_{i,j}\}_{i\in U_{1},j\in U_{2}}$, where $U_{1},U_{2}\subseteq\mathcal{V}$ and $T_{i,j}=I_{n_{i}\times n_{i}}$ if $i=j$ and $0_{n_{i}\times n_{j}}$ otherwise. Note that $\boldsymbol{z}_{-U}$ is supposed to represent the solution information that is complementary to $U$. The full complementary solution information includes the solution on $V\setminus U$. However, the variables and constraints in $U$ are coupled only with $N_{V}(U)$, so it suffices to incorporate information only for the coupled complementary region $N_{V}(U)$. Therefore, we will abuse the term complementary solution to represent the coupled complementary solution $\boldsymbol{z}_{-U}$.

The core part of the algorithm (line 3) consists of three steps: subproblem solution, solution restriction, and primal-dual exchange. In the first step, one formulates the subproblem for the $k$th subdomain as $P_{W_{k}}(\boldsymbol{d}_{W_{k}}-\boldsymbol{H}_{-W_{k}}\boldsymbol{z}^{(\ell)}_{-W_{k}})$ (this formulation will be justified later in Lemma 6). The subproblem incorporates complementary solution information $\boldsymbol{H}_{-W_{k}}\boldsymbol{z}^{(\ell)}_{-W_{k}}$, which is obtained during the third inner step of the previous iteration step. The subproblem is solved to obtain its solution $\boldsymbol{z}^{*}_{W_{k}}(\boldsymbol{d}_{W_{k}}-\boldsymbol{H}_{-W_{k}}\boldsymbol{z}^{(\ell)}_{-W_{k}})$. Here, we observe that solution multiplicity exists at the overlapping region. To remove such multiplicity, we restrict the solution in the second step. In particular, we abandon the primal-dual solutions associated with $W_{k}\setminus V_{k}$ (subdomain region acquired by expansion) and take only those solutions associated with $V_{k}$ (the original subdomain region). This procedure is represented by the restriction operator $T_{V_{k}\leftarrow W_{k}}$. After restriction, the solutions are assembled over $k\in\mathbb{I}_{1:K}$ to make the next guess of the solution $\boldsymbol{z}^{(\ell+1)}_{V}$. In the third step, the primal-dual solutions are exchanged across the subdomains to update the complementary information $\boldsymbol{H}_{-W_{k}}\boldsymbol{z}^{(\ell)}_{-W_{k}}$ for each subproblem. The schematic of the algorithm is shown in Fig. 1. The algorithm can be implemented in a fully decentralized manner, and different updating schemes can be used (e.g., Gauss-Seidel or asynchronous) [18].

In this section we analyze the convergence of the Schwarz algorithm. We will see that parametric sensitivity plays a central role in convergence behavior because, intuitively, primal-dual solutions of the neighbors of a subdomain enter as parametric perturbations. We first analyze the parametric sensitivity of $P_{U}(\cdot)$ in a way that the result can be generally applied to any node subset $U\subseteq\mathcal{V}$. We then exploit this sensitivity result to establish convergence.