Consider applying first-order methods to solve the smooth convex constrained optimization problem of the form $\min_{x \in X} F(x).$ For a simple closed convex set $X \subseteq \bbr^n$ which is easy to project onto, Nesterov proposed the Accelerated Gradient Descent (AGD) method to solve the constrained problem as efficiently as an unconstrained problem in terms of the number of gradient computations of $F$ (i.e., oracle complexity). For a more complicated $\mathcal{X}$ described by function constraints, i.e., $\mathcal{X} = \{x \in X: g(x) \leq 0\}$, where the projection onto $\mathcal{X}$ is not possible,  it is an open question whether the function constrained problem can be solved as efficiently as an unconstrained problem in terms of the number of gradient computations for $F$ and $g$. In this paper, we provide an affirmative answer to the question by proposing a single-loop Accelerated Constrained Gradient Descent (ACGD) method. The ACGD method modifies the AGD method by changing the descent step to a constrained descent step, which adds only a few linear constraints to the prox mapping. It enjoys almost the same oracle complexity as the optimal one for minimizing the optimal Lagrangian function, i.e., the Lagrangian multiplier $\lambda$ being fixed to the optimal multiplier $\lambda^*$. These upper oracle complexity bounds are shown to be unimprovable under a certain optimality regime with new lower oracle complexity bounds. To enhance its efficiency for large-scale problems with many function constraints, we introduce an ACGD with Sliding (ACGD-S) method which replaces the possibly computationally demanding constrained descent step with a sequence of basic matrix-vector multiplications. The ACGD-S method shares the same oracle complexity as the ACGD method, and its computation complexity, measured by the number of matrix-vector multiplications, is also unimprovable. So the method provides an almost complete characterization of the hardness of smooth function-constrained optimization from both the oracle complexity and the computation complexity perspectives.

Keywords: constrained optimization,nonlinear optimization, first-order algorithm, convex optimization, lower complexity.

§ INTRODUCTION

Consider convex smooth function-constrained optimization of the form

\begin{equation}\label{eq:prob}
\begin{split}
\min_{x \in X}\ & \{F(x) :=f(x) + u(x)\}\\
s.t.\ & g(x) \leq 0,
\end{split}
\end{equation}

where both $f:\R^n \rightarrow \R$ and $g := [g_1, g_2, \ldots, g_m]^\top:\R^n \rightarrow \R^m$ are convex and differentiable with Lipschitz continuous gradients, and the domain $X$ is convex and closed. The regularization function $u(x)$ is assumed to be both simple [14] and $\alpha$-strongly convex for some $\alpha\geq 0$. Problem of this type finds a wide range of applications in, for example, the Neyman-Pearson classification problem [22], the fairness-constrained classification [25], and the risk-constrained portfolio optimization [5]. Since the dimensions $n$ and $m$ are large in many applications, we focus on first-order methods to find an approximation solution. Specifically, $(f, g)$ is assumed to be accessible only via a black-box first-order oracle, which returns $(f(x), g(x); \grad f(x), \grad g(x))$ when queried at some $x \in \R^n$, and the goal is to find an $(\ep; \ep/c)$-optimal solution (or, in short, an $\ep$-solution):

\begin{equation}\label{eq:opt-criteria}
F(x^N) - F(\xstar) \leq \ep \text{ and} \norm{[g(x^N)]_+} \leq \ep/c,\end{equation}

where the scaling constant $c \geq 1$ represents the modeler's preference for constraint violation relative to sub-optimality.

This paper intends to develop fast methods and to understand the requisite computation cost. To ensure practical efficiency, special attention is paid to simple methods involving only oracle evaluations, projections onto $X$, and basic vector operations like matrix-vector multiplication.
In particular, we consider
consider two kinds of cost involved in solving (<ref>): a) the oracle complexity, i.e., the total number of queries to the first-order oracle, and b) the computation complexity, i.e., the total number of matrix-vector multiplications (each costing at most $O(mn)$ FLOPs). During implementation, these complexities translate to different computation burdens, so the dominating cost depends on the context. For instance,
if the constraint $g$ is complicated, say a large finite-sum function in the fairness-constrained problem, the gradient evaluation could be the bottleneck. On the other hand, if both $n$ and $m$ are large and $g$ is an affine function, the matrix-vector multiplication might be the bottleneck. Thus the ideal optimization method should excel in both directions.

Ideal Complexity for Solving (<ref>)

Case     Oracle Complexity     Computation Complexity

Non-strongly Convex $\alpha=0$     $\bigO(1/\sqrt{\ep})$     $\bigO(1/\ep)$

Strongly Convex $\alpha > 0$     $\bigO(\log(1/\ep))$     $\bigO(1/\sqrt{\ep})$

Since convex optimization has been studied extensively, e.g. [19], we can conjecture the best possible complexities for (<ref>) by reducing it to simpler cases for which the optimal results are available. Recall from [16, 17], the number of oracle evaluations required to find an $\ep$-optimal solution to a smooth problem without function constraints is $\Theta(1/\sqrt{\ep})$ in general, and $\Theta(\sqrt{\kappa}\log(1/\ep))$ if the problem is also strongly convex, where $\kappa$ denotes the condition number.
One way to reduce (<ref>)  is to consider the Lagrange dual formulation:

\begin{equation}\label{eq:La}
\min_{x\in X} \max_{\lambda \in \R^m_+} \{\La(x; \lambda):=f(x) + u(x) + \lambda^\top (g(x))\}.\end{equation}

If the multiplier $\lam$ is fixed to the optimal dual multiplier $\lamstar$, the optimal Lagrangian function $\La(x; \lamstar)$ is smooth and does not have any function constraint. Since some minimizer to $\La(x; \lamstar)$ is the solution to (<ref>), in the ideal case, optimizing (<ref>) should have the same oracle complexity as optimizing $\La(x; \lamstar)$, i.e., $\bigO(1/\sqrt{\ep})$ for the non-strongly case ($\alpha =0$) and   $\bigO(\sqrt{\kappa}\log(1/\ep))$ for the strongly convex case ($\alpha > 0$), where $\kappa$ is the condition number associated with $\La(x; \lamstar)$. On the other hand, inspired by [19] and [15], Xu and Ouyang constructed novel large-scale linearly constrained quadratic programs in [21] to show the tight computation complexity for the linearly constrained smooth problem is $\Theta(1/\ep)$ if $\alpha =0 $, and $\Theta(1/\sqrt{\ep})$ if  $\alpha > 0$.  Since the nonlinear constraint function includes the affine function as a special case,  the computation complexity for (<ref>) should be at least as expensive as the affine case. These conjectured complexities, summarized in Table <ref>, lead naturally to the research question:

Can we solve function-constrained problems with the same oracle complexities as those without function constraints, while maintaining the same computation complexities as solving linearly constrained problems?

However, despite much research effort on the subject from different directions, the question remains open. Broadly speaking, the current results can be grouped according to reformulations (see the summary in Table <ref>). The methods based on the Lagrangian formulation (<ref>) are usually simple to implement and have a low per-iteration cost. For instance, the ConEx method in [4] and the APD method in [6] are single-loop algorithms with optimal computation complexities, but their oracle complexities are worse than the ideal ones by an order of magnitude. The current best method based on the augmented Lagrangian reformulation [23], a three-loop algorithm constructed from the inexact augmented Lagrangian method, the ellipsoid method, and the accelerated gradient descent (AGD) method, is closer to the ideal oracle complexities. However, the proposed method in [23] scales poorly with the number of constraints. In fact, the method is advantageous to the APD method only when $m \leq 5$ in the numerical experiments [23].

Another line of research [13, 20], called the level-set method, reformulates (<ref>) to a root finding problem associated with a certain mini-max problem. These methods typically require all constraint functions to share the same strong convexity modulus as the objective function. Such an assumption can be quite restrictive because it is violated when there exists one affine constraint. Assuming the uniform strong-convexity condition holds, the method proposed by Nesterov in Section 2.3.5 of [20] can achieve an $\bigO(\log(1/\ep))$ oracle complexity when $\alpha >0$. However, the method might be too computationally demanding to implement for the large-scale setting because it requires the exact solution to a quadratic program (QP) in each iteration, and the exact solution to a quadratic-constrained quadratic program (QCQP) from time to time. Lin et al. relaxes the expensive computation oracle assumption in [13], but the oracle complexity also becomes worse by an order of magnitude.  To sum up, there exist two major deficiencies: a) current methods fail to match the ideal oracle complexity, and b) methods that are close to the ideal oracle complexity are impractical in the large-scale setting.

Complexities for Smooth Constrained Optimization (<ref>)

2c|Strongly Convex $\alpha > 0$     2c|Convex $\alpha =0$     Strong Oracle &

Method     Oracle     Computation     Oracle      Computation      Assumption

Level Set1 [20]             $\bigO(\tfrac{1}{\sqrt{\ep}}\log(\tfrac{1}{\ep}))$          QCQP & QP Oracle

Level Set1 [13]             $\bigO(1/\ep)$     $\bigO(1/\ep)$

iALM [23]     $\bigO(\tfrac{m}{\sqrt{\ep}}\log^3(\tfrac{1}{\ep}))$
$\bigO(\tfrac{m}{\sqrt{\ep}}\log^3(\tfrac{1}{\ep}))$     $\bigO(m\log^3(\tfrac{1}{\ep}))$     $\bigO(m\log^3(\tfrac{1}{\ep}))$     $m=\bigO(1)$ constraints

APD [6]     $\bigO({1}/{\sqrt{\ep}})$     $\bigO({1}/{\sqrt{\ep}})$     $\bigO({1}/{\ep})$     $\bigO({1}/{\ep})$

ConEx [4]    $\bigO({1}/{\sqrt{\ep}})$     $\bigO({1}/{\sqrt{\ep}})$     $\bigO({1}/{\ep})$     $\bigO({1}/{\ep})$

ACGD [*]     $\bigO(\log(1/\ep))$          $\bigO(1/\sqrt{\ep})$         QP Oracle

ACGD-S [*]      $\bigO(\log(1/\ep))$     $\bigO(1/\sqrt{\ep})$     $\bigO(1/\sqrt{\ep})$     $\bigO(1/\ep)$

1 If, in addition to $f$, all constraint functions, $g_1,\ldots,g_m$, are $\alpha$-strongly convex, the oracle complexity is $\bigO(\log(1/\ep))$ for the level-set method in [20], and is  $\bigO(1/\sqrt{\ep})$ for [13].

In this paper, we provide an affirmative answer to the research question by proposing efficient algorithms to achieve the ideal complexities. An important observation for our development is that the ideal oracle complexity and the ideal computation complexity, shown in Table <ref>, have different orders of magnitudes. This explains why the single loop algorithms which carry out $\bigO(1)$ matrix-vector operations per oracle evaluation,  for example, the APD method [6] and the ConEx method [4], can only achieve the ideal computation complexity, but not the ideal oracle complexity. This also explains why the more complicated algorithms with better oracle complexities, for example, the level set method in [20] and the iALM method in [23], require stronger computation oracles than basic matrix-vector operations. Following this observation, we first ask if the ideal oracle complexity in Table <ref> is attainable with any strong computation oracle. The question leads us to the oracle-efficient Accelerated Constrained Gradient Descent (ACGD).

The ACGD method is based on the Lagrangian formulation (<ref>). The method can be motivated by how Nesterov adapts the AGD method to solve a minimax problem of the form $\min_{x}\max_{i}\{f_1(x), ..., f_m(x)\}$ in Section 2.3 of [20].  Since the Lagrangian function is also a minimax problem, it also satisfies the following max-type smoothness condition (see Lemma 2.3.1 in [20]):

\begin{equation}
\max_{\lam \in \R^m_+} \La(x, \lam) - \max_{\lam \in \R^m_+} \{l_f(x;\xbar) + \sumi \lami l_g(x; \xbar)\} \leq \tfrac{L(\R^m_+)}{2} \normsq{x -\xbar},
\end{equation}

where $l_f$ and $l_g$ denote the linearization of $f$ and $g$ at $\xbar$, and the aggregate smoothness constant $L(\R^m_+)$ denotes the maximum, over $\lambda \in \R^m_+$, of the Lipschitz smoothness constants of $\La(x; \lambda)$ for a fixed $\lambda$. To obtain the desired oracle complexity, Nesterov proposes to modify the AGD method by replacing the descent step  with a max-type descent step (see (2.3.12) in [20] or the DRAO method in [12]), which in our case becomes

