[1]Zachary Frangella These authors contributed equally to this work.

These authors contributed equally to this work.

1]Management Science and Engineering, Stanford University, 475 Via Ortega,
Stanford, 94305, CA

2]Systems Engineering, Cornell University, 136 Hoy Rd, Ithaca, 14850, NY

3]Electrical Engineering and Computer Science, Massachusetts Institute of
Technology, 50 Vassar St, Cambridge, Cambridge, 02139, MA

4]Operations Research and Financial Engineering, Princeton University, 98
Charlton St, Princeton, 08540, NJ

# On the (linear) convergence of Generalized Newton Inexact ADMM

<EMAIL_ADDRESS>Shipu Zhao<EMAIL_ADDRESS>Theo Diamandis
<EMAIL_ADDRESS>Bartolomeo Stellato<EMAIL_ADDRESS>Madeleine Udell
<EMAIL_ADDRESS>[  [  [  [

###### Abstract

This paper presents GeNI-ADMM, a framework for large-scale composite convex
optimization that facilitates theoretical analysis of both existing and new
approximate ADMM schemes. GeNI-ADMM encompasses any ADMM algorithm that solves
a first- or second-order approximation to the ADMM subproblem inexactly. GeNI-
ADMM exhibits the usual $\mathcal{O}(1/t)$-convergence rate under standard
hypotheses and converges linearly under additional hypotheses such as strong
convexity. Further, the GeNI-ADMM framework provides explicit convergence
rates for ADMM variants accelerated with randomized linear algebra, such as
NysADMM and sketch-and-solve ADMM, resolving an important open question on the
convergence of these methods. This analysis quantifies the benefit of improved
approximations and can inspire the design of new ADMM variants with faster
convergence.

###### keywords:

ADMM, Inexactness, Linearization, Randomized Numerical Linear Algebra

##  1 Introduction

The Alternating Direction Method of Multipliers (ADMM) is one of the most
popular methods for solving composite optimization problems, as it provides a
general template for a wide swath of problems and converges to an acceptable
solution within a moderate number of iterations [1]. Indeed, [1] implicitly
promulgates the vision that ADMM provides a unified solver for various convex
machine learning problems. Unfortunately, for the large-scale problem
instances routinely encountered in the era of Big Data, ADMM scales poorly and
cannot provide a unified machine learning solver for problems of all scales.
The scaling issue arises as ADMM requires solving two subproblems at each
iteration, whose cost can increase superlinearly with the problem size. As a
concrete example, in the case of $\ell_{1}$-logistic regression with an
$n\times d$ data matrix, ADMM requires solving an $\ell_{2}$-regularized
logistic regression problem at each iteration [1]. With a fast-gradient
method, the total complexity of solving the subproblem is
$\tilde{\mathcal{O}}(nd\sqrt{\kappa})$ [2], where $\kappa$ is the condition
number of the problem. When $n$ and $d$ are in the tens of thousands or
larger—a moderate problem size by contemporary machine learning standards—and
$\kappa$ is large, such a high per-iteration cost becomes unacceptable. Worse,
ill-conditioning is ubiquitous in machine learning problems; often
$\kappa=\Omega(n)$, in which case the cost of the subproblem solve becomes
superlinear in the problem size.

Randomized numerical linear algebra (RandNLA) offers promising tools to scale
ADMM to larger problem sizes. Recently [3] proposed the algorithm NysADMM,
which uses a randomized fast linear system solver to scale ADMM up to problems
with tens of thousands of samples and hundreds of thousands of features. The
results in [3] show that ADMM combined with the RandNLA primitive runs 3 to 15
times faster than state-of-the-art solvers on machine learning problems from
LASSO to SVM to logistic regression. Unfortunately, the convergence of
randomized or approximate ADMM solvers like NysADMM is not well understood.
NysADMM approximates the $x$-subproblem using a linearization based on a
second-order Taylor expansion, which transforms the $x$-subproblem into a
Newton-step, i.e., a linear system solve. It then solves this system
approximately (and quickly) using a randomized linear system solver. The
convergence of this scheme, which combines linearization and inexactness, is
not covered by prior theory for approximate ADMM; prior theory covers either
linearization [4] or inexact solves [5] but not both.

In this work, we bridge the gap between theory and practice to explain the
excellent performance of a large class of approximate linearized ADMM schemes,
including NysADMM [3] and many methods not previously proposed. We introduce a
framework called Generalized Newton Inexact ADMM, which we refer to as GeNI-
ADMM (pronounced _genie-ADMM_). GeNI-ADMM includes NysADMM and many other
approximate ADMM schemes as special cases. The name is inspired by viewing the
linearized $x$-subproblem in GeNI-ADMM as a generalized Newton-step. GeNI-ADMM
allows for inexactness in both the $x$-subproblem and the $z$-subproblem. We
show GeNI-ADMM exhibits the usual $\mathcal{O}(1/t)$-convergence rate under
standard assumptions, with linear convergence under additional assumptions.
Our analysis also clarifies the value of using curvature in the generalized
Newton step: approximate ADMM schemes that take advantage of curvature
converge faster than those that do not at a rate that depends on the
conditioning of the subproblems. As the GeNI-ADMM framework covers any
approximate ADMM scheme that replaces the $x$-subproblem by a linear system
solve, our convergence theory covers any ADMM scheme that uses fast linear
system solvers. Given the recent flurry of activity on fast linear system
solvers within the (randomized) numerical linear algebra community [6, 7, 8],
our results will help realize these benefits for optimization problems as
well. To demonstrate the power of the GeNI-ADMM framework, we establish
convergence of NysADMM and another RandNLA-inspired scheme, sketch-and-solve
ADMM, whose convergence was left as an open problem in [9].

###  1.1 Contributions

Our contributions may be summarized concisely as follows:

  1. 1.

We provide a general ADMM framework GeNI-ADMM, that encompasses prior
approximate ADMM schemes as well as new ones. It can take advantage of second-
order information and allows for inexact subproblem solves.

  2. 2.

We show that GeNI-ADMM converges ergodically with the usual ergodic
$\mathcal{O}(1/t)$ rate, despite all the approximations it makes. Further, it
converges at a linear rate under strong convexity.

  3. 3.

We apply our framework to show some RandNLA-based approximate ADMM schemes
converge at the same rate as vanilla ADMM, answering some open questions
regarding their convergence [9, 3]. In the case of sketch-and-solve ADMM, we
show modifications to the naive scheme are required to ensure convergence.

###  1.2 Roadmap

In Section 2, we formally state the optimization problem that we focus on in
this paper and briefly introduce ADMM. In Section 3 we introduce the
Generalized Newton Inexact ADMM framework and review ADMM and its variants.
Section 4 gives various technical backgrounds and assumptions needed for our
analysis. Section 5 establishes that GeNI-ADMM converges at an
$\mathcal{O}(1/t)$-rate in the convex setting. Section 6 then shows that, when
the objective is strongly convex, GeNI-ADMM converges linearly. In Section 7,
we apply our theory to establish convergence rates for two methods that
naturally fit into our framework, and we illustrate these results numerically
in Section 8.

###  1.3 Notation and preliminaries

We call a matrix psd if it is positive semidefinite. We denote the convex cone
of $n\times n$ real symmetric psd matrices by $\mathbb{S}^{+}_{n}$. We denote
the Loewner ordering on $\mathbb{S}^{+}_{n}$ by $\preceq$, that is $A\preceq
B$ if and only if $B-A$ is psd. Given a matrix $H$, we denote its spectral
norm by $\|H\|$. If $f$ is a smooth function we denote its smoothness constant
by $L_{f}$. We say a positive sequence $\\{\varepsilon^{k}\\}_{k\geq 1}$ is
summable if $\sum_{k=1}^{\infty}\varepsilon^{k}<\infty$.

##  2 Problem statement and ADMM

Let $\mathcal{X},\mathcal{Z},\mathcal{H}$ be finite-dimensional inner-product
spaces with inner-product $\langle\cdot,\cdot\rangle$ and norm $\|\cdot\|$. We
wish to solve the convex constrained optimization problem

$\begin{array}[]{ll}\mbox{minimize}&f(x)+g(z)\\\ \mbox{subject to}&Mx+Nz=0\\\
&x\in X,~{}z\in Z,\end{array}$ (1)

with variables $x$ and $z$, where $X\subset\mathcal{X}$ and
$Z\subset\mathcal{Z}$ are closed convex sets, $f$ is a smooth convex function,
$g$ is a convex proper lower-semicontinuous (lsc) function, and
$M:\mathcal{X}\mapsto\mathcal{H}$ and $N:\mathcal{Z}\mapsto\mathcal{H}$ are
bounded linear operators.

######  Remark 1.

Often, the constraint is presented as $Mx+Nz=c$ for some non-zero vector $c$
however, by increasing the dimension of $\mathcal{Z}$ by $1$ and replacing $N$
by $[N~{}c]$, we can make $c=0$. Thus, our setting of $c=0$ is without loss of
generality.

We can write problem (1) as the saddle point problem

$\underset{{x\in X,z\in Z}}{\mbox{minimize}}\;\underset{y\in
Y}{\mbox{maximize}}\;f(x)+g(z)+\langle y,Mx+Nz\rangle,$ (2)

where $Y\subset\mathcal{Z}$ is a closed convex set. The saddle-point
formulation will play an important role in our analysis. Perform the change of
variables $u=y/\rho$ and define the Lagrangian

$L_{\rho}(x,z,u)\coloneqq f(x)+g(z)+\langle\rho u,Mx+Nz\rangle.$

Then (2) may be written concisely as

$\underset{x\in X,z\in Z}{\mbox{minimize}}\;\underset{{u\in
U}}{\mbox{maximize}}\;L_{\rho}(x,z,u).$ (3)

The ADMM algorithm (Algorithm 1) is a popular method for solving (1). Our
presentation uses the scaled form of ADMM [1, 10], which uses the change of
variables $u=y/\rho$, this simplifies the algorithms and analysis, so we
maintain this convention throughout the paper.

Algorithm 1 ADMM

convex proper lsc functions $f$ and $g$, constraint matrix $M$, stepsize
$\rho$

repeat

$x^{k+1}=\underset{x\in X}{\mathop{\rm
argmin}}\\{f(x)+\frac{\rho}{2}\|Mx+Nz^{k}+u^{k}\|^{2}\\}$

$z^{k+1}=\underset{z\in Z}{\mathop{\rm
argmin}}\\{g(z)+\frac{\rho}{2}\|Mx^{k+1}+Nz+u^{k}\|^{2}\\}$

$u^{k+1}=u^{k}+Mx^{k+1}+Nz^{k+1}$

until convergence

solution $(x^{\star},z^{\star})$ of problem (1)

##  3 Generalized Newton Inexact ADMM

As shown in Algorithm 1, at each iteration of ADMM, two subproblems are solved
sequentially to update variables $x$ and $z$. ADMM is often the method of
choice when the $z$-subproblem has a closed-form solution. For example, if
$g(x)=\|x\|_{1}$, the $z$-subproblem is soft thresholding, and if
$g(x)=1_{\mathcal{S}}(x)$ is the indicator function of a convex set
$\mathcal{S}$, the $z$-subproblem is projection onto $\mathcal{S}$ [11, §6].
However, it may be expensive to compute even with a closed-form solution. For
example, when $g(x)=1_{\mathcal{S}}(x)$ and $\mathcal{S}$ is the psd cone, the
$z$-subproblem requires an eigendecomposition to compute the projection, which
is prohibitively expensive for large problems [12].

Let us consider the $x$-subproblem

$x^{k+1}=\underset{x\in X}{\mathop{\rm
argmin}}\left\\{f(x)+\frac{\rho}{2}\|Mx+Nz^{k}+u^{k}\|^{2}\right\\}.$ (4)

In contrast to the $z$-subproblem, there is usually no closed-form solution
for the $x$-subproblem. Instead, an iterative scheme is often used to solve it
inaccurately, especially for large-scale applications. This solve can be very
expensive when the problem is large. To reduce computational effort, many
authors have suggested replacing this problem with a simplified subproblem
that is easier to solve. We highlight several strategies to do so below.

##### Augmented Lagrangian linearization.

One strategy is to linearize the augmented Lagrangian term
$\frac{\rho}{2}\|Mx+Nz^{k}+u^{k}\|^{2}$ in the ADMM subproblem and replace it
by the quadratic penalty $\frac{1}{2}\|x-x^{k}\|^{2}_{P}$ for some (carefully
chosen) positive definite matrix $P$. More formally, the strategy adds a
quadratic term to form a new subproblem

$x^{k+1}=\underset{x\in X}{\mathop{\rm
argmin}}\left\\{f(x)+\frac{\rho}{2}\|Mx+Nz^{k}+u^{k}\|^{2}+\frac{1}{2}\|x-x^{k}\|^{2}_{P}\right\\},$

which is substantially easier to solve for an appropriate choice of $P$. One
canonical choice is $P=\eta I-\rho M^{T}M$, where $\eta>0$ is a constant. For
this choice, the quadratic terms involving $M$ cancel, and we may omit
constants with respect to $x$, resulting in the subproblem

$x^{k+1}=\underset{x\in X}{\mathop{\rm argmin}}\left\\{f(x)+\rho\langle
Mx,Mx^{k}-z^{k}+u^{k}\rangle+\frac{\eta}{2}\|x-x^{k}\|^{2}\right\\}.$

Here, an isotropic quadratic penalty replaces the augmented Lagrangian term in
(4). This strategy allows the subproblem solve to be replaced by a proximal
operator with a (possibly) closed-form solution. This strategy has been
variously called _preconditioned ADMM_ , proximal ADMM, and (confusingly)
linearized ADMM [13, 14, 4].

##### Function approximation.

The second strategy to simplify the $x$-subproblem, is to approximate the
function $f$ by a first- or second-order approximation [4, 10, 3], forming the
new subproblem

$x^{k+1}=\underset{x\in X}{\mathop{\rm argmin}}\left\\{f(x^{k})+\langle\nabla
f(x^{k}),x-x^{k}\rangle+\frac{1}{2}\|x-x^{k}\|_{H}^{2}+\frac{\rho}{2}\|Mx+Nz^{k}+u^{k}\|^{2}\right\\}$
(5)

where $H$ is the Hessian of $f$ at $x^{k}$. The resulting subproblem is
quadratic and may be solved by solving a linear system or (for $M=I$)
performing a linear update, as detailed below in Section 3.1.1.

##### Inexact subproblem solve.

The third strategy is to solve the ADMM subproblems inexactly to achieve some
target accuracy, either in absolute error or relative error. An absolute-error
criterion chooses the subproblem error a priori [5], while a relative error
criterion requires the subproblem error to decrease as the algorithm nears
convergence, for example, by setting the error target at each iteration
proportional to $\|u^{k+1}-u^{k}-\rho(z^{k+1}-z^{k})\|$ [15].

##### Approximations used by GeNI-ADMM.

The GeNI-ADMM framework allows for any combination of the three strategies:
augmented Lagrangian linearization, function approximation, and inexact
subproblem solve. Consider the generalized second-order approximation to $f$

$\displaystyle f(x)\approx f_{1}(x)+f_{2}(\tilde{x}^{k})+\langle\nabla
f_{2}(\tilde{x}^{k}),x-\tilde{x}^{k}\rangle+\frac{1}{2}\|x-\tilde{x}^{k}\|^{2}_{\Theta^{k}},$
(6)

where $f_{1}+f_{2}=f$ and $\\{\Theta^{k}\\}_{k\geq 1}$ is a sequence of psd
matrices that approximate the Hessian of $f$. GeNI-ADMM uses this
approximation in the $x$-subproblem, resulting in the new subproblem

$\displaystyle\tilde{x}^{k+1}=\underset{x\in X}{\textrm{argmin}}$
$\displaystyle\\{f_{1}(x)+\langle\nabla
f_{2}(\tilde{x}^{k}),x-\tilde{x}^{k}\rangle+\frac{1}{2}\|x-\tilde{x}^{k}\|^{2}_{\Theta^{k}}$
(7) $\displaystyle+\frac{\rho}{2}\|Mx+N\tilde{z}^{k}+\tilde{u}^{k}\|^{2}\\}.$

We refer to (7) as a generalized Newton step. The intuition for the name is
made plain when $X=\mathbb{R}^{d}$, $f_{1}=0,f_{2}=f$ in which case the update
becomes

$\displaystyle\left(\Theta^{k}+\rho
M^{T}M\right)\tilde{x}^{k+1}=\Theta^{k}\tilde{x}^{k}-\nabla
f(\tilde{x}^{k})-\rho M^{T}(\tilde{u}^{k}-\tilde{z}^{k}).$ (8)

Equation (8) shows that the $x$-subproblem reduces to a linear system solve,
just like the Newton update. We can also interpret GeNI-ADMM as a linearized
proximal augmented Lagrangian (P-ALM) method [16, 17]. From this point of
view, GeNI-ADMM replaces $f$ in the P-ALM step by its linearization and adds a
specialized penalty defined by the $\Theta^{k}$-norm.

GeNI-ADMM also replaces the $z$-subproblem of ADMM with the following problem

$\tilde{z}^{k+1}=\underset{z\in Z}{\mathop{\rm
argmin}}\\{g(z)+\frac{\rho}{2}\|M\tilde{x}^{k+1}+Nz+\tilde{u}^{k}\|^{2}+\frac{1}{2}\|z-\tilde{z}^{k}\|^{2}_{\Psi^{k}}\\}.$
(9)

Incorporating the quadratic $\Psi^{k}$ term allows us to linearize the
augmented Lagrangian term, which is useful when $N$ is very complicated.

However, even with the allowed approximations, it is unreasonable to assume
that (7), (9) are solved exactly at each iteration. Indeed, the hallmark of
the NysADMM scheme from [3] is that it solves (7) inexactly (but efficiently)
using a randomized linear system solver. Thus, the GeNI-ADMM framework allows
for inexactness in the $x$ and $z$-subproblems, as seen in Algorithm 2.

Algorithm 2 Generalized Newton Inexact ADMM (GeNI-ADMM)

penalty parameter $\rho$, sequence of psd matrices $\\{\Theta^{k}\\}_{k\geq
1},\\{\Psi^{k}\\}_{k\geq 1}$, forcing sequences
$\\{\varepsilon_{x}^{k}\\}_{k\geq 1},\\{\varepsilon_{z}^{k}\\}_{k\geq 1},$

repeat

$\tilde{x}^{k+1}\overset{\varepsilon_{x}^{k}}{\approx}\underset{x\in
X}{\mathop{\rm argmin}}\\{f_{1}(x)+\langle\nabla
f_{2}(\tilde{x}^{k}),x-\tilde{x}^{k}\rangle+\frac{1}{2}\|x-\tilde{x}^{k}\|^{2}_{\Theta^{k}}+\frac{\rho}{2}\|Mx+N\tilde{z}^{k}+\tilde{u}^{k}\|^{2}\\}$

$\tilde{z}^{k+1}\overset{\varepsilon_{z}^{k}}{\approxeq}\underset{z\in
Z}{\mathop{\rm
argmin}}\\{g(z)+\frac{\rho}{2}\|M\tilde{x}^{k+1}+Nz+\tilde{u}^{k}\|^{2}+\frac{1}{2}\|z-\tilde{z}^{k}\|^{2}_{\Psi^{k}}\\}$

$\tilde{u}^{k+1}=\tilde{u}^{k}+M\tilde{x}^{k+1}+N\tilde{z}^{k+1}$

until convergence

solution $(x^{\star},z^{\star})$ of problem (1)

Algorithm 2 differs from ADMM (Algorithm 1) in that a) the $x$ and
$z$-subproblems are now given by (7) and (9) and b) both subproblems may be
solved inexactly. Given the inexactness and the use of the generalized Newton
step in place of the original $x$-subproblem, we refer to Algorithm 2 as
Generalized Newton Inexact ADMM (GeNI-ADMM). The inexactness schedule is
controlled by the forcing sequences $\\{\varepsilon_{x}^{k}\\}_{k\geq
1},\\{\varepsilon_{z}^{k}\\}_{k\geq 1}$ that specify how accurately the $x$
and $z$-subproblems are solved at each iteration. These subproblems have
different structures that require different notions of accuracy. To
distinguish between them, we make the following definition.

######  Definition 1 ($\varepsilon$-minimizer and $\varepsilon$-minimum).

Let $h:\mathcal{T}\mapsto\mathbb{R}$ be strongly-convex and let
$t^{\star}=\mathop{\rm argmin}_{t^{\prime}\in\mathcal{T}}h(t^{\prime})$.

  * •

($\varepsilon$-minimizer) Given $t\in\mathcal{T}$, we say $t$ is an
$\varepsilon$-_minimizer_ of $\textrm{minimize}_{t\in\mathcal{T}}h(t)$ and
write

$t\overset{\varepsilon}{\approx}\mathop{\rm
argmin}_{t^{\prime}\in\mathcal{T}}h(t^{\prime})\quad\textrm{if and only
if}\quad\|t-t^{\star}\|\leq\varepsilon.$

In words, $t$ is nearly equal to the argmin of $h(t)$ in set $\mathcal{T}$.

  * •

($\varepsilon$-minimum) Given $t\in\mathcal{T}$, we say $t$ gives an
$\varepsilon$-_minimum_ of $\textrm{minimize}_{t\in\mathcal{T}}h(t)$ and write

$t\overset{\varepsilon}{\approxeq}\mathop{\rm
argmin}_{t^{\prime}\in\mathcal{T}}h(t^{\prime})\quad\textrm{if and only
if}\quad h(t)-h(t^{\star})\leq\varepsilon.$

In words, $t$ produces nearly the same objective value as the argmin of $h(t)$
in set $\mathcal{T}$.

Thus, from Definition 1 and Algorithm 2, we see for each iteration $k$ that
$\tilde{x}^{k+1}$ is an $\varepsilon_{x}^{k}$-minimizer of the $x$-subproblem,
while $\tilde{z}^{k+1}$ gives an $\varepsilon_{z}^{k}$-minimum of the
$z$-subproblem.

###  3.1 Related work

The literature on the convergence of ADMM and its variants is vast, so we
focus on prior work most relevant to our setting. Section 3.1 lists some of
the prior work that developed and analyzed the approximation strategies
described in Section 3. GeNI-ADMM differs from all prior work in Section 3.1
by allowing (almost) all these approximations and more. It also provides
explicit rates of convergence to support choices between algorithms. Moreover,
many of these algorithms can be recovered from GeNI-ADMM.

Reference |  Convergence |  |  |  |
---|---|---|---|---|---
rate |  Problem |  |  |  |
class |  Augmented Lagrangian Linearization |  Function |  |  |
approximation |  Subproblem |  |  |  |
inexactness |  |  |  |  |
Eckstein and Bertsekas [5] |  ✗ |  Convex |  ✗ |  ✗ |  $x,z$
absolute error |  |  |  |  |
He and Yuan [14] |  Ergodic $\mathcal{O}(1/t)$ |  Convex |  $x$ |  ✗ |  ✗
Monteiro and Svaiter [18] |  Ergodic $\mathcal{O}(1/t)$ |  Convex |  ✗ |  ✗ |  ✗
Ouyang et al. [19] |  $\mathcal{O}(1/\sqrt{t})$ |  Stochastic |  |  |
convex |  ✗ |  $f$ |  |  |
stochastic first-order |  ✗ |  |  |  |
Ouyang et al. [4] |  Ergodic $\mathcal{O}(1/t)$ |  Convex |  $x$ |  $f$ |
first-order |  ✗ |  |  |  |
Deng and Yin [13] |  Linear |  Strongly convex |  $x,z$ |  ✗ |  ✗
Eckstein and Yao [15] |  ✗ |  Convex |  ✗ |  ✗ |  $x$
relative error |  |  |  |  |
Hager and Zhang [20] |  Ergodic $\mathcal{O}(1/t),$ |  |  |  |
$\mathcal{O}(1/t^{2})$ |  Convex, |  |  |  |
Strongly convex |  $x$ |  ✗ |  $x$ |  |
relative error |  |  |  |  |
Yuan et al. [21] |  Locally |  |  |  |
linear |  Convex |  $x,z$ |  ✗ |  ✗ |
Ryu and Yin [10] |  ✗ |  Convex |  $x,z$ |  $f,g$ |
partial first-order |  ✗ |  |  |  |
This work |  Ergodic $\mathcal{O}(1/t),$ |  |  |  |
Linear |  Convex, |  |  |  |
Strongly convex |  $x,z$ |  $f$ |  |  |
partial generalized second-order |  $x,z$ |  |  |  |
absolute error |  |  |  |  |

Table 1: A structured comparison of related work on the convergence of ADMM
and its variants. The “$x$” (“$z$”) in the table denotes that a paper uses the
corresponding strategy of the column to simplify the $x$-subproblem
($z$-subproblem). In the “Function approximation” column, the “$f$” (“$g$”)
indicates that a paper approximates function “$f$” (“$g$”) in the
$x$-subproblem ($z$-subproblem). “stochastic first-order” means a paper uses
first-order function approximation but replace the gradient term with a
stochastic gradient. “partial generalized second-order” means a paper uses
second-order function approximation as (7), but $H$ is not necessarily the
Hessian.

####  3.1.1 Algorithms recovered from GeNI-ADMM

Various ADMM schemes in the literature can be recovered by appropriately
selecting the parameters in Algorithm 2. Let us consider a few special cases
to highlight important prior work on approximate ADMM, and provide concrete
intuition for the general framework provided by Algorithm 2.

##### NysADMM

The NysADMM scheme [3] assumes the problem is unconstrained $M=I,N=-I$. We
maintain this assumption for simplicity of exposition, although it is
unnecessary. GeNI-ADMM specializes to NysADMM taking
$\Theta^{k}=\eta(H_{f}^{k}+\sigma I)$, and $\Psi^{k}=0$, where $H_{f}^{k}$
denotes the Hessian of $f$ at the $k$th iteration. Unlike the original NysADMM
scheme of [3], we include the regularization term $\sigma I$, where
$\sigma>0$. This inclusion is required for theoretical analysis but seems
unnecessary in practice (see Section 7.1 for a detailed discussion).
Substituting this information into (8) leads to the following update for the
$x$-subproblem.

$\tilde{x}^{k+1}=\tilde{x}^{k}-\left(\eta
H_{f}^{k}+(\rho+\eta\sigma)I\right)^{-1}\left(\nabla
f(\tilde{x}^{k})+\rho(\tilde{x}^{k}-\tilde{z}^{k}+\tilde{u}^{k})\right).$ (10)

##### Sketch-and-solve ADMM

If $M=I,N=-I$, $\Psi^{k}=0$, and $\Theta^{k}$ is chosen to be a matrix such
that $\Theta^{k}+\rho I$ is easy to factor, we call the resulting scheme
_sketch-and-solve ADMM_. The name _sketch-and-solve ADMM_ is motivated by the
fact that such a $\Theta^{k}$ can often be obtained via sketching techniques.
However, the method works for any $\Theta^{k}$, not only ones constructed by
sketching methods. The update is given by

$\tilde{x}^{k+1}=\tilde{x}^{k}-\left(\Theta^{k}+\rho I\right)^{-1}\left(\nabla
f(\tilde{x}^{k})+\rho(\tilde{x}^{k}-\tilde{z}^{k}+\tilde{u}^{k})\right).$

We provide the full details of sketch-and-solve ADMM in Section 7.2. To our
knowledge, sketch-and-solve ADMM has not been previously proposed or analyzed.

##### Proximal/Generalized ADMM

Set $f_{1}=f,f_{2}=0$, and fix $\Theta^{k}=\Theta$ and $\Psi^{k}=\Psi$, where
$\Theta$ and $\Psi$ are symmetric positive definite matrices. Then GeNI-ADMM
simplifies to

$\displaystyle\tilde{x}^{k+1}=\mathop{\rm argmin}_{x\in
X}\left\\{f(x)+\frac{\rho}{2}\|Mx+N\tilde{z}^{k}+\tilde{u}^{k}\|^{2}+\frac{1}{2}\|x-\tilde{x}^{k}\|_{\Theta}^{2}\right\\},$
$\displaystyle\tilde{z}^{k+1}=\mathop{\rm argmin}_{z\in
Z}\left\\{g(z)+\frac{\rho}{2}\|M\tilde{x}^{k+1}+Nz+\tilde{u}^{k}\|^{2}+\frac{1}{2}\|z-\tilde{z}^{k}\|_{\Psi}^{2}\right\\},$

which is the Proximal/Generalized ADMM of Deng and Yin [13].

##### Linearized ADMM

Set $\Theta^{k}=\alpha^{-1}I-\rho M^{T}M$ and $\Psi=\beta^{-1}I-\rho N^{T}N$,
for all $k\geq 1$. Then GeNI-ADMM simplifies to

$\displaystyle\tilde{x}^{k+1}=\textup{prox}_{\alpha
f}\left(\tilde{x}^{k}-\alpha\rho
M^{T}(M\tilde{x}^{k}+N\tilde{z}^{k}+\tilde{u}^{k})\right),$
$\displaystyle\tilde{z}^{k+1}=\textup{prox}_{\beta
g}\left(\tilde{z}^{k}-\beta\rho
N^{T}(M\tilde{x}^{k+1}+N\tilde{z}^{k}+\tilde{u}^{k})\right),$

which is exactly Linearized ADMM [11, §4.4.2].

##### Primal Dual Hybrid Gradient

The Primal Dual Hybrid Gradient (PDHG) or Chambolle-Pock algorithm of
Chambolle and Pock [22] is a special case of GeNI-ADMM since PDHG is a special
case of Linearized ADMM (see section 3.5 of [10]).

##### Gradient descent ADMM

Consider two linearization schemes from Ouyang et al. [4]. The first scheme we
call _gradient descent ADMM_ (GD-ADMM) is useful when it is cheap to solve a
least-squares problem with $M$. GD-ADMM is obtained from GeNI-ADMM by setting
$f_{1}=0,f_{2}=f$, and $\Theta^{k}=\eta I$ for all $k$. The $\tilde{x}^{k+1}$
update (8) for GD-ADMM simplifies to

$\tilde{x}^{k+1}=\tilde{x}^{k}-(\rho M^{T}M+\eta I)^{-1}\left(\nabla
f(\tilde{x}^{k})+\rho
M^{T}(M\tilde{x}^{k}+N\tilde{z}^{k}+\tilde{u}^{k})\right).$ (11)

The second scheme, _linearized-gradient descent ADMM_ (LGD-ADMM), is useful
when $M$ is not simple, so that the update (11) is no longer cheap. To make
the $x$-subproblem update cheaper, it linearizes the augmented Lagrangian term
by setting $\Theta^{k}=\alpha^{-1}I-\rho M^{T}M$ for all $k$ in addition to
linearizing $f$. In this case, (8) yields the $\tilde{x}^{k+1}$ update

$\tilde{x}^{k+1}=\tilde{x}^{k}-\alpha\left(\nabla f(\tilde{x}^{k})+\rho
M^{T}(M\tilde{x}^{k}+N\tilde{z}^{k}+\tilde{u}^{k})\right).$ (12)

Observe in the unconstrained case, when $M=I$, the updates (11) and (12) are
equivalent and generate the same iterate sequences when initialized at the
same point [23]. Indeed, they are both generated by performing a gradient step
on the augmented Lagrangian (7), for suitable choices of the parameters.
Notably, this terminology differs from [4], who refer to (11) as “linearized
ADMM” (L-ADMM) and (12) as “linearized preconditioned ADMM” (LP-ADMM). We
choose our terminology to emphasize that GD-ADMM accesses $f$ via its
gradient, as in the literature the term “linearized ADMM” is usually reserved
for methods that access $f$ through its prox operator [14, 11, 13].

In the remainder of this paper, we will prove convergence of Algorithm 2 under
appropriate hypotheses on the sequences $\\{\Theta^{k}\\}_{k\geq
1},\\{\Psi^{k}\\}_{k\geq 1},\\{\varepsilon_{x}^{k}\\}_{k\geq 1},$ and
$\\{\varepsilon_{z}^{k}\\}_{k\geq 1}$.

##  4 Technical preliminaries and assumptions

We introduce some important concepts that will be central to our analysis. The
first is $\Theta$-relative smoothness, which is crucial to establish that
GeNI-ADMM benefits from curvature information provided by the Hessian.

######  Definition 2 ($\Theta$-relative smoothness).

We say $f:\mathcal{D}\to\mathbb{R}$ is $\Theta$-relatively smooth with respect
to the bounded function $\Theta:\mathcal{D}\to\mathbb{S}_{+}^{n}$ if there
exists $\hat{L}_{\Theta}>0$ such that for all $x,y\in\mathcal{D}$

$f(x)\leq f(y)+\langle\nabla
f(y),x-y\rangle+\frac{\hat{L}_{\Theta}}{2}\|x-y\|^{2}_{\Theta(y)}.$ (13)

That is, the function $f$ is smooth with respect to the $\Theta$-norm.
Definition 2 generalizes _relative smoothness_ , introduced in [24] to analyze
Newton’s method. The definition in [24] takes $\Theta$ to be the Hessian of
$f$, $H_{f}$. When $f$ belongs to the popular family of generalized linear
models, then (13) holds with a value of $\hat{L}_{H_{f}}$ that is independent
of the conditioning of the problem [24]. For instance, if $f$ is quadratic and
$\Theta(y)=H_{f}$, then (13) holds with equality for $\hat{L}_{H_{f}}=1$.
Conversely, if we take $\Theta=I$, which corresponds to GD-ADMM, then
$\hat{L}_{\Theta}=L_{f}$, the smoothness constant of $f$. Our theory uses the
fact that $\hat{L}_{\Theta}$ is much smaller than the smoothness constant
$L_{f}$ for methods that take advantage of curvature, and relies on
$\hat{L}_{\Theta}$ to characterize the faster convergence speed of these
methods.

The other important idea we need is the notion of an $\varepsilon$-subgradient
[25, 26].

######  Definition 3 ($\varepsilon$-subgradient).

Let $r:\mathcal{D}\to\mathbb{R}$ be a convex function and $\varepsilon>0$. We
say that $s\in\mathcal{D}^{*}$ is an $\varepsilon$-subgradient for $r$ at
$z\in\mathcal{D}$ if, for every $z^{\prime}\in\mathcal{D}$, we have

$r(z^{\prime})-r(z)\geq\langle s,z^{\prime}-z\rangle-\varepsilon.$

We denote the set of $\varepsilon$-subgradients for $r$ at $z$ by
$\partial_{\varepsilon}r(z)$.

Clearly, any subgradient is an $\varepsilon$-subgradient, so Definition 3
provides a natural weakening of a subgradient. The $\varepsilon$-subgradient
is critical for analyzing $z$-subproblem inexactness, and our usage in this
context is inspired by the convergence analysis of inexact proximal gradient
methods [27]. We shall need the following proposition whose proof may be found
in [25, Proposition 4.3.1].

######  Proposition 1.

Let $r$, $r_{1}$, and $r_{2}$ be convex functions. Then for any $z$, the
following holds:

  1. 1.

$0\in\partial_{\varepsilon}r(z)$ if and only if
$r(z)\overset{\varepsilon}{\approxeq}\mathop{\rm
argmin}_{z^{\prime}}r(z^{\prime})$, that is $z$ gives an $\varepsilon$-minimum
of $\textup{minimize}_{z^{\prime}}r(z^{\prime})$.

  2. 2.

$\partial_{\varepsilon}(r_{1}+r_{2})(z)\subset\partial_{\varepsilon}r_{1}(z)+\partial_{\varepsilon}r_{2}(z).$

With Proposition 1 recorded, we prove the following lemma in Section A.2 of
the appendix, which will play a critical role in establishing the convergence
of GeNI-ADMM.

######  Lemma 1.

Let $\tilde{z}^{k+1}$ give an $\varepsilon^{k}_{z}$-minimum of the
$z$-subproblem. Then there exists an $\tilde{s}$ with $\|\tilde{s}\|\leq
C\sqrt{\varepsilon_{z}^{k}}$ such that

$-\rho
N^{T}\tilde{u}^{k+1}+\Psi_{k}(\tilde{z}^{k}-\tilde{z}^{k+1})+\tilde{s}\in\partial_{\varepsilon_{z}^{k}}g(\tilde{z}^{k+1}).$

###  4.1 Assumptions

In this section, we present the main assumptions required by our analysis.

######  Assumption 1 (Existence of saddle point).

There exists an optimal primal solution $(x^{\star},z^{\star})\in X\times Z$
for (1) and an optimal dual solution $u^{\star}\in U$ such that
$(x^{\star},z^{\star},u^{\star})$ is a saddle point of (2). Here,
$U\subset\mathcal{Z}$ is a closed convex set and $\rho U=Y$. We denote the
optimal objective value of (1) as $p^{\star}$.

Assumption 1 is standard and merely assumes that (1) has a solution.

######  Assumption 2 (Regularity of $f$ and $g$).

The function $f$ is twice-continuously differentiable and is $1$-relatively
smooth with respect to $\Theta$. The function $g$ is finite-valued, convex,
and lower semi-continuous.

Assumption 2 is also standard. It ensures that it makes sense to talk about
the Hessian of $f$ and that $f$ is relatively smooth. Note assuming
$\hat{L}_{\Theta}=1$ is without loss of generality, for we can always redefine
$\Theta^{\prime}=\hat{L}_{\Theta}\Theta$, and $f$ will be $1$-relatively
smooth with respect to $\Theta^{\prime}$.

######  Assumption 3 (Forcing sequence summability and approximate subproblem
oracles).

Let $\\{\varepsilon^{k}_{x}\\}_{k\geq 1}$ and
$\\{\varepsilon^{k}_{z}\\}_{k\geq 1}$ be given forcing sequences, we assume
they satisfy

$\mathcal{E}_{x}=\sum_{k=1}^{\infty}\varepsilon_{x}^{k}<\infty,\quad\mathcal{E}_{z}=\sum^{\infty}_{k=1}\sqrt{\varepsilon_{z}^{k}}<\infty.$

Further we define the constants $K_{\varepsilon_{x}}\coloneqq\sup_{k\geq
1}\varepsilon_{x}^{k}$, $K_{\varepsilon_{z}}\coloneqq\sup_{k\geq
1}\sqrt{\varepsilon_{z}^{k}}$. Observe $K_{\varepsilon_{x}}$ and
$K_{\varepsilon_{z}}$ are finite owing to the summability hypotheses.
Moreover, we assume Algorithm 2 is equipped with oracles for solving the $x$
and $z$-subproblems, which at each iteration produce approximate solutions
$\tilde{x}^{k+1}$, $\tilde{z}^{k+1}$ satisfying:

$\displaystyle\tilde{x}^{k+1}\overset{\varepsilon_{x}^{k}}{\approx}\underset{x\in
X}{\mathop{\rm argmin}}\\{f_{1}(x)+\langle\nabla
f_{2}(\tilde{x}^{k}),x-\tilde{x}^{k}\rangle+\frac{1}{2}\|x-\tilde{x}^{k}\|^{2}_{\Theta^{k}}+\frac{\rho}{2}\|Mx+N\tilde{z}^{k}+\tilde{u}^{k}\|^{2}\\},$
$\displaystyle\tilde{z}^{k+1}\overset{\varepsilon_{z}^{k}}{\approxeq}\underset{z\in
Z}{\mathop{\rm
argmin}}\\{g(z)+\frac{\rho}{2}\|M\tilde{x}^{k+1}+Nz+\tilde{u}^{k}\|^{2}+\frac{1}{2}\|z-\tilde{z}^{k}\|^{2}_{\Psi^{k}}\\}.$

The conditions on the $x$ and $z$ subproblem oracles are consistent with those
of [5], which requires the sum of the errors in the subproblems to be
summable. The approximate solution criteria of Assumption 3 are also easily
met in practical applications, as the subproblems are either simple enough to
solve exactly, or can be efficiently solved approximately via iterative
algorithms. For instance, with LGD-ADMM the $x$-subproblem is simple to solve
enough exactly, while for NysADMM the $x$-subproblem is efficiently solved via
conjugate gradient with a randomized preconditioner.

######  Assumption 4 (Regularity of $\\{\Theta^{k}\\}_{k\geq 1}$ and
$\\{\Psi^{k}\\}_{k\geq 1}$).