\begin{align}
\xt \leftarrow \argmin_{x \in X} \max_{\lam \in \R^m_+}  l_f(x; \xundert) + \sum_{i=1}^m \lami (l_{g_i}(x; \xundert)) + u(x) + \tfrac{\etat}{2} \normsq{x - \xtt}, \label{eq:prox-linear}
% &\Leftrightarrow \xt \leftarrow \min_{x \in X} \{f(\xbar) + \inner{\grad f(\xbar)}{x - \xbar}  + u(x) + \tfrac{\etat}{2} \normsq{x - \xtt}, \text{s.t.}\   \lami (g_i(\xbar) + \inner{\grad g_i(\xbar)}{x - \xbar} \leq 0\ \forall i \in [m]\}, \notag
\end{align}

where the linearization center $\xundert$ is a certain convex combination of $\{\xt[0],\ldots, \xtt\}$. Since the maximization of $\lambda$ is over $\R^m_+$, (<ref>) is equivalent to  the following quadratic program,

$$\xt \leftarrow \min_{x \in X} \{l_f(x; \xundert)  + u(x) + \tfrac{\etat}{2} \normsq{x - \xtt} \text{s.t.}\   l_{g_i}(x; \xundert) \leq 0\ \forall i \in [m]\}.

We call it the constrained descent step and the so-modified AGD method the Accelerated Constrained Gradient Descent method.
Just like the AGD method, the stepsize parameter $\etat$ should be selected to be proportional to the aggregate smoothness constant $L(\R^m_+).$ Note, however, that one difficulty with this approach is $L(\R^m_+) = \infty$, so $\etat$ is infinite and the convergence rate is $\bigO(\infty/N^2)$, i.e., the method for  may not even converge. By taking a primal-dual perspective, we show $\etat$ only needs to be proportional to  the smoothness constant of $\La(x, \lamstar)$, denoted by $L(\lamstar)$. As a result, $\etat$ remains finite and the oracle complexity of the proposed ACGD method almost matches the optimal oracle complexity for optimizing $\La(x;\lamstar)$.

We also investigate the optimality of the proposed ACGD method by constructing new lower oracle complexity bounds. In the strongly convex case with $\alpha> 0$, the lower bound shows the oracle complexity bound of the ACGD method to be unimprovable for all problem parameters. In the non-strongly case with $\alpha =0$, the lower bound shows the oracle complexity bound of the ACGD method to be tight up to a factor of $\bigO(\sqrt{\norm{\lamstar} /c}+1)$, i.e., the ACGD method has tight oracle complexity with respect to all problem parameters if $c \geq \norm{\lamstar}$. If  $\lamstar$ is interpreted as the shadow price of constraint violation around $\xstar$, such a choice of $c$ ensures that the increase in objective value incurred by moving an $\ep/c$-feasible solution to feasibility is roughly $\bigO(\ep).$

To enhance its efficiency for the large-scale setting with many constraints, we use the sliding technique [7, 12] to extend the ACGD method to the ACGD with Sliding  (ACGD-S) method. Given a linearization center $\xundert$, instead of solving (<ref>) to optimality in each iteration, the ACGD-S method solves the bilinear saddle point inexactly by performing only a finite number of  $\lambda$-prox mappings and $x$-prox mappings (see (<ref>)), the most expensive operation during which is matrix-vector multiplication. Particularly, after the $t$th oracle evaluation, i.e., $\{\grad g(\xundert), \grad f(\xundert); f(\xundert), g(\xundert)\}$, the $x$-prox mapping and $\lambda$-prox mapping are repeated only $\bigO(t)$ times if $\alpha = 0$, and $\bigO(\sqrt{\theta^t})$ times for some $\theta > 1$ if $\alpha > 0$. The repetitive updating of prox-mappings and judicious computation of gradients allow the ACGD-S to achieve the optimal computation complexity (matching the lower bound for linearly constrained problem [21]), while maintaining the same oracle complexity as the ACGD method. Therefore the ACGD-S method provides an almost complete characterization of both the computation and the oracle complexity for solving a smooth function-constrained problem. Moreover, the intricate step-size choice to achieve both the optimal  $\bigO(\log(1/\ep))$ oracle evaluations (outer loops) and the optimal $\bigO(1/\sqrt{ \ep})$ matrix-vector multiplications (inner iterations) appears to be new to the sliding technique [7, 11, 12], so it is of independent interest.

One potential concern for implementing ACGD and ACGD-S is that their stepsize selection depends on the Lipschitz smoothness constant of $\La(x; \lamstar)$ (see (<ref>)), which could be hard to estimate in practice. To address this issue, we proposed a binary guess-and-check scheme to search for the elusive problem parameter.  Then, given the strong convexity parameter $\alpha$, under the assumption $X$ is bounded and its radius $\DX$ is known, both the ACGD and ACGD-S methods can be implemented without the knowledge of any other problem parameters to find an $(\ep; \ep/c)$-solution with the same oracle and computation complexities as those listed in Table <ref>.

The rest of the is organized as follows. Section <ref> proposes the ACGD method, and Section <ref> develops matching lower bounds. Section <ref> extends it to a computationally efficient ACGD-S method, and some concluding remarks are made in Section <ref>.

§.§ Notations & Assumptions

The following assumptions and notations will be used throughout the paper.

* The set of optimal solutions to (<ref>), $X^*$, is nonempty, and $\xstar$ is an arbitrary optimal solution. $F_*$ denotes the optimal objective, $F(\xstar)$.

* The Fenchel conjugate (see [1]) of a convex function g(x) is defined as $g^*(\pi) := \max_{x \in \R^n} \inner{x}{\pi} - g(x).$

* $U_h$ denotes the Bregman distance function (see [1]) generated by a convex function $h$, i.e., $U_h(\pi; \pibar):= h(\pi) - h(\pibar) - \inner{h'(\pibar)}{\pi - \pibar}$ where $h'(\pibar)$ is some fixed subgradient in $\partial h(\pibar).$ If $g$ is vector-valued, i.e.,  $g(\nu):=[g_1(\nu), g_2(\nu),\ldots,g_m(\nu)]$, $U_g$ is vector-valued with its $i$th component being the Bregman distance function generated by $g_i$, namely, $U_{g_i}$.

* We refer to the following computation as either a prox-mapping or a projection:

\begin{equation}\label{eq:prox-def}
\hat w \leftarrow \argmin_{w \in W} \inner{y}{w} + h(w) + \tau U(w;\bar w),\end{equation}

where the vector $y$ represents some “descent direction” (the gradient for example), and
$h(w)$ is a simple convex function [8]. $U$ is the Bregman distance function, $\bar w$ is a prox center, and $\tau$ is a stepsize parameter. Together they ensure the output $\hat w$ is close to $\bar w$. In particular, the following will be referred to as the $X$-projection:

\begin{equation}\label{eq:X-proj-def}
\hat x \leftarrow \argmin_{x \in X} \inner{y}{x} + u(x) + \tfrac{\tau}{2} \normsq{x - \xbar}.\end{equation}

§ THE ACCELERATED CONSTRAINED GRADIENT DESCENT METHOD

We present in this section the ACGD method. Specifically, Subsection <ref> introduces a novel primal-dual perspective to motivate the ACGD method, and Subsection <ref> presents the convergence results. Subsection <ref> proposes the guess-and-check scheme to look for the problem parameter $L(\lambda^*)$, and Subsection <ref> contains the detailed proofs of the convergence results.

§.§ The ACGD method

This subsection presents the primal-dual perspective to motivate the ACGD method. Such a perspective is important for understanding the finite stepsize parameter $\etat$ discussed in Section <ref>. Specifically, we first introduce a novel nested Lagrangian function which reformulates (<ref>) as a $\min-\max-\max$ trilinear saddle point problem. To search for the saddle point, we propose a primal-dual type method similar to [12, 26]. The ACGD method then simply follows from rewriting the proposed method in the primal form.

First, we need to assume the existence of a KKT point to the Lagrangian function defined in (<ref>) throughout this paper:

There exists a $\lamstar \in \R^m_+$ and $\xstar \in X$ such that $\grad f(\xstar) + \sumi \lamstar_i\grad g_i(\xstar)+ u' \in N_X(\xstar)$ for some $u' \in \partial u(\xstar)$, $g(\xstar) \leq 0$, and $\lami^* g_i(\xstar) = 0\ \forall i$, where $N_X$ denotes the normal cone to $X.$

Note  $\xstar$ and $\lamstar$ could be interpreted as an arbitrary element in $X^*$ and $\Lam^*:= \argmax_{\lam \in \R^m_+} \min_{x\in X} \La(x; \lambda)$ for the rest of the paper,
because any $(\xbar, \bar{\lam})$ with  $\xbar \in X^*$ and $\lambar \in \Lam^* $ also constitutes a KKT point (see Proposition 3.4.1 in [2]).

To motivate the nested Lagrangian function, consider the simplified case where the optimal dual multiplier $\lamstar$ is known. Fix $\lambda$ to $\lamstar$, an optimal solution $\xstar$ can be found by solving the following simplified problem under certain regularity conditions:

\begin{equation}\label{eq:unconstrained}
\min_{x \in X} f(x) + \inner{\lamstar}{g(x)} +  u(x).\end{equation}

One useful framework for designing an optimal algorithm for the constrained problem in (<ref>) is to consider a bilinear saddle-point reformulation

\begin{equation}\label{eq:fenchel}
\min_{x\in X} \max_{  \nu \in V, \pi \in \Pi} \{ \inner{x}{\pi} - f^*(\pi) + u(x) + \inner{\lamstar}{\nu x - g^*(\nu)} \},\end{equation}

where $f^*$ is the Fenchel conjugate function to $f$ and $\Pi$ is its domain, namely,  $\{\pi \in \R^n: f^*(\pi) < \infty \}$, and $g^*:=[g_1^*, \ldots, g_m^*]$ is the vector-valued Fenchel conjugate to $g$ and $V$ is its domain, namely, $\{v \in \R^{m\times n}: g_i^*(v_i) < \infty\ \forall i \in [m] \}.$

Moreover, since the dual variables $\pi$ and $\nu$ are associated with a common primal variable $x$, it is sometimes helpful consider the following joint domain $[V, \Pi]$:

\begin{equation}\label{eq:joint-domain}
[V, \Pi] = \{(\grad g(x), \grad f(x)): x\in \R^n\}.\end{equation}

However, $\lamstar$ is unknown in practice. So we propose to consider the following nested Lagrangian reformulation which combines (<ref>) and (<ref>):

\begin{equation}\label{eq:nLa}
\min_{x\in X} \max_{\lambda \in \R^m_+, \pi \in \Pi, \nu \in V} \{\La(x; \lambda, \pi, \nu):= \inner{x}{\pi} - f^*(\pi) + u(x) + \inner{\lambda}{\nu x - g^*(\nu)}\},\end{equation}

where $f^*$, $g^*$, $\Pi$ and $V$ are defined in the same way as (<ref>). Notice a common notation $\La$ is used for the nested Lagrangian and the ordinary Lagrangian, but the exact meaning should be clear from the context.
Let $Z$ denote the joint domain, $Z:= X \times \R^m_+ \times \Pi \times V$ and $\zstar:= (\xstar; \lamstar, \nu^* :=\grad g(\xstar), \pistar:=\grad f(\xstar))$. It is not hard to see that $\zstar$ is saddle point to (<ref>), and a useful criterion to measure the optimality of an iterate $z^t=(\xt; \lamt, \vt, \pit ) \in Z$ is to compare it to some reference point $z \in Z$ in the following gap function:

\begin{equation}\label{eq:Q_def}
Q(\zt; z):= \La(\xt; \lam, \nu, \pi) - \La(x; \lamt, \vt, \pit).
\end{equation}

Indeed, the saddle point $\zstar$ satisfies $Q(\zstar; z) \leq 0\ \forall z \in Z$.

A crucial observation for our development is that for the convergence to an $\ep$-solution, it is sufficient to consider only the reference $\lam$'s inside a certain bounded set rather than the positive orthant.

The next lemma shows the $Q$ function still provides upper bounds for both the feasibility violation and the optimality gap.

Let $\zt=(\xt; \lamt, \nut, \pit) \in Z$ be given and let $\Lamr$ denote a certain set of reference $\lambda$'s,

\begin{equation}\label{eq:ref_lambs}
\Lamr =\{0\} \cup \{\lambda^* + \lambda: \lambda \in B_{r}(0) \cap \R^m_+\}.\end{equation}

If $\max_{\lam \in \Lamr, (\nu,\pi) \in [V, \Pi]}Q(\zt; (\xstar; \lam, \nu, \pi)) \leq \ep$, we have $F(\xt) - F(\xstar) \leq \ep$ and $ \norm{[g(\xt)]_+} \leq \ep/r$.

Consider Part a). Fixing $\hpi = \grad f(\xt)$ and $\hnu = \grad g(\xt)$, the given condition implies that $Q(\zt; (\xstar; \lam, \hnu, \hpi)) \leq \ep\ \forall \lam \in \Lamr$. It then follows from the conjugate duality relationship (see Section 4.1 and 4.2 of [1]) that

$$f(\xt) + \inner{\lam}{g(\xt)} + u(\xt) - [f(\xstar) + \inner{\lamt}{g(\xstar)} + u(\xstar)] \leq \ep\ \forall \lam \in \Lamr.$$

Since $\xstar$ is feasible, i.e., $g(\xstar)\leq 0$, we have $\inner{\lamt}{g(\xstar)} \leq 0$. Then taking $\lam = 0$ leads to

$$f(\xt) + u(\xt) - f(\xstar) - u(\xstar)\leq f(\xt) + u(\xt)  - [f(\xstar) + \inner{\lamt}{g(\xstar)} + u(\xstar)] \leq \ep.$$

Next since $(\xstar; \lamstar)$ is a saddle point to (<ref>), we have
$0 \leq F(\xt) + \inner{\lamstar}{g(\xt)} - [F(\xstar) + \inner{\lamt}{g(\xstar)}]$. Setting $\hlam = \lamstar+ r [g(\xt)]_+/ \norm{[g(\xt)]_+} \in \Lamr$, we get

\begin{align*}
r\norm{[g(\xt)]_+} &\leq F(\xt) +  \inner{\lamstar+ r[g(\xt)]_+/ \norm{[g(\xt)]_+}}{g(\xt)} - [F(\xstar) + \inner{\lamt}{g(\xstar)}] \\
&\leq Q(\zt; (\xstar; \hlam, \hnu, \hpi)) \leq \ep.
\end{align*}

We remark that a lower bound to $F(\xt) - \Fstar$ can also be derived in similar fashion, see [9, 24].

Now, let us move on to consider minimizing $Q$ to find a saddle point to (<ref>). An essential feature of the nested Lagrangian function (<ref>) is the trilinear term $\inner{\lam}{\nu x - g^*(\nu)}$. The problem cannot be simplified to a $\min-\max$ saddle point problem by combining the $\nu$ and $\lam$ into a single dual block, because their joint maximization is difficult to compute. Similar problems have been studied by Lan and Zhang in [12, 26, 27] and the key to handle the $\min-\max-\max$ trilinear structure is to decompose the $Q$ function into the sub-gap functions and optimize each sub-gap function sequentially. Specifically, the following decomposition of (<ref>) into sub-gap functions associated with $x$, $\lam$, $\nu$ and $\pi$ is useful.

\begin{equation}\label{eq:Q_decomp}
Q(\zt; z) = Q_x(\zt; z) + Q_\lam(\zt; z) + Q_\pi(\zt; z) + Q_\nu(\zt; z)
\end{equation}

\begin{align}
% \begin{split}
Q_\pi(\zt; z)&:= \La(\xt; \lam, \nu, \pi) - \La(\xt; \lam, \nu, \pit) = \inner{\xt}{\pi} - f^*(\pi) \boxed{- [\inner{\xt}{\pit} - f^*(\pit)]}, \label{eq:Qpi}\\
Q_\nu(\zt; z)&:= \La(\xt; \lam, \nu, \pit) - \La(\xt; \lam, \nut, \pit) = \sumi \lami (\inner{\nui}{\xt} - g^*_i(\nui)] \boxed{- \sumi \lami[i] [\inner{\nuit}{\xt} - g_i^*(\nuit]}), \label{eq:Qnu}\\
Q_\lambda(\zt; z)&:= \La(\xt; \lam, \nut, \pit) - \La(\xt; \lamt, \nut, \pit) = \inner{\lam}{\nut \xt - g^*(\nut)} \boxed{- \inner{\lamt}{\nut \xt - g^*(\nut)}}, \label{eq:Qlam}\\
Q_x(\zt; z)&:= \La(\xt; \lamt, \nut, \pit) - \La(x; \lamt, \nut, \pit) = \boxed{\inner{\pit + \sumi \lamti \nuti}{ \xt} + u(\xt)} - \inner{\pit + \sumi \lamti \nuti}{ x} - u(x). \label{eq:Qx}
% \end{split}
\end{align}

Similar to the DRAO method in [12], we propose to reduce the boxed terms associated with $Q_\nu$ and $Q_\pi$, and $Q_x$ and $Q_\lam$ using the following sequence of prox-mappings in each iteration:

\begin{align}\label{ls:anla}
\begin{split}
\xtilt  &\leftarrow \xtt + \thetat (\xtt - \xt[t-2]);\\
\pit &\leftarrow \argmax_{\pi  \in \Pi} \inner{\pi}{\xtilt} - f^*(\pi) - \taut U_{f^*}(\pi; \pit[t-1]);\\
\nuti &\leftarrow  \argmax_{\nui \in V_i} \inner{\nui}{\xtilt} - g_i^*(\nui) - \taut U_{g_i^*}(\nui; \nuit[t-1])\ \forall i \in [m];\\
(\xt, \lamt) &\leftarrow \text{arg}\min_{x \in X} \max_{\lam \in \R^m_+} \inner{\pit}{x} + \inner{\lam}{\nut x - g^*(\nut)} + u(x) + \etat \normsq{x - \xtt}/2.\end{split}
\end{align}

In the above listing, $\thetat$, $\taut$ and $\etat$ are non-negative stepsize parameters, and $U_{g_i^*}$ and $U_{f^*}$ are Bregman distance functions generated by $g_i^*$ and $f^*$ respectively. Particularly, with the yet available $\xt$ replaced by some proxy $\xtilt$, the above $\pit$ update corresponds to the minimization of the variable $\pit$ in (<ref>) subject to a prox term $\taut U_{f^*}(\pi; \pit[t-1])$.
Similarly, the $\nut$ update corresponds to the minimization of the variable $\nut$ in (<ref>) subject to a prox term $\taut \sumi \lami U_{g_i^*}(\nui; \nuit[t-1])$. Since the summation weights $\{\lami\}$ are non-negative and the maximizations are separable, the $\nut$ update is written equivalently as individual $\vit$ updates in (<ref>). Finally, since both the maximization of $\lamt$ in (<ref>) and the minimization of $\xt$ in (<ref>) do not require any oracle information related to either $f$ or $g$, we evaluate them simultaneously subject to a prox term $\etat \normsq{x - \xtt}/2$. This leads to the joint $(\xt, \lamt)$ update in (<ref>).

We remark  $\xtilt$, the proxy for $\xt$, is chosen as the momentum extrapolation from $\xtt$. Such a choice helps us to achieve acceleration ( e.g. see Section 3.4 in [8] for the connection).

The implementable version of (<ref>) is shown in Algorithm <ref>. It employs two additional simplifications. First, we initialize the dual variables to some gradients, i.e.,  $\pit[0] = \grad f(\xt[0]) $ and $\nut[0] = \grad g(\xt[0])$. With $U_{g^*}$ and $U_{f^*}$ selected as prox-functions, we can show recursively that  $\pit$ and $\vt$ in (<ref>) are same as the gradients at some averaged point (see Lemma 2 in [27]). Thus the $\nuit$ and $\pit$ computation in (<ref>) simplifies to gradient evaluations in Line 4 of Algorithm <ref>. Second,  the $(\xt, \lamt)$-saddle point problem in (<ref>) is formulated as a linearly constrained quadratic program in Line 5 of Algorithm <ref>. Since  $\nut = \grad g(\xundert)$ in Algorithm <ref> implies the relation $g_i(\xundert) + g_i^*(\nuit) = \inner{\xundert}{\nuit}$ (see Theorem 4.20 in [1]), we have $g^*(\nut)=\nut\xundert - g(\xundert)$. Interestingly, other than the additional linear constraint associated with $g$ in the descent step (Line 5), Algorithm <ref> is the same as Nesterov's AGD method [18]. Therefore we name it the Accelerated Constrained Gradient Descent method.
Accelerated Cconstrained Gradient Descent Method

[1]
$\xt[-1]= \xundert[0] = \xt[0]\in X$,  stepsizes $\{\thetat\}$, $\{\etat\}$, $\{\taut\}$, and  weights $\{\wt\}$.
Set $\pit[0] = \grad f(\xt[0]) $, $\nut = \grad g(\xt[0])$.
$t =  1,2,3 ... N$
Set $\xundert \leftarrow (\taut \xundert[t-1] + \xtilt)/(1 + \taut)$ where  $\xtilt  = \xtt + \thetat (\xtt - \xt[t-2])$.
Set $\pit \leftarrow \grad f(\xundert)$ and $\nut \leftarrow \grad g(\xundert)$.
Solve $\xt \leftarrow \argmin_{x \in X} \{\inner{\pit}{x} +u(x) + \etat \normsq{x - \xtt}/2\ s.t.\ \nut (x - \xundert) + g(\xundert) \leq 0\} .$
$\xbart[N]:= \tsum_{t=1}^N \wt \xt / (\sumt \wt). $

§.§ The Convergence Results

Next, we present convergence results for the ACGD method. The detailed analysis is deferred to Subsection <ref>. The next proposition states some conditions required for the $Q$ function in (<ref>) to converge.
Since these conditions depend on the set of reference multiplier $\lambda$'s under consideration (see Lemma <ref>), it is useful to define an aggregate Lipschitz smoothness constant as a function of $\Lambda$:

\begin{equation}\label{eq:agg-L}
L(\Lambda):= \max_{\lam \in \Lambda} L_\lambda, \text{ where } L_\lambda \norm{\xbar - \hx} \geq \norm{\grad_x \La(\xbar; \lambda) - \grad_x \La(\hx; \lambda)} \forall \xbar, \hx \in \R^n .
\end{equation}

Let a set of reference multipliers $\Lambda \subset \R^m_+$ be given and let the aggregate smoothness constant $\Llam$ be defined in (<ref>). Consider the iterates ${\zt:=(\xt; \lamt, \nut, \pit)}$ generated by Algorithm <ref>, where $\lamt$ is defined in (<ref>). Suppose the following conditions are satisfied by the stepsizes together with some non-negative weights $\wt \geq 0$ for all $t \geq 2$:

\begin{align}
% \begin{split}
& \wt \etat \leq \wtt (\etatt + \alpha), \label{stp-con:eta-recur}\\
&\wt \taut \leq \wtt (\tautt + 1), \label{stp-con:Q-recur}\\
&\etat[N](\taut[N] + 1) \geq \Llam, \etatt \taut \geq \theta_t\Llam \text{ with }  \thetat :=\wtt /\wt\label{stp-con:Q-cancel}.
% \end{split}
\end{align}

Consider the ergodic iterate $\zbart[N]= (\xbarN; \lambart[N], \nubarti[N][], \pibart[N])$ specified according to

\begin{equation}\label{eq:ergo_solns}
\begin{split}
&\xbarN := \sumt \wt \xt / (\sumwt),\ \lambart[N]:= \sumwt \lamt/ (\sumwt),\ \pibart[N]:= \sumwt \pit / (\sumwt),\\
&\nutibar[N]:=  \begin{cases}\sumwt \lamti \nuit /(\sumwt \lamti) & \text{o.w.} \\
\grad g_i(x^0) &\text{if } \lamti = 0\ \forall t. \end{cases}
\end{split}
\end{equation}

Then the following convergence bound is valid for any reference point $z=(x; \lam,  \nu, \pi) \in X \times \Lam \times [V, \Pi]$ with $[V, \Pi]$ defined in (<ref>):

\begin{equation}\label{eq:Q-conv}
\begin{split}
(\sumwt)& Q(\zbart, z) + \wt[N](\etat[N] + \alpha) \normsq{\xt[N]-x}/2 \\
\leq& \wt[1]\etat[1] \normsq{\xt[0] - x}/2 + \wt[1]\taut[1](U_f^*(\pi; \pit[0]) + \inner{\lambda}{U_g^*(\nu; \nut[0])}).
\end{split}\end{equation}

Proposition <ref>, together with Lemma <ref>, show that it is possible to select a finite stepsize $\etat$ related to $L(\Lamr)$ in place of the infinite $L(\R^m_+)$.
To provide convergence guarantees for both the feasibility violation and the optimality gap associated with the function-constrained problem in  (<ref>),
Lemma <ref> states that we should only consider the $Q$ gap function defined with respect to $\Lamr$, some small neighborhood of reference $\lambda$'s around $\lamstar$.  Proposition <ref> shows such an $Q$ function converges when the stepsize choice satisfies certain conditions related to $L(\Lamr)$, rather than $L(\R^m_+)$. Thus a finite $\etat$ proportional to $L(\Lamr)$ is sufficient for our purpose, and the next theorem states the result more precisely.

Let a smooth constrained optimization problem (<ref>) be given and let $\Llamr$ be defined in  (<ref>) and (<ref>). Denote its condition number by $\kar= \Llamr/\alpha$ (We set $\kar = \infty$ if $\alpha=0$).  Suppose the solution iterates $\{\xt\}$ are generated by Algorithm <ref> with the following stepsizes for $t \geq 1$

\begin{equation}\label{stp-choice:ANLA}
\taut = \min\{\tfrac{t-1}{2}, \sqrt{\kar}\},\ \etat = \tfrac{\Llamr}{\taut[t+1]},\ \thetat = \tfrac{\taut}{\taut[t-1] + 1},\ \wt = \begin{cases} \wt[t-1] /\thetat & \text{if } t  \geq 2,\\
1 &\text{if } t = 1.
\end{cases}
\end{equation}

Then the ergodic average $\xbarN$ solution satisfies

\begin{equation}\label{eq:thm-non-strong}
\max\{F(\xbarN) - \Fstar, r\norm{[g(\xbarN)]_+}\} \leq \frac{2\Llamr}{N(N+1)} \normsq{\xt[0] - \xstar}.
\end{equation}

Moreover, in the strongly convex case with $\alpha >0 $, $\xbarN$ also satisfies

\begin{align}
&\max\{F(\xbarN) - \Fstar, r\norm{[g(\xbarN)]_+}\} \leq \frac{ \sqrt{\Llamr \alpha} \normsq{\xt[0] - \xstar}}{(1+ 1/\sqrt{\kar})^{N-4} - 1}, \label{eq:thm-strong}\\
&\normsq{\xbarN - \xstar} \leq \frac{2 \sqrt{\kar} \normsq{\xt[0] - \xstar}}{(1+ 1/\sqrt{\kar})^{N-4} - 1}. \label{eq:thm-strong-dist}
\end{align}

As a consequence of the preceding theorem, we can derive upper bounds on the number of iterations required of the ACGD method to find an $(\ep; \ep/c)$-solution (see (<ref>)). The next corollary focuses on the non-strongly convex case, i.e., $\alpha=0$.

Suppose $\{\xt\}$ are generated by Algorithm <ref> using the stepsizes choice in (<ref>), with $\Llamr=\Llamr[c]$ (see (<ref>)). The ergodic average solution $\xbarN$ is an $(\ep, \ep/c)$-solution if

$$N \geq \sqrt{\tfrac{2\Llamr[c]}{\ep}} \norm{\xt[0]- \xstar}.$$

For the strongly convex problem with $\alpha > 0$, it is advantageous to choose a small $r$ to ensure a small aggregate smoothness constant $\Llamr$ and hence a small condition number $\kar$. We set it to $\Llamr[1]$ in the next corollary.

Suppose $\{\xt\}$ are generated by Algorithm <ref> using the stepsizes choice in (<ref>), with $\Llamr=\Llamr[1]$ (see (<ref>)). Then the ergodic average solution $\xbarN$ is an $(\ep, \ep/c)$-solution if

$$N \geq  \min\{\sqrt{\tfrac{2\max\{c, 1\} \Llamr[1]}{\ep}} \norm{\xt[0] - \xstar}, [\sqrt{\tfrac{\Llamr[1]}{\alpha}} + 1]\log[\tfrac{\max\{c,1\}\sqrt{\Llamr[1]\alpha}\normsq{\xt[0] - \xstar}}{\ep} + 1] + 4\}.$$

Note that for a small $\lamstar$, i.e., $\norm{\lamstar}$ being significantly less than $1$, it may be worthwhile to choose $\Llamr = \Llamr[\norm{\lamstar}]$ to further reduce the condition number from $\kar[1] = \Llamr[1]/\alpha$ to $\kar[\norm{\lamstar}] = \Llamr[\norm{\lamstar}]/\alpha$.

Since each iteration of the ACGD method requires only one gradient evaluation, the preceding two corollaries establish the desired $\bigO(1/\sqrt{\ep})$ and $\bigO(\log(1/\ep))$ oracle complexities for the non-strongly convex and the strongly convex problems, respectively.

§.§ The Binary Search for ${\Llamr}$

A crucial limitation to the ACGD method is that the desired aggregate smoothness constants, $\Llamr[c]$ and $\Llamr[1]$, may be unavailable for the stepsize calculation in (<ref>) during implementation. To address this issue, we propose a binary search outer-loop method to look for the suitable constant in an online fashion. Starting with a guess $\Ltil$, we run the ACGD method for $N(\Ltil)$ iterations and check if the ergodic average solution $\xbart[N(\Ltil)]$ satisfies a certain termination criterion. The outer loop would be stopped if the termination criterion is met, otherwise the process would be repeated with the guess doubled, i.e., $2\Ltil$. To ensure the soundness of the search routine, we need to choose the iteration limit function $N(\Ltil)$ to guarantee the satisfaction of the termination criterion when $\Ltil$ were larger than the desired aggregate smoothness constant. Additionally, the termination criterion should also be both easily computable and valid, i.e., a test solution $\xbart[N]$ which satisfies the criterion should be $(\ep; \ep/c)$-optimal. Toward these ends, we need the following additional assumptions.

a) $X$ is bounded and an estimate $\DX$ of its radius, i.e., $\DX \geq \max_{x,y \in X} \norm{x - y},$ is known