We assume there exists constants
$\theta_{\textup{max}},~{}\theta_{\textup{min}},~{}\psi_{\textup{max}},$ and
$\nu>0$, such that

$\theta_{\textup{max}}I\succeq\Theta^{k}\succeq\theta_{\textup{min}}I,~{}\psi_{\textup{max}}I\succeq\Psi^{k},~{}\Psi^{k}+\rho
N^{T}N\succeq\nu I,\quad\text{for all}~{}k\geq 1.$ (14)

Moreover, we also assume that The sequences $\\{\Theta^{k}\\}_{k\geq
1},\\{\Psi^{k}\\}_{k\geq 1}$ satisfy

$\displaystyle\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k}}\leq(1+\zeta^{k-1})\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k-1}},$
(15)
$\displaystyle\|\tilde{z}^{k}-z^{\star}\|^{2}_{\Psi^{k}}\leq(1+\zeta^{k-1})\|\tilde{z}^{k}-z^{\star}\|^{2}_{\Psi^{k-1}},$

where $\\{\zeta^{k}\\}_{k\geq 1}$ is a non-negative summable sequence, that
is, $\sum^{\infty}_{k=1}\zeta^{k}<\infty$.

The first half of Assumption 4 is standard, and essentially requires that
$\Theta^{k}$ and $\Psi^{k}+\nu N^{T}N$ define norms. For most common choices
of $\Theta^{k}$ and $\Psi^{k}$ these assumptions are satisfied, they are also
readily enforced by adding on a small multiple of the identity. The addition
of a small regularization term is common practice in popular ADMM solvers,
such as OSQP [28], as it avoids issues with degeneracy and increases numerical
stability.

The second part of the assumption requires that $\Theta^{k}$ ($\Psi^{k})$ and
$\Theta^{k-1}$ ($\Psi^{k-1}$) eventually not differ much on the distance of
$\tilde{x}^{k}$ ($\tilde{z}^{k}$) to $x^{\star}$ ($z^{\star}$). Assumptions of
this form are common in analyses of optimization algorithms that use variable
metrics. For instance, He et al. [29] which develops an alternating directions
method for monotone variational inequalities, assumes their equivalents of the
sequences $\\{\Theta^{k}\\}$ and $\\{\Psi^{k}\\}$ satisfy

$(1-\zeta^{k-1})\Theta^{k-1}\preceq\Theta^{k}\preceq(1+\zeta^{k-1})\Theta^{k-1},~{}(1-\zeta^{k-1})\Psi^{k-1}\preceq\Psi^{k}\preceq(1+\zeta^{k-1})\Psi^{k-1}.$
(16)

More recently, Rockafellar [30] analyzed the proximal point method with
variable metrics, under the assumption the variable metrics satisfy (16).
Thus, Assumption 4 is consistent with the literature and is, in fact, weaker
than prior work, as (16) implies (15). Assumption 4 may always be enforced by
changing the $\Theta^{k}$’s and $\Psi^{k}$’s only finitely many times.

We also define the following constants which shall be useful in our analysis,

$\tau_{\zeta}\coloneqq\prod_{k\geq
2}(1+\zeta^{k}),\quad\mathcal{E}_{\zeta}=\tau_{\zeta}\left(\sum^{\infty}_{k=1}\zeta^{k}\right)<\infty.$

Note Assumption 4 implies $\tau_{\zeta}<\infty$. Moreover, for any $k\geq 1$
it holds that

$\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k}}\leq\tau_{\zeta}\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{1}},~{}\|\tilde{z}^{k}-z^{\star}\|^{2}_{\Psi^{k}}\leq\tau_{\zeta}\|\tilde{z}^{k}-z^{\star}\|^{2}_{\Psi^{1}}.$

##  5 Sublinear convergence of GeNI-ADMM

This section establishes our main theorem, Theorem 1, which shows that
Algorithm 2 enjoys the same $\mathcal{O}(1/t)$-convergence rate as standard
ADMM.

######  Theorem 1 (Ergodic convergence).

Define constants
$d_{x^{\star},\Theta^{1}}=\|\tilde{x}^{1}-x^{\star}\|_{\Theta^{1}}$,
$d_{u^{\star}}=\|\tilde{u}^{1}-u^{\star}\|$,
$d_{z^{\star},\Psi^{1}_{\rho,N}}=\|\tilde{z}^{1}-z^{\star}\|_{\Psi^{1}+\rho
N^{T}N}$. Let $p^{\star}$ denote the optimum of (1). For each $t\geq 1$,
denote $\bar{x}^{t+1}=\frac{1}{t}\sum^{t+1}_{k=2}\tilde{x}^{k}$, and
$\bar{z}^{t+1}=\frac{1}{t}\sum^{t+1}_{k=2}\tilde{z}^{k}$, where
$\\{\tilde{x}^{k}\\}_{k\geq 1}$ and $\\{\tilde{z}^{k}\\}_{k\geq 1}$ are the
iterates produced by Algorithm 2 with forcing sequences
$\\{\varepsilon^{k}_{x}\\}_{k\geq 1}$ and $\\{\varepsilon^{k}_{z}\\}_{k\geq
1}$ with $\tilde{u}^{1}=0$ and $M\tilde{x}^{1}=-N\tilde{z}^{1}$. Instate
Assumptions 1-4. Then, the suboptimality gap satisfies

Furthermore, the feasibility gap satisfies

$\|M\bar{x}^{t+1}+N\bar{z}^{t+1}\|\leq\frac{2}{t}\sqrt{\frac{\Gamma}{\rho}+d^{2}_{u^{\star}}},$

Consequently, after $\mathcal{O}(1/\epsilon)$ iterations,

$f(\bar{x}^{t+1})+g(\bar{z}^{t+1})-p^{\star}\leq\epsilon,\;\text{and}\;\;\;\|M\bar{x}^{t+1}+N\bar{z}^{t+1}\|\leq\epsilon.$

Theorem 1 shows that, with a constant value of $\eta$ and appropriate forcing
sequences, the suboptimality gap and the feasibility residuals both go to zero
at a rate of $\mathcal{O}(1/t)$. Hence, the overall convergence rate of ADMM
is preserved despite all the approximations involved in GeNI-ADMM. However,
some schemes yield better constants than others.

To see this difference, consider the case when $f$ is a convex quadratic
function and $M=I$. Let $H_{f}$ be the Hessian of $f$. Consider (a) NysADMM
$\Theta^{k}=H_{f}$ (10) and (b) GD-ADMM $\Theta^{k}=I$ (11). Further, suppose
that NysADMM and GD-ADMM are initialized at $0$. The rates of convergence
guaranteed by Theorem 1 for Algorithm 2 for both methods are outlined in Table
2.

Table 2: Convergence rate comparison of NysADMM and GD-ADMM when initialized
at $0$ for quadratic $f$.
Method |  NysADMM $\left(\Theta^{k}=H_{f}\right)$ |  GD-ADMM $\left(\Theta^{k}=I\right)$
---|---|---
Feasibility gap | $\mathcal{O}\left(\frac{1}{t}\left(\sqrt{\frac{2}{\rho}\|x^{\star}\|^{2}_{H_{f}+\sigma I}}\right)\right)$ | $\mathcal{O}\left(\frac{1}{t}\left(\sqrt{\frac{2}{\rho}\lambda_{1}\left(H_{f}\right)\|x^{\star}\|^{2}}\right)\right)$
Suboptimality gap | $\mathcal{O}\left(\frac{1}{2t}\|x^{\star}\|^{2}_{H_{f}+\sigma I}\right)$ | $\mathcal{O}\left(\frac{1}{2t}\lambda_{1}\left(H_{f}\right)\|x^{\star}\|^{2}\right)$

In the first case, the relative smoothness constant satisfies
$\hat{L}_{\Theta}=1$, while in the second,
$\hat{L}_{\Theta}=\lambda_{1}\left(H_{f}\right)$, which is the largest
eigenvalue of $H_{f}$. Comparing the rates in Table 2, we see NysADMM improves
over GD-ADMM whenever

$\left(\lambda_{1}\left(H_{f}\right)-\sigma\right)\|x^{\star}\|^{2}\geq\|x^{\star}\|^{2}_{H_{f}}.$

Hence NysADMM improves significantly over GD-ADMM when $H_{f}$ exhibits a
decaying spectrum, provided $x^{\star}$ is not concentrated on the top
eigenvector of $H_{f}$. We formalize the latter property in the following
definition.

######  Definition 4 ($\mu_{\mathbf{v}}$-incoherence).

Let $(\lambda_{1}(H_{f}),\mathbf{v})$ be the largest eigenpair of $H_{f}$. We
say $x^{\star}$ is $\mu_{\mathbf{v}}$-incoherent if there exists
$0\leq\mu_{\mathbf{v}}<1$ such that:

$|\langle\mathbf{v},x^{\star}\rangle|^{2}\leq\mu_{\mathbf{v}}\|x^{\star}\|^{2}.$
(17)

Definition 4 is a weak form of the incoherence condition from compressed
sensing and matrix completion, which plays a key role in signal and low-rank
matrix recovery [31, 32]. In words, $x^{\star}$ is
$\mu_{\mathbf{v}}$-incoherent if its energy is not solely concentrated on the
top eigenvector of $H_{f}$ and can be expected to hold generically. The
parameter $\mu_{\mathbf{v}}$ controls the allowable concentration. When
$\mu_{\mathbf{v}}=0$, $x^{\star}$ is orthogonal to $\mathbf{v}$, so its energy
is distributed amongst the other eigenvectors. Conversely, the closer
$\mu_{\mathbf{v}}$ is to $1$, the more $x^{\star}$ is allowed to concentrate
on $\mathbf{v}$.

Using $\mu_{\mathbf{v}}$-incoherence, we can say more about how NysADMM
improves on GD-ADMM.

######  Proposition 2.

Suppose $x^{\star}$ is $\mu_{\mathbf{v}}$-incoherent and
$\sigma=\tau\lambda_{1}(H_{f})$ where $\tau\in(0,1)$. Then, the following
bound holds

$\frac{\left(\lambda_{1}(H_{f})-\sigma\right)\|x^{\star}\|^{2}}{\|x^{\star}\|^{2}_{H_{f}}}\geq(1-\tau)\left(\mu_{\mathbf{v}}+\frac{\lambda_{2}(H_{f})}{\lambda_{1}(H_{f})}\right)^{-1}.$

Hence if
$\mu_{\mathbf{v}}+\lambda_{2}(H_{f})/\lambda_{1}(H_{f})\leq(1-\tau)\alpha^{-1}$,
where $\alpha\geq 1$, then

$\frac{\lambda_{1}(H_{f})\|x^{\star}\|^{2}}{\|x^{\star}\|^{2}_{H_{f}}}\geq\alpha.$

The proof of Proposition 2 is given in Section A.2. Proposition 2 shows when
$x^{\star}$ is $\mu_{\mathbf{v}}$-incoherent and $H_{f}$ has a decaying
spectrum, NysADMM yields a significant improvement over GD-ADMM. As a concrete
example, consider when $\alpha=2$, then Proposition 2 implies the ergodic
convergence of NysADMM is twice as fast as GD-ADMM. We observe this
performance improvement in practice; see Section 8 for corroborating numerical
evidence. Just as Newton’s method improves on gradient descent for ill-
conditioned problems, NysADMM is less sensitive to ill-conditioning than GD-
ADMM.

We now move to proving Theorem 1.

###  5.1 Our approach

To prove Theorem 1, we take the approach in Ouyang et al. [4], and analyze
GeNI-ADMM by viewing Eq. 1 through its formulation as saddle point problem Eq.
2. Let $W=X\times Z\times U$, $\hat{w}=(\hat{x},\hat{z},\hat{u})$, and
$w=(x,z,u)$, where $w,\hat{w}\in W$. Define the _gap function_

$Q(\hat{w},w)\coloneqq[f(x)+g(z)+\langle\rho\hat{u},Mx+Nz\rangle]-[f(\hat{x})+g(\hat{z})+\langle\rho
u,M\hat{x}+N\hat{z}\rangle].$ (18)

The gap function naturally arises from the saddle-point formulation Eq. 2. By
construction, it is concave in its first argument and convex in the second,
and satisfies the important inequality

$Q(w^{\star},w)=L_{\rho}(x,z,u^{\star})-L_{\rho}(x^{\star},z^{\star},u)\geq
0,$

which follows by definition of $(x^{\star},z^{\star},u^{\star})$ being a
saddle-point. Hence the gap function may be viewed as measuring the distance
to the saddle $w^{\star}$.

Further, given a closed set $U\subset\mathcal{Z}$ and $v\in\mathcal{Z}$ we
define

$\displaystyle\ell_{U}(v,w)$ $\displaystyle\coloneqq\sup_{\hat{u}\in
U}\\{Q(\hat{w},w)+\langle
v,\rho\hat{u}\rangle\mid\hat{w}=(\hat{x},\hat{z},\hat{u}),~{}\hat{x}=x^{\star},~{}\hat{z}=z^{\star}\\}$
(19) $\displaystyle=f(x)+g(z)-p^{\star}+\sup_{\hat{u}\in
U}\langle\rho\hat{u},v-(Mx+Nz)\rangle.$

The following lemma of Ouyang et al. [4] relates the gap function to the
suboptimality and feasibility gaps. As the proof is elementary, we provide it
in Section A.2 for completeness.

######  Lemma 2.

For any $U\subseteq\mathcal{Z}$, suppose
$\ell_{U}(Mx+Nz,w)\leq\epsilon<\infty$ and $\|Mx+Nz\|\leq\delta$, where
$w=(x,z,u)\in W$. Then

$f(x)+g(z)-p^{\star}\leq\epsilon.$ (20)

In other words, $(x,z)$ is an approximate solution of (1) with suboptimality
gap $\epsilon$ and feasibility gap $\delta$. Further, if $U=\mathcal{Z}$, for
any $v$ such that $\ell_{U}(v,w)\leq\epsilon<\infty$ and $\|v\|\leq\delta$, we
have $v=Mx-z$.

Lemma 2 shows that if we can find $w$ such that
$\ell_{U}(Mx+Nz,w)\leq\epsilon$ and $\|Mx+Nz\|\leq\delta$, then we have an
approximate optimal solution to (1) with gaps $\epsilon$ and $\delta$, that
is, $\ell_{U}(Mx+Nz,w)$ controls the suboptimality and feasibility gaps.

###  5.2 Controlling the gap function

Lemma 2 shows the key to establishing convergence of GeNI-ADMM is to achieve
appropriate control over the gap function. To accomplish this, we use the
optimality conditions of the $x$ and $z$ subproblems. However, as the
subproblems are only solved approximately, the inexact solutions satisfy
_perturbed_ optimality conditions. To be able to reason about the optimality
conditions under inexact solutions, the iterates must remain bounded. Indeed,
if the iterates are unbounded, they can fail to satisfy the subproblem
optimality conditions arbitrarily badly. Fortunately, the following
proposition shows the iterates remain bounded.

######  Proposition 3 (GeNI-ADMM iterates remain bound).

Let $\\{\varepsilon_{x}^{k}\\}_{k\geq 1},\\{\varepsilon_{z}^{k}\\}_{k\geq
1},\\{\Theta^{k}\\}_{k\geq 1},\\{\Psi^{k}\\}_{k\geq 1}$, and $\rho>0$ be
given. Instate Assumptions 1-4. Run Algorithm 2, then the output sequences
$\\{\tilde{x}^{k}\\}_{k\geq 1},\\{\tilde{z}^{k}\\}_{k\geq
1},\\{\tilde{u}^{k}\\}_{k\geq 1}$ are bounded. That is, there exists $R>0$,
such that

$\\{\tilde{x}^{k}\\}_{k\geq 1}\subset
B(x^{\star},R),\quad\\{\tilde{z}^{k}\\}_{k\geq 1}\subset
B(z^{\star},R),\quad\\{\tilde{u}^{k}\\}_{k\geq 1}\subset B(u^{\star},R).$

The proof of Proposition 3 is provided in Section A.2.

As the iterates remain bounded, we can show that the optimality conditions are
approximately satisfied at each iteration. The precise form of these perturbed
optimality conditions is given in Lemmas 3 and 4. Detailed proofs establishing
these lemmas are given in Section A.2.

######  Lemma 3 (Inexact $x$-optimality condition).

Instate Assumptions 1-4. Suppose $\tilde{x}^{k+1}$ is an
$\varepsilon^{k}_{x}$-minimizer of the $x$-subproblem under Assumption 3. Then
for some absolute constant $C_{x}>0$ we have

$\displaystyle\langle\nabla f_{1}(\tilde{x}^{k+1})+\nabla
f_{2}(\tilde{x}^{k}),x^{\star}-\tilde{x}^{k+1}\rangle\leq\langle\Theta^{k}(\tilde{x}^{k+1}-\tilde{x}^{k}),x^{\star}-\tilde{x}^{k+1}\rangle$
(21) $\displaystyle+\rho\langle
M\tilde{x}^{k+1}+N\tilde{z}^{k}+\tilde{u}^{k},M(x^{\star}-\tilde{x}^{k+1})\rangle+C_{x}\varepsilon_{x}^{k}.$

######  Lemma 4 (Inexact $z$-optimality condition).

Instate Assumptions 1-4. Suppose $\tilde{z}^{k+1}$ is an
$\varepsilon^{k}_{z}$-minimum of the $z$-subproblem under Definition 1. Then
for some absolute constant $C_{z}>0$ we have

$g(z^{\star})-g(\tilde{z}^{k+1})\geq\langle-\rho
N^{T}\tilde{u}^{k+1}+\Psi^{k}(\tilde{z}^{k}-\tilde{z}^{k+1}),z^{\star}-\tilde{z}^{k+1}\rangle-
C_{z}\sqrt{\varepsilon_{z}^{k}}.$ (22)

When $\varepsilon_{x}^{k}=\varepsilon_{z}^{k}=0$ the approximate optimality
conditions of Lemmas 3 and 4 collapse to the exact optimality conditions. We
also note while Lemma 3 (Lemma 4) necessarily holds when $\tilde{x}^{k+1}$
($\tilde{z}^{k+1})$ is an $\varepsilon_{x}^{k}$-approximate minimizer
($\varepsilon_{z}^{k}$-minimum), the converse does not hold. With Lemmas 3 and
4 in hand, we can establish control of the gap function for one iteration.

######  Lemma 5.

Instate Assumptions 1-4.. Let
$\tilde{w}^{k+1}=(\tilde{x}^{k+1},\tilde{z}^{k+1},\tilde{u}^{k+1})$ denote the
iterates generated by Algorithm 2 at iteration $k$. Set
$w=(x^{\star},z^{\star},u)$, then the gap function $Q$ satisfies

$\displaystyle Q(w,\tilde{w}^{k+1})$
$\displaystyle\leq\frac{1}{2}\left(\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k}}-\|\tilde{x}^{k+1}-x^{\star}\|^{2}_{\Theta^{k}}\right)+\frac{1}{2}\left(\|\tilde{z}^{k}-z^{\star}\|_{\Psi^{k}+\rho
N^{T}N}^{2}-\|\tilde{z}^{k+1}-z^{\star}\|_{\Psi^{k}+\rho N^{T}N}^{2}\right)$
$\displaystyle+\frac{\rho}{2}\left(\|\tilde{u}^{k}-u\|^{2}-\|\tilde{u}^{k+1}-u\|^{2}\right)+C_{x}\varepsilon_{x}^{k}+C_{z}\sqrt{\varepsilon_{z}^{k}}.$

###### Proof.

From the definition of $Q$,

$\displaystyle
Q(w,\tilde{w}^{k+1})=f(\tilde{x}^{k+1})-f(x^{\star})+g(\tilde{z}^{k+1})-g(z^{\star})+\langle\rho
u,M\tilde{x}^{k+1}+N\tilde{z}^{k+1}\rangle-\langle\rho\tilde{u}^{k+1},Mx^{\star}+Nz^{\star}\rangle.$

Our goal is to upper bound $Q(w,\tilde{w}^{k+1})$. We start by bounding
$f(\tilde{x}^{k+1})-f(x^{\star})$ as follows:

$\displaystyle
f(\tilde{x}^{k+1})-f(x^{\star})=\left(f_{1}(\tilde{x}^{k+1})-f_{1}(x^{\star})\right)+\left(f_{2}(\tilde{x}^{k+1})-f_{2}(\tilde{x}^{k})\right)+\left(f_{2}(\tilde{x}^{k})-f_{2}(x^{\star})\right)$
$\displaystyle\overset{(1)}{\leq}\langle\nabla
f_{1}(\tilde{x}^{k+1}),\tilde{x}^{k+1}-\tilde{x}^{\star}\rangle+\langle\nabla
f_{1}(\tilde{x}^{k}),\tilde{x}^{k+1}-\tilde{x}^{k}\rangle+\frac{1}{2}\|\tilde{x}^{k+1}-\tilde{x}^{k}\|^{2}_{\Theta_{k}}+f_{2}(\tilde{x}^{k})-f_{2}(x^{\star})$
$\displaystyle\overset{(2)}{\leq}\langle\nabla
f_{1}(\tilde{x}^{k+1}),\tilde{x}^{k+1}-\tilde{x}^{\star}\rangle+\langle\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-\tilde{x}^{k}\rangle+\frac{1}{2}\|\tilde{x}^{k+1}-\tilde{x}^{k}\|^{2}_{\Theta_{k}}+\langle\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k}-x^{\star}\rangle$
$\displaystyle=\langle\nabla f_{1}(\tilde{x}^{k+1})+\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-x^{\star}\rangle+\frac{1}{2}\|\tilde{x}^{k+1}-\tilde{x}^{k}\|^{2}_{\Theta^{k}},$

where $(1)$ uses convexity of $f_{1}$ and $1$-relative smoothness of $f_{2}$,
and $(2)$ uses convexity of $f_{2}$. Inserting the upper bound on
$f(\tilde{x}^{k+1})-f(x)$ into the expression for $Q(w,\tilde{w}^{k+1})$, we
find

$\displaystyle Q(w,\tilde{w}^{k+1})\leq$ $\displaystyle\langle\nabla
f_{1}(\tilde{x}^{k+1})+\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-x^{\star}\rangle+\frac{1}{2}\|\tilde{x}^{k+1}-\tilde{x}^{k}\|^{2}_{\Theta^{k}}+g(\tilde{z}^{k+1})-g(z^{\star})$
(23) $\displaystyle+\langle\rho
u,M\tilde{x}^{k+1}-\tilde{z}^{k+1}\rangle-\langle\rho\tilde{u}^{k+1},Mx^{\star}-z^{\star}\rangle.$

Now, using the inexact optimality condition for the $x$-subproblem (Lemma 3),
the above display becomes

$\displaystyle
Q(w,\tilde{w}^{k+1})\leq\langle\Theta^{k}(\tilde{x}^{k+1}-\tilde{x}^{k}),x^{\star}-\tilde{x}^{k+1}\rangle+\rho\langle
M\tilde{x}^{k+1}+N\tilde{z}^{k}+\tilde{u}^{k},M(x^{\star}-\tilde{x}^{k+1})\rangle+C_{x}\varepsilon_{x}^{k}$
(24)
$\displaystyle+\frac{1}{2}\|\tilde{x}^{k+1}-\tilde{x}^{k}\|^{2}_{\Theta^{k}}+g(\tilde{z}^{k+1})-g(z^{\star})+\langle\rho
u,M\tilde{x}^{k+1}+N\tilde{z}^{k+1}\rangle-\langle\rho\tilde{u}^{k+1},Mx^{\star}+Nz^{\star}\rangle.$
(25)

Similarly, applying the inexact optimality for the $z$-subproblem (Lemma 4),
we further obtain