b) The dual multiplier to the quadratic program in Line 5 of Algorithm <ref>, $\lamt$, is available.

Given iterates $\{(\xt, \xundert), \lamt,    \nut, \pit\}_{t=1}^N$ generated by the ACGD method, our termination criterion is the two conditions:

\begin{equation}\label{cr:check-ANLA}
\norm{[g(\xbarN)]_+} \leq c \ep \text{ and } F(\xbarN) - \Fundert \leq \ep,
\end{equation}

where $\xbarN$ is the ergodic average solution and $\fundert$ is given by

\begin{align}\label{eq:lower-obj}
\begin{split}
\Fundert:= \argmin_{x \in X}&\    \tfrac{1}{\sumwt} \sumt\wt[t][\inner{\pit}{x- \xundert} + f(\xundert) + u(x))]\\
s.t.&\  \tfrac{1}{\sumwt\lamti} \sumt\wt[t]\lamti[t][i] [\inner{\nuit}{x- \xundert} + g_i(\xundert)] \leq 0\ \forall i \in [m].
\end{split}
\end{align}

Notice  Assumption <ref>.a)  ensures $\Fundert$ to be finite, and the dual solutions in Assumption <ref>.b) are used for constructing the linearized constraints above.
Since its objective and constraints are lower linear approximations to $f$ and $g_1, g_2,\ldots,g_m$ respectively, (<ref>) is a relaxation to (<ref>), and the minimum value $\Fundert$ is a lower bound to $F_*$. Therefore the satisfaction of (<ref>) guarantees the solution $\xbarN$ to be $(\ep, \ep/c)$-optimal.

Binary Search for the Aggregate Smoothness Constant $L$.

[1]
$\xbar \in X$, $\Ltil> 0$ , and an iteration limit function  $N(L)$.

Set stepsizes to (<ref>) with $\Llamr$ replaced by $\Ltil$.
Run the ACGD method for $N(\Ltil)$ iterations starting from $\xbar$.
Set $\xbar \leftarrow \xbart[N(\Ltil)]$, where $\xbart[N(\Ltil)]$ is the ergodic average solution.
$\xbar$ satisfies the termination criterion (<ref>)
Set $\Ltil \leftarrow 2 \Ltil$.

Putting the pieces together, the binary search scheme is listed in Algorithm <ref>. When it terminates, the output solution is  $(\ep, \ep/c)$-optimal. Now we suggest some concrete iteration limit functions, $N(L)$, and derive upper bounds on the total numbers of iterations, across all ACGD runs, required for termination. Their proof is deferred to Subsection <ref>.

Consider a smooth constrained optimization problem of form (<ref>). Suppose its aggregate smoothness constants, $\Llamr[c]$ and $\Llamr[1]$, are defined according to (<ref>) and (<ref>), and Assumption <ref> is valid. Let $\Ltil$ be initialized to some $\Ltil_0$, then Algorithm <ref> terminate finitely for the following choices of the iteration limit function $N(L)$.

a) In the non-strongly convex case with $\alpha=0$, if
$N(L)= \ceil{\sqrt{\tfrac{2L}{\ep}}\DX},$
the total number of ACGD  iterations required for termination is upper bounded by

\begin{equation}\label{eq:bin-search-upper-complexity-non-strong}
\max\{7 \sqrt{\tfrac{\Llamr[c]}{\alpha}}\Dx+ \ceil{\log_2(\tfrac{\Llamr[c]}{\Ltil_0})}, \sqrt{\tfrac{2\Ltil_0}{\ep}}\DX + 1\}.\end{equation}

b) In the strongly convex case with $\alpha>0$, if

$$N(L)= \ceil{(\sqrt{\tfrac{L}{\alpha}} + 1)\log(\tfrac{\max\{c, 1\}\sqrt{L \alpha} \DX^2}{\ep} + 1)} + 4,$$

the total number of ACGD  iterations required for termination is upper bounded by

\begin{align}\label{eq:bin-search-upper-complexity-strong}\begin{split}
\max\left\{9 \sqrt{\tfrac{\Llamr[1]}{\alpha}}\log(\tfrac{\max\{c, 1\}\sqrt{2 \Llamr[1] \alpha} \DX^2}{\ep} + 1) + 5 \ceil{\log_2 \tfrac{\Llamr[1]}{\Ltil_0}}, \right.\\
\left. \ceil{(\sqrt{\tfrac{\Ltil_0}{\alpha}} + 1)\log(\tfrac{\max\{c, 1\}\sqrt{\Ltil_0 \alpha} \DX^2}{\ep} + 1)} + 4\right\}.\end{split}
\end{align}

We make two remarks regarding the result. First, the proposed $N(L)$'s depend only on the target accuracy $(\ep; \ep/c)$ and some easy-to-estimate problem parameters, that is, $\DX$ in the non-strongly convex case, and $\DX$ and $\alpha$ in the strongly convex case. So Algorithm <ref> is easy to implement. Second, the total numbers of ACGD iterations required by Algorithm <ref>, (<ref>) and (<ref>), are nearly independent of the initial guess $\Ltil_0$ for all $\Ltil_0 \leq \Llamr[0]$. The FO oracle complexities of Algorithm <ref> have the same orders as those in Corollary <ref> and <ref>, so the total cost of the multiple ACGD runs in Algorithm <ref> is on the same order as a single ACGD run with some carefully selected stepsizes.

§.§ Convergence Analysis

We present the detailed proofs to convergence results of Algorithm <ref> and Algorithm <ref>.

Proof of Proposition <ref>
Let $Q_x$, $Q_\lam$, $Q_\nu$, and $Q_\pi$ be defined in (<ref>). It is useful to view the updates in Algorithm <ref> from the perspective of prox-mappings in (<ref>). First, let's consider $Q_\nu$. We have from the definition of $\xtilt$ that

\begin{align*}
\inner{\nui -\nuit }{(\xtilt-\xt)}=& -\inner{\nui - \nuit}{(\xt -\xtt)} + \thetat \inner{\nui -\nuitt}{(\xtt -\xt[t-2])}\\
&+\thetat \inner{\nuitt -\nuit}{(\xtt -\xt[t-2])}.
\end{align*}

Since $g^*_i$ has a strong convexity modulus of $1$ with respect to $U_{g_i^*}$, the definition of $\nuit$ prox-mapping in (<ref>) implies a three-point inequality (see Lemma 3.5 in [8]):

\begin{align*}
\inner{\nui - \nuit}{ \xt}& + \gistar(\nuit) - \gistar(\nui) + \inner{\nui -\nuit }{\xtilt-\xt} \\
&\leq \taut U_{g_i^*}(\nui; \nuitt) - (\taut + 1) U_{g_i^*}(\nui; \nuit) - \taut U_{g_i^*}(\nuit; \nuitt).
\end{align*}

So, combining the above two relations,
taking the $\wt$ weighted sum of the resulting inequalities and using the conditions $\wtt=\wt\thetat$ and $\wt\taut \leq \wtt (\tautt  + \mu)$, we obtain

\begin{align*}
\sumt \wt &(\inner{\nui - \nuit}{ \xt} + \gistar(\nuit) - \gistar(\nui)) \\
\leq& - (\wt[N](\taut[N] + 1) U_{g_i^*}(\nui; \nuit[N]) -\wt[N]\inner{\nui - \nuit[N]}{ (\xt[N] -\xt[N-1])})\\
& -\tsum_{t=2}^{N} [\wt \taut U_{g_i^*}(\nuit; \nuitt) + \wtt \inner{\nuitt -\nuit}{(\xtt -\xt[t-2])} ]\\
&+ \wt[1]\taut[1] U_{g_i^*}(\nui; \nuit[0]).
\end{align*}

A $\lami$-weighted sum of the above inequality leads to the desired $Q_\nu$ convergence bound given by

\begin{align*}
\begin{split}
\sumt \wt& Q_{\nu}(\zt;z) \\
\leq& - (\wt[N](\taut[N] + 1) \sumi \lami U_{g_i^*}(\nui; \nuit[N]) -\wt[N]\inner{\sumi \lami (\nui - \nuit[N])}{\xt[N] -\xt[N-1]})\\
& -\tsum_{t=2}^N [\wt \taut\sumi \lami U_{g_i^*}(\nuit; \nuitt) + \wtt \inner{\sumi \lami (\nuitt -\nuit)}{\xtt -\xt[t-2]} ]\\
&+ \wt[1]\taut[1] \sumi \lami U_{g_i^*}(\nui; \nuit[0]) .
\end{split}
\end{align*}

A $Q_\pi$ bound can be derived similarly. Taken together, we get

\begin{align}\label{pr1:Qnu-Qpi}
\begin{split}
&\sumt \wt[t] [Q_{\nu}(\zt;z) + Q_\pi(\zt; z)] \\
&\ \leq - (\wt[N](\taut[N] + 1) [\sumi \lami U_{g_i^*}(\nui; \nuit) + \Ufstar(\pi; \pit)] -\wt[N]\inner{\pi - \pit[N] + \sumi \lami (\nui - \nuit[N])}{\xt[N] -\xt[N-1]})\\
&\quad \ -\tsum_{t=2}^N \{\wt \taut[t][\sumi \lami U_{g_i^*}(\nuit; \nuitt) + \Ufstar(\pit; \pitt)] + \wtt \inner{\pit[t-1] - \pit + \sumi \lami (\nuitt -\nuit)}{\xtt -\xt[t-2]} \}\\
&\quad \ + \wt[1]\taut[1] [\sumi \lami U_{g_i^*}(\nui; \nuit[0]) + \Ufstar(\pi; \pit[0])] .
\end{split}
\end{align}

Since $(\nut, \pit) = (\grad g(\xundert), \grad f(\xundert))\ \forall t$, applying Lemma <ref> with $\tilde{\lam}= [\lam; 1]$ and $\tilde{g} = [g; f]$ implies that

$$\sumi \lami U_{g_i^*}(\nuit; \nuitt) + \Ufstar(\pit; \pitt) \geq \tfrac{1}{2\Llam} \normsq{\pit - \pitt + \sumi \lami (\nuit - \nuitt)}.$$

Similarly, with $(\nu, \pi) \in [V, \Pi]$ (see (<ref>)), we get

$$\sumi \lami U_{g_i^*}(\nu; \nuit[N]) + \Ufstar(\pi; \pit[N]) \geq \tfrac{1}{2\Llam} \normsq{\pi - \pit[N] + \sumi \lami (\nui - \nuit[N])}.$$

Thus applying the Young's inequality to (<ref>) leads to

\begin{align}\label{eq:pr-pf-Qnu}
\begin{split}
\sumt& \wt[t] [Q_{\nu}(\zt;z) + Q_\pi(\zt; z)] \\
\leq&  \tfrac{\wt[N]\Llam }{2(\taut[N] + 1)} \normsq{\xt[N] - \xt[N-1]} + \tsum_{t=1}^{N-1} \tfrac{\wt\thetat \Llam }{2\taut} \normsq{\xt - \xtt}
+ \wt[1]\taut[1] [\sumi \lami U_{g_i^*}(\nui; \nuit[0]) + \Ufstar(\pi; \pit[0])] .
\end{split}
\end{align}