$\displaystyle
Q(w,\tilde{w}^{k+1})\leq\langle\Theta^{k}(\tilde{x}^{k+1}-\tilde{x}^{k}),x^{\star}-\tilde{x}^{k+1}\rangle+\rho\langle
M\tilde{x}^{k+1}+N\tilde{z}^{k}+\tilde{u}^{k},M(x^{\star}-\tilde{x}^{k+1})\rangle+C_{x}\varepsilon_{x}^{k}$
$\displaystyle+\frac{1}{2}\|\tilde{x}^{k+1}-\tilde{x}^{k}\|^{2}_{\Theta^{k}}-\rho\langle
N^{T}\tilde{u}^{k+1},\tilde{z}^{k+1}-z^{\star}\rangle+\langle\Psi^{k}(\tilde{z}^{k+1}-\tilde{z}^{k}),z^{\star}-\tilde{z}^{k+1}\rangle+C_{z}\sqrt{\varepsilon_{z}^{k}}$
(26) $\displaystyle+\langle\rho
u,M\tilde{x}^{k+1}+N\tilde{z}^{k+1}\rangle-\langle\rho\tilde{u}^{k+1},Mx^{\star}+Nz^{\star}\rangle.$

We now simplify (5.2) by combining terms. Some basic manipulations show the
terms on line 2 of (5.2) may be rewritten as

$\displaystyle-\rho\langle
N^{T}\tilde{u}^{k+1},\tilde{z}^{k+1}-z^{\star}\rangle+\rho\langle
u,M\tilde{x}^{k+1}+N\tilde{z}^{k+1}\rangle-\rho\langle\tilde{u}^{k+1},Mx^{\star}+Nz^{\star}\rangle$
$\displaystyle=\rho\langle
u-\tilde{u}^{k+1},\tilde{M}\tilde{x}^{k+1}+N\tilde{z}^{k+1}\rangle-\rho\langle\tilde{u}^{k+1},M(x^{\star}-\tilde{x}^{k+1})\rangle.$

We can combine the preceding display with the second term of line 1 in (5.2)
to reach

$\displaystyle\rho\langle
M\tilde{x}^{k+1}+N\tilde{z}^{k}+\tilde{u}^{k},M(x^{\star}-\tilde{x}^{k+1})\rangle-\rho\langle\tilde{u}^{k+1},M(x^{\star}-\tilde{x}^{k+1})\rangle+\rho\langle
u-\tilde{u}^{k+1},M\tilde{x}^{k+1}+N\tilde{z}^{k+1}\rangle$
$\displaystyle=\rho\langle\tilde{u}^{k+1}+N(\tilde{z}^{k}-\tilde{z}^{k+1}),M(x^{\star}-\tilde{x}^{k+1})\rangle-\rho\langle\tilde{u}^{k+1},M(x^{\star}-\tilde{x}^{k+1})\rangle+\rho\langle
u-\tilde{u}^{k+1},M\tilde{x}^{k+1}+N\tilde{z}^{k+1}\rangle$
$\displaystyle=\rho\langle
N(\tilde{z}^{k}-\tilde{z}^{k+1}),M(x^{\star}-\tilde{x}^{k+1})+\rho\langle
u-\tilde{u}^{k+1},M\tilde{x}^{k+1}+N\tilde{z}^{k+1}\rangle$
$\displaystyle=\rho\langle
N^{-}(\tilde{z}^{k+1}-\tilde{z}^{k}),M(x^{\star}-\tilde{x}^{k+1})\rangle+\rho\langle\tilde{u}^{k+1}-u,\tilde{u}^{k}-\tilde{u}^{k+1}\rangle,$

where $N^{-}=-N.$ Inserting the preceding simplification into (5.2), we reach

$\displaystyle Q(w,\tilde{w}^{k+1})$
$\displaystyle\leq\langle\Theta^{k}(\tilde{x}^{k+1}-\tilde{x}^{k}),x^{\star}-\tilde{x}^{k+1}\rangle+\rho\langle\tilde{u}^{k+1}-u,\tilde{u}^{k}-\tilde{u}^{k+1}\rangle$
(27)
$\displaystyle+\langle\Psi^{k}(\tilde{z}^{k+1}-\tilde{z}^{k}),z^{\star}-\tilde{z}^{k+1}\rangle+\rho\langle
N^{-}(\tilde{z}^{k+1}-\tilde{z}^{k}),M(x^{\star}-\tilde{x}^{k+1})\rangle$
$\displaystyle+\frac{1}{2}\|\tilde{x}^{k+1}-\tilde{x}^{k}\|_{\Theta^{k}}^{2}+C_{x}\varepsilon_{x}^{k}+C_{z}\sqrt{\varepsilon^{k}_{z}}.$

Now, we bound the first two leading terms in line 1 of (27) by invoking the
identity $\langle
a-b,\Upsilon(c-d)\rangle=1/2\left(\|a-d\|^{2}_{\Upsilon}-\|a-c\|^{2}_{\Upsilon}\right)+1/2\left(\|c-b\|^{2}_{\Upsilon}-\|d-b\|^{2}_{\Upsilon}\right)$
to obtain

$\displaystyle\langle\Theta^{k}(\tilde{x}^{k+1}-\tilde{x}^{k}),x^{\star}-\tilde{x}^{k+1}\rangle+\rho\langle\tilde{u}^{k+1}-u,\tilde{z}^{k+1}-M\tilde{x}^{k+1}\rangle=$
$\displaystyle\frac{1}{2}\left(\|x^{\star}-\tilde{x}^{k}\|^{2}_{\Theta^{k}}-\|x^{\star}-\tilde{x}^{k+1}\|^{2}_{\Theta^{k}}\right)-\frac{1}{2}\|\tilde{x}^{k+1}-\tilde{x}^{k}\|^{2}_{\Theta^{k}}+\frac{\rho}{2}\left(\|\tilde{u}^{k}-u\|^{2}-\|\tilde{u}^{k+1}-u\|^{2}-\|\tilde{u}^{k}-\tilde{u}^{k+1}\|^{2}\right).$

Similarly, to bound the third and fourth terms in (27), we again invoke
$(a-b)^{T}\Upsilon(c-d)=1/2\left(\|a-d\|^{2}_{\Upsilon}-\|a-c\|^{2}_{\Upsilon}\right)+1/2\left(\|c-b\|^{2}_{\Upsilon}-\|d-b\|^{2}_{\Upsilon}\right)$
which yields

$\displaystyle\langle\Psi^{k}(\tilde{z}^{k+1}-\tilde{z}^{k}),z^{\star}-\tilde{z}^{k+1}\rangle=\frac{1}{2}\left(\|\tilde{z}^{k}-z^{\star}\|^{2}_{\Psi^{k}}-\|\tilde{z}^{k+1}-z^{\star}\|^{2}_{\Psi^{k}}\right)-\frac{1}{2}\|\tilde{z}^{k+1}-\tilde{z}^{k}\|_{\Psi^{k}}^{2},$
$\displaystyle\rho\langle
N^{-}(\tilde{z}^{k}-\tilde{z}^{k+1}),M(\tilde{x}^{k+1}-x^{\star})\rangle$
$\displaystyle=\frac{\rho}{2}\left(\|N^{-}\tilde{z}^{k}-Mx^{\star}\|^{2}-\|\tilde{z}^{k+1}-Mx^{\star}\|^{2}+\|N^{-}\tilde{z}^{k+1}-M\tilde{x}^{k+1}\|^{2}-\|N^{-}\tilde{z}^{k}-M\tilde{x}^{k+1}\|^{2}\right)$
$\displaystyle=\frac{\rho}{2}\left(\|N\tilde{z}^{k}-Nz^{\star}\|^{2}-\|N\tilde{z}^{k+1}-Nz^{\star}\|^{2}-\|M\tilde{x}^{k+1}+N\tilde{z}^{k}\|^{2}\right)+\frac{\rho}{2}\|\tilde{u}^{k}-\tilde{u}^{k+1}\|^{2}.$

Putting everything together, we conclude

$\displaystyle Q(w,\tilde{w}^{k+1})$
$\displaystyle\leq\frac{1}{2}\left(\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k}}-\|\tilde{x}^{k+1}-x^{\star}\|^{2}_{\Theta^{k}}\right)+\frac{1}{2}\left(\|\tilde{z}^{k}-z^{\star}\|_{\Psi^{k}+\rho
N^{T}N}^{2}-\|\tilde{z}^{k+1}-z^{\star}\|_{\Psi^{k}+\rho N^{T}N}^{2}\right)$
$\displaystyle+\frac{\rho}{2}\left(\|\tilde{u}^{k}-u\|^{2}-\|\tilde{u}^{k+1}-u\|^{2}\right)+C_{x}\varepsilon_{x}^{k}+C_{z}\sqrt{\varepsilon_{z}^{k}}.$

as desired. ∎

###  5.3 Proof of Theorem 1

With Lemma 5 in hand, we are ready to prove Theorem 1.

###### Proof.

From Lemma 5, we have for each $k$ that

$\displaystyle Q(w,\tilde{w}^{k+1})$
$\displaystyle\leq\frac{1}{2}\left(\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k}}-\|\tilde{x}^{k+1}-x^{\star}\|^{2}_{\Theta^{k}}\right)+\frac{1}{2}\left(\|\tilde{z}^{k}-z^{\star}\|_{\Psi^{k}+\rho
N^{T}N}^{2}-\|\tilde{z}^{k+1}-z^{\star}\|_{\Psi^{k}+\rho N^{T}N}^{2}\right)$
$\displaystyle+\frac{\rho}{2}\left(\|\tilde{u}^{k}-u\|^{2}-\|\tilde{u}^{k+1}-u\|^{2}\right)+C_{x}\varepsilon_{x}^{k}+C_{z}\sqrt{\varepsilon_{z}^{k}}.$

Now, summing up the preceding display from $k=1$ to $t$ and using
$w=(x^{\star},z^{\star},u)$, we obtain

$\displaystyle\sum^{t+1}_{k=2}Q(x^{\star},z^{\star},u;\tilde{w}^{k})$
$\displaystyle\leq\frac{1}{2}\underbrace{\sum^{t}_{k=1}\left(\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k}}-\|\tilde{x}^{k+1}-x^{\star}\|^{2}_{\Theta^{k}}\right)}_{T_{1}}$
$\displaystyle+\underbrace{\frac{\rho}{2}\sum^{t}_{k=1}\left(\|\tilde{z}^{k}-z^{\star}\|_{\Psi^{k}+\rho
N^{T}N}^{2}-\|\tilde{z}^{k+1}-z^{\star}\|_{\Psi^{k}+\rho
N^{T}N}^{2}\right)}_{T_{2}}$
$\displaystyle+\underbrace{\frac{\rho}{2}\sum^{t}_{k=1}\left(\|\tilde{u}^{k}-u\|^{2}-\|\tilde{u}^{k+1}-u\|^{2}\right)}_{T_{3}}+C_{x}\mathcal{E}_{x}+C_{z}\mathcal{E}_{z}.$

We now turn to bounding $T_{1}$. Using the definition of $T_{1}$, we find

$\displaystyle T_{1}$
$\displaystyle=\frac{1}{2}\sum^{t}_{k=1}\left(\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k}}-\|\tilde{x}^{k+1}-x^{\star}\|^{2}_{\Theta^{k}}\right)$
$\displaystyle=\frac{1}{2}\left(\|\tilde{x}^{1}-x^{\star}\|^{2}_{\Theta^{1}}-\|\tilde{x}^{t+1}-x^{\star}\|^{2}_{\Theta^{t}}\right)+\frac{1}{2}\sum_{k=2}^{t}\left(\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k}}-\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k-1}}\right).$

Now using our hypotheses on the sequence $\\{\Theta^{k}\\}_{k}$ in (15), we
obtain

$\displaystyle\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k}}-\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k-1}}$
$\displaystyle=(\tilde{x}^{k}-x^{\star})^{T}(\Theta^{k}-\Theta^{k-1})(\tilde{x}^{k}-x^{\star})\leq\zeta^{k-1}\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k-1}}$
$\displaystyle\leq\tau_{\zeta}\zeta^{k-1}\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{1}}.$

Inserting the previous bound into $T_{1}$, we reach

$\displaystyle T_{1}$
$\displaystyle\leq\frac{1}{2}\|\tilde{x}^{1}-x^{\star}\|^{2}_{\Theta^{1}}+\frac{1}{2}\tau_{\zeta}\sum_{k=1}^{t}\zeta^{k-1}\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{1}}-\frac{1}{2}\|\tilde{x}^{t+1}-x^{\star}\|^{2}_{\Theta^{t}}$
$\displaystyle\leq\frac{1}{2}\left(\|\tilde{x}^{1}-x^{\star}\|^{2}_{\Theta^{1}}+C_{\zeta}\mathcal{E}_{\zeta}\right)-\frac{1}{2}\|\tilde{x}^{t+1}-x^{\star}\|^{2}_{\Theta^{t}}$
$\displaystyle\leq\frac{1}{2}\left(d^{2}_{x^{\star},\Theta^{1}}+C_{\zeta}\mathcal{E}_{\zeta}\right).$

Next, we bound $T_{2}$.

$\displaystyle T_{2}$
$\displaystyle=\frac{1}{2}\sum^{t}_{k=1}\left(\|\tilde{z}^{k}-z^{\star}\|^{2}-\|\tilde{z}^{k+1}-z^{\star}\|^{2}\right)+\frac{1}{2}\sum_{k=1}^{t}\left(\|\tilde{z}^{k}-z^{\star}\|_{\Psi^{k}+\rho
N^{T}N}^{2}-\|\tilde{z}^{k+1}-z^{\star}\|_{\Psi^{k}+\rho N^{T}N}^{2}\right)$
$\displaystyle\leq\frac{\rho}{2}\sum^{t}_{i=1}\left(\|\tilde{z}^{k}-z^{\star}\|^{2}-\|\tilde{z}^{k+1}-z^{\star}\|^{2}\right)+\frac{1}{2}\left(\|\tilde{z}^{1}-z^{\star}\|^{2}_{\Psi^{1}}+C_{\zeta}\mathcal{E}_{\zeta}\right)$
$\displaystyle=\frac{1}{2}\|\tilde{z}^{1}-z^{\star}\|_{\Psi^{1}+\rho
I}^{2}=\frac{1}{2}d_{z^{\star},\Psi^{1}_{\rho,N}}^{2}+\frac{C_{\zeta}}{2}\mathcal{E}_{\zeta}.$

Last, $T_{3}$ is a telescoping sum, hence

$\displaystyle T_{2}$
$\displaystyle=\frac{\rho}{2}\sum^{t}_{k=1}\left(\|\tilde{u}^{k}-u\|^{2}-\|\tilde{u}^{k+1}-u\|^{2}\right)=\frac{\rho}{2}\left(\|\tilde{u}^{1}-u\|^{2}-\|\tilde{u}^{t+1}-u\|^{2}\right).$

Using our bounds on $T_{1}$ through $T_{3}$, we find

$\displaystyle\sum^{t+1}_{k=2}Q(x^{\star},z^{\star},u;\tilde{w}^{k})$
$\displaystyle\leq\frac{1}{2}\left(d^{2}_{x^{\star},\Theta^{1}}+C_{\zeta}\mathcal{E}_{\zeta}\right)+\frac{1}{2}d_{z^{\star},\Psi^{1}_{\rho,N}}^{2}+\frac{1}{2}C_{\zeta}\mathcal{E}_{\zeta}+C_{x}\mathcal{E}_{x}+C_{z}\mathcal{E}_{z}$
$\displaystyle+\frac{\rho}{2}\left(\|\tilde{u}^{1}-u\|^{2}-\|\tilde{u}^{t+1}-u\|^{2}\right)$
$\displaystyle=\frac{1}{2}d^{2}_{x^{\star},\Theta^{1}}+\frac{1}{2}d_{z^{\star},\Psi^{1}_{\rho,N}}^{2}+C_{x}\mathcal{E}_{x}+C_{z}\mathcal{E}_{z}+C_{\zeta}\mathcal{E}_{\zeta}+\frac{\rho}{2}\left(\|\tilde{u}^{1}-u\|^{2}-\|\tilde{u}^{t+1}-u\|^{2}\right),$

Now, as $\bar{w}^{t+1}=\frac{1}{t}\sum^{t+1}_{k=2}\tilde{w}^{k}$, the
convexity of $Q$ in its second argument yields

$\displaystyle Q(x^{\star},z^{\star},u;\bar{w}^{t+1})$
$\displaystyle\leq\frac{1}{t}\sum_{k=2}^{t+1}Q(x^{\star},z^{\star},u;\tilde{w}^{k})$
(28)
$\displaystyle\leq\frac{1}{t}\left(\frac{1}{2}d^{2}_{x^{\star},\Theta^{1}}+\frac{\rho}{2}d_{z^{\star},\Psi_{\rho,N}^{1}}^{2}+C_{x}\mathcal{E}_{x}+C_{z}\mathcal{E}_{z}+C_{\zeta}\mathcal{E}_{\zeta}+\frac{\rho}{2}\left(\|\tilde{u}^{1}-u\|^{2}-\|\tilde{u}^{t+1}-u\|^{2}\right)\right).$

Define
$\Gamma\coloneqq\frac{1}{2}d^{2}_{x^{\star},\Theta^{1}}+\frac{1}{2}d^{2}_{z^{\star},\Psi^{1}_{\rho,N}}+C_{x}\mathcal{E}_{x}+C_{z}\mathcal{E}_{z}+C_{\zeta}\mathcal{E}_{\zeta}$.
Since $Q(w^{\star},\bar{w}^{t+1})\geq 0$, by (28) we reach

$\|\tilde{u}^{t+1}-u^{\star}\|^{2}\leq\frac{2}{\rho}\Gamma+d^{2}_{u^{\star}}.$

Let $\tilde{v}^{t+1}=\frac{1}{t}\left(\tilde{u}^{1}-\tilde{u}^{t+1}\right)$.
Then we can bound $\|\tilde{v}_{t+1}\|^{2}$ as

$\displaystyle\|\tilde{v}^{t+1}\|^{2}\leq$
$\displaystyle\frac{2}{t^{2}}\left(\|\tilde{u}^{1}-u^{\star}\|^{2}+\|\tilde{u}^{t+1}-u^{\star}\|^{2}\right)\leq\frac{4}{t^{2}}\left(\frac{\Gamma}{\rho}+d^{2}_{u^{\star}}\right).$

By (28), given the fact $\tilde{u}^{1}=0$, we also have

$\displaystyle Q(x^{\star},z^{\star},u;\bar{w}^{t+1})$
$\displaystyle\leq\frac{\Gamma-\rho\langle\tilde{u}^{1}-\tilde{u}^{t+1},u\rangle}{t}=\frac{\Gamma}{t}-\rho\langle\tilde{v}^{t+1},u\rangle,$

where the equality follows from the definition of $\tilde{v}^{t+1}$. Hence for
any $u$

$\displaystyle
Q(x^{\star},z^{\star},u;\bar{w}^{t+1})+\langle\tilde{v}^{t+1},\rho u\rangle$
$\displaystyle\leq\frac{\Gamma}{t},$

and therefore

$\ell_{U}(\tilde{v}^{t+1},\bar{w}^{t+1})\leq\frac{1}{t}\left(\frac{1}{2}d^{2}_{x^{\star},\Theta^{1}}+\frac{1}{2}d_{z^{\star},\Psi^{1}_{\rho,N}}^{2}+C_{x}\mathcal{E}_{x}+C_{z}\mathcal{E}_{z}+C_{\zeta}\mathcal{E}_{\zeta}\right).$

We finish the proof by invoking Lemma 2. $\Box$

##  6 Linear convergence of GeNI-ADMM

In this section, we seek to establish linear convergence results for Algorithm
2. In general, the linear convergence of ADMM relies on strong convexity of
the objective function [1, 33, 11]. Consistently, the linear convergence of
GeNI-ADMM also requires strong convexity. Many applications of GeNI-ADMM fit
into this setting, such as elastic net [34]. However, linear convergence is
not restricted to strongly convex problems. It has been shown that local
linear convergence of ADMM can be guaranteed even without strong convexity
[21]. Experiments in Section 8 show the same phenomenon for GeNI-ADMM: it
converges linearly after a couple of iterations when the iterates reach some
manifold containing the solution. The linear convergence theory of GeNI-ADMM
provides a way to understand this phenomenon. We first list the additional
assumptions required for linear convergence:

######  Assumption 5 (Optimization is over the whole space).

The sets $X$ and $Z$ in (1) satisfy

$X=\mathcal{X},~{}\text{and}~{}Z=\mathcal{Z}.$

Assumption 5 states that the optimization problem in (1) is over the entire
spaces $\mathcal{X}$ and $\mathcal{Z}$, not closed subsets. This assumption is
met in many practical optimization problems of interest, where
$\mathcal{X}=\mathcal{Z}=\mathcal{H}=\mathbb{R}^{d}$. Moreover, it is
consistent with prior analyses such as Deng and Yin [13], who specialize their
analysis to the setting of the last sentence.

######  Assumption 6 (Regularity of $f$).

The function $f$ is finite valued, strongly convex with parameter
$\sigma_{f}$, and smooth with parameter $L_{f}$.

Assumption 6 imposes standard regularity conditions on $f$, in addition to the
conditions of Assumption 2.

######  Assumption 7 (Non-degeneracy of constraint operators).

The linear operators $MM^{T}$ and $N^{T}N$ are invertible.

Assumption 7 is consistent with prior analyses of ADMM-type schemes under
strong convexity, such as Deng and Yin [13], who make this assumption in their
analysis of Generalized ADMM. Moreover, the assumption holds in many important
application problems, especially those that arise from machine learning where,
typically, $M=I$ and $N=-I$.

######  Assumption 8 (Geometric decay of the forcing sequences).

There exists a constant $q>0$ such that the forcing sequences
$\\{\varepsilon^{k}_{x}\\}_{k\geq 1}$, $\\{\varepsilon^{k}_{z}\\}_{k\geq 1}$
satisfy

$\varepsilon^{k+1}_{x}\leq\varepsilon^{k}_{x}/(1+q),\;\text{and}\;\;\varepsilon^{k+1}_{z}\leq\varepsilon^{k}_{z}/(1+q)^{2}.$
(29)

Moreover, we assume Algorithm 2 is equipped with oracles for solving the $x$
and $z$-subproblems, which at each iteration produce approximate solutions
$\tilde{x}^{k+1}$, $\tilde{z}^{k+1}$ satisfying:

$\displaystyle\tilde{x}^{k+1}\overset{\varepsilon_{x}^{k}}{\approx}\underset{x\in
X}{\mathop{\rm argmin}}\\{f_{1}(x)+\langle\nabla
f_{2}(\tilde{x}^{k}),x-\tilde{x}^{k}\rangle+\frac{1}{2}\|x-\tilde{x}^{k}\|^{2}_{\Theta^{k}}+\frac{\rho}{2}\|Mx+N\tilde{z}^{k}+\tilde{u}^{k}\|^{2}\\},$
$\displaystyle\tilde{z}^{k+1}\overset{\varepsilon_{z}^{k}}{\approxeq}\underset{z\in
Z}{\mathop{\rm
argmin}}\\{g(z)+\frac{\rho}{2}\|M\tilde{x}^{k+1}+Nz+\tilde{u}^{k}\|^{2}+\frac{1}{2}\|z-\tilde{z}^{k}\|^{2}_{\Psi^{k}}\\}.$

Assumption 8 replaces Assumption 3 and requires the inexactness sequences to
decay geometrically. Compared with the sublinear convergence result, linear
convergence requires more accurate solutions to the subproblems. Again, since
the $z$-subproblem inexactness is weaker than the $x$-subproblem inexactness,
$\\{\varepsilon^{k}_{z}\\}_{k\geq 1}$ should have a faster decay rate
$(1+q)^{2}$ than the decay rate ($1+q$) of $\\{\varepsilon^{k}_{x}\\}_{k\geq
1}$.

The requirement that the forcing sequences decay geometrically is somewhat
burdensome, as it leads to the subproblems needing to be solved to higher
accuracy sooner than if the forcing sequences were only summable. Fortunately,
this condition seems to be an artifact of the analysis; our numerical
experiments with strongly convex $f$ (Section 8.2) only use summable forcing
sequences but show linear convergence of GeNI-ADMM.

###  6.1 Our approach

Inspired by Deng and Yin [13], we take a Lyapunov function approach to proving
linear convergence. Let $\tilde{w}=(\tilde{x},\tilde{z},\tilde{u})$, and
$w^{\star}=(x^{\star},z^{\star},u^{\star})$. We define the Lyapunov function:

$\Phi^{k}=\frac{1}{\rho}\|\tilde{x}^{k}-x^{\star}\|_{\Theta^{k}}^{2}+\frac{1}{\rho}\|\tilde{z}^{k}-\tilde{z}^{\star}\|_{\Psi^{k}+\rho
N^{T}N}^{2}+\|\tilde{u}^{k}-u^{\star}\|^{2}=\|\tilde{w}^{k}-w^{\star}\|_{G^{k}}^{2},$

where

$G^{k}\coloneqq\begin{pmatrix}\frac{1}{\rho}\Theta^{k}&0&0\\\
0&\frac{1}{\rho}\Psi^{k}+N^{T}N&0\\\ 0&0&I\end{pmatrix}.$

Our main result in this section is the following theorem, which shows the
Lyapunov function converges linearly to $0$.

######  Theorem 2.

Instate Assumptions 1-2, and Assumptions 4-8. Moreover, suppose that
$\theta_{\textup{min}}$ in (14) satisfies
$\theta_{\textup{min}}>L^{2}_{f}/\sigma_{f}$. Then there exist constants
$\delta$ and $S>0$ such that if $q>\delta$,

$(1+\delta)^{k}\Phi^{k}\leq\tau_{\zeta}\Phi^{1}+S.$ (30)

Hence after
$k=\mathcal{O}\left(\frac{1}{\delta}\log\left(\frac{\tau_{\zeta}\Phi^{1}}{\epsilon}\right)\right)$
iterations,

$\Phi^{k}\leq\epsilon.$

The proof of Theorem 2 is deferred to Section 6.3. As $G^{k}\succ 0$ for all
$k$, Theorem 2 implies the iterates
$(\tilde{x}^{k},\tilde{z}^{k},\tilde{u}^{k})$ converge linearly to optimum
$(x^{\star},z^{\star},u^{\star})$. Thus, despite inexactly solving
approximations of the original ADMM subproblems, GeNI-ADMM still enjoys linear
convergence when the objective is strongly convex. An unattractive aspect of
Theorem 2 is the requirement that $\theta_{\textup{min}}$ satisfy
$\theta_{\textup{min}}>L_{f}^{2}/\sigma_{f}$. Although this condition can be
enforced by adding a damping term $\sigma I$ to $\Theta^{k}$, this is
undesirable as a large regularization term can lead GeNI-ADMM to converge
slower. We have found that this condition is unnecessary—our numerical
experiments do not enforce this condition, yet GeNI-ADMM achieves linear
convergence. We believe this condition to be an artifact of the analysis,
which stems from lower bounding the term $-\langle\nabla
f_{2}(\tilde{x}^{k+1})-\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-x^{\star}\rangle$ in Lemma 7.

###  6.2 Sufficient descent

From Theorem 2, to establish linear convergence of GeNI-ADMM, it suffices to
show $\Phi^{k}$ decreases geometrically. We take two steps to achieve this.
First, we show that
$\|\tilde{w}^{k}-w^{\star}\|_{G^{k}}^{2}-\|\tilde{w}^{k+1}-w^{\star}\|_{G^{k}}^{2}$
decreases for every iteration $k$ (Lemma 7). Second, we show that $\Phi^{k+1}$
decreases geometrically by a factor of $1/(1+\delta)$ with respect to
$\Phi^{k}$ and some small error terms that stem from inexactness (Lemma 8).

As in the convex case, the optimality conditions of the subproblems play a
vital role in the analysis. Since the subproblems are only solved
approximately, we must again consider the inexactness of the solutions in
these two steps. For the first step, we use strong convexity of $f$ and
convexity of $g$ with appropriate perturbations to account for the
inexactness. We call these conditions _perturbed convexity_ conditions, as
outlined in Lemma 6.

######  Lemma 6 (Perturbed convexity).

Instate the assumptions of Theorem 2. Let $\tilde{x}^{k+1}$ and
$\tilde{z}^{k+1}$ be the inexact solutions of $x$ and $z$-subproblems under
Definition 1. Recall $(x^{\star},z^{\star},u^{\star})$ is a saddle point of
(1). Then for some constant $C\geq 0$, the following inequalities are
satisfied:

  1. 1.

(Semi-inexact $f$-strong convexity)

$\displaystyle\langle\tilde{x}^{k}-\tilde{x}^{k+1},\tilde{x}^{k+1}-x^{\star}\rangle_{\Theta^{k}}+\rho\langle
N(\tilde{z}^{k+1}-\tilde{z}^{k})+u^{\star}-\tilde{u}^{k+1},M(\tilde{x}^{k+1}-x^{\star})\rangle+C\varepsilon_{x}^{k}$
(31)
$\displaystyle\geq\sigma_{f}\|\tilde{x}^{k+1}-x^{\star}\|^{2}-\langle\nabla
f_{2}(\tilde{x}^{k+1})-\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-x^{\star}\rangle,$

  2. 2.

(Semi-inexact $g$-convexity)

$\langle\tilde{z}^{k+1}-z^{\star},\rho
N^{T}(u^{\star}-\tilde{u}^{k+1})+\Psi^{k}(\tilde{z}^{k}-\tilde{z}^{k+1})\rangle\geq-C\sqrt{\varepsilon_{z}^{k}}.$
(32)

A detailed proof of Lemma 6 is presented in Section A.3. In Lemma 6, we call
(31) and (32) _semi-inexact_ (strong) convexity because $x^{\star}$
($z^{\star}$) is part of the exact saddle point but $\tilde{x}^{k+1}$
($\tilde{z}^{k+1}$) is an inexact subproblem solution. With Lemma 6, we
establish a descent-type inequality, which takes into inexactness.

######  Lemma 7 (Inexact sufficient descent).

Define $\varepsilon_{w}^{k}=\varepsilon_{x}^{k}+\sqrt{\varepsilon_{z}^{k}}$,
then the following descent condition holds.

$\displaystyle\|\tilde{w}^{k}-w^{\star}\|_{G^{k}}^{2}-\|\tilde{w}^{k+1}-w^{\star}\|^{2}_{G^{k}}+C\varepsilon_{w}^{k}\geq\frac{1}{2}\|\tilde{w}^{k}-\tilde{w}^{k+1}\|_{G^{k}}^{2}+\frac{2\sigma_{f}}{\rho}\|\tilde{x}^{k+1}-x^{\star}\|^{2}$
(33) $\displaystyle-\frac{2}{\rho}\langle\nabla f_{2}(\tilde{x}^{k+1})-\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-x^{\star}\rangle.$

###### Proof.

Adding the inequalities (31) and (32) together, and using the relation
$M(\tilde{x}^{k+1}-x^{\star})=\tilde{u}^{k+1}-\tilde{u}^{k}+N(z^{\star}-\tilde{z}^{k+1})$,
we reach

$\displaystyle\langle\tilde{x}^{k}-\tilde{x}^{k+1},\tilde{x}^{k+1}-x^{\star}\rangle_{\Theta^{k}}+\langle\tilde{z}^{k}-\tilde{z}^{k+1},\tilde{z}^{k+1}-z^{\star}\rangle_{\Psi^{k}+\rho
N^{T}N}+\rho\langle\tilde{u}^{k}-\tilde{u}^{k+1},\tilde{u}^{k+1}-u^{\star}\rangle+C\varepsilon_{x}^{k}$
$\displaystyle\geq\sigma_{f}\|\tilde{x}^{k+1}-x^{\star}\|^{2}-\langle\nabla
f_{2}(\tilde{x}^{k+1})-\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-x^{\star}\rangle-C\sqrt{\varepsilon_{z}^{k}}+\rho\langle\tilde{u}^{k}-\tilde{u}^{k+1},N(\tilde{z}^{k+1}-\tilde{z}^{k})\rangle$

Recalling the definitions of $\tilde{w}$, $\varepsilon_{w}^{k}$, and $G^{k}$,
and using the identity $\langle
a-b,\Upsilon(c-d)\rangle=1/2\left(\|a-d\|^{2}_{\Upsilon}-\|a-c\|^{2}_{\Upsilon}\right)+1/2\left(\|c-b\|^{2}_{\Upsilon}-\|d-b\|^{2}_{\Upsilon}\right)$,
we arrive at

$\displaystyle\|\tilde{w}^{k}-w^{\star}\|_{G^{k}}^{2}-\|\tilde{w}^{k+1}-w^{\star}\|_{G^{k}}^{2}+\varepsilon_{w}^{k}\geq\frac{2\sigma_{f}}{\rho}\|\tilde{x}^{k+1}-x^{\star}\|^{2}+\|\tilde{w}^{k+1}-\tilde{w}^{k}\|^{2}_{G^{k}}$
$\displaystyle+\langle\tilde{u}^{k}-\tilde{u}^{k+1},N(\tilde{z}^{k+1}-\tilde{z}^{k})\rangle-\frac{2}{\rho}\langle\nabla
f_{2}(\tilde{x}^{k+1})-\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-x^{\star}\rangle.$

Now, for the term
$\|\tilde{w}^{k+1}-\tilde{w}^{k}\|^{2}_{G^{k}}+\langle\tilde{u}^{k}-\tilde{u}^{k+1},N(\tilde{z}^{k+1}-\tilde{z}^{k})\rangle$,
Cauchy-Schwarz implies

$\displaystyle\|\tilde{w}^{k+1}-\tilde{w}^{k}\|_{G^{k}}^{2}+\langle\tilde{u}^{k}-\tilde{u}^{k+1},N(\tilde{z}^{k+1}-\tilde{z}^{k})\rangle$
$\displaystyle\geq\|\tilde{w}^{k+1}-\tilde{w}^{k}\|_{G^{k}}^{2}-\frac{1}{2}\|\tilde{u}^{k+1}-\tilde{u}^{k}\|^{2}-\frac{1}{2}\|\tilde{z}^{k+1}-\tilde{z}^{k}\|^{2}_{N^{T}N}$
$\displaystyle\geq\|\tilde{w}^{k+1}-\tilde{w}^{k}\|_{G^{k}}^{2}-\frac{1}{2}\|\tilde{u}^{k+1}-\tilde{u}^{k}\|^{2}-\frac{1}{2}\|\tilde{z}^{k+1}-\tilde{z}^{k}\|^{2}_{1/\rho\Psi^{k}+N^{T}N}$
$\displaystyle=\frac{1}{2}\|\tilde{w}^{k+1}-\tilde{w}^{k}\|_{G^{k}}^{2}.$

Hence we obtain

$\displaystyle\|\tilde{w}^{k}-w^{\star}\|_{G^{k}}^{2}-\|\tilde{w}^{k+1}-w^{\star}\|_{G^{k}}^{2}+C\varepsilon_{w}^{k}\geq\frac{2\sigma_{f}}{\rho}\|\tilde{x}^{k+1}-x^{\star}\|^{2}+\frac{1}{2}\|\tilde{w}^{k+1}-\tilde{w}^{k}\|^{2}_{G^{k}}$
$\displaystyle-\frac{2}{\rho}\langle\nabla f_{2}(\tilde{x}^{k+1})-\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-x^{\star}\rangle,$

as desired. $\Box$ ∎

Given the inexact sufficient descent condition (33), the next step in proving
linear convergence is to show (33) leads to a contraction relation between
$\Phi^{k+1}$ and $\Phi^{k}$.

######  Lemma 8 (Inexact Contraction Lemma).

Under the assumptions of Theorem 2, there exists constants $\delta>0$, and
$C\geq 0$ such that

$(1+\delta)\Phi^{k+1}\leq(1+\zeta^{k})\left(\Phi^{k}+C\varepsilon_{w}^{k}\right).$
(34)

The proof of Lemma 8, and an explicit expression for $\delta$, appears in
Section A.3. As in Theorem 2, the constant $\delta$ gives the rate of linear
convergence and depends on the conditioning of $f$ and the constraint matrices
$M$ and $N$. The better the conditioning, the faster the convergence, with the
opposite holding true as the conditioning worsens. With Lemma 8 in hand, we
now prove Theorem 2.

###  6.3 Proof of Theorem 2

###### Proof.

By induction on (34), we have

$\displaystyle(1+\delta)^{k}\Phi^{k}$
$\displaystyle\leq\left(\prod^{k}_{j=1}(1+\zeta^{j})\right)\Phi^{1}+C\sum^{k}_{j=1}\left(\prod^{k}_{i=j}(1+\zeta^{i})\right)(1+\delta)^{j-1}\varepsilon^{j}_{w}$
$\displaystyle\leq\tau_{\zeta}\Phi^{1}+C\tau_{\zeta}\sum_{j=1}^{k}(1+\delta)^{j-1}\varepsilon^{j}_{w}\leq\tau_{\zeta}\Phi^{1}+C\tau_{\zeta}\sum_{j=1}^{k}\left(\frac{1+\delta}{1+q}\right)^{j-1}$
$\displaystyle\leq\tau_{\zeta}\Phi^{1}+\frac{C\tau_{\zeta}}{1-\frac{1+\delta}{1+q}}=\tau_{\zeta}\Phi^{1}+S,$

where the second inequality uses Assumption 8 to reach
$\varepsilon^{j}_{w}\leq\frac{C}{(1+q)^{j}}$, and the third inequality uses
$q>\delta$ to bound the sum by the sum of the geometric series. Hence, we have
shown the first claim. The second claim follows immediately from the first via
a routine calculation. $\Box$

##  7 Applications

This section applies our theory to establish convergence rates for NysADMM and
sketch-and-solve ADMM.

###  7.1 Convergence of NysADMM

Algorithm 3 NysADMM

penalty parameter $\rho$, step-size $\eta$, regularization $\sigma\geq 0$,
forcing sequences $\\{\varepsilon_{x}^{k}\\}_{k\geq 1}$,
$\\{\varepsilon_{z}^{k}\\}_{k\geq 1}$

repeat

Find $\varepsilon_{x}^{k}$-approximate solution $\tilde{x}^{k+1}$ of  $(\eta
H_{f}^{k}+(\rho+\eta\sigma)I)x=\eta(H^{k}_{f}+\sigma I)\tilde{x}^{k}-\nabla
f(\tilde{x}^{k})+\rho(\tilde{z}^{k}-\tilde{u}^{k})$

$\tilde{z}^{k+1}\overset{\varepsilon_{z}^{k}}{\approxeq}\underset{z\in
Z}{\mathop{\rm
argmin}}\\{g(z)+\frac{\rho}{2}\|M\tilde{x}^{k+1}-z+\tilde{u}^{k}\|^{2}\\}$

$\tilde{u}^{k+1}=\tilde{u}^{k}+M\tilde{x}^{k+1}-\tilde{z}^{k+1}$

until convergence

solution $(x^{\star},z^{\star})$ of problem (1)

We begin with the NysADMM scheme from [3]. Recall NysADMM is obtained from
Algorithm 2 by setting $f_{1}=0,f_{2}=f$, using the exact Hessian
$\Theta^{k}=H_{f}(\tilde{x}^{k})=\eta(H_{f}^{k}+\sigma I)$, and setting
$\Psi^{k}=0$. Instantiating these selections into Algorithm 2, we obtain
NysADMM, presented as Algorithm 3. Compared to the original NysADMM, Algorithm
3 adds a regularization term $\sigma I$ to the Hessian. In theory, when $f$ is
only convex, this regularization term is required to ensure the condition
$\Theta^{k}\succeq\theta_{\textup{min}}$ of Assumption 2, as the Hessian along
the optimization path may fail to be uniformly bounded below. The addition of
the $\sigma I$ term removes this issue. However, the need for this term seems
to be an artifact of the proof, as Zhao et al. [3] runs Algorithm 3 on non-
strongly convex objectives with $\sigma=0$, and convergence is still obtained.
Similarly, we set $\sigma=0$ in all our experiments, and convergence
consistent with our theory is obtained. Hence, in practice, we recommend
setting $\sigma=0$, or some small value, say $\sigma=10^{-8}$, so that
convergence isn’t slowed by unneeded regularization.

The $x$-subproblem for NysADMM is a linear system. NysADMM solves this system
using the Nyström PCG method from [7]. The general convergence of NysADMM was
left open in [3], which established convergence only for quadratic $f$.
Moreover, the result in [3] does not provide an explicit convergence rate. We
shall now rectify this state of affairs using Theorem 1. We obtain the
following convergence guarantee by substituting the parameters defining
NysADMM (with the added $\sigma I$ term) into Theorem 1.

######  Corollary 1 (Convergence of NysADMM).

Instate the assumptions of Theorem 1. Let $\sigma>0$. Set $f_{1}=0,f_{2}=f$,
$\eta=\hat{L}_{f}$, $\Theta^{k}=\eta(H^{k}_{f}+\sigma I)$, and $\Psi^{k}=0$ in
Algorithm 2. Then

$f(\bar{x}^{t+1})+g(\bar{z}^{t+1})-p^{\star}\leq\frac{\Gamma}{t}$

and

$\|M\bar{x}^{t+1}+N\bar{z}^{t+1}\|\leq\frac{2}{t}\sqrt{\frac{\Gamma}{\rho}+d^{2}_{u^{\star}}}.$

Here $\Gamma$ and $d_{u^{\star}}^{2}$ are the same as in Theorem 1.

NysADMM converges at the same $\mathcal{O}(1/t)$-rate as standard ADMM,
despite all the approximations it makes. Thus, NysADMM offers the same level
of performance as ADMM, but is much faster due to its use of inexactness. This
result is empirically verified in Section 8, where NysADMM converges almost
identically to ADMM. Corollary 1 supports the empirical choice of a constant
step-size $\eta=1$, which was shown to have excellent performance uniformly
across tasks in [3]: the theorem sets $\eta=\hat{L}_{f}$ and $\hat{L}_{f}=1$
for quadratic functions, and satisfies $\hat{L}_{f}=\mathcal{O}(1)$ for loss
functions such as the logistic loss. We recommend setting $\eta=1$ as the
default value for GeNI-ADMM. Given NysADMM’s superb empirical performance in
[3] and the firm theoretical grounding given by Corollary 1, we conclude that
NysADMM provides a reliable framework for solving large-scale machine learning
problems.

###  7.2 Convergence of sketch-and-solve ADMM

Algorithm 4 Sketch-and-solve ADMM

penalty parameter $\rho$, step-size $\eta$, $\\{\varepsilon_{z}^{k}\\}_{k\geq
1}$

repeat

Find solution $\tilde{x}^{k+1}$ of
$(\eta\hat{H}^{k}+(\rho+\eta\gamma^{k})I)x=\eta\left(\hat{H}^{k}+\gamma^{k}I\right)\tilde{x}^{k}-\nabla
f(\tilde{x}^{k})+\rho(\tilde{z}^{k}-\tilde{u}^{k})$