Now let us move onto $Q_x$ and $Q_\lam$. Fix $\lamt$. The $\xt$-prox mapping implies a three-point inequality (Lemma 3.5 in [8]):

$$\inner{\pit + \sumi \lamti\nuti}{\xt - x} + u(\xt) - u(x) +  \tfrac{\etat + \alpha}{2} \normsq{\xt - \hx} \leq \tfrac{\etat}{2} \normsq{\xtt - \hx} - \tfrac{\etat}{2} \normsq{\xt - \xtt}.$$

Fix $\xt$. The optimality of $\lamt$ implies that:

$$\inner{\lam - \lamt}{\nut \xt - g^*(\nut)} \leq 0.$$

So we have

$$Q_x(\zt; z) + Q_\lam(\zt; z) +  \tfrac{\etat + \alpha}{2} \normsq{\xt - x} \leq \tfrac{\etat}{2} \normsq{\xtt - x} - \tfrac{\etat}{2} \normsq{\xt - \xtt}.$$

Summing across iterations with weight $\wt$ leads to

\begin{equation}\label{eq:pr-pf-Qx}
\sumwt[t] [Q_x(\zt; z) + Q_\lam(\zt; z)] +  \tfrac{\wt[N](\etat[N] + \alpha)}{2} \normsq{\xt - x} \leq \tfrac{\wt[1]\etat[1]}{2} \normsq{\xt[0] - x} - \sumt \tfrac{\wt\etat}{2} \normsq{\xt - \xtt}.
\end{equation}

Utilizing the stepsize assumption (<ref>), we can add it to (<ref>) to obtain a convergence bound for the $Q$ function

\begin{equation}\label{eq:pr-pf-Q}
\sumwt Q(\zt; z) + \tfrac{\wt[N](\etat[N] + \alpha)}{2} \normsq{\xt - x} \leq \tfrac{\wt[1]\etat[1]}{2} \normsq{\xt[0] - x} + \wt[1]\taut[1] [\sumi \lami U_{g_i^*}(\nui; \nuit[0]) + \Ufstar(\pi; \pit[0])].
\end{equation}

Moreover, the Jensen's inequality implies that

\begin{equation}\label{pr1:jensen}
(\sumwt)\La(\xbart[N]; \lam,\nu, \pi) \leq \sumwt \La(\xt; \lam, \nu, \pi),\end{equation}

\begin{align*}
\sumwt& \La(x; \lamt, \nut, \pit) \geq (\sumwt)[\inner{x}{\pibart[N]} - f^*(\pibart[N])+ u(x)] + \sumi(\sumwt\lami)[ \inner{\nubarti}{x}-g_i^*(\nubarti[N])] \\
\geq& \sumwt (\La(x; \lambart[N], \nubarti[N][],\pibart[N] )).
\end{align*}

Thus, we get $(\sumwt )Q(\zbart[N]; z) \leq \sumwt Q(\zt; z)$, and the desired inequality in (<ref>) follows from (<ref>).
Next, the proof of Theorem <ref> is a direct application of Proposition <ref>. The analysis is complicated by the switch from a diminishing stepsize to a constant stepsize in (<ref>).

Proof of Theorem <ref>
We apply Proposition <ref> to obtain the results. First, we verify that the requirements in (<ref>) and (<ref>) are satisfied by the stepsize choice in (<ref>). Since $\thetat=\taut/(\taut[t-1]+1)$ and $\etat = \Llamr/\taut[t+1]$, all requirements other than (<ref>) hold automatically.  Now
let $T = \ceil{2 \sqrt{\kar}} + 1$ denote the first iteration at which we switch to $\taut=\sqrt{\kar}.$ Since the other iterations and the case with $\alpha = 0$ are straightforward to check, we focus on iteration $T-1$ and $T$, and assume $\alpha >0$. For $t=T$, we have

\begin{align*}
\wt[T]\etat[T] = \wt[T-1] \etat[T] \tfrac{\taut[T-1] + 1}{\taut[T]} = \wt[T-1] \tfrac{\Llamr}{\sqrt{\kar}} \tfrac{\ceil{2\sqrt{\kar} } + 1}{2\sqrt{\kar}} = \wt[T-1] \alpha \tfrac{2\sqrt{\kar} + 2}{2} \leq \wt[T-1] (\alpha\sqrt{\kar} + \alpha) \leq \wt[T-1] (\etat[T-1] + \alpha).
\end{align*}

For $t=T-1$, we have

\begin{align*}
\wt[T-1]\etat[T-1] = \wt[T-2] \etat[T-1] \tfrac{\taut[T-2] + 1}{\taut[T-1]} = \wt[T-2] \tfrac{2\Llamr}{2\sqrt{\kar}} \tfrac{\ceil{2\sqrt{\kar}} -1}{\ceil{2\sqrt{\kar}} -2} \leq \wt[T-2] \tfrac{2 \Llamr}{\ceil{2\sqrt{\kar}} -2} = \wt[T-2] \etat[T-2] \leq \wt[T-2](\etat[T-2] + \alpha).
\end{align*}

Thus the requirements in Proposition <ref> are satisfied, and we have

\begin{equation}\label{thm-pf:ANLA-Q}
Q(\zbart, z)  \leq \Llamr \normsq{\xt[0] - x}/(\sumwt)\ \forall z=(x; \lam, \nu, \pi)\in X \times \Lamr \times [V, \Pi].
\end{equation}

We provide a lower bound to $(\sumwt)$. It is useful to show

$$\wt \geq \max\{t, (1 + 1/\sqrt{\kar})^{t -5}\}.$$

Since $\wt= 1/(\prod_{t=2}^{N} \thetat)$, the fact $\wt \geq t\ \forall t \geq 1$ is straightforward. Regarding the second lower bound, let us  first consider $t \leq T = \ceil{2 \sqrt{\kar}} +1.$ The algebraic fact in Lemma <ref> implies that
$\wt[t+2] \geq t \geq (1 + 1/\sqrt{\kar})^{t-3} \forall t \leq 2 \sqrt{\kar}$, so
$\wt \geq (1 + 1/\sqrt{\kar})^{t-5} \forall t \leq T = \ceil{2 \sqrt{\kar}} + 1$. For $t\geq T+1$, the relation is also valid since $1/\thetat = (1 + 1/\sqrt{\kar}).$ Therefore we get $\wt \geq \max\{t, (1 + 1/\sqrt{\kar})^{t -5}\} \forall t \geq 1$. Using these two lower bounds, it is easy to derive

\begin{equation}\label{eq:wt_lower_bound}
\sumwt \geq \max\{N(N+1)/2, \sqrt{\kar}[(1+ 1/\sqrt{\kar})^{N-4} -1]\}.\end{equation}

Substituting the preceding lower bound into (<ref>) and applying Lemma <ref> lead us to the convergence results in (<ref>) and (<ref>).

Now we deduce the convergence bound for $\norm{\xt - \xstar}$. Choosing $\hz = (\xstar; \lamstar, \grad g(\xbarN), \grad f(\xbart[N]))$ to be the reference point $z$ leads to

\begin{equation}\label{thm-pf:ANLA-Q-dist}
Q(\zbart; \hz) \leq \Llamr \normsq{\xt[0] - \xstar}/(\sumwt).\end{equation}