$\tilde{z}^{k+1}\overset{\varepsilon_{z}^{k}}{\approxeq}\underset{z\in
Z}{\mathop{\rm
argmin}}\\{g(z)+\frac{\rho}{2}\|M\tilde{x}^{k+1}-z+\tilde{u}^{k}\|^{2}\\}$

$\tilde{u}^{k+1}=\tilde{u}^{k}+M\tilde{x}^{k+1}-\tilde{z}^{k+1}$

until convergence

solution $(x^{\star},z^{\star})$ of problem (1)

Sketch-and-solve ADMM is obtained from GeNI-ADMM by setting $\Theta^{k}$ to be
an approximate Hessian computed by a sketching procedure. The two most popular
sketch-and-solve methods are the Newton sketch [35, 36, 37] and Nyström
sketch-and-solve [38, 39, 7]. Both methods use an approximate Hessian that is
cheap to compute and yields a linear system that is fast to solve. ADMM
together with sketching techniques, provides a compelling means to handle
massive problem instances. The recent survey [9] suggests using sketching to
solve a quadratic ADMM subproblem for efficiently solving large-scale inverse
problems. However, it was previously unknown whether the resulting method
would converge. Here, we formally describe the sketch-and-solve ADMM method
and prove convergence.

Sketch-and-solve ADMM is obtained from Algorithm 2 by setting

$\Theta^{k}=\eta\left(\hat{H}^{k}+\gamma^{k}I\right),$ (35)

where $\hat{H}^{k}$ is an approximation to the Hessian $H_{f}(\tilde{x}^{k})$
at the $k$th iteration, and $\gamma^{k}\geq 0$ is a constant chosen to ensure
convergence. The term $\gamma^{k}I$ ensures that the approximate linearization
satisfies the $\Theta$-relative smoothness condition when $\gamma^{k}$ is
chosen appropriately, as in the following lemma:

######  Lemma 9.

Suppose $f$ is $\hat{L}_{f}$-relatively smooth with respect to its Hessian
$H_{f}$. Construct $\\{\Theta^{k}\\}_{k\geq 1}$ as in (35) and select
$\gamma^{k}>0$ such that
$\gamma^{k}\geq\|E^{k}\|=\|H_{f}(\tilde{x}^{k})-\hat{H}^{k}\|$ for every $k$.
Then

$f(x)\leq f(\tilde{x}^{k})+\langle\nabla
f(\tilde{x}^{k}),x-\tilde{x}^{k}\rangle+\frac{\hat{L}_{f}}{2}\|x-\tilde{x}^{k}\|^{2}_{\Theta^{k}}.$

###### Proof.

By $\hat{L}_{f}$ relative smoothness we have,

$f(x)\leq f(\tilde{x}^{k})+\langle\nabla
f(\tilde{x}^{k}),x-\tilde{x}^{k}\rangle+\frac{\hat{L}_{f}}{2}\|x-\tilde{x}^{k}\|^{2}_{H_{f}(\tilde{x}^{k})}.$

Now,

$\displaystyle\|x-\tilde{x}^{k}\|^{2}_{H_{f}(\tilde{x}^{k})}$
$\displaystyle=\|x-\tilde{x}^{k}\|^{2}_{\hat{H}^{k}}+(x-\tilde{x}^{k})^{T}(H_{f}(\tilde{x}^{k})-\hat{H}^{k})(x-\tilde{x}^{k})$
$\displaystyle\leq\|x-\tilde{x}^{k}\|^{2}_{\hat{H}^{k}}+\|E^{k}\|\|x-\tilde{x}^{k}\|^{2}\leq\|x-\tilde{x}^{k}\|^{2}_{\hat{H}^{k}}+\gamma^{k}\|x-\tilde{x}^{k}\|^{2}$
$\displaystyle=\|x-\tilde{x}^{k}\|^{2}_{\Theta^{k}}.$

The desired claim now immediately follows. $\Box$

Lemma 9 shows we can ensure relative smoothness by selecting $\gamma^{k}>0$
appropriately. Assuming relative smoothness, we may invoke Theorem 1 to
guarantee that sketch-and-solve ADMM converges. Unlike with NysADMM, we find
it is necessary to select $\gamma^{k}$ carefully (such as in Lemma 9) to
ensure the relative smoothness condition holds, otherwise sketch-and-solve
ADMM will diverge, see Section 8.3 for numerical demonstration. This condition
is somewhat different from those in prior sketching schemes in optimization,
such as the Newton Sketch [35, 36], where convergence is guaranteed as long as
the Hessian approximation is invertible.

######  Corollary 2.

Suppose $f$ is $\hat{L}_{f}$-relatively smooth with respect to its Hessian
$H_{f}$ and instate the assumptions of Theorem 1. In Algorithm 2, Set
$\Theta^{k}=\eta(\hat{H}^{k}+\gamma^{k}I)$ with $\eta=\hat{L}_{f}$,
$\gamma^{k}=\|E^{k}\|=\|\hat{H}^{k}-H_{f}(\tilde{x}^{k})\|$, and $\Psi^{k}=0$.
Then

$f(\bar{x}^{t+1})+g(\bar{z}^{t+1})-p^{\star}\leq\frac{1}{t}\left(\hat{L}_{f}\left(d^{2}_{x^{\star},H_{f}(\tilde{x}^{1})}+\|E^{1}\|d^{2}_{x^{\star}}\right)+\frac{\rho}{2}d^{2}_{z^{\star}}+C_{x}\mathcal{E}_{x}+C_{z}\mathcal{E}_{z}+C_{\zeta}\mathcal{E}_{\zeta}\right)=:\frac{1}{t}\hat{\Gamma}$

and

$\|M\bar{x}^{t+1}-\bar{z}^{t+1}\|\leq\frac{2}{t}\sqrt{\frac{\hat{\Gamma}}{\rho}+d^{2}_{u^{\star}}}.$

Here variables $\bar{x}^{t+1}$ and $\bar{z}^{t+1}$, diameters
$d^{2}_{x^{\star},H_{f}(\tilde{x}^{1})}$, $d^{2}_{z^{\star}}$, and
$d^{2}_{u^{\star}}$, and constants $C_{x}$, $C_{z}$, $\mathcal{E}_{x}$,
$\mathcal{E}_{z}$, $\mathcal{E}_{\zeta}$, $R_{\star}$, and $p^{\star}$ are all
defined as in Theorem 1.

######  Remark 2.

It is not hard to see Corollary 2 also holds for any $\gamma>0$ if the
$\\{\Theta^{k}\\}_{k\geq 1}$ satisfy

$(1-\tau)(H_{f}^{k}+\gamma I)\preceq\Theta^{k}\preceq(1+\tau)(H_{f}^{k}+\gamma
I),$

for some $\tau\in(0,1)$, and $\eta\geq\hat{L}_{f}/(1-\tau)$. Many sketching
methods can ensure this relative error condition with high probability [40] by
choosing the sketch size proportional to the effective dimension of
$H_{f}+\gamma I$ [7, 36].

Corollary 2 shows sketch-and-solve ADMM obtains an
$\mathcal{O}(1/t)$-convergence rate. The main difference between NysADMM and
sketch-and-solve is the additional error term due to the use of the
approximate Hessian, which results in a slightly slower convergence rate. In
this sense, sketch-and-solve ADMM can be regarded as a compromise between
NysADMM ($\Theta^{k}=\eta(H_{f}(\tilde{x}^{k})+\sigma I)$) and gradient
descent ADMM ($\Theta^{k}=\eta I$) — with the convergence rate improving as
the accuracy of the Hessian approximation increases.

##  8 Numerical experiments

In this section, we numerically illustrate the convergence results developed
in Section 5 and Section 6 for several methods highlighted in Section 3.1.1
that fit into the GeNI-ADMM framework: sketch-and-solve ADMM (Algorithm 4),
NysADMM (Algorithm 3), and “gradient descent” ADMM (GD-ADMM) (11). As a
baseline, we also compare to exact ADMM (Algorithm 1) to see how various
approximations or inexactness impact the convergence rate. We conduct three
sets of experiments that verify different aspects of the theory:

  * •

Section 8.1 verifies that NysADMM, GD-ADMM and sketch-and-solve ADMM converge
sublinearly for convex problems. Moreover, we observe a fast transition to
linear convergence, after which all methods but GD-ADMM converge quickly to
high accuracy.

  * •

Section 8.2 verifies that NysADMM, GD-ADMM, and sketch-and-solve ADMM converge
linearly for strongly convex problems.

  * •

Section 8.3 verifies that, without the correction term, sketch-and-solve ADMM
diverges, showing the necessity of the correction term in Section 7.2.

We consider three common problems in machine learning and statistics in our
experiments: lasso [41], elastic net regression [42], and $\ell_{1}$-logistic
regression[43]. All experiments use the realsim dataset from LIBSVM [44],
accessed through OpenML [45], which has $72,309$ samples and $20,958$
features. Our experiments use a subsample of realsim, consisting of $10,000$
random samples, which ensures the objective is not strongly convex for lasso
and $\ell_{1}$-logistic regression.

All methods solve the $z$-subproblem exactly, and every method, except
NysADMM, solves their $x$-subproblem exactly. NysADMM solves the linear system
in Algorithm 3 inexactly using Nyström PCG. At each iteration, the following
absolute tolerance on the PCG residual is used for determining the termination
of Nyström PCG:

$\varepsilon^{k}_{x}=\min\left\\{\frac{\sqrt{r_{p}^{k-1}r_{d}^{k-1}}}{k^{1.5}},1\right\\},$

where $r_{p}^{k-1}$ and $r_{d}^{k-1}$ are the ADMM primal and dual residuals
at the previous iteration. We multiply a summable sequence by the geometric
mean of the primal and dual residuals at the previous iteration. The geometric
mean of the primal and dual residuals is a commonly used tolerance for PCG in
practical ADMM implementations such as OSQP [46]—despite this sequence not
being apriori summable.

For NysADMM, a sketch size of $50$ is used to construct the Nyström
preconditioner, and for sketch-and-solve ADMM, we use a sketch size $500$ to
form the Hessian approximation. For sketch-and-solve ADMM, the parameter
$\gamma^{k}$ in (35) is chosen by estimating the error of the Nyström sketch
using power iteration.

We use a bound on the duality gap as the stopping criterion for each
algorithm. The details of how the bound on the duality gap bound is computed
are presented in Appendix B. All our experiments use a tolerance of
$\epsilon=10^{-4}$. An algorithm is terminated when the bound on the duality
falls below $\epsilon$ or has completed 500 iterations.

All experiments are performed in the Julia programming language [47] on a
MacBook Pro with a M1 Max processor and 64GB of RAM. To compute the “true”
optimal values, we use the commercial solver Mosek [48] (with tolerances set
low for high accuracy and presolve turned off to preserve the problem scaling)
and the modeling language JuMP [49, 50].

###  8.1 GeNI-ADMM converges sublinearly and locally linearly on convex
problems

To illustrate the global sublinear convergence of GeNI-ADMM methods on convex
objectives, we look at the performance of ADMM, NysADMM, GD-ADMM, and sketch-
and-solve ADMM on solving a lasso and $\ell_{1}$-logistic regression problem
with the realsim dataset. Note, as the number of samples is smaller than the
number of features, the corresponding optimization problems are convex but not
strongly convex.

##### Lasso regression

(a) Primal residual convergence

(b) Dual residual convergence

(c) Dual gap convergence

(d) Objective value convergence

Figure 1: Convergence of lasso regression for NysADMM, sketch-and-solve ADMM,
and gradient descent ADMM.

The lasso regression problem is to minimize the $\ell_{2}$ error of a linear
model with an $\ell_{1}$ penalty on the weights:

minimize $\displaystyle({1}/{2})\|Ax-b\|_{2}^{2}+\gamma\|x\|_{1}.$

This can be easily transformed into the form (1) by taking
$f(x)=(1/2)\|Ax-b\|_{2}^{2}$, $g(z)=\gamma\|z\|_{1}$, $M=I$, $N=-I$, and
$X=Z=\mathbb{R}^{n}$. We set $\gamma=0.05\cdot\gamma_{\mathrm{max}}$, where
$\gamma_{\mathrm{max}}=\|A^{T}b\|_{\infty}$ is the value above which the all
zeros vector is optimal111This can be seen from the optimality conditions,
which we state in Appendix B of the appendix.. We stop the algorithm when the
gap is less than $10^{-4}$, or after $500$ iterations.

The results of lasso regression are illustrated in Figure 1. Figure 1 shows
all methods initially converge at sublinear rate, but quickly transition to a
linear rate of convergence after reaching some manifold containing the
solution. ADMM, NysADMM, and sketch-and-solve ADMM (which use curvature
information) all converge converge much faster then GD-ADMM, confirming the
predictions of Section 5 that methods which use curvature information will
converge faster than methods that do not. Moreover, the difference in
convergence between NysADMM and ADMM is negligible, despite the former having
a much cheaper iteration complexity due to the use of inexact linear system
solves.

##### Logistic Regression

(a) Primal residual convergence

(b) Dual residual convergence

(c) Dual gap convergence

(d) Objective value convergence

Figure 2: Convergence of logistic regression for NysADMM, sketch-and-solve
ADMM, and gradient descent ADMM.

For the logistic regression problem, with data
$(\tilde{a}_{i},\tilde{b}_{i})$, where $\tilde{b}_{i}\in\\{\pm 1\\}$. We
define

$\mathbb{P}\left({\tilde{b}_{i}\mid\tilde{a}_{i}}\right)=\frac{1}{1+\exp\left(\tilde{b}_{i}(\tilde{a}_{i}^{T}x)\right)}=\frac{1}{1+\exp(a_{i}^{T}x)},$

where we define $a_{i}=\tilde{b}_{i}\tilde{a}_{i}$. The loss is the negative
of the log likelihood,
$-\log\left(\prod_{i=1}^{m}\mathbb{P}(\tilde{b}_{i}\mid\tilde{a}_{i})\right)$.
Thus, the optimization problem is

minimize
$\displaystyle\sum_{i=1}^{m}\log\left(1+\exp(a_{i}^{T}x)\right)+\gamma\|x\|_{1}.$

We set $\gamma=0.05\cdot\gamma_{\mathrm{max}}$, where
$\gamma_{\mathrm{max}}=(1/2)\|A^{T}\mathbf{1}\|_{\infty}$ is the value above
which the all zeros vector is optimal. For NysADMM, the preconditioner is re-
constructed after every 20 iterations. For sketch-and-solve ADMM, we re-
construct the approximate Hessian at every iteration. Since the $x$-subproblem
is not a quadratic program, we use the L-BFGS [51] implementation from the
Optim.jl package [52] to solve the $x$-subproblem.

Figure 2 presents the results for logistic regression. The results are
consistent with the lasso experiment—all methods initially converge
sublinearly before quickly transitioning to linear convergence, and methods
using better curvature information, converge faster. In particular, although
sketch-and-solve-ADMM converges slightly faster than GD-ADMM, its convergence
is much slower than NysADMM, which more accurately captures the curvature due
to using the exact Hessian. The convergence of NysADMM and ADMM is essentially
identical, despite the former having a much cheaper iteration cost due to
approximating the $x$-subproblem and using inexact linear system solves.

###  8.2 GeNI-ADMM converges linearly on strongly convex problems

(a) Primal residual convergence

(b) Dual residual convergence

(c) Objective value convergence

Figure 3: Convergence of elastic net regression for NysADMM, sketch-and-solve
ADMM, and gradient descent ADMM.

To verify the linear convergence of GeNI-ADMM methods in the presence of
strong convexity, we experiment with the elastic net problem:

minimize
$\displaystyle({1}/{2})\|Ax-b\|_{2}^{2}+\gamma\|x\|_{1}+(\mu/2)\|x\|_{2}^{2}.$

We set $\mu=1$, and use the same problem data and value of $\gamma$ as the
lasso experiment. The results of the elastic-net experiment are presented in
Figure 3. Comparing Figures 1 and 3, we clearly observe the linear convergence
guaranteed by the theory in Section 6. Although in Figure 1, ADMM and NysADMM
quickly exhibit linear convergence, Figure 3 clearly shows strong convexity
leads to an improvement in the number of iterations required to converge.
Moreover, we see methods that make better use of curvature information
converge faster than methods that do not, consistent with the results of the
lasso and logistic regression experiments.

(a) Lasso regression

(b) Logistic regression

Figure 4: We show the convergence or lack thereof for sketch-and-solve ADMM
with (solid lines) and without (dashed lines) the correction term required by
the theoretical results (see section 7.2). When this term is not included, the
algorithm does not converge.

###  8.3 Sketch-and-solve ADMM fails to converge without the correction term

To demonstrate the necessity of the correction term in Section 7.2, we run
sketch-and-solve ADMM on lasso and $\ell_{1}$-logistic regression without the
correction term. Figure 4 presents the results of these simulations. Without
the correction term in (35), sketch-and-solve ADMM quickly diverges on the
lasso problem. For the logistic regression problem, it oscillates and fails to
converge. These results highlight the importance of selecting this term
appropriately as discussed in Section 7.2.

##  9 Conclusion and future work

In this paper, we have developed a novel framework, GeNI-ADMM, that
facilitates efficient theoretical analysis of approximate ADMM schemes and
inspires the design of new, practical approximate ADMM methods. GeNI-ADMM
generalizes prior approximate ADMM schemes as special cases by allowing
various approximations to the $x$ and $z$-subproblems, which can then be
solved inexactly. We have established the usual ergodic
$\mathcal{O}(1/t)$-convergence rate for GeNI-ADMM under standard hypotheses,
and linear convergence under strong convexity. We have shown how to derive
explicit rates of convergence for ADMM variants that exploit randomized
numerical linear algebra using the GeNI-ADMM framework, We have provided
convergence results for NysADMM and sketch-and-solve ADMM as special cases,
resolving whether these schemes are globally convergent. Numerical experiments
on real-world data generally show an initial sublinear phase followed by
linear convergence, validating the theory we developed.

There are several interesting directions for future work. Our analysis shows
that in approximate ADMM schemes, a better approximation usually leads to
faster convergence. However, a better approximation usually makes
$\\{\Theta^{k}\\}_{k\geq 1}$ more complex and the resulting subproblem harder
to solve. A comprehensive study of the trade-off between theoretical
convergence speed and subproblem complexity depending on the choice of
$\\{\Theta^{k}\\}_{k\geq 1}$ will be an important future research direction.
Our numerical experiments show that GeNI-ADMM can achieve linear convergence
locally in the neighborhood of the optimal solution, even for problems that
are not strongly convex. Existing studies have proved that exact ADMM exhibits
local linear convergence without strong convexity [21]. It would be
interesting to establish local linear convergence for GeNI-ADMM schemes
without strong convexity, a phenomenon observed in all our experiments
(Section 8). Extending GeNI-ADMM to the multiblock setting is another exciting
avenue for future work. Such an extension will likely be nuanced, as standard
multi-block ADMM can diverge on convex objectives [53]. Our theoretical
analysis shows that curvature information provides an advantage in approximate
ADMM schemes. This observation can help the design of efficient ADMM-type
methods for a broader range of applications, including stochastic and
nonconvex optimization.

##  Appendix A Proofs not appearing in the main paper

###  A.1 Proofs for Section 4

We begin with the following lemma, which plays a key role in the proof of
Lemma 1. For a proof, see Example 1.2.2. in Chapter XI of [26].

######  Lemma 10 ($\varepsilon$-subdifferential of a quadratic function).

Let $h(z)=\langle b,z\rangle+\frac{1}{2}\|z\|_{A}^{2}$, where $A$ is a
symmetric positive definite linear operator. Then

$\displaystyle\partial_{\varepsilon}h(z)$
$\displaystyle=\left\\{b+A(z+v)\hskip 2.0pt\middle|\hskip
2.0pt\frac{\|v\|_{A}^{2}}{2}\leq\varepsilon\right\\}.$

With Lemma 10 in hand, we now prove Lemma 1.

##### Proof of Lemma 1

###### Proof.

Observe the function defining the $z$-subproblem may be decomposed as
$G(z)=g(z)+h(z)$ with

$h(z)=\frac{\rho}{2}\|M\tilde{x}^{k+1}+Nz+\tilde{u}^{k}\|^{2}_{2}+\frac{1}{2}\|z-\tilde{z}^{k}\|_{\Psi^{k}}^{2}.$

Now, by hypothesis $\tilde{z}^{k+1}$ gives an $\varepsilon_{z}^{k}$-minimum of
$G(z)$. Hence by Proposition 1,

$0\in\partial_{\varepsilon_{z}^{k}}G(\tilde{z}^{k+1})~{}\text{and}~{}\partial_{\varepsilon^{k}_{z}}G(\tilde{z}^{k+1})\subset\partial_{\varepsilon^{k}_{z}}g(\tilde{z}^{k+1})+\partial_{\varepsilon^{k}_{z}}h(\tilde{z}^{k+1}).$

Thus, we have $0=s_{g}+s_{h}$ where
$s\in\partial_{\varepsilon_{z}^{k}}g(\tilde{z}^{k+1})$ and
$s_{h}\in\partial_{\varepsilon_{z}^{k}}h(\tilde{z}^{k+1})$. Applying Lemma 10,
with $A=\Psi^{k}+\rho N^{T}N$ and $b=\rho
N^{T}(M\tilde{x}^{k+1}+\tilde{u}^{k})-\Psi^{k}\tilde{z}^{k}$, we reach

$s_{h}=\rho
N^{T}(M\tilde{x}^{k+1}+N\tilde{z}^{k+1}+\tilde{u}^{k})+\Psi^{k}(\tilde{z}^{k+1}-\tilde{z}^{k})+v.$

The desired claim now immediately follows from using $s=-s_{h}$. $\Box$

###  A.2 Proofs for Section 5

##### Proof of Proposition 2

###### Proof.

Using the eigendecomposition, we may decompose $H_{f}$ as

$H_{f}=\lambda_{1}(H_{f})\mathbf{v}\mathbf{v}^{T}+\sum^{d}_{i=2}\lambda_{i}(H_{f})\mathbf{v}_{i}\mathbf{v}^{T}_{i}.$

So

$\displaystyle\|x^{\star}\|^{2}_{H_{f}}$
$\displaystyle=\lambda_{1}(H_{f})\langle\mathbf{v},x^{\star}\rangle^{2}+\sum^{d}_{i=2}\lambda_{i}(H_{f})\langle\mathbf{v}_{i},x_{\star}\rangle^{2}$
$\displaystyle\leq\lambda_{1}(H_{f})\mu_{\mathbf{v}}\|x^{\star}\|^{2}+\lambda_{2}(H_{f})\sum^{d}_{i=2}\langle\mathbf{v}_{i},x^{\star}\rangle^{2}$
$\displaystyle\leq\left(\mu_{\mathbf{v}}+\frac{\lambda_{2}(H_{f})}{\lambda_{1}(H_{f})}\right)\lambda_{1}(H_{f})\|x^{\star}\|^{2}.$

Recalling that $\sigma=\tau\lambda_{1}(H_{f})$ with $\tau\in(0,1)$, we may put
everything together to conclude

$\displaystyle\frac{(\lambda_{1}(H_{f})-\sigma)\|x^{\star}\|^{2}}{\|x^{\star}\|^{2}_{H_{f}}}\geq(1-\tau)\frac{\lambda_{1}(H_{f})\|x^{\star}\|^{2}}{\|x^{\star}\|^{2}_{H_{f}}}\geq(1-\tau)\left(\mu_{\mathbf{v}}+\frac{\lambda_{2}(H_{f})}{\lambda_{1}(H_{f})}\right)^{-1},$

proving the first claim. The second claim follows immediately from the first.
$\Box$

##### Proof of Lemma 2

###### Proof.

By definition of $\ell_{U}(v,w)$

$\ell_{U}(v,w)=f(x)+g(z)-p^{\star}+\sup_{\hat{u}\in
U}\langle\rho\hat{u},v-(Mx+Nz)\rangle.$ (36)

Hence if $v=Mx+Nz$ and $\ell_{U}(v,w)\leq\epsilon$, (36) yields

$f(x)+g(z)-p^{\star}\leq\epsilon.$

Thus, we obtain the first statement.

Now, suppose $\ell_{U}(v,w)\leq\epsilon$, $\|v\|\leq\delta$ and
$U=\mathcal{Z}$. Then, combining these hypotheses with (36) and using
$p^{\star}$ that is the optimum, we reach

$\sup_{\hat{u}\in\mathcal{Z}}\langle\rho\hat{u},v-(Mx+Nz)\rangle\leq\epsilon.$

Observe if $v\neq Mx+Nz$ then
$\sup_{\hat{u}\in\mathcal{Z}}\langle\rho\hat{u},v-(Mx+Nz)\rangle=\infty$,
contradicting the preceding inequality. Thus $v=Mx+Nz$, which yields the
second statement, and completes the proof. $\Box$

##### Proof of Proposition 3

######  Lemma 11 (Exact Gap function 1-step bound ).

Let $w^{k+1}=(x^{k+1},z^{k+1},u^{k+1})$ denote the iterates generated by
Algorithm 2 at iteration $k+1$, when we solve the subproblems exactly. Then,
the following inequality holds

$\displaystyle Q(w,w^{k+1})$
$\displaystyle\leq\frac{1}{2}\left(\|z-\tilde{x}^{k}\|^{2}_{\Theta^{k}}-\|x-x^{k+1}\|^{2}_{\Theta^{k}}\right)+\frac{1}{2}\left(\|z-\tilde{z}^{k}\|^{2}_{\Psi^{k}}-\|z-z^{k+1}\|^{2}_{\Psi^{k}}\right)$
$\displaystyle\frac{\rho}{2}\left(\|\tilde{z}^{k}-Mx\|^{2}-\|z^{k+1}-Mx\|^{2}\right)+\frac{\rho}{2}\left(\|\tilde{u}^{k}-u\|^{2}-\|u^{k+1}-u\|^{2}\right).$

###### Proof.

The argument is identical to the proof of Lemma 5. Indeed, the only difference
is $\varepsilon_{x}^{k}=\varepsilon_{z}^{k}=0$, as the subproblems are solved
exactly. So the same line of argumentation may be applied, except the exact
subroblem optimality conditions are used. $\Box$

We now prove Proposition 3.

###### Proof.

First, observe that item 2. is an immediate consequence of item 1, so it
suffices to show item 1. To this end, plugging in $w=w^{\star}$, using
$Q(w^{\star},w^{k+1})\geq 0$, and rearranging, we reach

$\displaystyle\|x^{k+1}-x^{\star}\|^{2}_{\Theta^{k}}+\|z^{k+1}-z^{\star}\|_{\Psi^{k}+\rho
N^{T}N}^{2}+\rho\|u^{k+1}-u^{\star}\|^{2}$
$\displaystyle\leq\|\tilde{x}^{k}-x^{\star}\|^{2}_{\Theta^{k}}+\|\tilde{z}^{k}-z^{\star}\|_{\Psi^{k}+\rho
N^{T}N}^{2}+\rho\|\tilde{u}^{k}-u^{\star}\|^{2}.$

Defining the norm
$\|w\|_{W,k}=\sqrt{\|x\|^{2}_{\Theta^{k}}+\|z\|^{2}_{\Psi^{k}+\rho
N^{T}N}+\rho\|u\|^{2}}$, the preceding inequality may be rewritten as

$\|w^{k+1}-w^{\star}\|_{W,k}\leq\|\tilde{w}^{k}-w^{\star}\|_{W,k}.$

Now, our inexactness hypothesis (Assumption 3) along with $\nu$-strong
convexity of the $z$-subproblem (which follows by Assumption 2) implies

$\|\tilde{x}^{k+1}-x^{k+1}\|\leq\varepsilon_{x}^{k},\quad\|\tilde{z}^{k+1}-z^{k+1}\|\leq\sqrt{\frac{\varepsilon_{z}^{k}}{\nu}}.$

Using $\tilde{u}^{k+1}=M\tilde{x}^{k+1}+N\tilde{z}^{k+1}+\tilde{u}^{k}$ and
$u^{k+1}=Mx^{k+1}+Nz^{k+1}+\tilde{u}^{k}$, we find

$\rho\|\tilde{u}^{k+1}-u^{k+1}\|\leq\varepsilon_{x}^{k}+C\sqrt{\varepsilon_{z}^{k}}.$

Hence we have,

$\|\tilde{w}^{k+1}-w^{\star}\|_{W,k}\leq\|w^{k+1}-w^{\star}\|_{W,k}+\|\tilde{w}^{k+1}-w^{k+1}\|_{W,k}\leq\|\tilde{w}^{k}-w^{\star}\|_{W,k}+C\left(\varepsilon_{x}^{k}+\sqrt{\varepsilon_{z}^{k}}\right).$

Defining the summable sequence
$\varepsilon_{w}^{k}=C(\varepsilon_{x}^{k}+\sqrt{\varepsilon_{z}^{k}})$, the
preceding display may be written as

$\|\tilde{w}^{k+1}-w^{\star}\|_{W,k}\leq\|\tilde{w}^{k}-w^{\star}\|_{W,k}+\varepsilon_{w}^{k}\quad(\triangle),$

Now, Assumption 4 implies that
$\|\tilde{w}^{k}-w^{\star}\|_{W,k}\leq\sqrt{1+\zeta^{k}}\|\tilde{w}^{k}-w^{\star}\|_{W,k-1}$.
Combining this inequality with induction on $(\triangle)$, we find

$\|\tilde{w}^{k+1}-w^{\star}\|_{W,k}\leq\left(\prod^{k}_{j=2}\sqrt{1+\zeta^{j}}\right)\|\tilde{w}^{1}-w^{\star}\|_{W,1}+\left(\prod^{k}_{j=2}\sqrt{1+\zeta^{j}}\right)\sum_{j=1}^{k}\varepsilon_{w}^{j}.$

Hence,

$\sqrt{\min\\{\theta_{\textup{min}},\nu\\}}\|\tilde{w}^{k+1}-w^{\star}\|\leq\|\tilde{w}^{k+1}-w^{\star}\|_{W,k}\leq\sqrt{\tau_{\zeta}}\left(\|\tilde{w}^{1}-w^{\star}\|_{W,1}+\mathcal{E}_{w}\right).$

It follows immediately that:

$\sup_{k\geq
1}\|\tilde{w}^{k+1}-w^{\star}\|\leq\max\left\\{\frac{\sqrt{\tau_{\zeta}}\left(\|\tilde{w}^{1}-w^{\star}\|_{W,1}+\mathcal{E}_{w}\right)}{\sqrt{\min\\{\theta_{\textup{min}},\nu\\}}},\|\tilde{w}^{0}-w^{\star}\|\right\\}<\infty,$

and so the sequence $\\{\tilde{w}^{k}\\}_{k\geq 0}$ is bounded, which in turn
implies $\\{\tilde{x}^{k}\\}_{k\geq 0},\\{\tilde{z}^{k}\\}_{k\geq 0},$ and
$\\{\tilde{u}^{k}\\}_{k\geq 0}$ are bounded. $\Box$

##### Proof of Lemma 3

###### Proof.

Throughout the proof, we shall denote by $S_{x}^{k}(x)$, the function defining
the $x$-subproblem at iteration $k$. The exact solution of the $x$-subproblem
shall be denoted by $x^{k+1}$. To begin, observe that:

$\nabla S^{k}_{x}(x)=\nabla f_{1}(x)+\nabla
f_{2}(\tilde{x}^{k})+(\Theta^{k}+\rho M^{T}M)(x-\tilde{x}^{k})+\rho
M^{T}(M\tilde{x}^{k}+N\tilde{z}^{k}+\tilde{u}^{k}).$

Now,

$\displaystyle\langle\nabla
S^{k}_{x}(\tilde{x}^{k+1}),x^{\star}-\tilde{x}^{k+1}\rangle$
$\displaystyle=\langle\nabla S^{k}_{x}(\tilde{x}^{k+1})-\nabla
S^{k}_{x}(x^{k+1})+\nabla
S^{k}_{x}(x^{k+1}),x^{\star}-x^{k+1}+x^{k+1}-\tilde{x}^{k+1}\rangle$
$\displaystyle=\langle\nabla S^{k}_{x}(\tilde{x}^{k+1})-\nabla
S^{k}_{x}(x^{k+1}),x^{\star}-x^{k+1}\rangle$ $\displaystyle+\langle\nabla
S^{k}_{x}(\tilde{x}^{k+1})-\nabla
S^{k}_{x}(x^{k+1}),x^{k+1}-\tilde{x}^{k+1}\rangle$
$\displaystyle+\langle\nabla
S^{k}_{x}(x^{k+1}),x^{\star}-x^{k+1}\rangle+\langle\nabla
S^{k}_{x}(x^{k+1}),x^{k+1}-\tilde{x}^{k+1}\rangle$
$\displaystyle\overset{(1)}{\geq}\langle\nabla f_{1}(\tilde{x}^{k+1})-\nabla
f_{1}(x^{k+1}),x^{\star}-\tilde{x}^{k+1}\rangle+\langle\tilde{x}^{k+1}-x^{k+1},x^{\star}-x^{k+1}\rangle_{\Theta^{k}+\rho
M^{T}M}$ $\displaystyle+\langle\nabla f_{1}(\tilde{x}^{k+1})-\nabla
f_{1}(x^{k+1}),\tilde{x}^{k+1}-x^{k+1}\rangle+\langle\tilde{x}^{k+1}-x^{k+1},\tilde{x}^{k+1}-x^{k+1}\rangle_{\Theta^{k}+\rho
M^{T}M}$ $\displaystyle+\langle\nabla
S^{k}_{x}(x^{k+1}),x^{k+1}-\tilde{x}^{k+1}\rangle$
$\displaystyle\overset{(2)}{\geq}-C_{1}\varepsilon_{x}^{k}-C_{2}(\varepsilon^{k}_{x})^{2}+\langle\nabla
S^{k}_{x}(x^{k+1}),x^{k+1}-\tilde{x}^{k+1}\rangle.$

Here $(1)$ uses $\nabla S_{x}^{k}(x)-\nabla S_{x}^{k}(y)=\nabla
f_{1}(x)-\nabla f_{1}(y)+(\Theta^{k}+\rho M^{T}M)(x-y)$, and that $x^{k+1}$ is
the exact solution of the $x$-subproblem, while $(2)$ uses that
$\tilde{x}^{k+1}$ is an $\varepsilon_{x}^{k}$ minimizer of the $x$-subproblem
and that $f_{1}$ has a Lipschitz gradient. Now, as the iterates all belong to
a compact set, it follows that $\sup_{k}\|\nabla S^{k}_{x}(x^{k+1})\|\leq
C_{3}$. So,

$\displaystyle\langle\nabla
S^{k}_{x}(\tilde{x}^{k+1}),x^{\star}-\tilde{x}^{k+1}\rangle\geq-
C_{1}\varepsilon_{x}^{k}-(\varepsilon^{k}_{x})^{2}-C_{3}\varepsilon_{x}^{k}\geq-
C_{x}\varepsilon_{x}^{k}.$

This establishes the first inequality, and the second inequality follows from
the first by plugging $\tilde{x}^{k+1}$ into the expression for $\nabla
S_{x}^{k}(x)$. $\Box$

##### Proof of Lemma 4

###### Proof.

Proof. By hypothesis, $\tilde{z}^{k+1}$ is an
$\varepsilon_{z}^{k}$-approximate minimizer of the $z$-subproblem, so Lemma 1
shows there exists a vector $\tilde{s}$ with $\|\tilde{s}\|\leq
C\sqrt{\varepsilon_{z}^{k}}$, such that

$g(z^{\star})-g(\tilde{z}^{k+1})\geq\rho\langle\tilde{u}^{k+1}+\Psi^{k}(\tilde{z}^{k}-\tilde{z}^{k+1})+\tilde{s},z^{\star}-\tilde{z}^{k+1}\rangle-\varepsilon_{z}^{k}.$

Rearranging, using Cauchy-Schwarz, along with the boundedness of the iterates,
the preceding display becomes

$g(z^{\star})-g(\tilde{z}^{k+1})-\rho\langle\tilde{u}^{k+1}+\Psi^{k}(\tilde{z}^{k}-\tilde{z}^{k+1}),z^{\star}-\tilde{z}^{k+1}\rangle\geq-C\sqrt{\varepsilon_{z}^{k}}-\varepsilon_{z}^{k}\geq-
C_{z}\sqrt{\varepsilon_{z}^{k}},$

where the last display uses $\varepsilon_{z}^{k}\leq K_{\varepsilon_{z}}$, and
absorbs constant terms to get the constant $C_{z}$. $\Box$

###  A.3 Proofs for Section 6

##### Proof of Lemma 6

###### Proof.

  1. 1.

We begin by observing that strong convexity of $f$ implies

$\displaystyle\langle\nabla f(\tilde{x}^{k+1})-\nabla
f(x^{\star}),\tilde{x}^{k+1}-x^{\star}\rangle\geq\sigma_{f}\|\tilde{x}^{k+1}-x^{\star}\|^{2}.$

Decomposing $f$ as $f=f_{1}+f_{2}$, the preceding display may be rewritten as

$\displaystyle\langle\nabla f_{1}(\tilde{x}^{k+1})+\nabla
f_{2}(\tilde{x}^{k})-\nabla
f(x^{\star}),\tilde{x}^{k+1}-x^{\star}\rangle+\langle\nabla
f_{2}(\tilde{x}^{k+1})-\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-x^{\star}\rangle$
$\displaystyle\geq\sigma_{f}\|\tilde{x}^{k+1}-x^{\star}\|^{2}.$

Now, the exact solution $x^{k+1}$ of the $x$-subproblem at iteration $k$
satisfies $\nabla S_{x}^{k}(x^{k+1})=0$. Moreover, $\tilde{x}^{k+1}$ is an
$\varepsilon_{x}^{k}$ minimizer of the $x$-subproblem, and $\nabla
S_{x}^{k}(x)$ is Lipschitz continuous. Combining these three properties, it
follows there exists a vector $\textrm{Err}^{k}_{\nabla}$ satisfying:

$\nabla S^{k}_{x}(\tilde{x}^{k+1})=\nabla
S_{x}^{k}(x^{k+1})+\textrm{Err}^{k}_{\nabla}=\textrm{Err}^{k}_{\nabla},\quad\text{where}~{}\|\textrm{Err}^{k}_{\nabla}\|\leq
C\varepsilon_{x}^{k}.$

Consequently, we have

$\nabla f_{1}(\tilde{x}^{k+1})+\nabla
f_{2}(\tilde{x}^{k})=\Theta^{k}(\tilde{x}^{k}-\tilde{x}^{k+1})+\rho
M^{T}N(\tilde{z}^{k+1}-\tilde{z}^{k})-\rho
M^{T}\tilde{u}^{k+1}+\textrm{Err}^{k}_{\nabla}.$

Utilizing the preceding relation, along with the fact that the iterates are
bounded, we reach

$\displaystyle\langle\tilde{x}^{k}-\tilde{x}^{k+1},\tilde{x}^{k+1}-x^{\star}\rangle_{\Theta^{k}}+\rho\langle
N(\tilde{z}^{k+1}-\tilde{z}^{k})+u^{\star}-\tilde{u}^{k+1},M(\tilde{x}^{k+1}-x^{\star})\rangle+C\varepsilon_{x}^{k}$
$\displaystyle\geq\sigma_{f}\|\tilde{x}^{k+1}-x^{\star}\|^{2}-\langle\nabla
f_{2}(\tilde{x}^{k+1})-\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-x^{\star}\rangle.~{}\Box$

  2. 2.

As $-\rho N^{T}u^{\star}\in\partial g(z^{\star})$, and $-\rho
N^{T}\tilde{u}^{k+1}+\Psi^{k}(\tilde{z}^{k}-\tilde{z}^{k+1})+\tilde{s}\in\partial_{\varepsilon_{z}^{k}}g(\tilde{z}^{k+1})$,
we have

$\displaystyle g(\tilde{z}^{k+1})-g(z^{\star})\geq\langle-\rho
N^{T}u^{\star},\tilde{z}^{k+1}-z^{\star}\rangle,$ $\displaystyle
g(z^{\star})-g(\tilde{z}^{k+1})\geq\langle-\rho
N^{T}\tilde{u}^{k+1}+\Psi^{k}(\tilde{z}^{k}-\tilde{z}^{k+1}),z^{\star}-\tilde{z}^{k+1}\rangle-\varepsilon_{z}^{k}.$

So, adding together the two inequalities and rearranging, we reach

$\langle\tilde{z}^{k+1}-z^{\star},\rho
N^{T}(u^{\star}-\tilde{u}^{k+1})+\Psi^{k}(\tilde{z}^{k}-\tilde{z}^{k+1})+\tilde{s}\rangle\geq-\varepsilon_{z}^{k}.$

Now using boundedness of the iterates and $\|\tilde{s}\|\leq
C\sqrt{\varepsilon_{z}^{k}}$, Cauchy-Schwarz yields

$\langle\tilde{z}^{k+1}-z^{\star},\rho
N^{T}(\tilde{u}^{k+1}-u^{\star})+\Psi^{k}(\tilde{z}^{k}-\tilde{z}^{k+1})\rangle\geq-C\sqrt{\varepsilon_{z}^{k}}.$

$\Box$

##### Proof of Lemma 8

We wish to show for some $\delta>0$, that the following inequality holds

$\|w^{k}-w^{\star}\|_{G^{k}}^{2}-\|w^{k+1}-w^{\star}\|_{G^{k}}^{2}+\varepsilon_{w}^{k}\geq\delta\|w^{k+1}-w^{\star}\|_{G^{k}}^{2}.$

To establish this result, it suffices to show the inequality

$\|w^{k+1}-w^{k}\|^{2}_{G^{k}}+\frac{2\sigma_{f}}{\rho}\|x^{k+1}-x^{\star}\|^{2}-\frac{2}{\rho}\langle\nabla
f_{2}(x^{k+1})-\nabla
f_{2}(x^{k}),x^{k+1}-x^{\star}\rangle\geq\delta\|w^{k+1}-w_{\star}\|_{G^{k}}^{2}.$

We accomplish this by upper bounding each term that appears in
$\|w^{k+1}-w_{\star}\|_{G^{k}}^{2}$, by terms that appear in the left-hand
side of the preceding equality. To that end, we have the following result.

######  Lemma 12 (Coupling Lemma).

Under the assumptions of Theorem 2, the following statements hold:

  1. 1.

Let $\mu_{1}\geq 2$,
$c_{1}=\frac{12L^{2}_{f}}{\rho^{2}\lambda_{\textup{min}}(MM^{T})},c_{2}=\frac{4\left(2L^{2}_{f}+\theta_{\textup{max}}^{2}\right)}{\rho^{2}\lambda_{\textup{min}}(MM^{T})},c_{3}=8\|M\|^{2}$.
Then for some constant $C\geq 0$, we have