Moreover, since $\xstar$ minimizes the optimal Lagrangian, we get  $\inner{\grad f(\xstar) + (\lamstar)^\top \grad g(\xstar) + u'(\xstar)}{\xbarN - \xstar} \geq 0$ for all $u'(\xstar) \in \partial u(\xstar).$ So we have

\begin{align}\label{thmpf:Q-strong}
\begin{split}
Q(\zbart; \hz) &\geq f(\xbarN) + (\lamstar)^\top g(\xbarN) + u(\xbarN) - [f(\xstar) + (\lamstar)^\top g(\xstar) + u(\xstar)]\\
&\quad - \inner{\grad f(\xstar) + (\lamstar)^\top \grad g(\xstar) + u'(\xstar)}{\xbarN - \xstar} \\
&= f(\xbarN) - f(\xstar) - \inner{\grad f(\xstar)}{\xbarN - \xstar} + \sumi \lamstar_i [g_i(\xbarN) - g_i(\xstar)- \inner{\grad g_i(\xstar)}{\xbarN - \xstar}] \\
&\quad + [u(\xbarN) - u(\xstar) - \inner{u'(\xstar)}{\xbarN - \xstar}]\\
&\geq \alpha \normsq{\xbarN - \xstar}/2,\end{split}
\end{align}

where the last inequality follows from the fact that the regularization function $u(x)$ is $\alpha$-strongly convex.
Combining the above two relations, we obtain the desired inequality in (<ref>):

$$\normsq{\xbarN - \xstar} \leq 2 \kar \normsq{\xt[0] - x}/(\sumwt)\leq \tfrac{2 \sqrt{\kar} \normsq{\xt[0] - \xstar}}{(1+ 1/\sqrt{\kar})^{N-4} - 1}.$$

Next, we move on to provide proofs to the binary search outer-loop in Algorithm <ref>.

Proof of Theorem <ref>:
First, let's consider the non-strongly problem. We show the soundness of the  Algorithm <ref>, that is, the output ergodic average solution from the ACGD method satisfies termination criterion (<ref>) if the estimate $\Ltil$ is larger than the desired aggregate smoothness constant $\Llamr[c]$. Let $\Ntil:= N(\Ltil)= \ceil{\Dx\sqrt{2\Ltil/\ep}}$ such that  $\xbart[\tilde N]$ denotes the output solution. Since $\DX \geq \norm{\xt[0] - \xstar}$, $\|[g(\xbart[\Ntil])]_+\| \leq  \ep/c$ is a direct consequence of Corollary <ref>. Now consider the objective criterion. Since $\Ltil \geq \Llamr[c]$, it follows from (<ref>) that

$$\sumwt Q(\zt; z)  \leq {\wt[1]\etat[1]} \normsq{\xt[0] - x}/2 \leq \wt[1]\etat[1]\DX^2/2, \forall z = (x; \lam, \nu, \pi) \in X \times \Llamr[c] \times [V , \Pi].$$

Substituting the stepsize choice in (<ref>) then leads to

$$\tfrac{\sumT[\Ntil] \wt Q(\zt, z)}{\sumT[\Ntil] \wt} \leq \tfrac{\Ltil \DX^2}{\sumT[\Ntil] \wt } \leq \ep, \forall z = (x; \lam, \nu, \pi) \in X \times \Llamr[c] \times [V , \Pi].$$

In particular, fixing $\hpi = \grad f(\xbart[\Ntil])$ and $\hnu = \grad g(\xbart[\Ntil])$, we have

$$F(\xbart[\Ntil]) \leq (\sumT[\Ntil] \wt \La(\xt; 0, \hnu, \hpi))(\sumT[\Ntil] \wt).$$

Moreover, since

\begin{align*}
\sumT[\Ntil] \tfrac{\wt}{\sumT[\Ntil] \wt} \La(x; \lamt, \nut, \pit) =&
\tfrac{1}{\sumwNtil} \sumt[\Ntil] \wt[t] [\inner{\pit}{x- \xundert} + f(\xundert) + u(\xt)]\\
&+ \sum_{i=1}^m[\frac{\sumwNtil \lamti}{\sumwNtil}] \frac{\sumwNtil \lamti[t][i] [\inner{\nuit}{x- \xundert} + g_i(\xundert)]}{\sumwNtil \lamti},
\end{align*}

the function is a Lagrangian relaxation to the linearly constrained problem in (<ref>), so

$$\min_{x \in X} \tfrac{\wt}{\sumT[\Ntil] \wt} \La(x; \lamt, \nut, \pit) \leq \fundert[\Ntil].$$

Putting these facts together, we get

$$F(\xbart[\Ntil]) - \Fundert[\Ntil] \leq \max_{z \in X \times \Llamr[c] \times [V, \Pi]}\tfrac{\sumT[\Ntil] \wt Q(\zt, z)}{\sumT[\Ntil] \wt} \leq \ep.$$

Thus $\xbart[\Ntil]$ must satisfy the termination criterion in (<ref>) if $\Ltil \geq \Llamr[c]$, and
the total number of ACGD iterations follows from a straightforward algebraic calculation.

Next, the finite termination of Algorithm <ref> for the strongly convex case can be deduced similarly. With the given choice of $N(L)$, we need to show the termination criterion (<ref>) is satisfied for the ergodic solution $\xbart[\Ntil]$ if the estimate $\Ltil$ were larger than $\Llamr[1]$.

§ LOWER ORACLE COMPLEXITY BOUND

In this section, we present the lower oracle complexity bounds, that is, the minimum number of queries to the FO oracle required to find an $(\ep; \ep/c)$-optimal solution. These results illustrate the optimality of the ACGD method under a certain optimality regime. We assume, for the sake of simplicity, that  $u(x)= \alpha\normsq{x}/2$ (see (<ref>)) and that $X$ is radially invariant, for example, the $\R^n$.

Similar to Nesterov's lower complexity computation model in [19], we consider the class of all first-order methods, $\calF$, verifying a linear-span update requirement. Given a (finite) memory of reachable points $\Mt[t-1]$ after the $t-1$th query to FO, the updated memory $\Mt$ after evaluating the FO oracle at some $y \in \Span(\Mt[t-1])$ needs to satisfy

\begin{equation}\label{req:fo-model}
\Mt \subset \{x + \etat \grad f(y) + \sumi \taut[i] \grad g_i(y): x \in \Span(\Mt[t-1]), \etat, \taut[i,t] \in \R \}.
\end{equation}

The freedom to choose arbitrary elements from the linear span allows $\calF$ to cover many first-order algorithms. For example, the ACGD method is a special case of $\calF$ because the generated points in the ACGD method, $\Mt =\{\xt[0], \xt[1], \ldots, \xt[t]\}$, satisfy the requirement (<ref>). Specifically, when the memory $\Mt[t-1] = \{\xt[0], \xt[1], \ldots, \xt[t-1]\}$, the evaluation point $\xundert$ in Line 3 is inside $\Span(\Mt[t-1])$. Moreover, the $\xt$-update in Line 5 of Algorithm <ref> can be expressed as

\begin{align*}
&\xt \leftarrow \argmin_{x \in X} \{\inner{\grad f(\xundert)}{x} + u(x) + \etat \normsq{x - \xtt}/2\ s.t.\ \grad g(\xundert) (x - \xundert) + g(\xundert) \leq 0\} \\
\stackrel{(a)}{\Leftrightarrow}\ &\xt  \leftarrow \argmin_{x \in X} \{\inner{\grad f(\xundert) + \sumi \lamti \grad g_i(\xundert)}{x}  + \alpha\norm{x}^2/2 + \etat \normsq{x - \xtt}/2  \}\\
\Leftrightarrow\ & \xt \leftarrow \argmin_{x \in X} \normsq{x - \tfrac{1}{\etat + \alpha}[\xt[t-1] - \tfrac{1}{\etat}(\grad f(\xundert) + \sumi \lamti \grad g_i(\xundert))]}\\
\stackrel{(b)}{\Leftrightarrow}\ & \xt =  \tfrac{1}{\gamma (\etat + \alpha)}[\xt[t-1] - \tfrac{1}{\etat}(\grad f(\xundert) + \sumi \lamti \grad g_i(\xundert))]\  \text{for some $\gamma>0$,}
\end{align*}

where the multiplier $\lamt$ in (a) is the optimal dual solution to the quadratic program, and (b) follows from $X$ being radially invariant. Therefore, $\xt$ is a member of the right-hand side of (<ref>), $\Mt := \Mt[t-1] \cup \{\xt\}$ satisfies the update requirement (<ref>), and the ACGD method is a member of $\calF$. In fact, similar arguments can be used to show that $\calF$ covers both the primal methods [13, 19],  and the primal-dual methods [3, 6, 23] in the literature.

Since the dependence on parameters of the smooth objective function is well-established, e.g. [16, 19], we will investigate the dependence of the lower complexity bounds on the parameters of the constraint function. Toward that end, we consider an affine $f$ in the objective, i.e. $\Lf = 0$, and study the dependence of the lower complexity bound on the norm of the optimal Lagrange multiplier $\norm{\lamstar}$ and the Lipschitz smoothness constant of the vector-valued constraint function $\Lgbar$, i.e., $\norm{\grad g(x) - \grad g(y) } \leq \Lgbar \norm{x - y} \forall x, y\in \R^n$. These parameters are used more often in the literature, but we will relate them to the aggregate smoothness constant $\Llamr$ (<ref>) for the ACGD method in the forthcoming discussion.

§.§ Strongly Convex Case

First, we study the strongly convex problem with $\alpha>0$. Since linear convergence is expected, different optimality criteria have little impact on the lower complexity bound, that is, they only have different constants inside the “$\log$”. So, without the loss of generality, we choose to focus on the convergence of $\norm{\xt - \xstar}$ in the next theorem.

Let problem parameters $\Lgbar > 0$, $l \geq  1$, and $\Lgbar l\geq  \alpha > 0 $ be given. There exists an infinite dimensional hard problem of the form (<ref>) with $\norm{\lamstar} = l$, $g$ being $\Lgbar$-Lipschitz smooth, and $f$ being affine, i.e., $L_f =0$, such that every first-order method in $\calF$ requires at least $\Omega(\sqrt{{\Lgbar l}/{\alpha}}\log(1/\ep))$ queries to the FO oracle to find an $\xt[N]$ with $\normsq{\xt[N] - \xstar} \leq \ep$ for all $\ep >0.$

Let $\gam = {\alpha}/({\Lgbar l})$, $\Delta = (1 - \sqrt{\gam})/(1 + \sqrt{\gam})$, $\beta = \Lgbar$ and $\xbar =  [\Delta, \Delta^2, \ldots, \Delta^i, \ldots] \in \R^{\infty}.$ Consider the following hard problem:

\begin{equation}
\begin{split}
\min_{x \in \R^\infty}& -\tfrac{l \beta - \alpha}{4} x_1 + \tfrac{\alpha}{2} \normsq{x}\\
s.t.&\ g(x):=h(x)  - h(\xbar) \leq 0 \text{ with } h(x)= \tfrac{\beta - \alpha /l}{8} [x_1^2 + \sumi[\infty] (x_i - x_{i+1})^2].
\end{split}
\end{equation}

Clearly, the objective has a strong convexity modulus of $\alpha$ and the constraint function $g(x)$ has a smoothness constant of $\Lgbar$. It is straightforward to verify that $\lamstar = l$ and $\xstar = \xbar$ satisfies the KKT condition:

\begin{align*}
g(\xstar)  = 0  \text{ and } \left(\tfrac{\lamstar}{l}\tfrac{l \beta - \alpha}{4} \begin{bmatrix} 2 & -1 \\ -1 & 2 & -1 \\ & -1 & 2 & -1& \dots\\& & & & \ddots \end{bmatrix} + \alpha \begin{bmatrix} 1\\ 1 \\ 1\\ \vdots\end{bmatrix} \right) \begin{bmatrix}x^*_1\\ x^*_2 \\ \xstar_3\\ \vdots\end{bmatrix}  = \tfrac{l \beta - \alpha}{4}  \begin{bmatrix} 1 \\ 0 \\ 0 \\ \vdots \end{bmatrix}.
\end{align*}

Starting from $\calM_0 = \{0\}$, at iteration $t$, the solution $x_t$ generated by any first-order method satisfying (<ref>) have non-zeros in only the first $t$ coordinates. Thus we have

\begin{equation}
\normsq{\xt - \xstar} \geq \tsum_{i=t+1}^\infty (\xstar_i)^2 = (\tfrac{1 - \sqrt{\gam}}{1 + \sqrt{\gam}})^{2t} \sumi[\infty] (\tfrac{1 - \sqrt{\gam}}{1 + \sqrt{\gam}})^{2i} = (\tfrac{1 - \sqrt{\gam}}{1 + \sqrt{\gam}})^{2t} \normsq{\xt[0] - \xstar} \geq (1 - \sqrt{\gamma})^{2t} \normsq{\xt[0] - \xstar},
\end{equation}

where the last inequality follows from $1\geq \gam \geq 0$. Therefore we require at least  $t = \Omega(\sqrt{{\Lgbar l}/{\alpha}}\log(1/\ep))$ iterations to find $\xt$ with $\normsq{\xt - \xstar} \leq \ep$.

Since, in this section, $f$ is assumed to be affine such that  $L_f = 0$, we have $\Llamr[1] \leq (1 + \norm{\lamstar})\Lgbar = \bigO(\norm{\lamstar}\Lgbar)$. The preceding theorem then shows the $\bigO(\sqrt{\Llamr[1]/\alpha}\log(1/\ep))$ upper complexity bound in Corollary <ref> to be unimprovable when $\norm{\lamstar} \geq 1$, that is, the ACGD method has a tight oracle complexity.

§.§ Non-strongly Convex Case

Now we move on to consider the non-strongly convex problem with $\alpha =0$. The next theorem states the lower oracle complexity bound to find an $(\ep; \ep/c)$-optimal solution.

Let problem parameters $\Lgbar > 0$, $\Rt[0] \geq 1$, $l > 0$, $c \geq 1$ and $\ep > 0$ be given. For a large enough problem dimension, $n > 2 \ceil{\Rt[0]\sqrt{\Lgbar c/ \ep} },$ there exists a hard problem of form (<ref>) with $\norm{\lamstar} = l$, $\norm{\xt[0] - \xstar} \leq \Rt[0]$, $g$ being $\Lgbar$-Lipschitz smooth, and $f$ being affine, i.e., $L_f =0$, such that every first order method in $\calF$ requires at least $\Omega(\sqrt{\Lgbar c}\Rt[0]/\sqrt{\ep})$ queries to the first order oracle to find an $(\ep;\ep/c)$-optimal solution.

Consider the following function-constrained problem parameterized by $\gamma> 0$, $\beta>0$,  and $k\in \mathbb{N}_+$:

\begin{align}\label{ls:lower-eg}
\begin{split}
\min_{x \in \R^{2k+1}}& - 2 l \gamma \beta x_1 \\
s.t.\ & g_1(x) :=\beta[x_1^2 + \sumi[2k] (x_i - x_{i+1})^2 + x_{2k+1}^2] - (\tfrac{2k+1}{2k+2}) \gam^2 \beta \leq 0 \\
\ & g_2(x) := \beta[-2x_1 \gam + x_1^2 + \sumi[2k] (x_i - x_{i+1})^2 + x_{2k+1}^2] + (\tfrac{2k+1}{2k+2}) \gam^2 \beta\leq 0,
\end{split}\end{align}

where $l$ is the given parameter in the theorem statement. Without loss of generality, we take $\xt[0] = 0$ and $\Mt[0] = \{0\}$. Let $\calK_i$ denote the subspace with non-zeros in only the first $i$th coordinates, i.e., $\{x \in \R^{2k+1}: x_j = 0\  \forall j > i\}$. Given a first-order method satisfying (<ref>), it is easy to show inductively that $\Mt \subset \calK_t\ \forall t \in [2k+1]$. This is because $\mathcal{G}_i = \{\grad f(\xbar), \grad g_1(\xbar), \grad g_2(\xbar): \xbar \in \calK_i\}$ are non-zero only in the first $i+1$ coordinates, i.e., $\mathcal{G}_i \subset \calK_{i+1}.$ Thus in $k$ iterations, we have the following lower bound on feasibility violation:

\begin{equation}\label{pf:feas-lowB}
\norm{[g(x^k)]_+} \geq  \min_{x \in \calK_k} g_2(x) \geq (\tfrac{1}{2k+2}) \beta \gam^2.\end{equation}

Now we calculate the problem parameters associated with (<ref>). It is straightforward to verify via the KKT condition that the optimal solution and the optimal dual multiplier are respectively:

$$\lamstar = [l, 0], \ \xstar_i = \gam [1 - \tfrac{i}{2k+2}]\ \forall i \in [2k+1].$$

So $\norm{\xt[0] -\xstar} \leq \gam \sqrt{k+1}$, $\norm{\lamstar} = l$, and the constraint function $g$ have a smoothness constant of $12\beta$. By selecting $k = \floor{\tfrac{1}{5} \sqrt{\tfrac{\Lgbar c}{\ep}} \Rt[0]} - 1$, $\beta = \Lgbar/12 $, $\gam = \Rt[0]/\sqrt{k+1}$, the problem satisfies the requirements in the theorem statement. Moreover, (<ref>) implies that in $k = \floor{\tfrac{1}{5} \sqrt{\tfrac{\Lgbar c}{\ep}} \Rt[0]} -1 = \Omega(\sqrt{\tfrac{\Lgbar c}{\ep}} \Rt[0])$ iterations, the feasibility violation are lower bounded by

$$\norm{[g(x^k)]_+} \geq (\tfrac{1}{2k+2}) \beta \gam^2 = (\tfrac{1}{24})(12 \beta) [\gam^2 (k+1)] [\tfrac{1}{(k+1)^2}] \geq (\tfrac{1}{24})(\Lgbar) \Rt[0]^2 \tfrac{25 \ep}{\Lgbar c \Rt[0]^2} \geq \tfrac{\ep}{c}.$$

This shows that $\Omega(\sqrt{\tfrac{\Lgbar c}{\ep}} \Rt[0])$ iterations are necessary for finding an $(\ep, \ep/c)$-optimal solution.
Since the choice among $\calF$ is arbitrary, the lower complexity bound is valid for all first-order methods in $\calF$.

Comparing the above lower bound of $\Omega\{\sqrt{\tfrac{\Lgbar c}{\ep}} \norm{\xstar - \xt[0]}\}$ to the upper bound of $\bigO\{\sqrt{\tfrac{\Llamr[c]}{\ep}} \norm{\xstar - \xt[0]}\}$ from Corollary <ref>, we see that the dependences of the oracle complexity of the ACGD method on $\ep$ and $\norm{\xstar - \xt[0]}$ are not improvable. Only the dependence on the smoothness constant $\Lgbar$ might be sub-optimal. Specifically, the (big-O) factor of sub-optimality can be characterized by the following function of $c$:

$$H(c):= \bigO(\sqrt{\tfrac{\Llamr[c]}{\Lgbar c}}).$$

Since, in this section, $\Lf = 0$ such that $\Llamr[c]\leq \Lgbar(\norm{\lamstar} + c)$ (see (<ref>)),  the big-O relationship between $\Llamr$ of the upper complexity bound and the $c \Lgbar$ of the lower complexity bound is shown in Figure <ref>.

[scale=0.5, font=]
[axis lines=middle,

xticklabels=1,$\norm{\lamstar}$, ]
[name path=F,blue,dashed,domain=0.5: 4] sqrt(x) node[pos=.8, below]$\sqrt{\Lgbar c}$;

[name path=G,red,domain=0.5:4] max(sqrt(2), sqrt()) + 0.1node[pos=0.7, above left]$\sqrt{\Lgbar(\norm{\lamstar} + c)}$;

[pattern=north west lines, pattern color=brown!50]fill between[of=F and G, soft clip=domain=0.5:4]
at (axis cs:1, 1.3) $\sqrt{\Llamr[c]}$;

(Big-O) Dependence of Complexities on $c$

Clearly $H(c)$ has two distinct regions. When  $c \in [1, \norm{\lamstar}]$, we have $1 \leq H(c) \leq \sqrt{\norm{\lamstar}/c}$, i.e., the oracle complexity of the ACGD method can be sub-optimal up to a  factor of $\sqrt{\norm{\lamstar}/c}$.  In practice, the factor $H(c)$ can be smaller. In fact, when the constraint functions have imbalanced Lipschitz-smoothness constants, $H(c)$ could be $1$ so that the ACGD method is optimal. For instance, consider a slightly modified version of (<ref>):

\begin{equation*}
\min_{x \in \R^{2k+1}} \{ - 2 l \gamma \beta x_1\  \text{s.t.}\ g_1(x) \leq 0, l g_2(x) \leq 0\},\\
\end{equation*}

where the second constraint $g_2$ is multiplied by $l$. Since the constraint $l g_2(x) \leq 0$ is not active,
the above modified problem still has the same $\xstar$ and $\lamstar = [l, 0]$ as (<ref>). Since the smoothness constants $\Lgbar$ is increased from $\bigO(\beta)$ to  $\bigO(  l \beta)$, $c\Lgbar = \bigO(l c \beta )$ and $\Llamr[c]= \bigO[(1 + c) l\beta]$ are now of the same order such that $H(c) = 1$.

When $c \geq \norm{\lamstar}$, we have $H(c) = 1$ such that the ACGD method is always optimal.  Since $\lamstar$ presents the shadow price of the constraints close to $\xstar$, the cost of changing an $(\ep/c)$-feasible solution, i.e. $g_i(\xt) \leq \tilde{\ep}_i$, to feasibility is roughly

$$\max_{\norm{\tilde{\ep}} \leq \ep /c}\sumi \lamstar \tilde{\ep}_i = \norm{\lamstar} \ep / c. $$

This shows that a scaling constant $c \geq \norm{\lamstar}$ should be selected to ensure that the $(\ep/c)$-feasibility requirement is comparable to the $\ep$-optimality requirement.  In this case, the proposed ACGD method is optimal.

§ THE ACGD-S METHOD

We extend the ACGD method to the ACGD with sliding (ACGD-S) method to handle the large-scale problem where both the problem dimension $n$ and the number of constraints $m$ are large. This section follows the same structure as Section <ref>. We first discuss the computation bottleneck in the large-scale setting. Then Subsection <ref> introduces the ACGD-S method and presents the convergence results, Subsection <ref> proposes a guess-and-check search scheme to look for some elusive problem parameters important for stepsizes selection, and Subsection <ref> contains the detailed proofs to the convergence results.

Despite its optimal oracle complexities, the ACGD method may be lacking in computation efficiency for the large-scale problem. The  bottleneck to Algorithm <ref> lies in  Line 5:

\begin{align}\label{eq:ANLA-QP}\begin{split}
\xt \leftarrow \argmin_{x \in X}& \inner{\pit}{x} +u(x) + \etat \normsq{x - \xtt}/2\\
s.t.&\ \nuit (x - \xundert) + g_i(\xundert) \leq 0\  \forall i \in [m].
\end{split}\end{align}

It amounts to a large-scale quadratic program (QP) if $X$ is linearly constrained, say a box, and a large-scale quadratic constrained quadratic program (QCQP) if $X$ is a Euclidean ball.

In this section, we address the bottleneck by replacing the large-scale QP with a sequence of basic matrix-vector operations, each requiring at most $\bigO(mn)$ FLOPs. The proposed ACGD-S method requires only a similar number of matrix-vector operations as solving a single linear constrained problem, i.e., $g(x)$ is affine, and maintains the same optimal oracle complexity as the ACGD method.
Towards that end, we need to assume the projection onto $X$ is easy, i.e., the following operation can be computed in $\bigO(n)$ FLOPs for any $\pi, \xbar \in \R^n$ and $\etat \geq 0.$

\begin{equation}\label{eq:x-projection}
x^t \leftarrow \targmin_{x \in X} \inner{\pi}{x} + u(x) + \etat\normsq{x - \xbar}/2.\end{equation}

For instance, if $u(x) = \alpha\normsq{x}/2$, the computation simplifies to component-wise thresholding if $X$ is a box,and to vector scaling if $X$ is a Euclidean ball. If $X$ is more challenging, we can model the complicated part using the function constraints.

The ACGD-S Method

[1]
$\xt[-1]= \xundert[0] = \yst[0][1] =\xt[0]\in X$,  stepsizes $\{\thetat\}$, $\{\etat\}$, $\{\taut\}$, and  weights $\{\wt\}$.
Set $\pit[0] = \grad f(\xt[0]) $, $\nut[0] = \grad g(\xt[0])$, $\lamst[-1][1] = \lamst[0][1] = 0 $.
$t =  1,2,3 ... N$
Set $\xundert \leftarrow (\taut \xundert[t-1] + \xtilt)/(1 + \taut)$ where  $\xtilt  = \xtt + \thetat (\xtt - \xt[t-2])$.
Set $\pit \leftarrow \grad f(\xundert)$ and $\nut \leftarrow \grad g(\xundert)$.
Calculate inner loop iteration limit $\{\St\}$ stepsizes $\{\betast\}$ and $\{\gamst\}$,  and weights $\{\delst\}$.
$s=1,2,..., \St$
Set $\htilst = \begin{cases} (\nut)^\top\lamst[0] + \rhost[1] (\nut[t-1])^\top (\lamst[0] - \lamst[-1]) &\text{if } s =1 \\
(\nut)^\top\lamst[s-1] + \rhost (\nut[t])^\top (\lamst[s-1] - \lamst[s-2]) &\text{ o.w.}\end{cases}$
Set $\yst \leftarrow \argmin_{y \in X}\inner{\htilst + \pit}{y} + u(y) + \etat\normsq{y - \xtt}/2 +  \betast \normsq{y - \yst[s-1]}/2. $
Set $\lamst \leftarrow \argmax_{\lam \in \R^m_+}\inner{\lam}{\nut[t] (\yst - \xundert) + g(\xundert)} + \gamst \normsq{\lam - \lamst[s-1]}/2.$
Set $\lamst[0][t+1] = \lamst[\St][t]$, $\lamst[-1][t+1]= \lamst[\St-1][t]$, $\yst[0][t+1] = \yst[\St][t]$
Set $\xt = \sums \delst \yst /(\sums \delst)$ and $\lamtilt = \sums \delst \lamst /(\sums \delst)$.

$\xbart[N]:= \tsum_{t=1}^N \wt \xt / (\sumt \wt). $

§.§ The ACGD-S Method and its Convergence Results

The ACGD-S method, listed in Algorithm <ref>, consisted of two loops. For clarity, we will call an outer iteration a phase, and an inner iteration an iteration for the rest of rest section. In phase $t$, the outer loop updates happen in Lines 3, 4, and 12; they are identical to the ACGD method except that the exact solution $\xt$ to the QP in (<ref>) is replaced by some average of inner iterates.
The others steps, Lines 5-11,  constitute the sliding subroutine. Its goal is to solve the Lagrangian reformulation to (<ref>), or the $( \xt, \lamt)$ saddle point problem in (<ref>), inexactly:

\begin{equation}\label{eq:slide-saddle}
(\xt, \lamt) \leftarrow \argmin_{y \in X} \argmax_{\lam \in \R^m_+} \inner{\pit}{y} + \inner{\lam}{\nut y - \xundert} + u(y) + \etat \normsq{y - \xtt}/2.\end{equation}

To avoid confusion, we use the dummy variable $y$ to emphasize it being used only in the inner loop. Specifically, Line 5 calculates the stepsize parameters and iteration number $\St$. Lines 7-9 carry out primal-dual type updates for $\St$ iterations. Line 7 computes a momentum extrapolation term $\htilst$ as a proxy for $(\nut)^\top  \lamst$. In Line 8, with the variable $(\nut)^\top \lambda$ being fixed to $\htilst$,   $\yst$ is generated by minimizing the variable $y$ in (<ref>) subject to a prox-function $\betast \normsq{y - \yst[s-1]}$. Then in Line 9,  with the variable $y$ fixed to $\yst$, $\lamst$ is generated by  maximizing (<ref>) with subject to a prox-function $\gamst \normsq{\lambda - \lamst[s-1]}$. After that, Line 11 prepares the initialization points for the inner loop in the next phase.

We highlight three features that are essential for achieving the desired computation efficiency. First, rather than being pre-specified, the inner loop step-size parameters and iteration limit $\St$ are calculated in an online fashion in Line 5. This allows the method to adjust dynamically to the varying difficulty of the saddle point problems (<ref>) from different phases. Second, the last operator $\nut[t-1]$, rather than $\nut$, is used for calculating the momentum extrapolation term at the first iteration $s=1$ in Line 7. This is characteristic of the sequential dual type algorithms [12, 27, 26] for solving the trilinear saddle point problem in (<ref>). Third, two primal iterates, $\xt$ and $\yst[\St]$, are stored after each inner loop to kick-start the next one. This is common to sliding-type algorithms [7, 11, 12].

Now we suggest certain stepsize choices to obtain concrete convergence rates for Algorithm <ref>.
The non-strongly convex case and the strongly convex case are presented in separate theorems.

Consider a non-strongly problem of form (<ref>) with $\alpha =0$. Let the aggregate smoothness constant $\Llamr[c]$ and the reference multiplier set $\Lamr[c]$ be defined in (<ref>) and (<ref>) respectively. Suppose Algorithm <ref> is performed with the following stepsizes. The outer-loop stepsizes are

\begin{equation}\label{stp:S-non-outer}
\taut = \tfrac{t-1}{2}, \etat = \tfrac{\Llamr[c]}{\taut[t+1]}, \wt = t, \thetat[t+1] = \wt[t+1] /\wt\  \forall t \geq 1.
\end{equation}

With $\MtU = \norm{\nut}$, $\rlamr[r] = \norm{\lamstar} + r$, and some $\Delta > 0$, the inner loop parameters in the phase $t$ are calculated according to

\begin{equation}\label{stp:S-non-inner}
\St = \ceil{\MtU \Delta t}, \Mtilt = \tfrac{\St}{\Delta t},\ \rhost = \begin{cases} \Mtilt/\Mtilt[t-1] &\text{if } s= 1\\
1 & \text{if } s \geq 2, \\
\end{cases}
\  \betast = \tfrac{\Mtilt \rlamr[c]}{\norm{\xt[0] - \xstar}},\ \gamst = \tfrac{\Mtilt^2}{\betast},\ \delst = 1\ \forall s \geq 1,
\end{equation}

Then $\MtU \leq \Mbar\ \forall t$, where $\Mbar$ is an upper bound of $\norm{\grad g(x)}$ for all $x$ in some bounded neighbourhood around $\xstar$. Moreover, the ergodic average solution $\xbarN$ satisfies

\begin{equation}\label{eq:S-non-strong-conv}
\max\{F(\xbarN) - \Fstar, c \norm{[g(\xbarN)]_+}\} \leq \tfrac{1}{N(N+1)} (\tfrac{2\rlamr[c]\norm{\xt[0] - \xstar}}{\Delta} + \Llamr[c] \normsq{\xt[0] - \xstar}).
\end{equation}

Under the setting of Theorem <ref>, if $\Delta = \tfrac{\rlamr[c]}{\norm{\xt[0] - \xstar} \Llamr[c]}$, the total numbers of operations required by the ACGD-S method to find an $(\ep; \ep/c)$-optimal solution are bounded by:

* $N_\ep = \ceil{\sqrt{\tfrac{3 \Llamr[c]}{\ep}}\norm{\xt[0] - \xstar}}$ FO-oracle evaluations.

* $C_\ep = \bigO\{(\sqrt{\tfrac{\Llamr[c]}{\ep}}\norm{\xt[0] - \xstar} + \tfrac{ \rlamr[c] \Mbar\norm{\xt[0] - \xstar}}{\ep})\}$ matrix-vector multiplications.

Three remarks are in order regarding the above results. First, the oracle complexity of the ACGD-S method matches that of the ACGD method, while its computation complexity, measured by the number of matrix-vector multiplications, matches the lower bound for solving a single linearly constrained problem [21]. Second, the stepsize choices in (<ref>) and (<ref>) only require an upper bound to $\Llamr[c]$. The misspecification of $\norm{\xt[0] - \xstar}$, $d(\Lamr[c])$ or $\Delta$ would still lead to an $\bigO(1/\sqrt{\ep})$ oracle complexity and an $\bigO(1/\ep)$ computation complexity. Third, the iteration limit function $\St$ in (<ref>)  adapts to the varying difficulty of the saddle point sub-problem (<ref>) from different phases. Specifically, $\St$ scales in proportion both to  $\norm{\nut}$, which characterizes the hardness of the saddle point sub-problem in  (<ref>), and to $t$, which captures the degree of accuracy required by the outer loop.

Consider a strongly convex problem of form (<ref>) with $\alpha >0$. Let the aggregate smoothness constant $\Llamr[r]$ and the reference multiplier set $\Lamr[r]$ be defined in (<ref>) and (<ref>) respectively, and let $\kar:=\Llamr/\alpha$ be the condition number. Suppose Algorithm <ref> is performed with the following stepsizes. The outer-loop stepsizes are

\begin{equation}\label{stp:S-str-outer}
\taut = \min\{\tfrac{t-1}{2}, \sqrt{2\kar}\}, \etat = \tfrac{\Llamr[c]}{\taut[t+1]}, \text{ and } \thetat = \tfrac{\taut}{\taut[t-1] + 1} \forall t \geq 1. \quad  \wt = \begin{cases} \wt[t-1] /\thetat & \text{if } t  \geq 2,\\
1 &\text{if } t = 1.\end{cases}
\end{equation}

For the inner loops, let $\MtU = \norm{\nut}$, and let a parameter $\Delta> 0$ be given. At the beginning, the iteration limit is set to $\St[1] = \min \{S\in \mathbb{N}_+: \tsum_{j=1}^{S} s \geq \wt[1] \MtU[1]^2 \Delta \}$, and the step-sizes for all $s \in [\St[1]]$ are set to

\begin{equation}\label{stp:S-str-init}
\delst[s][1] = \tfrac{s}{\Gamt[1]},\ \betast[s][1] = \tfrac{\alpha}{4} (s-1),\ \gamst[s][1] = \tfrac{4}{\alpha}\tfrac{\MtU[1]^2 \Gamt[1]}{\delst[s][1]},\text{ and } \Wt[2] = \tfrac{\delst[\St[1]][1]}{ \MtU[1]},
\end{equation}

where $\Gamt[1]$ is the non-negative root to $\sums[\St[1]] s = \Gamt[1]^2 (\wt[1] \MtU[1]^2 \Delta).$ Then for phase $t \geq 2$, the iteration limit $\St$ and the parameter $\Gamt \geq 0$ are specified to satisfy

\begin{equation}\label{stp:St-strong}
\St = \min\{S\in\mathbb{N}_+: \sums[\St] \Wt \MtU + (s-1) \geq \wt \MtU^2 \Delta\},\ \sums \Gamt (\Wt \MtU) + (s-1) = \Gamt^2 \wt \MtU^2 \Delta,
\end{equation}

and the stepsizes are chosen according to

\begin{equation}\label{stp:S-str-inner}
\delst = \Wt \MtU + \tfrac{1}{\Gamt} (s-1),\ \gamst[s] = \tfrac{4}{\alpha}\tfrac{\MtU^2 \Gamt}{\delst[s]},\ \betast = \begin{cases}\tfrac{\alpha}{4}(\Wt\Gamt\MtU) & \text{ if } s=1, \\ \tfrac{\alpha}{4}[\Wt\Gamt\MtU + (s-2)] &\text{otherwise,}\end{cases}\ \text{ and } \Wt[t+1] = \tfrac{\delst[\St]}{ \MtU}.
\end{equation}

Then we have $\MtU \leq \Mbar\ \forall t$, where $\Mbar$ is an upper bound for $\norm{\grad g(x)}$ for $x$ in some bounded neighbourhood around $\xstar$, and the ergodic average solution $\xbarN$ satisfies

\begin{equation}\label{eq:S-str-conv}
\begin{split}
&\max\{F(\xbarN) - \Fstar, r \norm{[g(\xbarN)]_+}\} \leq \tfrac{1}{\calWt} (\tfrac{2}{\alpha\Delta} [\rlamr[r]]^2 + \Llamr[r] \normsq{\xt[0] - \xstar}), \\
&\normsq{\xbarN - \xstar} \leq \tfrac{2}{\alpha\calWt} (\tfrac{2}{\alpha\Delta} [\rlamr[r]]^2 + \Llamr[r] \normsq{\xt[0] - \xstar}),
\end{split}
\end{equation}

where $\rlamr=\norm{\lamstar} + r$ and the denominator satisfies $\calWt \geq \max\{N(N+1)/2, \sqrt{2\kar}[(1+ 1/\sqrt{2\kar})^{N-4} -1]\}.$

Under the setting of Theorem <ref>, if $\Delta = \tfrac{[\rlamr]^2}{\Llamr \normsq{\xt[0] - \xstar} \alpha}$, the numbers of FO-oracle evaluations $N_\ep$ and of matrix-vector multiplications $C_\ep$ required by the ACGD-S method to find an $(\ep; \ep/c)$-optimal solution are bounded by:

* $N_\ep \leq  \min\{\sqrt{\tfrac{6 \Llamr \max\{c/r, 1\}}{\ep}} \normsq{\xt[0] - \xstar}, [\sqrt{\tfrac{2\Llamr}{\alpha}} + 1]\log(\tfrac{3\max\{c/r,1\}\sqrt{\Llamr[r]\alpha}\normsq{\xt[0] - \xstar}}{\ep} + 1) + 4\} + 1\}$.

* $C_\ep = \bigO\{\sqrt{\tfrac{\max\{c/r, 1\} [\rlamr]^2 \Mbar^2}{\alpha\ep}} + N_\ep\}$.

Under the setting of Theorem <ref>, if the goal is to find an $\epsilon$-close solution, i.e., $\normsq{\xbart[N] - \xstar} \leq \ep$,  we can choose $c =r =1$ and $\Delta = \tfrac{[\rlamr[1]]^2}{\Llamr[1] \normsq{\xt[0] - \xstar} \alpha}$ such that the numbers of required operations can be bounded by

* $N_\ep = \bigO\{\sqrt{\tfrac{2\Llamr[1]}{\alpha}} \log(\tfrac{\sqrt{\kar[1]}\normsq{\xt[0] - \xstar}}{\ep})\} $

* $C_\ep = \bigO\{{\tfrac{ \rlamr[1] \Mbar}{\alpha\sqrt{\ep}}} + \sqrt{\tfrac{2\Llamr[1]}{\alpha}} \log(\tfrac{\sqrt{\kar[1]}\normsq{\xt[0] - \xstar}}{\ep})\}.$

Again, we make a few remarks regarding the results. First, to find an $\ep$-close solution, Corollary (<ref>) implies that the ACGD-S method has the same oracle complexity as the ACGD method, and has the same computation complexity as that of the lower computation complexity bound for solving a single strongly-convex linearly constrained problem [21]. Second, the iteration limit function $\St$ is again adaptive to the varying difficulty of the saddle-point subproblem (<ref>) from different phases. Third, the rather complicated inner-loop stepsize choice in (<ref>) and (<ref>) is the first among sliding algorithms, e.g. [7, 10, 11, 12],  to achieve both the optimal inner loop complexity of $O(1/\sqrt{\epsilon})$ and the optimal outer loop complexity of $\bigO(\sqrt{\kappa}\log(1/\ep))$  without restarting. It is unclear if the same effect is achievable with simpler stepsize choices.  Notwithstanding that, the stepsize is easy to implement in practice because only conservative estimates of $\alpha$ and $\Llamr$ are required to obtain the $\bigO(\sqrt{\kar}\log(1/\ep))$ oracle complexity and the $\bigO(1/\sqrt{\ep})$ computation complexity.

§.§ The Binary Search for $\Llamr$ and $\rlamr$

To enhance its implementability, we propose a binary search scheme to find (upper) approximations to these hard-to-estimate problem parameters required by the ACGD-S method. Since both quantities are proportional to $\Lamr$, it is convenient to look for a common upper bound:

$$\Htil \geq \max\{\Llamr, \rlamr\}.$$

The search scheme has the same logic as Algorithm <ref>: starting with a guess $\Htil$, we first run the ACGD-S method for $N(\Htil)$ phases, with the stepsizes specified according to $\Htil$, to generate an ergodic solution $\xbart[N(\Htil)]$, then we check it against a certain termination criterion and update our guess $\Htil$ accordingly.
The following assumption is required for constructing an easy-to-verify termination criterion.

The feasible region $X$ is bounded and we know an estimate $\DX$ of its radius, i.e., $\DX \geq \max_{x,y \in X} \norm{x - y}.$

Compared to Assumption <ref>, we no longer require access to the dual solution $\lamt$ (Line 5 of Algorithm <ref>). This is because the sliding subroutine (Line 12 of Algorithm <ref>) provides a feasible dual iterate $\lamtilt$ as a by-product. So the new termination criterion, termed (<ref>)-(<ref>), requires only the following modification to (<ref>):

\begin{equation}\label{cr:check-S}
\boxed{\text{Change $\lamt_i$ in \eqref{eq:lower-obj} to $\lamtilt_i$.}}\end{equation}

Since $\lamtilt$ from Line 12 of Algorithm <ref> is dual feasible,  the minimum value $\Fundert$ to the modified linear program in (<ref>) (with $\lamtilt$ in place of $\lamt$) still bounds $\Fstar$ from below. The satisfaction of new termination criterion (<ref>)-(<ref>) then implies that $\xbart[N(\Htil)]$ is $(\ep; \ep/c)$-optimal. Accordingly, we modify Algorithm <ref> by switching to the ACGD-S method and to the termination criterion (<ref>)-(<ref>) to arrive at Algorithm <ref>.

Binary Search Outer Loop to ACGD-S.

[1]
$\xbar \in X$, $\Htil> 0$ , and a phase-limit function  $N(H)$.

With $\Dx$, $\Htil$, and $\Htil$ in place of $\norm{\xt[0] - \xstar}$, $\Llamr$, and $\rlamr$,  set stepsizes to (<ref>) and (<ref>) if $\alpha=0$, and to (<ref>), (<ref>) and (<ref>) if $\alpha > 0$.
Run the ACGD-S method for $N(\Htil)$ phases starting from $\xbar$.
Set $\xbar \leftarrow \xbart[N(\Htil)]$, where $\xbart[N(\Htil)]$ is the ergodic average solution.
$\xbar$ satisfies the condition (<ref>)-(<ref>)
Set $\Htil \leftarrow 2 \Htil$.

Clearly, the above algorithm terminates with an $(\ep, \ep/c)$-optimal solution. We only need to specify the phase-limit function $N(\Htil)$ to make it implementable. The next theorem provides some choices of $N(\Htil)$ and the corresponding complexities. Its derivation is similar to Theorem <ref>

Consider a problem of form (<ref>). Let Assumption <ref> hold and $\Mbar:= \max_{x \in X} \norm{\grad g(x)}$. Suppose the reference set $\Lamr$ is defined according to (<ref>). Let $\Llamr$  denote the aggregate smoothness constant (<ref>) and $\rlamr=\norm{\lamstar} + r$.  If $\Htil$ is initialized to some $\Htilt[0]$,  Algorithm <ref> then terminates finitely with an $(\ep; \ep/c)$-optimal solution under the following choices of phase limit function $N(\Htil)$.

a) In the non-strongly convex case with $\alpha=0$, if
$N(\Htil)= \ceil{\sqrt{\tfrac{3\Htil}{\ep}}\Dx},$ the total numbers of oracle evaluations (phases) and matrix-vector multiplications required for termination are  bounded  respectively by

$$\bigO(\sqrt{\tfrac{\max\{\Llamr[c], \rlamr[c], \Htilt[0] \}}{\ep}}\Dx) \text{, and } \bigO(\tfrac{\max\{\Llamr[c], \rlamr[c], \Htilt[0] \}\Dx \Mbar}{\ep}).$$

b) In the strongly convex case with $\alpha>0$, if
$N(\Htil)= \ceil{(\sqrt{\tfrac{2 \Htil}{\alpha}} + 1)\log(\tfrac{3\max\{c/r, 1\}\sqrt{\Htil \alpha} \DX^2}{\ep} + 1)} + 4,$
the total numbers of oracle evaluations (phases) and matrix-vector multiplications required for termination are bounded  respectively by

$$\bigO(\sqrt{\tfrac{\max\{\Htilt[0], \Llamr[1], \rlamr[1] \}}{\alpha}} \log \tfrac{c \max\{\Htilt[0], \Llamr[1], \rlamr[1] \}}{\ep}), \text{ and } \bigO(\tfrac{\sqrt{c} \max\{\Htilt[0], \Llamr[1], \rlamr[1] \} \Mbar \Dx}{\sqrt{\alpha \ep}}).$$

We remark that both the oracle complexity and the computation complexity across the multiple ACGD-S runs in Algorithm <ref> have the optimal order for all $\alpha \geq 0$. If the Lipschitz smoothness constants $L_f$ and $L_g$ for $f$ and $g$ are known, we can search instead for $\rlamr = \norm{\lamstar} + r$ and utilize $L_f + \rlamr L_g$ as an upper bound to $\Llamr$. The resulting search scheme would have the same complexities, with matching constant dependences, as those of Corollary <ref> and <ref>.

§.§ The Convergence Analysis

We first prove a generic result for the $Q$-function (<ref>) useful for both the non-strongly convex and the strongly convex cases.

Consider an $\alpha$-strongly convex problem of form (<ref>). Let a set of reference multipliers $\Lambda \in \R_+^m$ be given and let the aggregate smoothness constant $\Llam$ be defined in (<ref>). Let iterates $\zt:=\{\xt; \lamtilt, \nut, \pit\}$ be generated by Algorithm <ref>.  Suppose the following stepsize requirements are met.  For all $t \geq 1$, the outer-loop stepsize requirements are

\begin{align}\label{req:ANLA-s-outer}
\begin{split}
& \wt \etat \leq \wtt (\etatt + \alpha/2), \\
&\wt \taut \leq \wtt (\tautt + 1), \\
& \etatt \taut \geq \theta_t\Llam \text{ with }  \thetat :=\wtt /\wt,\\
&\etat[N](\taut[N] + 1) \geq \Llam.
\end{split}
\end{align}

For all $t \geq 1, s \geq 1$, the intra-phase stepsize requirements are

\begin{align}\label{req:intra-phase}
\begin{split}
&\delst (\betast + \alpha/2) \geq \delst[s+1] \betast[s+1],   \\
&\delst \gamst \geq \delst[s+1] \gamst[s+1], \\
&  \gamst \betast[s+1] \geq \rhost[s+1] \norm{\nut}^2,\  \rhost[s+1] = \delst / \delst[s+1],\\
&\gamst[\St[N]][N] (\betast[\St[N]][N] + \alpha/2) \geq  \normsq{\nut[N]}.
\end{split}
\end{align}

For all $t \geq 1$, the inter-phase requirements are

\begin{align}\label{req:inter-phase}
\begin{split}
&\wst[\St] (\betast[\St] + \alpha/2) \geq \wst[1][t+1] \betast[1][t+1],   \\
&\wst[\St] \gamst[\St] \geq \wst[1][t+1] \gamst[1][t+1], \\
& \gamst[\St][t] \betast[1][t+1] \geq \rhost[1][t+1] \norm{\nut[t-1]}^2,\  \rhost[1][t+1] = \wst[\St] / \wst[1][t+1],
\end{split}
\end{align}

where $\wst:= \wt \delst / (\sums \delst) $ denotes the aggregate weights. Then for any reference point $z=(x; \lam,  \nu, \pi) \in X \times \Lam \times [V, \Pi]$ with $[V, \Pi]$ being defined in (<ref>), we have

\begin{align}\label{eq:ANLA-S-Q-conv}
\begin{split}
\sumwt Q(\zt; z) + \tfrac{\wt[N]}{2} (\etat[N]  + \tfrac{\alpha}{2}) \normsq{\xt[N] - x} &\leq \tfrac{\wst[1][1] \betast[1][1] + \wt[1]\etat[1]}{2} \normsq{\xt[0] - x}  + \tfrac{\wst[1][1]\gamst[1][1]}{2} \normsq{\lamst[0][1] - \lam} \\
&+ \wt[1] \taut[1] [U_{f^*}(\pi; \pit[0]) + \lam^\top U_{g^*} (\nu; \nut[0])].
\end{split}
\end{align}

We first establish a convergence bound the inner loop within a phase. Fix $t \geq 1$. Consider the convergence of $\yst$. Since $u(y) + \etat \normsq{y-x^{t-1}}/2$ has a strong convexity modulus of $\alpha + \etat$, the $y$-prox mapping in Line 8 of Algorithm <ref> leads to a three point inequality (see Lemma 3.1 of [8]):

\begin{align*}
&\inner{\yst - x}{\htilst} + u(\yst) - u(x) +  \tfrac{\etat}{2}(\normsq{\yst - \xtt} - \normsq{x - \xtt})\notag\\
&\quad \tfrac{1}{2}[(\betast+\alpha+\etat) \normsq{x - \yst[s]} + \betast \normsq{\yst - \yst[s-1]} - \betast \normsq{\yst[s-1] -x }] \leq 0. \notag
\end{align*}

Equivalently, we have

\begin{align}
&\inner{\yst - x}{\htilst[s]}  +
\tfrac{1}{2}[(\betast + \alpha/2) \normsq{x - \yst[s]} + \betast \normsq{\yst - \yst[s-1]} - \betast \normsq{\yst[s-1] -x }]\notag\\
&\quad  + u(\yst) - u(x) +\tfrac{1}{2}[\etat \normsq{\yst - \xtt} +  (\etat + \alpha/2)\normsq{\yst - x} - \etat\normsq{x - \xtt} ] \leq 0.\label{p2:Q_X}
\end{align}

In particular, the definition of $\htilst[s]$ in Line 7 of Algorithm <ref> implies

\begin{aligned}
\inner{\yst - x}{\htilst[s]} = &\inner{\yst - x}{\sumi \lamist[s] \nuit} - \inner{\yst - x}{\sumi (\lamist[s] - \lamist[s-1]) \nuit} \\
& + \rhost[s] \inner{\yst[s] - \yst[s-1]}{\sumi (\lamist[s-1] - \lamist[s-2]) \nuit} \\
& + \rhost[s] \inner{\yst[s-1] - x}{\sumi (\lamist[s-1] - \lamist[s-2]) \nuit} ,\ \forall s \geq 2,\\
\end{aligned}

\begin{aligned}
\inner{\yst[1]- x}{\htilst[1]} = &\inner{\yst[1]- x}{\sumi \lamist[1] \nuit} - \inner{\yst[1]- x}{\sumi (\lamist[1] - \lamist[0]) \nuit} \\
& + \rhost[1] \inner{\yst[1] - \yst[0]}{\sumi (\lamist[0] - \lamist[-1]) \nuitt}\\
&+  \rhost[1] \inner{\yst[0] - x}{\sumi (\lamist[0] - \lamist[-1]) \nuitt}.
\end{aligned}

So, substituting them into (<ref>), summing up the resulting inequality with weight $\delst$, noting the stepsizes conditions in (<ref>), and utilizing Young's inequality,  we get

\begin{align}\label{p2:y_conv}
\begin{split}
\sum_{s=1}^{\St} &\delst \left(\La(\yst; \lamst, \nut, \pit) - \La(x; \lamst, \nut, \pit) + \tfrac{[\etat \normsq{\yst - \xtt} +  (\etat + \tfrac{\alpha}{2})\normsq{\yst - x} - \etat\normsq{x - \xtt}]}{2}\right)\\
&+ \delst[1] \rhost[1] \inner{\yst[0] - x}{\sumi (\lamist[0] - \lamist[-1]) \nuitt}-\delst[\St]  \inner{\yst[\St]- x}{\sumi (\lamist[\St] - \lamist[\St-1]) \nuit}\\
\leq& \tsum_{s=2}^{\St} \tfrac{\delst[s-1]\gamst[s-1]}{2} \normsq{\lamst[s-1] - \lamst[s-2]} + \tfrac{\delst[1] \rhost[1] \normsq{\nut[t-1]} }{2\betast[1]} \normsq{\lamst[0] - \lamst[-1]} \\
&-\tfrac{1}{2}[\delst[\St] (\betast[\St] + \alpha/2) \normsq{\yst[\St]- x}  - \delst[1] \betast[1] \norm{\yst[0] - x}^2].
\end{split} \end{align}

Next, consider the convergence of $\lamst$. The $\lam$-proximal mapping in Line 9 of Algorithm <ref> implies

\begin{equation*}
\La(\yst; \lam, \nut,  \pi^t) - \La(\yst; \lamst, \nut, \pi^t) + \tfrac{\gamst}{2} [\normsq{\lam- \lamst} + \normsq{\lamst- \lamst[s-1]}- \normsq{\lam- \lamst[s-1]}] \leq 0.
\end{equation*}

Due to the stepsize conditions in (<ref>), the $\delst$ weighted sum satisfies

\begin{equation*}%\label{p2:\lam_conv}
\begin{split}
&\sums \delst[s][t][\La(\yst; \lam, \nut, \pi^t) - \La(\yst; \lamst, \nut, \pi^t)] + \tfrac{\gamst[\St] \delst[\St]}{2} \normsq{\lam- \lamst[\St]}\\
&\quad + \sums \tfrac{\delst \gamst}{2}  \normsq{\lamst[s] - \lamst[s-1]} \leq \delst[1] \gamst[1]\normsq{\lam- \lamst[0]}.
\end{split}
\end{equation*}

Then, combining it with the $y$ convergence bound in (<ref>), we get

\begin{align*}
\begin{split}
\sums& \delst \left(\La(\yst; \lam, \nut, \pi^t) - \La(x; \lamst, \nut, \pi^t)  + \tfrac{{\eta_t} \normsq{\yst - \xtt} + (\eta_t + \alpha/2) \normsq{\yst - x} - \etat \normsq{\xtt -x}}{2}\right)\\
&+\delst[1] \rhost[1] \inner{\yst[0] - x}{\sumi \nuitt (\lam_i^0 - \lamist[-1])} - \delst[\St]\inner{\yst[\St] - x}{\sumi \nuit (\lamist[\St] - \lamist[\St - 1])}\\
\leq& \tfrac{\delst[1]\gamst[1]}{2} \normsq{\lam- \lamst[0]} - \tfrac{\delst[\St] \gamst[\St]}{2} [\normsq{\lam-\lamst[\St]} + \normsq{\lamst[\St]- \lamst[\St-1]}]  +\tfrac{\delst[1] \rhost[1] \normsq{\nut[t-1]}}{2\betast[1]} \normsq{\lam^{0} - \lamst[-1]} \\
&-\tfrac{1}{2}[\delst[\St] (\betast[\St] + \alpha/2) \normsq{\yst[\St]- x}  - \delst[1] \betast[1] \norm{\yst[0] - x}^2].
\end{split}
\end{align*}

Moreover, since $\La(\yst; \lam, \nut, \pi^t)$, $\normsq{\yst - \xtt}$ and $\normsq{\yst - x}$ are convex with respect to $\yst$ and
$\La(x; \lamst, \nut, \pi^t)$ is linear with respect to $\lamst$, multiplying both sides by $\wt/(\sums \delst)$ and applying the Jensen's inequality leads to

\begin{align*}
\begin{split}
\wt& \left(\La(\xt; \lam, \nut, \pi^t) - \La(x; \lamtilt, \nut, \pi^t)  + \tfrac{{\eta_t} \normsq{\xt - \xtt} + (\eta_t + \alpha/2) \normsq{\xt - x} - \etat \normsq{\xtt -x}}{2}\right)\\
&+\wst[1] \rhost[1] \inner{\yst[0] - x}{\sumi \nuitt (\lam_i^0 - \lamist[-1])} - \wst[\St]\inner{\yst[\St] - x}{\sumi \nuit (\lamist[\St] - \lamist[\St - 1])}\\
&\leq \tfrac{\wst[1]\gamst[1]}{2} \normsq{\lam- \lamst[0]} - \tfrac{\wst[\St] \gamst[\St]}{2} [\normsq{\lam-\lamst[\St]} + \normsq{\lamst[\St]- \lamst[\St-1]}]  +\tfrac{\wst \rhost[1] \normsq{\nut[t-1]}}{2\betast[1]} \normsq{\lam^{0} - \lamst[-1]} \\
&\ -\tfrac{1}{2}[\wst[\St] (\betast[\St] + \alpha/2) \normsq{\yst[\St]- x}  - \wst[1] \betast[1] \norm{\yst[0] - x}^2],
\end{split}
\end{align*}

where $\wst = \wt \delst/ (\sums \delst)$ represents the aggregate weight for the inner iterates.

Next, we consider the inner loops from different phases. The inter-phase stepsize condition in (<ref>) implies the sum of preceding inequality across $t$ satisfies

\begin{equation}\label{eq:sliding-QxQp}
\begin{split}
\sumt &{\wt} [Q_x(\zt; z)  + Q_\lam(\zt; z)] + \sumt \tfrac{\wt}{2} \etat\normsq{\xt - \xtt} + \wt[N](\etat[N] + \alpha/2) \normsq{\xt[N] - x} \\
&\leq \tfrac{\wst[1][1]}{2}(\gamst[1][1]\normsq{\lamst[1][0] - \lam} + \betast[1][1] \normsq{\yst[0][1] - x}) + \tfrac{\wt[1] \etat[1]}{2} \normsq{\xt[0] - x}.
\end{split}
\end{equation}

Observe that (<ref>) is almost identical to the $Q_x$ and $Q_\lam$ inequality in (<ref>). Thus a similar argument to Proposition <ref> and the outer-loop stepsize requirements in (<ref>) leads to the desired convergence result in <ref>.

Now we leverage the preceding proposition to prove the convergence of the ACGD-S method under the non-strongly convex setting.

Proof to Theorem <ref> and Corollary <ref>:
It is straightforward to verify that the outer loop stepsize in (<ref>) satisfies the condition (<ref>), and the adaptive inner loop stepsize in (<ref>) satisfies both the intra-phase condition (<ref>) and the inter-phase condition (<ref>). So it follows from (<ref>) that

\begin{align}\label{pfeq:thm-S-non-1}
\begin{split}
\sumwt Q(\zt; z) + \tfrac{\wt[N]\etat[N]}{2} \normsq{\xt[N] - x} &\leq \tfrac{\wst[1][1] \betast[1][1] + \wt[1]\etat[1]}{2} \normsq{\xt[0] - x}  + \tfrac{\wst[1][1]\gamst[1][1]}{2} \normsq{\lamst[0][1] - \lam}.
\end{split}
\end{align}

Now consider setting the reference point to $z^* = (\xstar; \lamstar, \nustar = \grad g(\xstar), \pistar= \grad f(\xstar))$ such that $Q(\zt, \zstar) \geq 0\ \forall t$. The preceding inequality implies that
$$\normsq{\xt[N] - \xstar} \leq \tfrac{1}{\Llamr[c]}[({\wst[1][1] \betast[1][1] + \wt[1]\etat[1]}) \normsq{\xt[0] - x}  + {\wst[1][1]\gamst[1][1]} \normsq{\lamst[0][1] - \lam}]\ \forall N \geq 1.$$
So $\xundert$, being the convex combination of $\{\xt\}$s, remains in a bounded ball around $\xstar$, and $\nut = \grad g(\xundert)$ is bounded for all $t \geq 1.$

Next, setting the ergodic average solution as $\zbart:= (\xbarN; \sumt \wt \lamtilt/ (\sumt\wt),\nubarti[N][], \pibart[N])$  with

\begin{equation}\label{eq:ergo_sol}
\pibart[N]:= \sumwt \pit / (\sumwt),\nutibar[N]:=  \begin{cases}\sumwt \lamti \nuit /(\sumwt \lamti) & \text{o.w.} \\
\grad g_i(x^0) &\text{if } \lamti = 0\ \forall t, \end{cases}\end{equation}

a similar application of the Jensen's inequality as (<ref>) and the stepsize choice in (<ref>) and (<ref>) lead to

\begin{equation}\label{thm4:Q-ave}
Q(\zbart[N]; (\xstar; \lam, \nu, \pi))  \leq \tfrac{1}{N(N+1)} (\tfrac{2\rlamr[c]\norm{\xt[0] - \xstar}}{\Delta} + \Llamr[c] \normsq{\xt[0] - \xstar}), \forall \lam \in  \Lamr[c], (\nu,\pi) \in [V, \Pi].  \end{equation}

The convergence in both the optimality gap and the feasibility violation in (<ref>) then follows from Lemma <ref>.

Consider now the given choice of $\Delta$ in Corollary <ref>. It follows from (<ref>) that at most $N_\ep = \ceil{\sqrt{\tfrac{3 \Llamr[c]}{\ep}}\norm{\xt[0] - \xstar}}$ phases are required to find an $(\ep; \ep/c)$-optimal solution.
Since each phase  requires one gradient evaluation for $f$ and $g$, the oracle complexity of the ACGD-S method is $N_\ep$. Moreover, since each inner iteration requires less than three matrix-vector multiplication,  the total number of matrix-vector multiplication across $N_\ep$ phases can bounded as

$$C_\ep = 3 \tsum_{t=1}^{N_\ep}S_t \leq N_\ep + N_\ep^2 \Mbar \Delta =  \bigO\{\sqrt{\tfrac{\Llamr[c]}{\ep}}\norm{\xt[0] - \xstar} + \tfrac{ \rlamr[c] \Mbar\norm{\xt[0] - \xstar}}{\ep}\}.$$

The next proof considers the strongly convex case.

Proof to Theorem <ref>, Corollary <ref> and Corollary <ref>:
It is straightforward to check that the outer-loop stepsize in (<ref>) satisfies the condition (<ref>), and the adaptive inner-loop stepsize in (<ref>) and (<ref>) satisfy the intra-phase condition (<ref>). Now we verify the inter-phase condition in (<ref>). Consider a fixed $t\geq 2$, we have

\begin{align*}
&\wst[1][t] = \tfrac{\wt\delst[1]}{\sums \delst} = \tfrac{\wt \delst[1]}{\wt \MtU^2 \Gamt \Delta } = \tfrac{\delst[1]}{\MtU^2\Gamt \Delta} = \tfrac{\Wt}{\MtU \Gamt \Delta},\\ &\wst[\St[t-1]][t-1]= \tfrac{\wt[t-1]\delst[\St[t-1]][t-1]}{\sums[\St[t-1]] \delst[s][t-1]} = \tfrac{\wt[t-1] \delst[\St[t-1]][t-1]}{\wtt \MtU[t-1]^2 \Gamt[t-1] \Delta } = \tfrac{\delst[\St[t-1]][t-1]}{\MtU[t-1]^2\Gamt[t-1] \Delta} = \tfrac{\Wt}{\MtU[t-1] \Gamt[t-1] \Delta}.
\end{align*}

\begin{align*}
\wst[\St[t-1]][t-1] (\betast[\St-1][t-1] + \alpha/2 ) &= \tfrac{\Wt}{\MtU[t-1] \Gamt[t-1] \Delta} \tfrac{\alpha}{4}(\Wt[t-1] \Gamt[t-1] \MtU[t-1] + \St) \geq  \tfrac{\Wt}{\MtU[t-1] \Gamt[t-1] \Delta} \tfrac{\alpha}{4}(\delst[\St[t-1]][t-1]\Gamt[t-1])  \\
& = \tfrac{\Wt}{\MtU[t-1] \Gamt[t-1] \Delta} \tfrac{\alpha}{4}(\MtU[t-1] \Wt \Gamt[t-1])= \tfrac{\alpha}{4} \tfrac{\Wt^2}{\Delta} = \tfrac{\Wt}{\MtU \Gamt \Delta} [\tfrac{\alpha}{4} \MtU \Gamt \Wt]= \wst[1]\betast[1].
\end{align*}

\begin{align*}
\wst[\St[t-1]][t-1] \gamst[\St[t-1]][t-1] = \tfrac{\delst[\St[t-1]][t-1]}{\MtU[t-1]^2\Gamt[t-1] \Delta} (\tfrac{4}{\alpha}) \tfrac{\MtU[t-1]^2 \Gamt[t-1] }{\delst[\St[t-1]][t-1]} = \tfrac{4}{\alpha \Delta} = \tfrac{\delst[1]}{\MtU^2\Gamt \Delta} \tfrac{\MtU[t]^2 \Gamt[t] }{\delst[1][t]} = \wst[1][t] \gamst[1].
\end{align*}

\begin{align*}
\betast[1] \gamst[\St[t-1]][t-1]  = (\tfrac{\alpha}{4} \MtU \Gamt \Wt) (\tfrac{4}{\alpha}) \tfrac{\MtU[t-1]^2 \Gamt[t-1] }{\MtU[t-1] \Wt} =  \MtU[t]\MtU[t-1] \Gamt \Gamt[t-1] \stackrel{(a)}{\geq} \tfrac{\MtU \Gamt}{\MtU[t-1] \Gamt[t-1]} \MtU[t-1]^2 = \tfrac{\wst[\St[t-1]][t-1]}{\wst[1]} \MtU[t-1]^2 = \rhost[1] \MtU[t-1]^2,
\end{align*}

where  the inequality in $(a)$ holds because we have $\Gamt[t-1] \geq 1$ as a consequence of its definition.
Thus all the requirements in Proposition <ref> are satisfied. We get from (<ref>) that

\begin{align}\label{pfeq:thm-S-str-1}
\begin{split}
\sumwt Q(\zt; z) + \tfrac{\wt[N]\etat[N]}{2} \normsq{\xt[N] - x} &\leq \tfrac{\wst[1][1] \betast[1][1] + \wt[1]\etat[1]}{2} \normsq{\xt[0] - x}  + \tfrac{\wst[1][1]\gamst[1][1]}{2} \normsq{\lamst[0][1] - \lam}.
\end{split}
\end{align}

Similar arguments as that of Theorem <ref> imply the boundedness of $\MtU$, and that the ergodic average solution $\zbart$ defined according to (<ref>) satisfies

\begin{equation}\label{thm5:Q-ave}
Q(\zbart[N]; (\xstar; \lam, \nu, \pi))  \leq \tfrac{1}{\sumwt} (\tfrac{2\rlamr\norm{\xt[0] - \xstar}}{\Delta} + \Llamr[r] \normsq{\xt[0] - \xstar}), \forall \lam \in  \Lamr, (\nu,\pi) \in [V, \Pi].  \end{equation}

Since $\sumwt \geq \max\{N(N+1)/2, \sqrt{2\kar}[(1+ 1/\sqrt{2\kar})^4 -1]\}$ (see (<ref>)), we get the optimality gap and feasibility violation convergence bound in (<ref>). Moreover, since since $\tfrac{\alpha}{2} \normsq{\xbart[N] - \xstar} \leq  Q(\zbart; (\xstar; \lamstar, \grad g(\xbart[N]), \grad f(\xbart[N])))$ (see (<ref>)), the convergence of $\xbart[N]$ to $\xstar$ in (<ref>) also follows from (<ref>).

Next we show Corollary <ref>. Fix a $r\geq 1$, and set $\Delta = \tfrac{[\rlamr]^2}{\Llamr \normsq{\xt[0] - \xstar} \alpha}$ . For any

$$N \geq \min\{\sqrt{\tfrac{6 \Llamr \max\{c/r, 1\}}{\ep}} \normsq{\xt[0] - \xstar}, [\sqrt{\tfrac{2\Llamr[r]}{\alpha}} + 1]\log(\tfrac{3\max\{c/r,1\}\sqrt{\Llamr[r]\alpha}\normsq{\xt[0] - \xstar}}{\ep} + 1) + 4\},$$

we get

\begin{align*}
Q(\zbart[N]; (\xstar; \lam, \nu, \pi)) \leq \min\{1, \tfrac{r}{c}\}\ep, \forall \lam \in  \Lamr[r], (\nu,\pi) \in [V, \Pi], \\
\Rightarrow \max\{F(\xbart[N]) - \Fstar, r \norm{[g(\xbart[N])]_+}\}  \leq \min\{1, \tfrac{r}{c}\}\ep,
\end{align*}

so that $\xbart[N]$ is an $(\ep; \ep/c)$ solution. Therefore, the least number of phases $\Nep$ required for such a solution admits the upper bounded in the corollary statement.

Now we consider the corresponding number of matrix-vector multiplications $\Cep$ in the $\Nep$ phases, i.e., $\Cep = \sumt[\Nep] \St.$ We can deduce from the preceding argument that

$$\sumt[\Nep] \wt \leq  \tfrac{6 \max\{1, \tfrac{c}{r}\} \Llamr[r]\normsq{\xt[0] - \xstar}}{\ep}.$$

Calculating the sum of $\St$ directly is challenging, so we consider an easier quantity $\Rt := [\St - 3]_+.$ An useful algebraic relation is

\begin{equation}\label{pfeq:Rt-lower-bound}
\tsum_{i=1}^{t-1} \Rt[i] \leq \Mbar \Wt.\end{equation}

The result can be deduced by  induction.
For $t=2$, the relation clearly holds since $\Wt[2] \Mbar \geq \Wt[2] \MtU[1] = \St[1] \geq \Rt[1]$. Assuming (<ref>) is valid up to $t\geq 2$, we have

\begin{align*}
\Wt[t+1] = \tfrac{\Wt \MtU + \tfrac{1}{\Gamt}(\St-1)}{\MtU} \stackrel{(a)}{\geq} \Wt[t] + \tfrac{[\St[t] - 3]_+}{\MtU} \stackrel{(b)}{\geq} \tfrac{1}{\Mbar}\{\tsum_{i=1}^{t-1} \Rt[i] + [\St[t] -3]_+\} = \tfrac{1}{\Mbar} \tsum_{i=1}^{t-1}[\Rt[i]],
\end{align*}

where (a) follows from the algebraic fact $\tfrac{\St}{\Gamt}\geq \St -2$ (see Lemma <ref>) and that $\St - 1\geq 0$ and $(b)$ follows from the induction hypothesis. Thus the principle of mathematical induction implies that (<ref>) is valid. Consequently, for $\Rt[\ep] = \sumt[\Nep] \Rt$,  we have

\begin{align*}
\Rt[\ep]^2/2 &\leq  \sums[\Rt[\ep]] s = \sumt[\Nep] \sums[\Rt] [(\tsum_{j=1}^{t-1} \Rt[j]) + s]
\leq \sumt[\Nep] \sums[\Rt] [\Mbar \Wt + s] = \sumt[\Nep] \sums[\St-3] [\Mbar \Wt +  s]\\
&= \sumt[\Nep] \tsum_{t=2}^{\St-2} [\Mbar \Wt +  (s-1)] \leq \sumt[\Nep] \tsum_{t=1}^{\St-1} [\Mbar \Wt +  (s-1)] \leq \Mbar \sumt[\Nep] \tsum_{t=1}^{\St-1} [\Wt +  \tfrac{(s-1)}{\MtU}]\\
&\stackrel{(a)}{\leq} \Mbar \sumt[\Nep] \MtU \wt \Delta  \leq \sumt[\Nep] \Mbar^2 \wt  \Delta \leq   \tfrac{6 \max\{1, \tfrac{c}{r}\} \Llamr[r]\normsq{\xt[0] - \xstar}}{\ep} \tfrac{[\rlamr]^2}{\Llamr \normsq{\xt[0] - \xstar} \alpha}\Mbar^2 \\
&\leq \tfrac{6 \max\{1, \tfrac{c}{r}\} \Mbar^2 [\rlamr]^2}{\ep \alpha},
\end{align*}

where (a) follows from the fact that $\sums[\St-1] \Wt \MtU + (s-1) < \wt \MtU^2 \Delta$ (see (<ref>)). Therefore we get $\Rt[\ep] \leq \sqrt{\tfrac{12}{\alpha\ep}} \Mbar \rlamr$. Since $\sumt[\Nep] \St \leq \Rt[\ep] + 3 \Nep$, the big-O bound on $\Cep$ follows immediately.

The complexity bounds in Corollary <ref> can be derived similarly.

§ CONCLUSION

To sum up, this paper proposes two efficient methods for large-scale function-constrained optimization. The simple ACGD method has the optimal oracle complexity, but it requires access to a QP solver. The more complicated ACGD-S method has both the optimal oracle complexity and the optimal computation complexity. Lower complexity bounds are provided to illustrate the oracle complexity of both ACGD and ACGD-S to be unimprovable for a general case of first-order methods. Together they provide a complete characterization of the difficulty of solving a smooth function-constrained optimization problem from both the oracle complexity and the computation complexity perspective.

[1]
A. Beck, First-order methods in optimization, vol. 25, SIAM, 2017.

[2]
D. Bertsekas, Convex optimization theory, vol. 1, Athena
Scientific, 2009.

[3]
D. Boob, Q. Deng, and G. Lan, Stochastic first-order methods for
convex and nonconvex functional constrained optimization, 2019.

[4]
D. Boob, Q. Deng, and G. Lan, Stochastic first-order methods for
convex and nonconvex functional constrained optimization, Mathematical
Programming,  (2022), pp. 1–65.

[5]
R. Gandy, Portfolio optimization with risk constraints, PhD thesis,
Universität Ulm, 2005.

[6]
E. Y. Hamedani and N. S. Aybat, A primal-dual algorithm with line
search for general convex-concave saddle point problems, SIAM Journal on
Optimization, 31 (2021), pp. 1299–1329.

[7]
G. Lan, Gradient sliding for composite optimization, Mathematical
Programming, 159 (2016), pp. 201–235.

[8]
height 2pt depth -1.6pt width 23pt, First-order and
stochastic Optimization Methods for Machine Learning, Springer-Nature, 2020.

[9]
G. Lan and R. D. Monteiro, Iteration-complexity of first-order
penalty methods for convex programming, Mathematical Programming, 138
(2013), pp. 115–139.

[10]
G. Lan and Y. Ouyang, Accelerated gradient sliding for structured
convex optimization, Computational Optimization and Applications,  (2020).
under revision.

[11]
G. Lan, Y. Ouyang, and Y. Zhou, Graph topology invariant gradient
and sampling complexity for decentralized and stochastic optimization, arXiv
preprint arXiv:2101.00143,  (2021).

[12]
G. Lan and Z. Zhang, Optimal methods for risk averse distributed
optimization, arXiv preprint arXiv:2203.05117,  (2022).

[13]
Q. Lin, S. Nadarajah, and N. Soheili, A level-set method for convex
optimization with a feasible solution path, SIAM Journal on Optimization, 28
(2018), pp. 3290–3311.

[14]
A. Nemirovski, Lectures on modern convex optimization, in Society
for Industrial and Applied Mathematics (SIAM, Citeseer, 2001.

[15]
A. S. Nemirovsky, On optimality of krylov's information when solving
linear operator equations, Journal of Complexity, 7 (1991), pp. 121–130.

[16]
A. S. Nemirovsky and D. B. Yudin, Problem complexity and method
efficiency in optimization., John Wiley UK/USA, 1983.

[17]
Y. Nesterov, A method for unconstrained convex minimization problem
with the rate of convergence o (1/k 2), in Doklady AN USSR, vol. 269,
1983, pp. 543–547.

[18]
height 2pt depth -1.6pt width 23pt, A method of solving
a convex programming problem with convergence rate o(k2), in Doklady
Akademii Nauk, vol. 269, Russian Academy of Sciences, 1983, pp. 543–547.

[19]
height 2pt depth -1.6pt width 23pt, Introductory
lectures on convex optimization: A basic course, vol. 87, Springer Science
& Business Media, 2003.

[20]
height 2pt depth -1.6pt width 23pt, Lectures on convex
optimization, vol. 137, Springer, 2018.

[21]
Y. Ouyang and Y. Xu, Lower complexity bounds of first-order methods
for convex-concave bilinear saddle-point problems, Mathematical Programming,
185 (2021), pp. 1–35.

[22]
P. Rigollet and X. Tong, Neyman-pearson classification, convexity
and stochastic constraints, Journal of Machine Learning Research,  (2011).

[23]
Y. Xu, First-order methods for problems with O(1) functional
constraints can have almost the same convergence rate as for unconstrained
problems, SIAM Journal on Optimization, 32 (2022), pp. 1759–1790.

[24]
S. Yang, X. Li, and G. Lan, Data-driven minimax optimization with
expectation constraints, arXiv preprint arXiv:2202.07868,  (2022).

[25]
M. B. Zafar, I. Valera, M. G. Rogriguez, and K. P. Gummadi, Fairness
constraints: Mechanisms for fair classification, in Artificial intelligence
and statistics, PMLR, 2017, pp. 962–970.

[26]
Z. Zhang, S. Ahmed, and G. Lan, Efficient algorithms for
distributionally robust stochastic optimization with discrete scenario
support, SIAM Journal on Optimization, 31 (2021), pp. 1690–1721.

[27]
Z. Zhang and G. Lan, Optimal algorithms for convex nested stochastic
composite optimization, arXiv preprint arXiv:2011.10076,  (2020).

§ APPENDIX

Let $\Lambda \subset \R^m_+$  and a convex vector-valued function $g:\R^n \rightarrow \R^m $ be given. If $\sumi \lam_i g_i$ is $L$-smooth for all $\lam \in \Lam$, i.e., $\norm{\sumi \lami \grad g_i(x) - g_i(\xbar) } \leq L \norm{x - \xbar} \forall x, \xbar \in R^n, \forall \lam \in \Lam$, the Bregman distance function generated by its (vector-valued) conjugate function $U_g^*$ satisfies

$$ \inner{\lam}{U_g^*(\nu, \nubar)} \geq \normsq{\sumi \lami (\nui - \nubar_i)}/(2L)\ \forall \nu, \nubar \in \{\grad g(x): x \in \R^n\}.$$

Let  $\lam \in \Lam$,  $\nu = \grad g(x),$ and $\nubar = \grad g(\xbar) $ be given. Consider the function $\barg := \sumi \lami g_i $. Clearly,  $\barg$ is $L$-smooth such that the Bregman distance function generated by its conjugate satisfies $U_{\barg^*}(v, \bar{v}) \geq \normsq{v - \bar{v}}/(2 L).$ Since $\grad \barg(x) = \lambda^\top \grad g(x)$, we have

\begin{align*}
\lambda^\top U_{g^*}(\nu; \nubar) &\stackrel{(a)}{=} \sumi \lami U_{g^*_i}(\nu_i; \nubar_i)=\sumi \lami U_{g_i}(\xbar; x) \\
&= \sumi \lami[i] [g_i(\xbar) - g_i(x)  - \inner{\grad g_i(x)}{\xbar - x}]\\
&= \barg(\xbar) - \barg(x) - \inner{\grad \barg(x)}{\xbar - x}\\
&\stackrel{(b)}{=} U_{\barg}(\xbar; x)
= U_{\barg^*}(\lam^\top \nu; \lam^\top \nubar) \geq \normsq{\sumi \lami (\nui - \nubar_i)}/(2L),
\end{align*}

where (a) and (b) follows from the algebraic identity between Bregman distance functions generated by Fenchel conjugates $(h, h^*)$, $U_h(y;\bar y) = U_{h^*}(\grad h(\bar y); \grad h(y)).$

The next two lemmas provide some basic algebraic identifies useful for deriving complexity bounds.

Given an $x > 0$, the following algebraic relation is valid:

\begin{equation}\label{al-fact:wt}
h(y):= (1 + 1/x)^{y -3} \leq y\  \forall\  2x \geq  y\geq 1.
\end{equation}

First, we show the relation for $y=2x$, i.e. $h(2x)= (1 + 1/x)^{2x -3} \leq 2x\  \forall x >0.$
Let's consider two cases. If $ x \geq 4$, we have

$$(1  + 1/x)^{2x -3} \leq [(1+ 1/x)^x]^2 \leq \exp(2) \leq 8 \leq 2x.$$

If $0<x < 4$, we have

$$ [(1 + 1/x)^{x -1.5}]^2 \leq [1/(1 - \tfrac{x-1.5}{x})]^2 = (x/1.5)^2 = (x/2.25) \times x \leq 2 x.$$

Thus $h(2x) \leq 2x$.
Since $h(1) \leq 1$, the relation in (<ref>) follows from the convexity of $h$ with respect to $y$.

Given non-negative parameters $\Delta_t$, $H$, and $h \in \{0,1\}$, suppose $\St = \min\{S\in \mathbb{N}_+: \sums[S] H + (s-h) \geq \Delta_t\}$, and $\Gamma$ is the non-negative root of $ \sums[\St]  [\Gamma H + (s-h)] = \Gamma^2\Delta_t$, then $\St$ satisfies

$$\St / \Gamma \geq \St -2.$$

Suppose for the sake of contradiction that $\St / \Gamma < \St -2$. On the one hand, the definition of $\St$ implies that

\begin{align*}
\Delta_t = H (\tfrac{\St}{\Gamma}) + \sums \tfrac{(s-h)}{\Gamma^2} = H (\tfrac{\St}{\Gamma}) +  \tfrac{1}{2}\tfrac{(\St-h+1)}{\Gamma} \tfrac{\St}{\Gamma} \geq H (\tfrac{\St}{\Gamma}) +  \tfrac{1}{2}(\tfrac{\St}{\Gamma})^2.
\end{align*}

On the other hand, the choice of $\St$ implies that $\sums[\St-1] [H + (s-h)] < \Delta_t$, thus

\begin{align*}
&H (\St - 1) + \tfrac{1}{2} (\St - 1)(\St-h -1) < \Delta_t\\
\Rightarrow&\ H (\tfrac{\St}{\Gamma}) + \tfrac{1}{2} (\tfrac{\St}{\Gamma})^2 < \Delta_t.
\end{align*}

These two relations leads to the desired contradiction.