$\|u^{k+1}-u^{\star}\|^{2}\leq
c_{1}\|\tilde{x}^{k+1}-x^{\star}\|^{2}+c_{2}\|\tilde{x}^{k}-\tilde{x}^{k+1}\|^{2}+c_{3}\|\tilde{z}^{k}-\tilde{z}^{k+1}\|_{1/\rho\Psi^{k}+N^{T}N}^{2}+C\varepsilon_{w}^{k}.$

  2. 2.

There exists constants
$c_{4}=2\left(\frac{\|N\|^{2}+\psi_{\textup{max}}/\rho}{\sigma^{2}_{\textrm{min}}(N)}\right)\|M\|^{2}$,
$c_{5}=c_{4}/\|M\|^{2}$, such that

$\|z^{k+1}-z^{\star}\|_{1/\rho\Psi^{k}+N^{T}N}^{2}\leq
c_{4}\|x^{k+1}-x^{\star}\|^{2}+c_{5}\|u^{k+1}-u^{k}\|^{2}.$

  3. 3.

For all $\mu>0$, we have

$\displaystyle\langle\nabla f_{2}(x^{k+1})-\nabla
f_{2}(x^{k}),x^{k+1}-x^{\star}\rangle\leq\mu/2\|x^{k+1}-x^{\star}\|^{2}+\frac{L^{2}_{f}}{2\mu}\|x^{k+1}-x^{k}\|^{2}.$

Taking Lemma 12 as given, let us prove Lemma 8.

###### Proof.

Observe Lemma 12 implies

$\displaystyle\frac{1}{2}\|w^{k+1}-w^{k}\|^{2}_{G^{k}}+\frac{2\sigma_{f}}{\rho}\|x^{k+1}-x^{\star}\|^{2}-\frac{2}{\rho}\langle\nabla
f_{2}(x^{k+1})-\nabla
f_{2}(x^{k}),x^{k+1}-x^{\star}\rangle-\delta\|w^{k+1}-w_{\star}\|_{G_{k}}^{2}$
$\displaystyle\geq\left(\frac{2\sigma_{f}-\mu_{2}}{\rho}-\delta(c_{1}+c_{4})\right)\|\tilde{x}^{k+1}-x^{\star}\|^{2}+\left(\frac{\theta_{\textup{min}}}{2\rho}-\frac{L_{f}^{2}}{2\rho\mu_{2}}-\delta
c_{2}\right)\|\tilde{x}^{k}-\tilde{x}^{k+1}\|^{2}$
$\displaystyle+\left(1/2-\delta
c_{3}\right)\|\tilde{z}^{k}-\tilde{z}^{k+1}\|_{1/\rho\Psi^{k}+N^{T}N}^{2}+\left(1/2-\delta
c_{5}\right)\|\tilde{u}^{k}-\tilde{u}^{k+1}\|^{2}-C\varepsilon_{w}^{k}.$

Setting $\mu_{2}=\sigma_{f}$, and using
$\theta_{\textup{min}}>L^{2}_{f}/\sigma_{f}=\kappa_{f}L_{f}$, we find by
setting

$\delta=\min\left\\{\frac{\sigma_{f}}{\rho(c_{1}+c_{4})},\frac{\theta_{\textup{min}}-\kappa_{f}L_{f}}{2\rho
c_{2}},\frac{1}{2c_{3}},\frac{1}{2c_{5}}\right\\}>0,$

that

$\|\tilde{w}^{k}-w^{\star}\|^{2}_{G^{k}}-\|\tilde{w}^{k+1}-w^{\star}\|^{2}_{G^{k}}+C\varepsilon_{w}^{k}\geq\delta\|\tilde{w}^{k+1}-w^{\star}\|^{2}_{G^{k}},$

as desired. $\Box$

We now turn to to the proof of Lemma 12.

##### Proof of Lemma 12

###### Proof.

  1. 1.

Observe by $L_{f}$-smoothness of $f$, that

$\displaystyle
2L^{2}_{f}\left(\|\tilde{x}^{k+1}-x^{\star}\|^{2}+\|\tilde{x}^{k}-x^{\star}\|^{2}\right)$
$\displaystyle\geq\|\nabla f_{1}(\tilde{x}^{k+1})-\nabla
f_{1}(x^{\star})+\nabla f_{2}(\tilde{x}^{k})-\nabla f_{2}(x^{\star})\|^{2}$
$\displaystyle=\|\Theta^{k}(\tilde{x}^{k}-\tilde{x}^{k+1})+\rho
M^{T}N(\tilde{z}^{k}-\tilde{z}^{k+1})+\rho
M^{T}(u^{\star}-\tilde{u}^{k+1})+\textrm{Err}^{k}_{\nabla}\|^{2}$
$\displaystyle\geq\frac{1}{2}\rho^{2}\|M^{T}(\tilde{u}^{k+1}-u^{\star})\|^{2}-\|\Theta^{k}(\tilde{x}^{k}-\tilde{x}^{k+1})+\rho
M^{T}N(\tilde{z}^{k}-\tilde{z}^{k+1})+\textrm{Err}^{k}_{\nabla}\|^{2}.$

Here the last inequality uses with $\mu=2$, the identity (valid for all
$\mu>1$)

$\|a+b\|^{2}\geq(1-\mu^{-1})\|a\|^{2}+(1-\mu)\|b\|^{2}.$

So, rearranging and using $\|a+b\|^{2}\leq 2\|a\|^{2}+2\|b\|^{2}$, we reach

$\displaystyle\frac{\rho^{2}}{2}\|M^{T}(\tilde{u}^{k+1}-u^{\star})\|^{2}$
$\displaystyle\leq
6L^{2}_{f}\|\tilde{x}^{k+1}-x^{\star}\|^{2}+2\left(2L^{2}_{f}+\|\Theta^{k}\|^{2}\right)\|\tilde{x}^{k}-\tilde{x}^{k+1}\|^{2}$
$\displaystyle+4\rho^{2}\|M\|^{2}\|N(\tilde{z}^{k}-\tilde{z}^{k+1})\|^{2}+4\|\textrm{Err}^{k}_{\nabla}\|^{2}.$

Now, using $\|\textrm{Err}^{k}_{\nabla}\|^{2}\leq C\varepsilon_{x}^{k}\leq
C\varepsilon_{w}^{k}$ and $\|v\|_{N^{T}N}\leq\|v\|_{1/\rho\Psi^{k}+N^{T}N}$,
we reach

$\displaystyle\|\tilde{u}^{k+1}-u^{\star}\|^{2}\leq\frac{12L^{2}_{f}}{\rho^{2}\lambda_{\textup{min}}(MM^{T})}\|\tilde{x}^{k+1}-x^{\star}\|^{2}+\frac{4\left(2L^{2}_{f}+\theta_{\textup{max}}^{2}\right)}{\rho^{2}\lambda_{\textup{min}}(MM^{T})}\|\tilde{x}^{k}-\tilde{x}^{k+1}\|^{2}$
$\displaystyle
8\|M\|^{2}\|\tilde{z}^{k}-\tilde{z}^{k+1}\|_{1/\rho\Psi^{k}+N^{T}N}^{2}+C\varepsilon_{w}^{k}.~{}\Box$

  2. 2.

This inequality is a straightforward consequence of the relation

$M(\tilde{x}^{k+1}-x^{\star})=\tilde{u}^{k+1}-\tilde{u}^{k}+N(z^{\star}-\tilde{z}^{k+1}).$

Indeed, using the identity $\|a+b\|^{2}\leq 2\|a\|^{2}+2\|b\|^{2}$, we reach

$\displaystyle\|N(\tilde{z}^{k+1}-z^{\star})\|^{2}\leq
2\|M\|^{2}\|\tilde{x}^{k+1}-x^{\star}\|+2\|\tilde{u}^{k}-\tilde{u}^{k+1}\|^{2}.$

Consequently

$\displaystyle\|\tilde{z}^{k+1}-z^{\star}\|_{1/\rho\Psi^{k}+N^{T}N}^{2}$
$\displaystyle\leq
2\left(\frac{\|N\|^{2}+\psi_{\textup{max}}/\rho}{\sigma^{2}_{\textrm{min}}(N)}\right)\|M\|^{2}\|\tilde{x}^{k+1}-x^{\star}\|$
$\displaystyle+2\left(\frac{\|N\|^{2}+\psi_{\textup{max}}/\rho}{\sigma^{2}_{\textrm{min}}(N)}\right)\|\tilde{u}^{k}-\tilde{u}^{k+1}\|^{2},$

which is precisely the desired claim. $\Box$

  3. 3.

Young’s inequality implies for all $\mu>0$, that

$\langle\nabla f_{2}(\tilde{x}^{k+1})-\nabla
f_{2}(\tilde{x}^{k}),\tilde{x}^{k+1}-x^{\star}\rangle\leq\frac{1}{2\mu}\|\nabla
f_{2}(\tilde{x}^{k+1})-\nabla
f_{2}(\tilde{x}^{k})\|^{2}+\frac{\mu}{2}\|\tilde{x}^{k+1}-x^{\star}\|^{2}.$

The desired claim now follows from $L_{f}$-smoothness of $f_{2}$. $\Box$

##  Appendix B Dual gap bounds for ADMM

In this section, we derive bounds for the duality gap for the machine learning
examples we consider in Section 8. The derivation is largely inspired by [54,
55]. We consider problems of the form

minimize
$\displaystyle\ell(x)=\sum_{i=1}^{m}f(a_{i}^{T}x-b_{i})+\gamma\|x\|_{1},$ (37)

where $f(\cdot)$ is some loss function, which we assume to be convex and
differentiable. (The latter is not necessary but is convenient for our
purposes.)

##### Dual problem.

We introduce new variable $y_{i}$ and reformulate the primal problem (37) as

minimize $\displaystyle l(x)=\sum_{i=1}^{m}f(y_{i})+\gamma\|x\|_{1}$ (38)
subject to $\displaystyle y=Ax-b.$

The Lagrangian is then

$L(x,y,\nu)=\sum_{i=1}^{m}f(y_{i})+\gamma\|x\|_{1}+\nu^{T}(Ax-b-y).$

Taking minimization over $x$ and $y$, we get the dual function

$g(\nu)=\begin{cases}-\sum_{i=1}^{m}f^{*}(\nu_{i})-b^{T}\nu,&\|A^{T}\nu\|_{\infty}\leq\gamma\\\
-\infty,&\text{otherwise}.\end{cases}$

Thus, the dual problem is

maximize $\displaystyle g(\nu)$ (39) subject to
$\displaystyle\|A^{T}\nu\|_{\infty}\leq\gamma.$

Importantly, for any dual feasible point $\nu$, we have that

$\ell(x)\geq\ell(x^{\star})=g(\nu^{\star})\geq g(\nu),$

where $x^{\star}$ and $\nu^{\star}$ are optimal solutions to (37) and (39)
respectively. It immediately follows that

$\frac{\ell(x)-g(\nu)}{\max\left(\ell(x),\lvert{g(\nu)}\rvert\right)}\geq\frac{\ell(x)-\ell(x^{\star})}{\max\left(\ell(x),\lvert{g(\nu)}\rvert\right)}.$
(40)

We will call this the _relative dual gap_.

##### Dual feasible points.

The optimality conditions of (38) are (after eliminating $y$)

$\nu_{i}^{\star}=f^{\prime}(a_{i}^{T}x^{\star}-b_{i})$ (41)

Inspired by the second condition, we construct a dual feasible point by taking
$\nu_{i}=f^{\prime}(a_{i}^{T}x-b_{i})$ and then projecting onto the norm ball
constraint. Specifically, we take

$\nu_{i}=\frac{\gamma}{\left\|\sum_{i=1}^{m}f^{\prime}(a_{i}^{T}x-b_{i})\cdot
a_{i}\right\|_{\infty}}\cdot f^{\prime}(a_{i}^{T}x-b_{i}).$

When $x\neq x^{\star}$, we can verify that $\nu$ is dual feasible:

$\|A^{T}\nu\|_{\infty}=\left\|\frac{\gamma}{\left\|\sum_{i=1}^{m}f^{\prime}(a_{i}^{T}x-b_{i})\cdot
a_{i}\right\|_{\infty}}\cdot\left(\sum_{i=1}^{m}f^{\prime}(a_{i}^{T}x-b_{i})\cdot
a_{i}\right)\right\|_{\infty}=\gamma.$

Assuming not all $x_{i}=0$, then this dual variable is optimal when
constructed with $x=x^{\star}$. This fact can be easily seen by using the
optimality conditions of (37), that $0$ is contained in the subdifferential,
i.e.,

$0\in\sum_{i=1}^{m}f^{\prime}(a_{i}^{T}x-b_{i})\cdot
a_{i}+\gamma\partial\|x\|_{1}.$

By rearranging, we have the conditions

$\left(\sum_{i=1}^{m}f^{\prime}(a_{i}^{T}x-b_{i})\cdot
a_{i}\right)_{i}\in\begin{cases}\\{+\gamma\\}&x_{i}<0\\\
\\{-\gamma\\}&x_{i}>0\\\ [-\gamma,\gamma]&x_{i}=0.\end{cases}$

Thus, when constructed with $x^{\star}$, $\nu$ is contained within the norm
ball and therefore optimal by (41). Note that these optimality conditions also
tell us that $x=0$ is a solution to (37) if and only if
$\|A^{T}f^{\prime}(-b)\|_{\infty}\leq\gamma$, where the function $f^{\prime}$
is applied elementwise.

This construction and the bound (40) suggest a natural stopping criterion. For
any primal feasible point, we construct a dual feasible point $\nu$. Using the
dual feasible point, we evaluate the duality gap. We then terminate the solver
if the relative gap is less than some tolerance $\epsilon$ as

$\frac{l(x)-g(\nu)}{\max\left(\lvert{g(\nu)\rvert},l(x)\right)}\leq\epsilon.$

If this condition is met, we are guaranteed that the true relative error is
also less than $\epsilon$ from (40).

## References

  *   * Boyd et al. [2011] Boyd, S., Parikh, N., Chu, E., Peleato, B., Eckstein, J.: Distributed optimization and statistical learning via the alternating direction method of multipliers. Foundations and Trends in Machine Learning 3, 1–122 (2011)
  * Bubeck [2015] Bubeck, S.: Convex optimization: Algorithms and complexity. Foundations and Trends® in Machine Learning 8(3-4), 231–357 (2015)
  * Zhao et al. [2022] Zhao, S., Frangella, Z., Udell, M.: NysADMM: faster composite convex optimization via low-rank approximation. In: International Conference on Machine Learning, vol. 162, pp. 26824–26840 (2022). PMLR
  * Ouyang et al. [2015] Ouyang, Y., Chen, Y., Lan, G., Pasiliao, E.: An accelerated linearized alternating direction method of multipliers. SIAM Journal on Imaging Sciences 8(1), 644–681 (2015)
  * Eckstein and Bertsekas [1992] Eckstein, J., Bertsekas, D.P.: On the Douglas—Rachford splitting method and the proximal point algorithm for maximal monotone operators. Mathematical Programming 55(1), 293–318 (1992)
  * Lacotte and Pilanci [2020] Lacotte, J., Pilanci, M.: Effective dimension adaptive sketching methods for faster regularized least-squares optimization. In: Advances in Neural Information Processing Systems (2020)
  * Frangella et al. [2021] Frangella, Z., Tropp, J.A., Udell, M.: Randomized Nyström preconditioning. arXiv preprint arXiv:2110.02820 (2021)
  * Meier and Nakatsukasa [2022] Meier, M., Nakatsukasa, Y.: Randomized algorithms for Tikhonov regularization in linear least squares. arXiv preprint arXiv:2203.07329 (2022)
  * Buluc et al. [2021] Buluc, A., Kolda, T.G., Wild, S.M., Anitescu, M., Degennaro, A., Jakeman, J.D., Kamath, C., Kannan, R.-R., Lopes, M.E., Martinsson, P.-G., et al.: Randomized algorithms for scientific computing (RASC). Technical report, Oak Ridge National Lab.(ORNL), Oak Ridge, TN (United States) (2021)
  * Ryu and Yin [2022] Ryu, E.K., Yin, W.: Large-scale Convex Optimization: Algorithms & Analyses Via Monotone Operators, (2022). Cambridge University Press
  * Parikh and Boyd [2014] Parikh, N., Boyd, S.: Proximal algorithms. Foundations and Trends® in Optimization 1(3), 127–239 (2014)
  * Rontsis et al. [2022] Rontsis, N., Goulart, P., Nakatsukasa, Y.: Efficient semidefinite programming with approximate admm. Journal of Optimization Theory and Applications 192(1), 292–320 (2022)
  * Deng and Yin [2016] Deng, W., Yin, W.: On the global and linear convergence of the generalized alternating direction method of multipliers. Journal of Scientific Computing 66(3), 889–916 (2016)
  * He and Yuan [2012] He, B., Yuan, X.: On the o(1/n) convergence rate of the Douglas–Rachford alternating direction method. SIAM Journal on Numerical Analysis 50(2), 700–709 (2012)
  * Eckstein and Yao [2018] Eckstein, J., Yao, W.: Relative-error approximate versions of douglas–rachford splitting and special cases of the admm. Mathematical Programming 170(2), 417–444 (2018)
  * Hermans et al. [2019] Hermans, B., Themelis, A., Patrinos, P.: QPALM: a Newton-type proximal augmented Lagrangian method for quadratic programs. In: 2019 IEEE 58th Conference on Decision and Control (CDC), pp. 4325–4330 (2019). IEEE
  * Hermans et al. [2022] Hermans, B., Themelis, A., Patrinos, P.: QPALM: a proximal augmented Lagrangian method for nonconvex quadratic programs. Mathematical Programming Computation, 1–45 (2022)
  * Monteiro and Svaiter [2013] Monteiro, R.D., Svaiter, B.F.: Iteration-complexity of block-decomposition algorithms and the alternating direction method of multipliers. SIAM Journal on Optimization 23(1), 475–507 (2013)
  * Ouyang et al. [2013] Ouyang, H., He, N., Tran, L., Gray, A.: Stochastic alternating direction method of multipliers. In: International Conference on Machine Learning, pp. 80–88 (2013). PMLR
  * Hager and Zhang [2020] Hager, W.W., Zhang, H.: Convergence rates for an inexact admm applied to separable convex optimization. Computational Optimization and Applications 77(3), 729–754 (2020)
  * Yuan et al. [2020] Yuan, X., Zeng, S., Zhang, J.: Discerning the linear convergence of ADMM for structured convex optimization through the lens of variational analysis. J. Mach. Learn. Res. 21, 83–1 (2020)
  * Chambolle and Pock [2011] Chambolle, A., Pock, T.: A first-order primal-dual algorithm for convex problems with applications to imaging. Journal of Mathematical Imaging and Vision 40, 120–145 (2011)
  * Zhao et al. [2021] Zhao, S., Lessard, L., Udell, M.: An automatic system to detect equivalence between iterative algorithms. arXiv preprint arXiv:2105.04684 (2021)
  * Gower et al. [2019] Gower, R.M., Kovalev, D., Lieder, F., Richtárik, P.: RSN: randomized subspace Newton. In: Advances in Neural Information Processing Systems (2019)
  * Bertsekas et al. [2003] Bertsekas, D., Nedic, A., Ozdaglar, A.: Convex Analysis and Optimization vol. 1, (2003). Athena Scientific
  * Hiriart-Urruty and Lemaréchal [1993] Hiriart-Urruty, J.-B., Lemaréchal, C.: Convex Analysis and Minimization Algorithms II vol. 305, (1993). Springer
  * Schmidt et al. [2011] Schmidt, M., Roux, N., Bach, F.: Convergence rates of inexact proximal-gradient methods for convex optimization. Advances in Neural Information Processing Systems 24 (2011)
  * Stellato et al. [2020] Stellato, B., Banjac, G., Goulart, P., Bemporad, A., Boyd, S.: OSQP: an operator splitting solver for quadratic programs. Mathematical Programming Computation 12(4), 637–672 (2020)
  * He et al. [2002] He, B., Liao, L.-Z., Han, D., Yang, H.: A new inexact alternating directions method for monotone variational inequalities. Mathematical Programming 92, 103–118 (2002)
  * Rockafellar [2023] Rockafellar, R.T.: Generic linear convergence through metric subregularity in a variable-metric extension of the proximal point algorithm. Computational Optimization and Applications, 1–20 (2023)
  * Candes and Romberg [2007] Candes, E., Romberg, J.: Sparsity and incoherence in compressive sampling. Inverse problems 23(3), 969 (2007)
  * Candes and Recht [2012] Candes, E., Recht, B.: Exact matrix completion via convex optimization. Communications of the ACM 55(6), 111–119 (2012)
  * Nishihara et al. [2015] Nishihara, R., Lessard, L., Recht, B., Packard, A., Jordan, M.: A general analysis of the convergence of ADMM. In: International Conference on Machine Learning, pp. 343–352 (2015). PMLR
  * Friedman et al. [2010] Friedman, J., Hastie, T., Tibshirani, R.: Regularization paths for generalized linear models via coordinate descent. Journal of Statistical Software 33(1), 1 (2010)
  * Pilanci and Wainwright [2017] Pilanci, M., Wainwright, M.J.: Newton sketch: A near linear-time optimization algorithm with linear-quadratic convergence. SIAM Journal on Optimization 27(1), 205–245 (2017)
  * Lacotte et al. [2021] Lacotte, J., Wang, Y., Pilanci, M.: Adaptive Newton sketch: Linear-time optimization with quadratic convergence and effective Hessian dimensionality. In: International Conference on Machine Learning, pp. 5926–5936 (2021). PMLR