# Implicit augmented Lagrangian and generalized optimization

Alberto De Marchi University of the Bundeswehr Munich, Department of Aerospace
Engineering, Institute of Applied Mathematics and Scientific Computing, 85577
Neubiberg, Germany<EMAIL_ADDRESS>orcid0000-0002-3545-6898.

###### Abstract

Generalized nonlinear programming is considered without any convexity
assumption, capturing a variety of problems that include nonsmooth objectives,
combinatorial structures, and set-membership nonlinear constraints. We extend
the augmented Lagrangian framework to this broad problem class, preserving an
implicit formulation and introducing slack variables merely as a formal
device. This, however, gives rise to a generalized augmented Lagrangian
function that lacks regularity, due to the marginalization with respect to
slack variables. Based on parametric optimization, we develop a tailored
stationarity concept to better qualify the iterates, generated as approximate
solutions to a sequence of subproblems. Using this variational
characterization and the lifted representation, a suitable multiplier update
rule is derived, and then asymptotic properties and convergence guarantees are
established for a safeguarded augmented Lagrangian scheme. An illustrative
numerical example showcases the modelling versatility gained by dropping
convexity assumptions and indicates the practical benefits of the advocated
implicit approach.

Keywords. Nonsmooth nonconvex optimization, generalized nonlinear programming,
multiplier methods, generalized augmented Lagrangians, proximal algorithms

AMS subject classifications. 65K10, 90C26, 65K05

##  1 Introduction

Optimization problems often bring plenty of mathematical structure to work
with. In this paper we investigate the setting of generalized programming
without any convexity, namely finite-dimensional optimization problems of the
form

$\operatorname*{minimize}_{x\in\mathbb{R}^{n}}{}\quad{}\varphi(x):=f(x)+g(c(x)),$
(P)

where $x$ is the decision variable, $f$ and $c$ are smooth mappings and $g$ is
merely lower semicontinuous. Problem data $f$ and $g$ are allowed to be
nonconvex mappings, the nonsmooth cost $g$ is not necessarily continuous nor
real-valued, and the mapping $c$ is potentially nonlinear. The composition of
a nonsmooth term $g$ with a smooth mapping $c$ results in great modeling
flexibility and encompasses a broad spectrum of problems [36]; both
regularizers and constraints can be encoded, as well as combinations thereof
[9, 7, 13]. Prominent special cases of (P) are nonlinear [4] and disjunctive
[2] programming, structured optimization [32, 40], and constrained structured
optimization [15]. The abstract problem (P) is also referred to as “extended
nonlinear programming” in [36], promoting better representation of structures
found in applications by featuring a composite term, beyond the conventional
format of nonlinear programming.

A fundamental technique for solving _constrained_ optimization problems is the
augmented Lagrangian (AL) framework [35, 4, 6]. Introducing an auxiliary
variable $z\in\mathbb{R}^{m}$, (P) can be equivalently rewritten in the form

$\operatorname*{minimize}_{x\in\mathbb{R}^{n},\
z\in\mathbb{R}^{m}}{}\quad{}f(x)+g(z){}\qquad{}\operatorname{subject~{}to}{}\quad{}c(x)-z=0,$
(P${}_{\text{S}}$)

which has a simpler, separable objective function but is subject to some
(explicit) constraints. It has been observed in [15, Lemma 3.1] that
incorporating slack variables in this manner does not affect minimizers and
stationary points when compared to the original (P) — somewhat remarkably in
light of [3]. Constrained minimization problems such as (P${}_{\text{S}}$) are
amenable to be addressed by means of AL methods, which we embrace as they can
effortlessly handle nonsmoothness [18, 19]. As the AL approach constitutes a
framework, rather than a single algorithm, several methods have been presented
in the past decades, expressing the foundational ideas in different flavors.

Since their inception [24, 33, 34], AL schemes have been deeply and
extensively investigated [10, 4, 20, 6, 21], also in infinite dimensions [26].
In the convex setting, it was recognized that the “method of multipliers” [24,
33] corresponds to solving an associated dual problem via the proximal point
algorithm (PPA) [34, 35], with duality theory playing a key role in this
parallel. Then, sufficient optimality conditions and local duality led to
analogous results in nonlinear programming by Bertsekas [4]. Pivoting on some
kind of local duality, in turn enabled by some local convexity, the AL scheme
can be traced back to an application of the PPA, thereby establishing
convergence.

This pattern was recently followed by Rockafellar [37, 38] to extend these
ideas to the much broader setting of (P) with $g$ convex. In contrast with
this line of work, which goes back at least to [36, Section 2], we are
interested here in investigating AL methods for generalized nonlinear
programming without any convexity assumption. With $g$ potentially nonconvex,
the abstract setting (P) not only covers a variety of problems, but it also
provides a more versatile and expressive modeling paradigm, with no need for
“suitable reformulations”. Dropping the convexity of $g$ comes at a price,
though, as it appears more difficult to leverage the perspective of PPA.
Nonetheless, the shifted-penalty approach underpinning the seminal “method of
multipliers” still applies, and can be exploited as such. Yet difficulties
arise as the nonconvexity of $g$ translates into potential set-valuedness of
its proximal mapping, which in turn leads to a lack of regularity. Here, we
seek a better understanding of the consequences on the analysis and design of
numerical methods for (P), concerning theoretical properties, convergence
guarantees, and practical aspects.

The following blanket assumptions are considered throughout, without further
mention. Technical definitions are given in Section 2.1.

######  Assumption 1.

The following hold in (P): (i) $f\colon\mathbb{R}^{n}\to\mathbb{R}$ and
$c\colon\mathbb{R}^{n}\to\mathbb{R}^{m}$ are continuously differentiable; (ii)
$g\colon\mathbb{R}^{m}\to\overline{\mathbb{R}}$ is proper, lower
semicontinuous, and prox-bounded; (iii) $\inf\varphi\in\mathbb{R}$.

Assumption 1(iii) is relevant to the well-posedness of the problem, as the
cost function $\varphi:=f+g\circ c$ has nonempty domain and is bounded from
below. Prox-boundedness in Assumption 1(ii) guarantees that the proximal
mapping of $g$ is well-defined, for some suitable stepsizes. Moreover, we work
under the practical assumption that (only) the following computational oracles
are available or simple to evaluate:

  * •

given $x\in\mathbb{R}^{n}$, cost function value $f(x)\in\mathbb{R}$ and
gradient $\nabla f(x)\in\mathbb{R}^{n}$;

  * •

given $v\in\mathbb{R}^{m}$ and $\gamma\in(0,\gamma_{g})$, with $\gamma_{g}$
being the prox-boundedness threshold of $g$, proximal point
$z\in\operatorname{prox}_{\gamma g}(v)$ and cost function value
$g(z)\in\mathbb{R}$ therein;

  * •

given $x\in\mathbb{R}^{n}$ and $v\in\mathbb{R}^{m}$, constraint function value
$c(x)\in\mathbb{R}^{m}$ and transposed Jacobian-vector product
$c^{\prime}(x)^{\top}v\in\mathbb{R}^{n}$, where
$c^{\prime}\colon\mathbb{R}^{n}\to\mathbb{R}^{m\times n}$ denotes the Jacobian
of $c$.

Relying only on these primitives, the numerical method presented in the
following is first-order and matrix-free by construction; as such, it involves
only simple operations and has low memory footprint.

#### Outline

Related works and topics are mentioned in Section 1.1, followed by our
contributions in Section 1.2. Notation, preliminary material, and useful tools
are given in Section 2.1, before introducing relevant optimality and
stationarity notions in Section 2.2. The proposed augmented Lagrangian scheme
is detailed and characterized in Section 3. A solver for the inner AL problems
is developed and analysed in Section 4. Finally, an illustrative numerical
example is reported in Section 5, before some concluding remarks in Section 6.

###  1.1 Related Work

The template (P) covers many classical problem classes. These include
structured two-term problems $f+g$ [32] and three-term problems $f+g+h$, as
well as settings with the nonsmooth term composed with a linear operator, with
two-term $g+h\circ L$, three-term problems $f+g+h\circ L$ [28], or with linear
constraints as for the alternating direction method of multipliers [7]. Most
of these problem classes have been extensively studied in the convex setting,
while nonconvexity has been mostly considered for structured two-term problems
$f+g$ only [40].

In conjunction with AL schemes, non-Lipschitzian or merely lower
semicontinuous cost functions have been recently considered in [8, 19, 15].
The issue of attentive convergence, defined in Section 2.1 and seemingly
necessary to perform limiting analysis, was overlooked in [8] and brought to
attention in [15]. Since it is difficult to attest a priori, as a sufficient
to guarantee attentive convergence one may additionally assume $g$ continuous
relative to its domain. While much weaker than convexity, this assumption is
satisfied by many practical problems, as pointed out in [16].

Another issue is that of well-definedness of the AL algorithm. In fact, the AL
subroblems could be ill-posed, as there is no a priori guarantee that the AL
function is bounded from below, due to the constraints relaxation [6]. As we
are concerned with the role of $g$, in this work, we neglect this possibility
and assume that all subproblems are well-posed. Several approaches for moving
beyond this assumption may be applicable, although possibly imposing other
requirements. Among these are localized minimizations by exploiting non-
relaxable constraints [1] or a trust-region included in the subproblems, as
well as different AL schemes [21, 19].

Constrained structured optimization [15] entails the minimization of the sum
$f(x)+g(x)$ over $x\in\mathbb{R}^{n}$ subject to set-membership constraints
$c(x)\in D$, where $c$ is a smooth mapping and $D$ a nonempty closed set. As
such, it is more general than nonlinear and disjunctive programming [4, 2],
but covered by the template (P). Conversely, possibly introducing some
auxiliary variables, (P) can always be reformulated as a constrained
structured program. Thus, the problem class investigated in [15] is as general
as (P) under Assumption 1, but our solution strategy here is essentially
different.

###  1.2 Approach and Contributions

The core idea of this work traces back to the author’s doctoral dissertation
[13, Chapter 1]. With the aim of improving the analysis and patching weak
reasoning, especially resolving the issue of a suitable termination criterion
for the subsolver, further developments took a different spirit and eventually
became [15]. Therein slack variables are introduced as in (P${}_{\text{S}}$),
making the constraints _explicit_ and then applying a traditional AL scheme.
To the contrary, here we maintain the _implicit_ approach advocated for in
[13], offering a rigorous convergence analysis. In particular, we use the
lifted reformulation (P${}_{\text{S}}$) only as a theoretical device to
develop our strategy for solving (P) and investigate its properties via
sharpened variational tools.

Marginalizing the AL function of (P${}_{\text{S}}$) with respect to the slack
variable $z$ yields the generalized (or, implicit) AL function for (P) [18,
13, 38]. However, by the nonconvexity of $g$, the generalized AL function
fails to be continuously differentiable (or prox-friendly) in general. For
overcoming this lack of regularity,

  * •

we devise the concept of $\Upsilon$-stationarity, which allows us to
characterize precisely the sought (approximate) solutions to the AL
subproblems,

thereby better capturing the intended philosophy. Rooted in parametric
optimization, $\Upsilon$-stationarity is strong enough for the outer AL scheme
but still computable and easy to handle numerically. Indeed, we show this
concept provides an intermediate qualification between M-stationarity for the
implicit and explicit formulations.

With the concept of $\Upsilon$-stationarity at hand, we focus on solving the
AL subproblems, namely minimizing the AL function, which constitutes the major
computational task and a key step in AL algorithms. We explore the possibility
to directly solve the AL subproblems, without slack variables. In particular,

  * •

we design and analyze a subsolver that exploits the proximal oracle of $g$ for
minimizing the generalized AL function and, by generating approximate
$\Upsilon$-stationary points, complies with the requirements of the AL scheme.

Reminiscent of _magical steps_ [11, Section 10.4.1], the rationale behind the
subsolver is to include slack variables only as algorithmic by-products, as
opposed to decision variables. These are then adopted within an iterative
scheme to generate updates that yield sufficient improvement and to terminate
with a suitable certificate of (approximate) stationarity.

Overall, the advantages of weaker assumptions and stronger stationarity
concepts would not be attractive without convergence results. Therefore,

  * •

we investigate the asymptotic properties of a safeguarded implicit AL scheme
for (P) under Assumption 1 and establish (global) convergence guarantees on
par with those in classical nonlinear programming.

Moreover, we illustrate how following an implicit approach without introducing
slack variables can result in improved modeling flexibility and practical
performance. Not only does sticking to a compact problem formulation avoid
increasing the problem size, but it also leads to better use of the available
computational oracles and tends to avoid spurious (local) minima and
stationary points.

##  2 Notation and Fundamentals

In this section, we comment on notation, preliminary definitions and useful
results.

###  2.1 Preliminaries

With $\mathbb{R}$ and $\overline{\mathbb{R}}:=\mathbb{R}\cup\\{\infty\\}$ we
denote the real and extended-real line, respectively. The vectors in
$\mathbb{R}^{n}$ with all elements equal to $0$ or $1$ are denoted as $0_{n}$
and $1_{n}$; whenever $n$ is clear from context we simply write $0$ and $1$,
respectively. With $\operatorname{\mathbb{B}}_{r}(x)$ we indicate the closed
ball centered at $x$ with radius $r$. The _effective domain_ of an extended-
real-valued function $h\colon\mathbb{R}^{n}\to\overline{\mathbb{R}}$ is
denoted by $\operatorname{dom}h:=\\{x\in\mathbb{R}^{n}\,|\,h(x)<\infty\\}$. We
say that $h$ is _proper_ if $\operatorname{dom}h\neq\emptyset$ and _lower
semicontinuous_ (lsc) if $h(\bar{x})\leq\liminf_{x\to\bar{x}}h(x)$ for all
$\bar{x}\in\mathbb{R}^{n}$. For some constant $\tau\in\mathbb{R}$,
$\operatorname{lev}_{\leq\tau}h:=\\{x\in\mathbb{R}^{n}\,|\,h(x)\leq\tau\\}$
denotes the $\tau$-_sublevel set_ associated with $h$. Function $h$ is (lower)
_level-bounded_ if for every $\tau\in\mathbb{R}$ its $\tau$-sublevel set is
bounded. Then, we introduce the following notion of _uniform level
boundedness_ , which plays an important role in parametric minimization [39,
Definition 1.16]. A function
$f\colon\mathbb{R}^{n}\times\mathbb{R}^{m}\to\overline{\mathbb{R}}$ with
values $f(x,u)$ is _level-bounded in $x$ locally uniformly in $u$_ if for each
$\alpha\in\mathbb{R}$ and $\bar{u}\in\mathbb{R}^{m}$ there exists
$\varepsilon>0$ such that the set
$\\{(x,u)\,|\,f(x,u)\leq\alpha,\|u-\bar{u}\|\leq\varepsilon\\}$ is bounded in
$\mathbb{R}^{n}\times\mathbb{R}^{m}$. We use the notation
$F\colon\mathbb{R}^{n}\rightrightarrows\mathbb{R}^{m}$ to indicate a point-to-
set function $F\colon\mathbb{R}^{n}\to\mathcal{P}(\mathbb{R}^{m})$, where
$\mathcal{P}(\mathbb{R}^{m})$ is the power set of $\mathbb{R}^{m}$ (the set of
all subsets of $\mathbb{R}^{m}$).

Given a proper and lsc function
$h\colon\mathbb{R}^{n}\to\overline{\mathbb{R}}$ and a point $\bar{x}$ with
$h(\bar{x})$ finite, we may avoid to assume $h$ continuous and instead appeal
to $h$-_attentive_ convergence of a sequence $\\{x^{k}\\}$, denote as
$x^{k}\overset{h}{\to}\bar{x}$ and given by $x^{k}\to\bar{x}$ with
$h(x^{k})\to h(\bar{x})$. Following [39, Definition 8.3], we denote by
$\hat{\partial}h\colon\mathbb{R}^{n}\rightrightarrows\mathbb{R}^{n}$ the
_regular subdifferential_ of $h$, where

$v\in\hat{\partial}h(\bar{x}){}\quad:\Leftrightarrow\quad{}\liminf_{\begin{subarray}{c}x\to\bar{x}\\\
x\neq\bar{x}\end{subarray}}\frac{h(x)-h(\bar{x})-\langle
v,x-\bar{x}\rangle}{\|x-\bar{x}\|}\geq 0.$

The (limiting) _subdifferential_ of $h$ is $\partial
h\colon\mathbb{R}^{n}\rightrightarrows\mathbb{R}^{n}$, where $v\in\partial
h(\bar{x})$ if and only if $\bar{x}\in\operatorname{dom}h$ and there exist
sequences $\\{x^{k}\\}$ and $\\{v^{k}\\}$ such that
$x^{k}\overset{h}{\to}\bar{x}$ and $v^{k}\in\hat{\partial}h(x^{k})$ with
$v^{k}\to v$. The subdifferential of $h$ at $\bar{x}$ satisfies
$\partial(h+h_{0})(\bar{x})=\partial h(\bar{x})+\nabla h_{0}(\bar{x})$ for any
$h_{0}\colon\mathbb{R}^{n}\to\overline{\mathbb{R}}$ continuously
differentiable around $\bar{x}$ [39, Exercise 8.8].

Given a parameter value $\gamma>0$, the _proximal_ mapping
$\operatorname{prox}_{\gamma h}$ is defined by

$\operatorname{prox}_{\gamma
h}(x){}:={}\operatorname*{\arg\min}_{z}\left\\{h(z)+\frac{1}{2\gamma}\|z-x\|^{2}\right\\}.$
(1)

We say that $h$ is _prox-bounded_ if it is proper and
$h+\|\cdot\|^{2}/(2\gamma)$ is bounded below on $\mathbb{R}^{n}$ for some
$\gamma>0$ [39, Definition 1.23]. The supremum of all such $\gamma$ is the
threshold $\gamma_{h}$ of prox-boundedness for $h$. In particular, if $h$ is
bounded below by an affine function, then $\gamma_{h}=\infty$. When $h$ is
lsc, for any $\gamma\in(0,\gamma_{h})$ the proximal mapping
$\operatorname{prox}_{\gamma h}$ is locally bounded, nonempty- and compact-
valued [39, Theorem 1.25]. The value function of the minimization problem
defining the proximal mapping is the _Moreau envelope_ with parameter
$\gamma>0$, denoted $h^{\gamma}\colon\mathbb{R}^{n}\to\mathbb{R}$, namely

$h^{\gamma}(x):=\inf_{z}\left\\{h(z)+\frac{1}{2\gamma}\|z-x\|^{2}\right\\}.$

Some tools of variational analysis will be exploited in order to describe the
geometry of some nonempty, closed, but not necessarily convex, set
$D\subset\mathbb{R}^{m}$. The _projection_ mapping
$\operatorname{\Pi}_{D}\colon\mathbb{R}^{m}\rightrightarrows\mathbb{R}^{m}$
and the _distance_ function
$\operatorname{dist}_{D}\colon\mathbb{R}^{m}\to\mathbb{R}$ are defined by

$\operatorname{\Pi}_{D}(v){}:={}\operatorname*{\arg\min}_{z\in
D}\|z-v\|{}\qquad\text{and}\qquad{}\operatorname{dist}_{D}(v){}\equiv{}\operatorname{dist}(v,D){}:={}\inf_{z\in
D}\|z-v\|.$

The former is a set-valued mapping whenever $D$ is nonconvex, whereas the
latter is always single-valued. The _indicator_ of $D$ is denoted
$\operatorname{\iota}_{D}\colon\mathbb{R}^{m}\to\overline{\mathbb{R}}$ and
defined by $\operatorname{\iota}_{D}(v)=0$ if $v\in D$, and
$\operatorname{\iota}_{D}(v)=\infty$ otherwise. If $D$ is nonempty and closed,
then $\operatorname{\iota}_{D}$ is proper and lsc. The proximal mapping of
$\operatorname{\iota}_{D}$ is the projection $\operatorname{\Pi}_{D}$; thus,
$\operatorname{\Pi}_{D}$ is locally bounded. Given $z\in D$, the _limiting
normal cone_ to $D$ at $z$ is the closed cone

$\mathcal{N}^{\text{lim}}_{D}(z):=\limsup_{v\to
z}\;\operatorname{cone}\left(v-\operatorname{\Pi}_{D}(v)\right).$

The limiting normal cone is robust in the sense that for all $z\in D$ it
satisfies $\mathcal{N}^{\text{lim}}_{D}(z)=\limsup_{v\to
z}\mathcal{N}^{\text{lim}}_{D}(v)$. For any proper and lsc function
$h\colon\mathbb{R}^{n}\to\overline{\mathbb{R}}$ and point
$\bar{x}\in\operatorname{dom}h$, we have

$\partial
h(\bar{x})=\left\\{v\in\mathbb{R}^{n}\,|\,(v,-1)\in\mathcal{N}^{\text{lim}}_{\operatorname{epi}h}(\bar{x},h(\bar{x}))\right\\},$
(2)

where
$\operatorname{epi}h:=\\{(x,\alpha)\in\mathbb{R}^{n}\times\mathbb{R}\,|\,h(x)\leq\alpha\\}$
denotes the _epigraph_ of $h$ [39, Theorem 8.9].

###  2.2 Stationarity Concepts

We now define some basic concepts and discuss optimality and stationarity
conditions for (P). As the cost function $\varphi$ is possibly extended-real-
valued, the concept of feasibility for a point refers to the domain of this
cost function. A point $x^{\ast}\in\mathbb{R}^{n}$ is called _feasible_ for
(P) if $x^{\ast}\in\operatorname{dom}\varphi$, that is, if
$c(x^{\ast})\in\operatorname{dom}g$. With (P), we seek points where the
minimum of $\varphi$ is attained, which are necessarily feasible. Notions of
optimality can have global or local scope. A feasible point
$x^{\ast}\in\mathbb{R}^{n}$ is said _locally optimal_ , or called a _local
minimizer_ , for (P) if there exists $\delta>0$ such that
$\varphi(x^{\ast})\leq\varphi(x)$ holds for all feasible
$x\in\operatorname{\mathbb{B}}_{\delta}(x^{\ast})$. Additionally, if the
inequality holds for all feasible $x\in\mathbb{R}^{n}$, then $x^{\ast}$ is
said _(globally) optimal_ for (P).

Focusing on the development of an “affordable” algorithm, our method is
designed to seek _candidate_ minimizers, guided by necessary optimality
conditions. Working without any convexity assumptions, a plausible
stationarity concept is that based on limiting subdifferentials; cf. [29,
Section 3] and [30, Theorem 5.48]. With respect to the minimization of a
proper lsc function $h\colon\mathbb{R}^{n}\to\overline{\mathbb{R}}$, we say
that $x^{\ast}\in\operatorname{dom}h$ is _M-stationary_ if $0\in\partial
h(x^{\ast})$, which constitutes a necessary condition for the optimality of
$x^{\ast}$ [39, Theorem 10.1]. When addressing (P), the composite structure
gives rise to first-order optimality conditions expressed in Lagrangian terms
[37, 38], as it intuitively follows from (P${}_{\text{S}}$). Although closely
related to M-stationarity, we prefer the nomenclature _KKT_ -stationarity to
emphasize the Lagrangian structure.

######  Definition 2.1 (KKT-stationarity).

Relative to (P), a point $x^{\ast}\in\mathbb{R}^{n}$ is called _KKT-
stationary_ if there exists a multiplier $y^{\ast}\in\mathbb{R}^{m}$ such that
$\displaystyle 0{}={}$ $\displaystyle\nabla
f(x^{\ast})+c^{\prime}(x^{\ast})^{\top}y^{\ast},$ (3a) $\displaystyle
y^{\ast}{}\in{}$ $\displaystyle\partial g(c(x^{\ast})).$ (3b)

Notice that (3b) implicitly requires the feasibility of $x^{\ast}$, namely
$c(x^{\ast})\in\operatorname{dom}g$, for otherwise the subdifferential
$\partial g(c(x^{\ast}))$ is empty. Moreover, note that Definition 2.1
coincides with the celebrated Karush-Kuhn-Tucker (KKT) conditions in nonlinear
programming if $g$ in (P) is the indicator of a nonempty closed convex set [4,
38, 15].

Although potentially ambiguous in the nonconvex setting, one could refer to
the Lagrange multiplier $y$ as to the _dual_ variable and to $(x,y)$ as to a
_primal-dual_ pair. In fact, a sort of local duality emerged while
investigating convergence of AL methods in nonconvex settings, as revealed by
Bertsekas [4] and recently by Rockafellar in a broader context [38, 37]—but
this seems premature for (P) under Assumption 1.

Given some tolerance $\varepsilon\geq 0$, an approximate M-stationarity
concept for the minimization of $h$ refers to $\operatorname{dist}(0,\partial
h(x^{\ast}))\leq\varepsilon$. We may also consider relaxed conditions for an
approximate counterpart of Definition 2.1.

######  Definition 2.2 (Approximate KKT-stationarity).

Let $\varepsilon\geq 0$ be given. Relative to (P), a point
$x^{\ast}\in\mathbb{R}^{n}$ is called _$\varepsilon$ -KKT-stationary_ if there
exist $z^{\ast}\in\operatorname{dom}g$ and a multiplier
$y^{\ast}\in\mathbb{R}^{m}$ such that $\displaystyle\|\nabla
f(x^{\ast})+c^{\prime}(x^{\ast})^{\top}y^{\ast}\|{}\leq{}$
$\displaystyle\varepsilon,$ (4a)
$\displaystyle\operatorname{dist}(y^{\ast},\partial g(z^{\ast})){}\leq{}$
$\displaystyle\varepsilon,$ (4b)
$\displaystyle\|c(x^{\ast})-z^{\ast}\|{}\leq{}$ $\displaystyle\varepsilon.$
(4c)

It should be noticed that $\varepsilon$-KKT-stationarity with $\varepsilon=0$
recovers the notion of KKT-stationarity, analogously to M-stationarity. This
fact follows from $\partial g(z)$ being closed at any point
$z\in\operatorname{dom}g$ [39, Theorem 8.6], which implies that, for any
$v\in\mathbb{R}^{m}$, $\operatorname{dist}(v,\partial g(z))=0$ is equivalent
to the inclusion $v\in\partial g(z)$.

The approximation allowed by Definition 2.2 entails not only stationarity (4a)
and complementarity (4b), but affects also feasibility: owing to (4c) it may
be that $c(x^{\ast})\notin\operatorname{dom}g$. Encompassing this inexactness
is particularly important in practice as it might be difficult to deal with
the composition $g\circ c$ as soon as $c$ is a nonlinear mapping. On the
contrary, it is dramatically easier to deal with (P) when $c$ is the identity
mapping. In this particular case, we apply a narrower definition compared to
Definition 2.2, as we deem feasibility of an approximate M-stationary point a
significant property. Therefore, relative to (P) with $c:=\operatorname{Id}$,
a point $x^{\ast}\in\mathbb{R}^{n}$ is called _$\varepsilon$ -M-stationary_ if
$\operatorname{dist}(-\nabla f(x^{\ast}),\partial
g(x^{\ast}))\leq\varepsilon$. Along with feasibility, namely
$x^{\ast}\in\operatorname{dom}g$, this demands a tighter approximation than
(4). As it suffices to invoke the proximal mapping of $g$ to reconcile the
“constraint” $x\in\operatorname{dom}g$, many feasible methods exist for such
two-term problems $f+g$; see e.g. [9, 32, 40, 17] and references therein.

We add to our theoretical toolbox also an asymptotic counterpart of KKT-
stationarity. The following closely patterns the _approximate_ KKT conditions
of [6, Definition 3.1], and so does the meaning of primal iterates $x^{k}$ and
Lagrange multipliers $y^{k}$.

######  Definition 2.3 (AKKT-stationarity).

Relative to (P), a feasible point $x^{\ast}\in\mathbb{R}^{n}$ is called
_asymptotically KKT-stationary_ if there exist sequences
$\\{x^{k}\\}_{k\in\mathbb{N}}\subseteq\mathbb{R}^{n}$ and
$\\{y^{k}\\}_{k\in\mathbb{N}},\\{z^{k}\\}_{k\in\mathbb{N}}\subseteq\mathbb{R}^{m}$
such that $z^{k}\overset{g}{\to}c(x^{\ast})$ and $\displaystyle\nabla
f(x^{k})+c^{\prime}(x^{k})^{\top}y^{k}{}\to{}$ $\displaystyle 0,$ (5a)
$\displaystyle y^{k}{}\in{}$ $\displaystyle\partial g(z^{k}),$ (5b)
$\displaystyle c(x)^{k}-z^{k}{}\to{}$ $\displaystyle 0.$ (5c)

Notice that Definition 2.3 does not require the sequence
$\\{y^{k}\\}_{k\in\mathbb{N}}$ to converge. A measure of (primal)
infeasibility is suggested by (5c), where $z^{k}$ is introduced to account in
(5b) for the fact that the condition $c(x^{k})\in\operatorname{dom}g$ can be
violated along the iterates, though it should hold asymptotically. Similarly,
with (5a) one can monitor the inexactness in the stationarity condition, or
_dual_ infeasibility. Finally, the attentive convergence
$z^{k}\overset{g}{\to}c(x^{\ast})$ plays a role when taking the limit in (5),
since we aim to recover the limiting subdifferential for KKT-stationarity; cf.
Definition 2.1.

A local minimizer for (P) is KKT-stationary under validity of a suitable
qualification condition, which, by non-Lipschitzness of $g$, will depend on
function $g$ as well, see [22]. However, we can show that each local minimizer
of (P) is always AKKT-stationary. Related results can be found in [27, Theorem
6.2], [29, Section 5.1] and [15, Proposition 2.5].

######  Proposition 2.4.

Let $x^{\ast}\in\mathbb{R}^{n}$ be a local minimizer for (P). Then, $x^{\ast}$
is AKKT-stationary.

###### Proof.

By local optimality of $x^{\ast}$ for (P) there exists $\delta>0$ such that
$\varphi(x^{\ast})\leq\varphi(x)$ is valid for all feasible
$x\in\operatorname{\mathbb{B}}_{\delta}(x^{\ast})$. Consequently, $x^{\ast}$
is the unique global minimizer of the localized problem

$\operatorname*{minimize}_{x\in\operatorname{\mathbb{B}}_{\delta}(x^{\ast})}{}\quad{}\varphi(x)+\frac{1}{2}\|x-x^{\ast}\|^{2}.$
(6)

Let us now consider the penalized surrogate problem

$\operatorname*{minimize}_{\begin{subarray}{c}x\in\operatorname{\mathbb{B}}_{\delta}(x^{\ast})\\\
z\in\operatorname{\mathbb{B}}_{1}(c(x^{\ast}))\end{subarray}}{}\quad{}f(x)+g(z)+\frac{\varrho_{k}}{2}\|c(x)-z\|^{2}+\frac{1}{2}\|x-x^{\ast}\|^{2}$
(P$(k)$)

where $k\in\mathbb{N}$ is arbitrary, $\varrho_{k}>0$, and the sequence
$\\{\varrho_{k}\\}_{k\in\mathbb{N}}$ satisfies $\varrho_{k}\to\infty$ as
$k\to\infty$. Noting that the objective function of this optimization problem
is lsc while its feasible set is nonempty (by feasibility of $x^{\ast}$) and
compact, it possesses a global minimizer $(x^{k},z^{k})$ for each
$k\in\mathbb{N}$. Therefore, without loss of generality, up to extracting a
subsequence, we assume $x^{k}\to\widehat{x}$ and $z^{k}\to\widehat{z}$ for
some $\widehat{x}\in\operatorname{\mathbb{B}}_{\delta}(x^{\ast})$ and
$\widehat{z}\in\operatorname{dom}g\cap\operatorname{\mathbb{B}}_{1}(c(x^{\ast}))$.

We now argue that $\widehat{x}=x^{\ast}$ and
$\widehat{z}=c(\widehat{x})=c(x^{\ast})$. To this end, we note that
$(x^{\ast},c(x^{\ast}))$ is feasible to (P$(k)$), which yields the estimate

$f(x^{k})+g(z^{k})+\frac{\varrho_{k}}{2}\|c(x^{k})-z^{k}\|^{2}+\frac{1}{2}\|x^{k}-x^{\ast}\|^{2}\leq\varphi(x^{\ast})$
(7)

for each $k\in\mathbb{N}$. Using lower semicontinuity of $g$, finiteness of
$\varphi(x^{\ast})$ as well as the convergences $c(x^{k})\to c(\widehat{x})$
and $z^{k}\to\widehat{z}$, taking the limit for $k\to\infty$ in (7) gives
$c(\widehat{x})=\widehat{z}\in\operatorname{dom}g$. Therefore, $\widehat{x}$
is feasible for (6) and local optimality of $x^{\ast}$ for (P) implies
$\varphi(x^{\ast})\leq\varphi(\widehat{x})$. Furthermore, exploiting (7) and
the optimality of each
$(x^{k},z^{k})\in\mathbb{R}^{n}\times\operatorname{dom}g$, we find that
$\widehat{x}=x^{\ast}$, since

$\varphi(\widehat{x})+\frac{1}{2}\|\widehat{x}-x^{\ast}\|^{2}\leq\liminf\limits_{k\to\infty}\left(f(x^{k})+g(z^{k})+\frac{\varrho_{k}}{2}\|c(x^{k})-z^{k}\|^{2}+\frac{1}{2}\|x^{k}-x^{\ast}\|^{2}\right)\\\
\leq\varphi(x^{\ast})\leq\varphi(\widehat{x}).$

Moreover, noting that (7) also gives $f(x^{k})+g(z^{k})\leq\varphi(x^{\ast})$
for each $k\in\mathbb{N}$, by

$\displaystyle\varphi(x^{\ast})\leq\liminf\limits_{k\to\infty}f(x^{k})+g(z^{k})\leq\limsup\limits_{k\to\infty}f(x^{k})+g(z^{k})\leq\varphi(x^{\ast})$

we have the attentive convergence of $\\{z^{k}\\}_{k\in\mathbb{N}}$ to
$\widehat{z}=c(x^{\ast})$, that is, $z^{k}\overset{g}{\to}c(x^{\ast})$.

Now, due to $x^{k}\to x^{\ast}$ and $z^{k}\to c(x^{\ast})$, we may assume
without loss of generality that $\\{x^{k}\\}_{k\in\mathbb{N}}$ and
$\\{z^{k}\\}_{k\in\mathbb{N}}$ are taken from the interior of
$\operatorname{\mathbb{B}}_{\delta}(x^{\ast})$ and
$\operatorname{\mathbb{B}}_{1}(c(x^{\ast}))$, respectively, as this is
eventually the case. Thus, for each $k\in\mathbb{N}$, $(x^{k},z^{k})$ is an
unconstrained (that is, with the constraints being inactive) minimizer of

$(x,z)\mapsto
f(x)+g(z)+\frac{\varrho_{k}}{2}\|c(x)-z\|^{2}+\frac{1}{2}\|x-x^{\ast}\|^{2},$

whose relevant necessary conditions for optimality read

$\displaystyle 0{}={}$ $\displaystyle\nabla
f(x^{k})+\varrho_{k}c^{\prime}(x^{k})^{\top}[c(x^{k})-z^{k}]+x^{k}-x^{\ast},$
$\displaystyle 0{}\in{}$ $\displaystyle\partial
g(z^{k})-\varrho_{k}[c(x^{k})-z^{k}].$

Thus, setting $y^{k}:=\varrho_{k}[c(x^{k})-z^{k}]$ for each $k\in\mathbb{N}$,
the conditions in (5) are a consequence of $x^{k}\to x^{\ast}$ and $z^{k}\to
c(x^{\ast})$. Overall, this shows that $x^{\ast}$ is AKKT-stationary for (P).
∎

In order to guarantee that local minimizers for (P) are not only AKKT- but
also KKT-stationary, the presence of a qualification condition is necessary.
Analogously to [15, Definition 2.6], the following regularity notion
generalizes the (comparatively weak) constraint qualification from [29,
Section 3.2] to the non-Lipschitzian setting of (P) and is closely related to
the _uniform qualification condition_ introduced in [27, Definition 6.8].

######  Definition 2.5 (AKKT-regularity).

Let $\bar{x}\in\mathbb{R}^{n}$ be a feasible point for (P), namely
$c(\bar{x})\in\operatorname{dom}g$, and define the set-valued mapping
$\mathcal{M}\colon\mathbb{R}^{n}\times\mathbb{R}^{m}\rightrightarrows\mathbb{R}^{n}$
by $\mathcal{M}(x,z):=c^{\prime}(x)^{\top}\partial g(z)$. Then, $\bar{x}$ is
called _asymptotically KKT-regular_ for (P) if
$\limsup_{\begin{subarray}{c}x\to\bar{x}\\\
z\overset{g}{\to}c(\bar{x})\end{subarray}}\mathcal{M}(x,z)\subset\mathcal{M}(\bar{x},c(\bar{x})).$

Sufficient conditions for the validity of the more general qualification
condition from Definition 2.5 can be captured following the arguments of [27].
A direct consequence of Proposition 2.4 is the next result, analogous in
spirit to [27, Proposition 6.9] and [15, Corollary 2.7].

######  Corollary 2.6.

Let $x^{\ast}\in\mathbb{R}^{n}$ be an AKK-regular AKKT-stationary point for
(P). Then, $x^{\ast}$ is a KKT-stationary point for (P). In particular, each
AKKT-regular local minimizer for (P) is KKT-stationary.

##  3 Augmented Lagrangian Method

Algorithms exhibiting a nested structure naturally arise in the AL framework,
which builds upon a sequence of subproblems. We develop our implicit approach
and derive the generalized AL function in Section 3.1. Then, in Section 3.2 we
discuss stationarity notions relevant for the AL subproblem, in order to
precisely qualify a subproblem “solution”, and thus the output expected from a
subsolver. The following Section 3.3 contains a detailed statement of our
algorithmic framework, whose convergence analysis is presented in Section 3.4.
Suitable termination criteria are discussed in Section 3.5. The subsequent
Section 4 is dedicated to the numerical solution of the AL subproblems.

###  3.1 Implicit Augmented Lagrangian

Introducing an auxiliary variable $z\in\mathbb{R}^{m}$, (P) can be
equivalently rewritten in the form (P${}_{\text{S}}$), which has a structured
objective function and some (explicit) equality constraints. For some penalty
parameter $\mu>0$, let us define the AL function
$\mathcal{L}^{\text{S}}_{\mu}\colon\mathbb{R}^{n}\times\mathbb{R}^{m}\times\mathbb{R}^{m}\to\overline{\mathbb{R}}$
associated to (P${}_{\text{S}}$) by

$\mathcal{L}^{\text{S}}_{\mu}(x,z,y){}:={}f(x)+g(z)+\langle
y,c(x)-z\rangle+\frac{1}{2\mu}\|c(x)-z\|^{2}\\\
=f(x)+g(z)+\frac{1}{2\mu}\|c(x)+\mu y-z\|^{2}-\frac{\mu}{2}\|y\|^{2}.$ (8)

Within AL schemes [4, 10, 6], it is required to minimize (at least, up to
approximate stationarity) the AL function. Entailing the sum of a smooth and a
prox-friendly function, the minimization of $\mathcal{L}^{\text{S}}_{\mu}$ can
be tackled by exploiting its structure and the computational oracles via,
e.g., proximal-gradient methods [9, 32, 40]. Owing to the convenient structure
of $\mathcal{L}^{\text{S}}_{\mu}$, the method developed in [15] relies on this
“explicit” approach. Here, instead of performing a joint minimization of
$\mathcal{L}^{\text{S}}_{\mu}(\cdot,\cdot,y)$ over both $(x,z)$, we exploit
the structure arising from the original problem (P), formally eliminating the
slack variable $z$, on the vein of the proximal AL approach [18, 13, 23, 14].
Given some $\mu>0$ and $y\in\mathbb{R}^{m}$, we define the oracle
$\operatorname{\mathbf{O}}_{\mu}(\cdot,y)$ as the mapping associated to the
pointwise minimization of $\mathcal{L}^{\text{S}}_{\mu}(\cdot,\cdot,y)$ with
respect to $z$. Owing to the definitions in (1) and (8), the set-valued
mapping
$\operatorname{\mathbf{O}}_{\mu}(\cdot,y)\colon\mathbb{R}^{n}\rightrightarrows\mathbb{R}^{m}$
is given by

$\operatorname{\mathbf{O}}_{\mu}(x,y){}:={}\operatorname*{\arg\min}_{z\in\mathbb{R}^{m}}\mathcal{L}^{\text{S}}_{\mu}(x,z,y){}={}\operatorname*{\arg\min}_{z\in\mathbb{R}^{m}}\left\\{g(z)+\frac{1}{2\mu}\|c(x)+\mu
y-z\|^{2}\right\\}\\\ {}={}\operatorname{prox}_{\mu g}\left(c(x)+\mu
y\right).$ (9)

Injecting back into $\mathcal{L}^{\text{S}}_{\mu}(x,\cdot,y)$ any arbitrary
element of $\operatorname{\mathbf{O}}_{\mu}(x,y)$, namely evaluating the AL on
the set corresponding to the minimization over $z$, we obtain the (single-
valued) AL function
$\mathcal{L}_{\mu}\colon\mathbb{R}^{n}\times\mathbb{R}^{m}\to\mathbb{R}$
associated to (P):

$\mathcal{L}_{\mu}(x,y){}:={}\inf_{z\in\mathbb{R}^{m}}\mathcal{L}^{\text{S}}_{\mu}(x,z,y){}={}\mathcal{L}^{\text{S}}_{\mu}(x,z,y)\Big{|}_{z\in\operatorname{\mathbf{O}}_{\mu}(x,y)}\\\
{}={}f(x)+g^{\mu}(c(x)+\mu y)-\frac{\mu}{2}\|y\|^{2}.$ (10)

We highlight that the Moreau envelope
$g^{\mu}\colon\mathbb{R}^{m}\to\mathbb{R}$ of $g$ is real-valued and strictly
continuous [31], [39, Example 10.32], but not continuously differentiable in
general, as the proximal mapping of $g$ is possibly set-valued. Regarded as
the AL for (P), $\mathcal{L}_{\mu}$ was called _generalized_ AL function in
[38] (for $g$ convex). Compared to $\mathcal{L}^{\text{S}}_{\mu}$ with the
_explicit_ slack variable $z$ for the lifted “constrained” formulation
(P${}_{\text{S}}$), we may also refer to $\mathcal{L}_{\mu}$ as to the
_implicit_ AL for (P).

The algorithm we are going to present in Section 3.3 fits in the AL framework
[4, 6, 38] in that it generates sequences of iterates $x^{k}$ and Lagrange
multipliers $y^{k}$ obeying the recursion

$x^{k+1}\mathrel{\mathchoice{\reflectbox{$\displaystyle\rightsquigarrow$}}{\reflectbox{$\textstyle\rightsquigarrow$}}{\reflectbox{$\scriptstyle\rightsquigarrow$}}{\reflectbox{$\scriptscriptstyle\rightsquigarrow$}}}\operatorname*{\arg\min}_{x\in\mathbb{R}^{n}}\mathcal{L}_{\mu_{k}}(x,\hat{y}^{k}),\qquad
y^{k+1}\mathrel{\mathchoice{\reflectbox{$\displaystyle\rightsquigarrow$}}{\reflectbox{$\textstyle\rightsquigarrow$}}{\reflectbox{$\scriptstyle\rightsquigarrow$}}{\reflectbox{$\scriptscriptstyle\rightsquigarrow$}}}\hat{y}^{k}+\frac{1}{\mu_{k}}\partial_{y}\mathcal{L}_{\mu_{k}}(x^{k+1},\hat{y}^{k}),$
(11)

with respect to suitable sequences of parameter values $\mu_{k}>0$ and
multiplier estimates $\hat{y}^{k}\in\mathbb{R}^{m}$. Here we intend to give
the reader only a prospect about the algorithm before looking into the
details; by using the symbol
“$\mathrel{\mathchoice{\reflectbox{$\displaystyle\rightsquigarrow$}}{\reflectbox{$\textstyle\rightsquigarrow$}}{\reflectbox{$\scriptstyle\rightsquigarrow$}}{\reflectbox{$\scriptscriptstyle\rightsquigarrow$}}}$”
in (11) we emphasize that the notation should not be taken too seriously. Most
especially, for the primal update in (11) we do not require a _global_ , nor
an _exact_ local minimizer of $\mathcal{L}_{\mu_{k}}(\cdot,\hat{y}^{k})$, but
only an approximate stationary point thereof (in the sense described in
Section 3.2). Regarding the multiplier, the update rule in (11) reduces to the
traditional one when (P) specializes to nonlinear programming [4, 6], but if
$g$ is nonconvex it may be necessary to select specific elements from the
relevant set-valued subdifferentials. This potential choice is also very much
related to the stationarity concept introduced in Section 3.2.

At the $k$-th iteration of the AL scheme, a subproblem in the form

$\operatorname*{minimize}_{x\in\mathbb{R}^{n}}\
\mathcal{L}_{\mu_{k}}(x,\hat{y}^{k}),$ (12)

with some given $\mu_{k}>0$ and $\hat{y}^{k}\in\mathbb{R}^{m}$, is to be
approximately solved. Being $g$ potentially nonconvex, the AL may lack
regularity. The following Section 3.2 introduces a stationary concept _weaker_
than M-stationarity for (12) but sharp enough to guarantee that the iterates
generated by the AL scheme detailed in Section 3.3 exhibit opportune
asymptotic behavior and properties, as demonstrated later in Section 3.4.

###  3.2 Stationarity for Parametric Optimization

In view of the marginalization step to derive $\mathcal{L}_{\mu}$ in (10), we
now study the AL subproblem (12) with the lenses of parametric optimization
[39, Sections 1.F, 10.C]. Let us consider an objective
$p\colon\mathbb{R}^{n}\to\mathbb{R}$ and an oracle
$\mathbf{O}\colon\mathbb{R}^{n}\rightrightarrows\mathbb{R}^{m}$ given by

$p(x):=\inf_{z\in\mathbb{R}^{m}}P(x,z)\qquad\text{and}\qquad\mathbf{O}(x):=\operatorname*{\arg\min}_{z\in\mathbb{R}^{m}}P(x,z)$
(13)

for a proper, lsc function
$P\colon\mathbb{R}^{n}\times\mathbb{R}^{m}\to\overline{\mathbb{R}}$. We are
interested in minimizing $p$, but first in devising suitable stationarity
concepts to qualify a “solution”. This setting is apparently mirroring
subproblem (12), with the explicit and implicit AL functions
$\mathcal{L}^{\text{S}}$ and $\mathcal{L}$ playing the role of $P$ and $p$,
respectively. Moreover, it highlights our intention to keep the solution
process implicit, effectively hiding the auxiliary variable $z$, by utilizing
the oracle $\mathbf{O}$.

Approximate M-stationarity of a point $\bar{x}\in\operatorname{dom}p$ for some
given $\varepsilon\geq 0$ amounts to requiring that
$\operatorname{dist}(0,\partial p(\bar{x}))\leq\varepsilon$. However, because
of the parametric nature of $p$, the subdifferential $\partial p$ is not a
simple object in general. Let us consider some results on subdifferentiation
in parametric optimization. Recalling from Section 2.1 that the notion of
uniform level-boundedness corresponds to a parametric extension of level-
boundedness [39, Definition 1.16], we suppose that $P(x,z)$ is level-bounded
in $z$ locally uniformly in $x$. Then, from [39, Theorem 10.13] we have for
every $\bar{x}\in\operatorname{dom}p$ the inclusion

$\partial
p(\bar{x}){}\subseteq{}\Upsilon(\bar{x}),\qquad\text{where}\qquad\Upsilon(\bar{x}){}:={}\bigcup_{\bar{z}\in\mathbf{O}(\bar{x})}\\{v\,|\,(v,0)\in\partial
P(\bar{x},\bar{z})\\},$ (14)

which implies the inequality

$\operatorname{dist}(0,\partial
p(\bar{x}))\geq\operatorname{dist}(0,\Upsilon(\bar{x})).$ (15)

This offers a lower bound on the M-stationarity measure, and it does not
provide in general an upper bound suitable for a termination condition.
Regardless, in the parametric optimization setting of (13), we may resort to
the following stationarity concept, here denominated _$\Upsilon$
-stationarity_.

######  Definition 3.1 ($\Upsilon$-stationarity).

Let $\varepsilon\geq 0$ be fixed and consider a proper, lsc function
$P\colon\mathbb{R}^{n}\times\mathbb{R}^{m}\to\overline{\mathbb{R}}$ with
$P(x,z)$ level-bounded in $z$ locally uniformly in $x$. Define
$p\colon\mathbb{R}^{n}\to\overline{\mathbb{R}}$ and
$\Upsilon\colon\mathbb{R}^{n}\rightrightarrows\mathbb{R}^{n}$ as in (13) and
(14), respectively. Then, relatively to $p$, a point
$x^{\ast}\in\operatorname{dom}p$ is called _$\varepsilon$
-$\Upsilon$-stationary_ if
$\operatorname{dist}(0,\Upsilon(x^{\ast}))\leq\varepsilon$. In the case
$\varepsilon=0$, such point $x^{\ast}$ is said _$\Upsilon$ -stationary_.

Notice that the inclusion $x^{\ast}\in\operatorname{dom}p$ is implicitly
required to have $\Upsilon(x^{\ast})$ nonempty. Furthermore, in the exact
case, namely $\varepsilon=0$, the requirement for $\Upsilon$-stationarity
reduces to $0\in\Upsilon(x^{\ast})$, by closedness of $\Upsilon$ [39, Theorem
10.13].

From the definition of $\Upsilon$ in (14), it appears that an
$\varepsilon$-$\Upsilon$-stationary point $x^{\ast}\in\mathbb{R}^{n}$ can be
associated to a (possibly non-unique) _certificate_
$z^{\ast}\in\mathbf{O}(x^{\ast})$ and a _pair_ $(x^{\ast},z^{\ast})$ that
satisfy

$\operatorname{dist}(0,\Upsilon(x^{\ast})){}={}\min_{v\in\Upsilon(x^{\ast})}\|v\|{}={}\min_{\begin{subarray}{c}v\in\mathbb{R}^{n}\\\
z\in\mathbf{O}(x^{\ast})\end{subarray}}\left\\{\|v\|\,|\,(v,0)\in\partial
P(x^{\ast},z)\right\\}\\\
{}\leq{}\min_{v\in\mathbb{R}^{n}}\left\\{\|v\|\,|\,(v,0)\in\partial
P(x^{\ast},z^{\ast})\right\\}\leq\varepsilon.$ (16)

If any such pair $(x^{\ast},z^{\ast})$ is found, the upper bound in (16)
certificates the $\varepsilon$-$\Upsilon$-stationarity of $x^{\ast}$.

The following result establishes that (approximate) M-stationarity provides a
stronger qualification than (approximate) $\Upsilon$-stationarity, meaning
that the former is a sufficient condition, not necessary in general, for the
latter. The sought implication readily follows from the upper bound in (15) on
the $\Upsilon$-stationarity measure.

######  Proposition 3.2.

Let $\varepsilon\geq 0$ be fixed. Then, relative to the minimization of $p$,
any $\varepsilon$-M-stationary point is also
$\varepsilon$-$\Upsilon$-stationary.

The converse implication does not hold in general, but it does under some
(possibly local) convexity assumptions; see, e.g., [39, Theorem 2.26]. Similar
considerations are valid with regard to the following result, which shows that
(approximate) $\Upsilon$-stationarity relative to $p$ is a stronger
qualification than (approximate) M-stationarity relative to $P$.

######  Proposition 3.3.

Let $\varepsilon\geq 0$ be fixed and $(x^{\ast},z^{\ast})$ be any
$\varepsilon$-$\Upsilon$-stationary pair associated to $p$. Then, relative to
the minimization of $P$, $(x^{\ast},z^{\ast})$ is $\varepsilon$-M-stationary.

###### Proof.

By $\varepsilon$-$\Upsilon$-stationarity relative to $p$, the pair
$(x^{\ast},z^{\ast})$ necessarily satisfies $x^{\ast}\in\operatorname{dom}p$
and $z^{\ast}\in\mathbf{O}(x^{\ast})$. Then, exploiting the upper bound in
(16), we can argue that

$\operatorname{dist}(0,\partial
P(x^{\ast},z^{\ast})){}={}\min_{\begin{subarray}{c}v\in\mathbb{R}^{n}\\\
u\in\mathbb{R}^{m}\end{subarray}}\left\\{\|(v,u)\|\,|\,(v,u)\in\partial
P(x^{\ast},z^{\ast})\right\\}\\\
{}\leq{}\min_{v\in\mathbb{R}^{n}}\left\\{\|v\|\,|\,(v,0)\in\partial
P(x^{\ast},z^{\ast})\right\\}{}\leq{}\varepsilon,$

where the first inequality is due to the restriction $u:=0$. Showing that
$(x^{\ast},z^{\ast})$ satisfies $\operatorname{dist}(0,\partial
P(x^{\ast},z^{\ast}))\leq\varepsilon$, the assertion is proved. ∎

Notice, however, that in general the statement of Proposition 3.3 is not valid
for all $z^{\ast}\in\mathbf{O}(x^{\ast})$, even if $x^{\ast}$ is
$\varepsilon$-$\Upsilon$-stationary for $p$. In fact, some pairs
$(x^{\ast},z^{\ast})$ may fail to satisfy (16), and those are not necessarily
$\varepsilon$-M-stationary for $P$.

At this point it appears clear that, in relation to (8) and (10), the oracle
$\mathbf{O}$ refers to $\operatorname{\mathbf{O}}_{\mu}(\cdot,\hat{y})$, for
some given $\mu>0$ and $\hat{y}\in\mathbb{R}^{m}$. Moreover, we shall focus
our attention on the special case of $P$ having the structure

$\displaystyle P(x,z){}:={}$ $\displaystyle F(x,z)+g(z),$ (17a) $\displaystyle
F(x,z){}:={}$ $\displaystyle f(x)+\frac{1}{2\mu}\|c(x)+\mu\hat{y}-z\|^{2},$
(17b)

where $F\colon\mathbb{R}^{n}\times\mathbb{R}^{m}\to\mathbb{R}$ is continuously
differentiable.

###  3.3 Algorithm

This section details an AL method for the solution of generalized programs of
the form (P). In this work we focus on a safeguarded AL scheme inspired by [6,
Algorithm 4.1] and investigate its convergence properties. Compared to the
classical AL or multiplier penalty approach for the solution of nonlinear
programs [10, 4], this variant uses a safeguarded update rule for the Lagrange
multipliers and has stronger global convergence properties [26]. Although we
restrict our analysis to this specific algorithm, analogous results can be
obtained for others with minor changes, e.g., [19, Algorithm 1].

The overall method is stated in Algorithm 3.1. Although a starting point is
not explicitly required, the subproblems at Algorithm 3.1 should be solved
warm-starting from the current primal estimate $x^{k-1}$, thus exploiting
initial guesses. The AL subproblem appears at Algorithm 3.1, involving the
implicit AL function with current penalty parameter $\mu_{k}$ and estimate
$\hat{y}^{k}$. Crucially, Algorithm 3.1 requires to compute an approximate
$\Upsilon$-stationary _pair_ $(x^{k},z^{k})$, not merely a point $x^{k}$, for
$\mathcal{L}_{\mu_{k}}(\cdot,\hat{y}^{k})$. As it may be employed to verify
approximate stationarity, it seems likely that a certificate $z^{k}$ for
$x^{k}$ can be obtained with no additional effort — this is the case for the
subsolver developed in Section 4.

Data: $\epsilon>0$

Result: $\epsilon$-KKT-stationary point $x^{\ast}$

1Select $\mu_{0}>0$, $\theta,\kappa\in(0,1)$, $Y\subset\mathbb{R}^{m}$
nonempty bounded

2 for  _$k=0,1,2\ldots$_ do

3 Select $\hat{y}^{k}\in Y$ and $\varepsilon_{k}\geq 0$

4 Compute an $\varepsilon_{k}$-$\Upsilon$-stationary pair $(x^{k},z^{k})$ for
$\mathcal{L}_{\mu_{k}}(\cdot,\hat{y}^{k})$

5 Set $y^{k}\leftarrow\hat{y}^{k}+\mu_{k}^{-1}[c(x^{k})-z^{k}]$ and
$V^{k}\leftarrow\|c(x^{k})-z^{k}\|$

6 if  _$\max\\{\varepsilon_{k},V^{k}\\}\leq\epsilon$_ then return
$x^{\ast}\leftarrow x^{k}$

7 if  _$k=0$ or $V^{k}\leq\theta V^{k-1}$_ then $\mu_{k+1}\leftarrow\mu_{k}$,
else $\mu_{k+1}\leftarrow\kappa\mu_{k}$

8

Algorithm 3.1 Safeguarded implicit AL method for (P)

The safeguarded Lagrange multiplier estimate $\hat{y}^{k}$ is drawn from a
bounded set $Y\subset\mathbb{R}^{m}$ at Algorithm 3.1. Although not necessary,
the choice of $\hat{y}^{k}$ should also depend on the current estimate
$y^{k-1}$. Moreover, the choice of $Y$ can take advantage of _a priori_
knowledge of $\partial g$ and its structure, in order to generate better
Lagrange multiplier estimates. In practice, it is advisable to choose the
safeguarded estimate $\hat{y}^{k}$ as the projection of the Lagrange
multiplier $y^{k-1}$ onto $Y$, thus effectively adopting the classical
approach as long as $y^{k-1}$ remains within $Y$; we refer to [6, Section 4.1]
and [26, Section 3] for a more articulated discussion.

Inspired by the classical first-order Lagrange multiplier estimate, the update
rule at Algorithm 3.1 is designed around (5a) and the identity

$\nabla_{x}\mathcal{L}^{\text{S}}_{\mu}(x,z,\hat{y})=\nabla
f(x)+\frac{1}{\mu}c^{\prime}(x)^{\top}[c(x)+\mu\hat{y}-z].$

Given any $(x^{k},z^{k})$ obtained at Algorithm 3.1, the specific choice of
$y^{k}$ leads by construction to

$\|\nabla f(x^{k})+c^{\prime}(x^{k})^{\top}y^{k}\|\leq\varepsilon_{k},$

effectively connecting inner and outer loops, and providing an upper bound on
the dual infeasibility at $(x^{k},y^{k})$. Algorithm 3.1 encapsulates suitable
termination conditions for returning an approximate KKT-stationary point;
these are discussed in Section 3.5. The monotonicity test at Algorithm 3.1 is
adopted to monitor primal infeasibility along the iterates and update the
penalty parameter accordingly. Aimed at driving $V^{k}:=\|c(x^{k})-z^{k}\|$ to
zero, the penalty parameter $\mu_{k}$ is reduced at Algorithm 3.1 in case of
insufficient decrease.

###  3.4 Convergence Analysis

Throughout our convergence analysis, we assume that Algorithm 3.1 is well-
defined, thus requiring that each subproblem at Algorithm 3.1 admits an
approximate stationary point. Moreover, the following statements assume the
existence of some accumulation point $x^{\ast}$ for a sequence
$\\{x^{k}\\}_{k\in\mathbb{N}}$ generated by Algorithm 3.1. In general,
coercivity or (level) boundedness arguments should be adopted to verify these
preconditions.

Due to their practical importance, we focus on _local_ , or affordable,
solvers, which return stationary points as candidate local minimizers, for the
subproblems at Algorithm 3.1. Instead, we do not investigate the case where
the subproblems are solved to _global_ optimality, whose analysis would
pattern the classical results in [6, Chapter 5] and [26].

As it is often the case for penalty-type methods in the nonconvex setting,
Algorithm 3.1 may generate accumulation points that are infeasible for (P).
Patterning standard arguments, the following result gives conditions that
guarantee feasibility of limit points; cf. [5, Problem 4.12], [15, Proposition
3.2].

######  Proposition 3.4.

Let Assumption 1 hold and consider a sequence $\\{x^{k}\\}_{k\in\mathbb{N}}$
of iterates generated by Algorithm 3.1 with $\epsilon=0$. Then, each
accumulation point $x^{\ast}$ of $\\{x^{k}\\}_{k\in\mathbb{N}}$ is feasible
for (P) if (at least) one of the following conditions holds: (i)
$x^{k}\overset{\varphi}{\to}_{K}x^{\ast}$ for some subsequence
$\\{x^{k}\\}_{k\in K}$; (ii) there exists some $B\in\mathbb{R}$ such that
$\mathcal{L}_{\mu_{k}}(x^{k},\hat{y}^{k})\leq B$ for all $k\in\mathbb{N}$;
(iii) $\\{\mu_{k}\\}_{k\in\mathbb{N}}$ is bounded away from zero and
$\operatorname{dom}g$ is closed; (iv) $\\{\mu_{k}\\}_{k\in\mathbb{N}}$ is
bounded away from zero and
$\mathcal{L}_{\mu_{k}}(x^{k},\hat{y}^{k})\leq\mathcal{L}_{\mu_{k}}(x^{k-1},\hat{y}^{k})$
for all $k\in\mathbb{N}$.

###### Proof.

Let $x^{\ast}\in\mathbb{R}^{n}$ be an arbitrary accumulation point of
$\\{x^{k}\\}$ and $\\{x^{k}\\}_{k\in K}$ a subsequence such that
$x^{k}\to_{K}x^{\ast}$. We shall argue that $x^{\ast}$ is feasible for (P),
that is, $c(x^{\ast})\in\operatorname{dom}g$. Before considering each
assertion, we point out that, if $\\{\mu_{k}\\}$ is bounded away from zero,
the conditions at Algorithm 3.1 of Algorithm 3.1 imply that
$V^{k}:=\|c(x^{k})-z^{k}\|\to 0$ for $k\to\infty$. Thus, since
$x^{k}\to_{K}x^{\ast}$, by continuity arguments we have also
$z^{k}\to_{K}c(x^{\ast})$. Finally, let us conveniently denote by
$\Delta_{Y}\in\mathbb{R}$ the finite size of the bounded set
$Y\subset\mathbb{R}^{m}$ initialized in Algorithm 3.1, such that
$\|\hat{y}\|\leq\Delta_{Y}$ for all $\hat{y}\in Y$.

  1. (i)

The $\varphi$-attentive convergence of $\\{x^{k}\\}_{k\in K}$ to $x^{\ast}$
requires that $x^{\ast}\in\operatorname{dom}\varphi$. Hence, by
$\varphi:=g\circ c$ and continuity of $c$, it must be
$c(x^{\ast})\in\operatorname{dom}g$.

  2. (ii)

In view of the preliminary observation regarding $\\{\mu_{k}\\}$ bounded away
from zero, it remains to consider the case $\mu_{k}\to 0$ for showing that
$z^{k}\to_{K}c(x^{\ast})$. Using (10) and
$z^{k}\in\operatorname{\mathbf{O}}_{\mu_{k}}(x^{k},\hat{y}^{k})$, the upper
bound yields

$B\geq\mathcal{L}_{\mu_{k}}(x^{k},\hat{y}^{k})=f(x^{k})+g(z^{k})+\frac{1}{2\mu_{k}}\|c(x^{k})+\mu_{k}\hat{y}^{k}-z^{k}\|^{2}-\frac{\mu_{k}}{2}\|\hat{y}^{k}\|^{2}$

for all $k\in\mathbb{N}$. Rearranging terms and considering $\mu_{k}\to 0$, we
get

$\|c(x^{k})+\mu_{k}\hat{y}^{k}-z^{k}\|^{2}\leq
2\mu_{k}\left[B+\frac{\mu_{k}}{2}\|\hat{y}^{k}\|^{2}-f(x^{k})-g(z^{k})\right]\to_{K}0,$

where the limit follows from boundedness of
$\\{\hat{y}^{k}\\}_{k\in\mathbb{N}}$ and $\\{z^{k}\\}_{k\in K}$, the latter
remaining bounded by $x^{k}\to_{K}x^{\ast}$ and local boundedness of the
proximal mapping. Using once again that $\mu_{k}\hat{y}^{k}\to 0$, the limit
implies $c(x^{k})-z^{k}\to_{K}0$. Thus, since $x^{k}\to_{K}x^{\ast}$, by
continuity arguments it is $z^{k}\to_{K}c(x^{\ast})$. Now, owing to
$x^{k}\to_{K}x^{\ast}$, boundedness of $\\{\hat{y}^{k}\\}_{k\in\mathbb{N}}$,
and $\mu_{k}\leq\mu_{0}$, from the upper bound condition we also obtain that
$\\{g(z^{k})\\}_{k\in K}$ remains bounded above:

$\limsup_{k\to_{K}\infty}g(z^{k})\leq
B+\limsup_{k\to_{K}\infty}\frac{\mu_{k}}{2}\|\hat{y}^{k}\|^{2}-f(x^{k})\leq
B+\frac{\mu_{0}}{2}\Delta_{Y}^{2}-f(x^{\ast})<\infty.$

Combined with $z^{k}\to_{K}c(x^{\ast})$ and lower semicontinuity of $g$, this
implies that $c(x^{\ast})\in\operatorname{dom}g$.

  3. (iii)

As $\\{\mu_{k}\\}_{k\in\mathbb{N}}$ remains bounded away from zero, we have
that $z^{k}\to_{K}c(x^{\ast})$. Then, since
$z^{k}\in\operatorname{\mathbf{O}}_{\mu_{k}}(x^{k},\hat{y}^{k})\subseteq\operatorname{dom}g$
for all $k\in\mathbb{N}$ by Algorithm 3.1, it follows from
$\operatorname{dom}g$ being closed that $c(x^{\ast})\in\operatorname{dom}g$.

  4. (iv)

As $\\{\mu_{k}\\}_{k\in\mathbb{N}}$ remains bounded away from zero, we have
that $z^{k}\to_{K}c(x^{\ast})$ and, owing to the conditions at Algorithm 3.1,
$\\{\mu_{k}\\}_{k\in\mathbb{N}}$ is asymptotically constant, namely for some
finite index $k_{\ast}\in\mathbb{N}$ it is $\mu_{k}=\mu_{\ast}>0$ for all
$k\geq k_{\ast}$. In turn, the primal infeasibility measure $V^{k}$ decreases
linearly, as it must be $V^{k}\leq\theta V^{k-1}$ for all $k\geq k_{\ast}$.
Now, rearranging (10) and exploiting the explicit minimization over $z$, we
obtain the upper bound

$\displaystyle\mathcal{L}_{\mu_{\ast}}(x^{k-1},\hat{y}^{k}){}={}$
$\displaystyle
f(x^{k-1})+\min_{z}\left\\{g(z)+\langle\hat{y}^{k},c(x^{k-1})-z\rangle+\frac{1}{2\mu_{\ast}}\|c(x^{k-1})-z\|^{2}\right\\}$
$\displaystyle{}\leq{}$ $\displaystyle
f(x^{k-1})+g(z^{k-1})+\langle\hat{y}^{k},c(x^{k-1})-z^{k-1}\rangle+\frac{1}{2\mu_{\ast}}\|c(x^{k-1})-z^{k-1}\|^{2}$
$\displaystyle{}={}$
$\displaystyle\mathcal{L}_{\mu_{\ast}}(x^{k-1},\hat{y}^{k-1})+\langle\hat{y}^{k}-\hat{y}^{k-1},c(x^{k-1})-z^{k-1}\rangle$

for all $k>k_{\ast}$. Invoking the Cauchy-Schwarz inequality, the descent
condition gives

$\displaystyle\mathcal{L}_{\mu_{\ast}}(x^{k},\hat{y}^{k}){}\leq{}\mathcal{L}_{\mu_{\ast}}(x^{k-1},\hat{y}^{k}){}\leq{}$
$\displaystyle\mathcal{L}_{\mu_{\ast}}(x^{k-1},\hat{y}^{k-1})+2\Delta_{Y}V^{k-1}$

for all $k>k_{\ast}$, where we used that $\|\hat{y}^{k}-\hat{y}^{k-1}\|\leq
2\Delta_{Y}$. Then, telescoping this argument, the linear decrease of $V^{k}$
yields

$\mathcal{L}_{\mu_{\ast}}(x^{k+\ell},\hat{y}^{k+\ell})-\mathcal{L}_{\mu_{\ast}}(x^{k-1},\hat{y}^{k-1}){}\leq{}2\Delta_{Y}\sum_{m=0}^{\ell}V^{k+m-1}{}\leq{}2\Delta_{Y}V^{k-1}\sum_{m=0}^{\ell}\theta^{\ell}$

for all $k>k_{\ast}$ and $\ell\in\mathbb{N}$. Hence, since $\theta\in(0,1)$,
for any arbitrary but fixed $k>k_{\ast}$ taking the limit results in

$\limsup_{\ell\to\infty}\mathcal{L}_{\mu_{\ast}}(x^{k+\ell},\hat{y}^{k+\ell})\leq\mathcal{L}_{\mu_{\ast}}(x^{k-1},\hat{y}^{k-1})+\frac{2\Delta_{Y}}{1-\theta}V^{k-1}<\infty,$

showing that, up to discarding early iterates, there is an upper bound to
$\mathcal{L}_{\mu_{k}}(x^{k},\hat{y}^{k})$. Then, feasibility of $x^{\ast}$
follows as in (ii). ∎

The following result concerns the Lagrange multiplier estimates. Considering
the particular case of $g$ being the indicator of some set $D$ with nonempty
interior, Proposition 3.5 shows that the Lagrange multipliers vanish for
constraints that are (strictly) inactive in the limit. This indicates that the
update rule at Algorithm 3.1 of Algorithm 3.1 is plausible with regard to
their expected behavior as shifts [6, Section 4].

######  Proposition 3.5.

Suppose $g\equiv\operatorname{\iota}_{D}$ for some nonempty and closed set
$D\subseteq\mathbb{R}^{m}$. Let Assumption 1 hold and consider a sequence
$\\{x^{k}\\}_{k\in\mathbb{N}}$ of iterates generated by Algorithm 3.1 with
$\epsilon=0$. Let $x^{\ast}$ be an accumulation point of
$\\{x^{k}\\}_{k\in\mathbb{N}}$ and $\\{x^{k}\\}_{k\in K}$ a subsequence such
that $x^{k}\to_{K}x^{\ast}$. If $c(x^{\ast})\in\operatorname{int}D$, then
$y^{k}=0$ for all $k\in K$ large enough.

###### Proof.

Let $c(x^{\ast})\in\operatorname{int}D\neq\emptyset$ and notice that the
proximal mapping of $\operatorname{\iota}_{D}$ coincides with the projection
onto $D$, namely
$\operatorname{prox}_{\gamma\operatorname{\iota}_{D}}=\operatorname{\Pi}_{D}$
for all $\gamma>0$. Consider the two cases.

  * •

If $\mu_{k}\to 0$, by boundedness of $\\{\hat{y}^{k}\\}_{k\in\mathbb{N}}$ and
continuity of $c$, there exists $k_{1}\in K$ such that
$c(x^{k})+\mu_{k}\hat{y}^{k}\in\operatorname{int}D$ for all $k\in K$, $k\geq
k_{1}$.

  * •

If $\\{\mu_{k}\\}_{k\in\mathbb{N}}$ is bounded away from zero, it follows from
Algorithm 3.1 that $c(x^{k})-z^{k}\to_{K}0$. Then, by continuity of $c$ and
$x^{k}\to_{K}x^{\ast}$, $z^{k}\to_{K}c(x^{\ast})\in\operatorname{int}D$
follows. Hence, there exists $k_{2}\in K$ such that
$z^{k}\in\operatorname{int}D$ for all $k\in K$, $k\geq k_{2}$.

In both cases, since
$z^{k}\in\operatorname{\mathbf{O}}_{\mu_{k}}(x^{k},\hat{y}^{k})=\operatorname{\Pi}_{D}(c(x^{k})+\mu_{k}\hat{y}^{k})$
and $z^{k}\in\operatorname{int}D$ for all $k\in K$ large enough, it must be
$z^{k}=c(x^{k})+\mu_{k}\hat{y}^{k}$ and, consequently, $y^{k}=0$ by Algorithm
3.1. ∎

Notice that the assertion can be easily refined by exploiting the separable
structure of $D$, if any. For instance, if $D$ is a hyperbox, every step in
the proof can be applied componentwise, recovering [6, Theorem 4.1].

The following convergence result provides fundamental theoretical support to
Algorithm 3.1. It shows that, under subsequential attentive convergence,
accumulation points are AKKT-stationary for (P). For any accumulation point
$x^{\ast}$ and for any subsequence $\\{x^{k}\\}_{k\in K}$ such that
$x^{k}\to_{K}x^{\ast}$, Proposition 3.4(i) assures the feasibility of
$x^{\ast}$ under the assumption that
$x^{k}\overset{\varphi}{\to}_{K}x^{\ast}$. However, the $\varphi$-attentive
convergence also entails $g(c^{k})\to_{K}g(c(x^{\ast}))$ as
$x^{k}\to_{K}x^{\ast}$, which is trivially satisfied, e.g., if ($x^{\ast}$ is
feasible and) $g$ is continuous relative to $\operatorname{dom}g$, in the
sense that, for all $x\in\operatorname{dom}g$ and
$\\{x^{\ell}\\}\subseteq\operatorname{dom}g$, whenever $x^{\ell}\to x$ it
holds that $g(x^{\ell})\to g(x)$. We refer to [15, Example 3.4] for an
illustration of the importance of $\varphi$-attentive convergence.

######  Theorem 3.6.

Let Assumption 1 hold and consider a sequence $\\{x^{k}\\}_{k\in\mathbb{N}}$
of iterates generated by Algorithm 3.1 with $\epsilon=0$ and
$\varepsilon_{k}\to 0$. Let $x^{\ast}\in\mathbb{R}^{n}$ be an accumulation
point of $\\{x^{k}\\}_{k\in\mathbb{N}}$ and $\\{x^{k}\\}_{k\in K}$ a
subsequence such that $x^{k}\overset{\varphi}{\to}_{K}x^{\ast}$. Then,
$x^{\ast}$ is an AKKT-stationary point for (P).

###### Proof.

We claim that the subsequences $\\{x^{k}\\}_{k\in K}$, $\\{y^{k}\\}_{k\in K}$
and $\\{z^{k}\\}_{k\in K}$ generated by Algorithm 3.1 satisfy the properties
in Definition 2.3 and therefore show that $x^{\ast}$ is an AKKT-stationary
point for (P). First of all, by Proposition 3.4(i), the $\varphi$-attentive
convergence of $\\{x^{k}\\}_{k\in K}$ implies that $x^{\ast}$ is feasible for
(P). From Algorithms 3.1 and 3.1 of Algorithm 3.1, we have by
$\varepsilon_{k}$-$\Upsilon$-stationarity of $x^{k}$ with certificate $z^{k}$
that

$\eta^{k}:=\nabla
f(x^{k})+c^{\prime}(x^{k})^{\top}y^{k}=\nabla_{x}\mathcal{L}^{\text{S}}_{\mu_{k}}(x^{k},z^{k},\hat{y}^{k})$

satisfies $\|\eta^{k}\|\leq\varepsilon_{k}$. Then, by construction, we have
$\eta^{k}\to_{K}0$. Further, from Algorithm 3.1 we obtain that, for all
$k\in\mathbb{N}$,

$y^{k}=\hat{y}^{k}+\frac{c(x^{k})-z^{k}}{\mu_{k}}\in\partial g(z^{k}),$

where the inclusion follows from
$z^{k}\in\operatorname{\mathbf{O}}_{\mu_{k}}(x^{k},\hat{y}^{k})$ (in fact,
from the necessary condition for the optimality of $z^{k}$). It remains to
show that $c(x^{k})-z^{k}\to_{K}0$. To this end, we consider the two cases.

  * •

If $\\{\mu_{k}\\}_{k\in\mathbb{N}}$ is bounded away from zero, the conditions
at Algorithm 3.1 imply that $V^{k}:=\|c(x^{k})-z^{k}\|\to 0$, hence in
particular $c(x^{k})-z^{k}\to_{K}0$.

  * •

If instead $\mu_{k}\to 0$, we exploit continuity of $c$ and boundedness of
$\\{\hat{y}^{k}\\}_{k\in\mathbb{N}}$. These properties result in
$c(x^{k})+\mu_{k}\hat{y}^{k}\to_{K}c(x^{\ast})\in\operatorname{dom}g$, where
the inclusion follows from feasibility of $x^{\ast}$. Then, since
$z^{k}\in\operatorname{\mathbf{O}}_{\mu_{k}}(x^{k},\hat{y}^{k})$ for all
$k\in\mathbb{N}$, from (9) and [16, Lemma 3.3] we have that

$z^{k}\in\operatorname{prox}_{\mu_{k}g}(c(x^{k})+\mu_{k}\hat{y}^{k})\to_{K}c(x^{\ast})$

while $\mu_{k}\to 0$. Thus, the limit $c(x^{k})-z^{k}\to_{K}0$ follows.

Overall, this proves that $x^{\ast}$ is an AKKT-stationary point for (P). ∎

Constrained optimization algorithms aim at finding feasible points and
minimizing the objective function subject to constraints. Employing affordable
local optimization techniques, one cannot expect to find global minimizers of
any infeasibility measure. Nevertheless, the next result proves that Algorithm
3.1 with bounded $\\{\varepsilon_{k}\\}_{k\in\mathbb{N}}$ finds at least
stationary points of an infeasibility measure. Notice that this property does
not require $\varepsilon_{k}\to 0$, but only boundedness; cf. [6, Theorem
6.3].

######  Proposition 3.7.

Let Assumption 1 hold and suppose $\operatorname{dom}g$ is closed. Consider a
sequence $\\{x^{k}\\}_{k\in\mathbb{N}}$ of iterates generated by Algorithm 3.1
with $\epsilon=0$ and $\\{\varepsilon_{k}\\}_{k\in\mathbb{N}}$ bounded. Let
$x^{\ast}\in\mathbb{R}^{n}$ be an accumulation point of
$\\{x^{k}\\}_{k\in\mathbb{N}}$ and $\\{x^{k}\\}_{k\in K}$ a subsequence such
that $x^{k}\to_{K}x^{\ast}$. Then, $x^{\ast}$ is an M-stationary point of
$\operatorname*{minimize}_{x\in\mathbb{R}^{n}}\,\operatorname{dist}^{2}(c(x),\operatorname{dom}g).$

###### Proof.

By Proposition 3.4(iii) and $\operatorname{dom}g$ being closed, if
$\\{\mu_{k}\\}_{k\in\mathbb{N}}$ is bounded away from zero, then each
accumulation point $x^{\ast}$ is feasible for (P), namely it is a global
minimizer of the feasibility problem, and so it is an M-stationary point
thereof by strict continuity of the objective function [30, Proposition 5.3].
Hence, it remains to consider the case $\mu_{k}\to 0$.

Let us rewrite the feasibility problem with the equivalent lifted
representation

$\operatorname*{minimize}_{\begin{subarray}{c}x\in\mathbb{R}^{n}\\\
(z,\alpha)\in\operatorname{epi}g\end{subarray}}{}\quad{}\frac{1}{2}\|c(x)-z\|^{2}$
(18)

obtained by epigraphical reformulation. The inclusion
$(z,\alpha)\in\operatorname{epi}g$ requires that $g(z)\leq\alpha<\infty$, thus
encoding $z\in\operatorname{dom}g$. M-stationarity of $x^{\ast}$ for (18)
requires the existence of some $z^{\ast}\in\mathbb{R}^{m}$ such that

$\displaystyle c^{\prime}(x^{\ast})^{\top}[c(x^{\ast})-z^{\ast}]{}={}$
$\displaystyle 0,$ (19a) $\displaystyle(c(x^{\ast})-z^{\ast},0){}\in{}$
$\displaystyle\mathcal{N}^{\text{lim}}_{\operatorname{epi}g}(z^{\ast},g(z^{\ast})).$
(19b)

Owing to Algorithm 3.1 of Algorithm 3.1 it is
$\|\eta^{k}\|\leq\varepsilon_{k}$, where

$\eta^{k}:=\nabla
f(x^{k})+\frac{1}{\mu_{k}}c^{\prime}(x^{k})^{\top}\left[c(x^{k})+\mu_{k}\hat{y}^{k}-z^{k}\right],$

for all $k\in K$. Taking the limit of $\mu_{k}\eta^{k}$ for $k\to_{K}\infty$,
while using $\mu_{k}\to 0$ and boundedness of
$\\{\eta^{k}\\}_{k\in\mathbb{N}}$ and $\\{\nabla f(x^{k})\\}_{k\in K}$, it
follows that (19a) is satisfied. Furthermore, by the minimizing property of
$z^{k}{}\in{}\operatorname{\mathbf{O}}_{\mu_{k}}(x^{k},\hat{y}^{k})=\operatorname{prox}_{\mu_{k}g}(c(x^{k})+\mu_{k}\hat{y}^{k})$,
it is $z^{k}\in\operatorname{dom}g$ and

$\frac{1}{\mu_{k}}\left[c(x^{k})+\mu_{k}\hat{y}^{k}-z^{k}\right]\in\partial
g(z^{k})$

for all $k\in\mathbb{N}$. By boundedness of
$\\{c(x^{k})+\mu_{k}\hat{y}^{k}\\}_{k\in K}$ and local boundedness of the
proximal mapping, $\\{z^{k}\\}_{k\in K}$ remains bounded too. Hence, since
$z^{k}\in\operatorname{dom}g$ for all $k\in\mathbb{N}$ and
$\operatorname{dom}g$ is closed, we may assume without loss of generality that
$z^{k}\to_{K}z^{\ast}$ for some $z^{\ast}\in\operatorname{dom}g$. Using the
relationship stated in (2), multiplying by $\mu_{k}>0$ and exploiting that
$\mathcal{N}^{\text{lim}}_{\operatorname{epi}g}(z^{k},g(z^{k}))$ is a cone, we
have

$(c(x^{k})+\mu_{k}\hat{y}^{k}-z^{k},-\mu_{k})\in\mathcal{N}^{\text{lim}}_{\operatorname{epi}g}(z^{k},g(z^{k}))$

for all $k\in\mathbb{N}$. Taking the limit for $k\to_{K}\infty$, the
robustness of the limiting normal cone [30, Proposition 1.20] and the
convergence $z^{k}\to_{K}z^{\ast}\in\operatorname{dom}g$ yield (19b),
concluding the proof. ∎

###  3.5 Termination Criteria

Algorithm 3.1 of Algorithm 3.1 involves the minimization of the implicit AL
function, defined in (10). Then, the multiplier update at Algorithm 3.1 allows
one to draw conclusions with respect to the original problem (P), as shown by
Theorem 3.6. Moreover, recalling the AKKT-stationarity conditions (5), in view
of Algorithm 3.1 one may construct a sequence
$\\{\varepsilon_{k}\\}_{k\in\mathbb{N}}\subset\mathbb{R}_{++}$ such that
$\varepsilon_{k}\to 0$. Then, given some user-specified tolerances
$\varepsilon^{\text{dual}},\varepsilon^{\text{prim}}>0$, it is reasonable to
declare successful convergence when the conditions
$\varepsilon_{k}\leq\varepsilon^{\text{dual}}$ and
$V^{k}:=\|c(x^{k})-z^{k}\|\leq\varepsilon^{\text{prim}}$ are satisfied.
Theorem 3.6 demonstrates that these termination criteria (the latter, in
particular) are satisfied in finitely many iterations if any subsequence of
$\\{x^{k}\\}_{k\in K}$ accumulates at a feasible point
$x^{\ast}\in\mathbb{R}^{n}$. As this might not be the case, a mechanism for
(local) infeasibility detection is needed, and usually included in practical
implementations.

Given some tolerances, Algorithm 3.1 can be equipped with relaxed conditions
on decrease requirements at Algorithm 3.1 and optimality at Algorithm 3.1 [5].
At Algorithm 3.1 the inner tolerance $\varepsilon_{k}$ can stay bounded away
from zero, as long as $\varepsilon_{k}\leq\varepsilon^{\text{dual}}$ for large
$k\in\mathbb{N}$. Similarly, the condition at Algorithm 3.1 can be relaxed by
adding the (inclusive) possibility that $V^{k}\leq\varepsilon^{\text{prim}}$,
thus avoiding unnecessarily strong penalization.

##  4 Inner Problem and Solver

In this section we elaborate upon the minimization of the implicit AL
function, defined in (10), focusing on the numerical solution of subproblems
(12) for some fixed $\mu_{k}>0$ and $\hat{y}^{k}\in\mathbb{R}^{m}$. As
remarked in Section 3.1, this task could be performed by jointly minimizing
the explicit AL function
$\mathcal{L}^{\text{S}}_{\mu_{k}}(\cdot,\cdot,\hat{y}^{k})$ in (8), that is,
with respect to both primal and slack variables. Owing to its composite
structure, any suitable methods for composite optimization could be adopted
for this task; as pointed out in [15, Section 3.4], a prominent example are
proximal gradient algorithms.

We take a different approach, aimed at minimizing
$\mathcal{L}_{\mu_{k}}(\cdot,\hat{y}^{k})$ directly, in order to solve inner
problems without auxiliary variables. To do so, we juggle the implicit and
explicit formulations, linked by the oracle, as in (13). In Section 4.1, an
algorithm fo the minimization of $p$ is discussed, that exploits the available
oracles and the structure of $P$ in (17). Reminiscent of [13, Section 1.5.4],
the scheme compensates for the lack of regularity by pivoting around the
oracle’s optimality, as a kind of reliable _magical step_ [11, 6]. In fact, a
surrogate gradient of $p$ is built based on that of $F$, the smooth part of
$P$, evaluated at the oracle’s response. Then, with a bit of abstraction, it
turns out that any off-the-shelf algorithm for unconstrained smooth
optimization can be directly adopted for the task of minimizing the implicit
AL function, in a remarkably transparent way.

Throughout this section, let $\mu_{k}>0$ and $\hat{y}^{k}\in\mathbb{R}^{m}$ be
some fixed given parameters. Denote $p\colon\mathbb{R}^{n}\to\mathbb{R}$ and
$F\colon\mathbb{R}^{n}\times\mathbb{R}^{m}\to\mathbb{R}$ as the implicit AL
and its smooth part, respectively. Also, let
$P\colon\mathbb{R}^{n}\times\mathbb{R}^{m}\to\overline{\mathbb{R}}$ be the
explicit AL; see Section 3.2. Patterning the transition from
$\mathcal{L}^{\text{S}}_{\mu}(\cdot,\cdot,y)$ in (8) to
$\mathcal{L}_{\mu}(\cdot,y)$ in (10), we suppose that an oracle $\mathbf{O}$
is available that performs the explicit (exact, global) minimization of $P$
with respect to $z\in\mathbb{R}^{m}$, apparently alluding to (9). In relation
to the $k$th iteration our AL method, we write
$\mathbf{O}=\operatorname{\mathbf{O}}_{\mu_{k}}(\cdot,\hat{y}^{k})$. Notice
that, for all $x\in\mathbb{R}^{n}$, it is
$\mathbf{O}(x)\subseteq\operatorname{dom}g$ and thus $p(x)<\infty$. We collect
some basic properties under Assumption 1.

######  Proposition 4.1.

Consider (13), (17), and let Assumption 1 hold. Then, (i)
$F\colon\mathbb{R}^{n}\times\mathbb{R}^{m}\to\mathbb{R}$ is continuously
differentiable; (ii) $g\colon\mathbb{R}^{m}\to\overline{\mathbb{R}}$ is
proper, lower semicontinuous, and prox-bounded; (iii)
$P\colon\mathbb{R}^{n}\times\mathbb{R}^{m}\to\overline{\mathbb{R}}$, $P(x,z)$
is level-bounded in $z$ locally uniformly in $x$; (iv)
$\mathbf{O}\colon\mathbb{R}^{n}\rightrightarrows\mathbb{R}^{m}$ is locally
bounded, nonempty- and compact-valued; (v)
$p\colon\mathbb{R}^{n}\to\mathbb{R}$ is strictly continuous (i.e. locally
Lipschitz).

###### Proof.

The assertions readily follow from Assumption 1(i)–(ii) and the definitions in
(13) and (17). In particular, (iii) is due to prox-boundedness of $g$ and that
$\mu_{k}\leq\gamma_{g}$, along with $\hat{y}^{k}$ being fixed and smoothness
of $f$ and $c$. Assertion (iv) is a consequence of the smoothness of $c$ by
Assumption 1(i) and of the properties of the proximal mapping [39, Theorem
1.25]. Assertion (v) follows from the identity
$p(\cdot)=f(\cdot)+g^{\mu_{k}}(c(\cdot)+\mu_{k}\hat{y}^{k})$ and the fact that
the Moreau envelope is strictly continuous [39, Example 10.32]. ∎

Owing to Assumption 1 and Proposition 4.1(v), the objective function $p$ is
proper and lsc. However, as the existence of a minimizer for $p$ is not
guaranteed a priori, we also consider the following assumption. In general, it
is verified under some coercivity or (level) boundedness conditions.

######  Assumption 2.

Consider (13). Then, $\inf p\in\mathbb{R}$.

######  Remark 4.1.

The problem format (P) could incorporate also an abstract constraint $x\in X$,
with the set $X\subseteq\mathbb{R}^{n}$ being nonempty, closed, and _simple_.
Although not included for the sake of a simpler exposition, the results in
this paper could be readily extended to account for a _convex_ $X$. As noted
in [37], a possibility is to embed $x\in X$ in the template (P) by replacing
$c(x)$ by $\bar{c}(x):=(c(x),x)$ and $g(z)$ by
$\bar{g}(z,z^{\prime}):=g(z)+\operatorname{\iota}_{X}(z^{\prime})$. This
approach, however, introduces an extra slack component $z^{\prime}$, with an
extra multiplier component $y^{\prime}$, while relaxing the simple constraint
$x\in X$. Following [1] instead, a sharper strategy is to treat $x\in X$
directly as an unrelaxable constraint for the AL subproblem, and not
incorporating it into the implicit AL function. Exploiting the projection onto
$X$, the inner solver could generate iterates that satisfy the abstract
constraint by design. With potentially nonconvex $X$, though, it appears
nontrivial to formulate a suitable termination criterion for attesting
approximate stationarity.

The possibility to localize the AL subproblem within a compact set $X$ helps
in securing Assumption 2. Together with properness and lower semicontinuity of
$p$, compactness of $X\subset\mathbb{R}^{n}$ ensures that, once localized, the
minimization of $p$ is well-posed. In the spirit of trust region methods [11],
however, it would be necessary to adapt the admissible region $X$ to avoid
interfering with the outer AL loop convergence. These considerations are left
for a future effort.

###  4.1 Descent Methods with an Oracle

In this section we detail and discuss an algorithm for the minimization of
$p$, which proceeds by taking steps along directions of descent for $p$. The
oracle-equipped descent method is stated in Algorithm 4.1, which closely
patterns [6, Algorithm 8.1] when variable $z$ is neglected — recovering a
descent method for unconstrained smooth optimization. We prove well-
definedness of the algorithm and investigate its convergence properties in
Section 4.2.

Data: $x^{0}\in\mathbb{R}^{n}$, $\varepsilon>0$

Result: $\varepsilon$-$\Upsilon$-stationary point $x^{\ast}\in\mathbb{R}^{n}$
with certificate $z^{\ast}\in\mathbf{O}(x^{\ast})$

1 Select parameters $\alpha,\beta,\theta\in(0,1)$, $\omega,\nu\in(0,1]$ and
set $\Phi_{0}\leftarrow p(x^{0})$

2 for  _$j=0,1,\ldots$_ do

3 Select $z^{j}\in\mathbf{O}(x^{j})$ and set $\nabla
p^{j}\leftarrow\nabla_{x}F(x^{j},z^{j})$

4 if  _$\|\nabla p^{j}\|\leq\varepsilon$_ then return
_$(x^{\ast},z^{\ast})\leftarrow(x^{j},z^{j})$_

5 Compute $d^{j}\in\mathbb{R}^{n}$ such that  $\langle\nabla
p^{j},d^{j}\rangle{}\leq{}-\theta\|\nabla
p^{j}\|\|d^{j}\|\quad\text{and}\quad\|d^{j}\|{}\geq{}\omega\|\nabla p^{j}\|$

6 Compute $x^{j+1}\in\mathbb{R}^{n}$ and the largest
$\gamma_{j}=\beta^{\ell_{j}}$, $\ell_{j}\in\mathbb{N}$, such that
$p(x^{j+1})\leq
p(x^{j}+\gamma_{j}d^{j})\leq\Phi_{j}+\alpha\gamma_{j}\langle\nabla
p^{j},d^{j}\rangle$

7 Set $\Phi_{j+1}\leftarrow(1-\nu)\Phi_{j}+\nu p(x^{j+1})$

8

Algorithm 4.1 Nonmonotone Descent Method with an Oracle

Some comments are in order. First, the method requires a starting point
$x^{0}\in\mathbb{R}^{n}$ for variable $x$ only, not for $z$. Then, at
Algorithm 4.1, a surrogate gradient $\nabla p^{j}\in\mathbb{R}^{n}$ of $p$ at
$x^{j}$ is evaluated: the oracle returns an arbitrary element
$z^{j}\in\mathbf{O}(x^{j})$ and, in light of the structure (17), we let
$\nabla p^{j}$ be the partial derivative $\partial_{x}P$ at $(x^{j},z^{j})$.
At Algorithm 4.1 a suitable search direction $d^{j}\in\mathbb{R}^{n}$ has to
be selected, possibly among many. The simplest choice is the _steepest_
direction given by the negative gradient, i.e. $d^{j}:=-\nabla p^{j}$, which
is always valid with respect to the conditions at Algorithm 4.1, for any
choice of $\theta\in(0,1)$ and $\omega\in(0,1]$. Another possibility is to
consider Newton-like directions, based on local quadratic approximations.
Barzilai-Borwein stepsizes and (limited memory) BFGS directions are popular
strategies in this domain; see [40, Section 5.1.2] for an overview.
Afterwards, the inner loop at Algorithm 4.1 seeks an appropriate stepsize
$\gamma_{j}>0$ such that the update $x^{j}+\gamma_{j}d^{j}$ yields some
improvement. In particular, this backtracking linesearch procedure is herein
included as a globalization mechanism. Finally, with $\nu\in(0,1)$, the update
rule at Algorithm 4.1 allows to relax the requirement of a monotone decrease
for the objective values, thus reducing conservatism of the linesearch and
possibly encouraging longer steps and improved practical performance [6,
Section 8.2], [40]. However, a monotone behavior for
$\\{p(x^{j})\\}_{j\in\mathbb{N}}$ can be enforced with $\nu:=1$.

Crucially, the tentative update and the backtracking condition at Algorithm
4.1 for adapting the stepsize $\gamma_{k}$ are concerned with variable $x$
only. Also, the oracle evaluation at Algorithm 4.1 can be formally excluded
from the algorithm and incorporated into the _surrogate_ gradient evaluation
of $p$ at $x^{j}$. This is possible thanks to the arbitrariness of the element
$z^{j}$ chosen by the oracle $\mathbf{O}$, which allows us to keep variable
$z$ implicit. These observations show that Algorithm 4.1 is equivalent to a
descent method applied directly to the minimization of $p$ (with a surrogate
gradient), thus making the oracle-equipped algorithm a formal tool only,
transparent to the user.

###  4.2 Convergence Analysis

In this section we investigate the convergence properties of Algorithm 4.1,
exploiting continuous differentiability of $F$ and the oracle’s optimality to
prove sufficient decrease with respect to $x$. We are going to show first that
Algorithm 4.1 is well defined, namely that each iteration terminates in finite
time. This is an essential theoretical and practical condition that an
implementable algorithm must satisfy.

######  Lemma 4.2 (Well definedness).

Let Assumption 1 hold and consider the iterates generated by Algorithm 4.1.
Then, (i) at every iteration, the number of backtrackings at Algorithm 4.1 is
finite;  (ii) at the end of the $j$-th iteration, one has
$p(x^{j+1})\leq\Phi_{j}+\delta_{j}\qquad\text{and}\qquad
p(x^{j+1})\leq\Phi_{j+1}\leq\Phi_{j}+\nu\delta_{j}$ where
$\delta_{j}:=\alpha\gamma_{j}\langle\nabla p^{j},d^{j}\rangle<0$;  (iii) every
iterate $x^{j}$ remains within the sublevel set
$\\{x\in\mathbb{R}^{n}\,|\,p(x)\leq\Phi_{0}<\infty\\}$.

###### Proof.

Recall the basic properties collected in Proposition 4.1.

  * (i)

Let us discard the trivial case by assuming that $\|\nabla
p^{j}\|>\varepsilon$. Since $\theta\in(0,1)$ and $\omega\in(0,1]$, the choice
$d^{j}\leftarrow-\nabla p^{j}$ is always viable, hence Algorithm 4.1 is well-
defined. Denote $x^{j+}:=x^{j}+\gamma_{j}d^{j}$ and
$z^{j+}\in\mathbf{O}(x^{j+})$. Then,

$p(x^{j+})=F(x^{j+},z^{j+})+g(z^{j+})\leq F(x^{j+},z)+g(z)$

for all $z\in\mathbb{R}^{m}$. In particular, $p(x^{j+})\leq
F(x^{j+},z^{j})+g(z^{j})$. By continuous differentiability of $F$, combining
Algorithms 4.1 and 4.1 yields

$\lim\limits_{\gamma\to 0}\frac{F(x^{j}+\gamma
d^{j},z^{j})-F(x^{j},z^{j})}{\gamma}=\langle\nabla_{x}F(x^{j},z^{j}),d^{j}\rangle=\langle\nabla
p^{j},d^{j}\rangle<0.$

Thus,

$\lim\limits_{\gamma\to 0}\frac{F(x^{j}+\gamma
d^{j},z^{j})-F(x^{j},z^{j})}{\gamma\langle\nabla p^{j},d^{j}\rangle}=1.$

Therefore, since $\alpha\in(0,1)$, for some $\gamma_{j}=\beta^{\ell_{j}}>0$
small enough we have that

$\frac{F(x^{j}+\gamma_{j}d^{j},z^{j})-F(x^{j},z^{j})}{\gamma_{j}\langle\nabla
p^{j},d^{j}\rangle}\geq\alpha.$

As $\langle\nabla p^{j},d^{j}\rangle<0$ by Algorithm 4.1 and $\gamma_{j}>0$,
we deduce that

$F(x^{j+},z^{j})\leq F(x^{j},z^{j})+\alpha\gamma_{j}\langle\nabla
p^{j},d^{j}\rangle.$

Now, invoking the oracle to get $z^{j+}\in\mathbf{O}(x^{j+})$ we obtain

$p(x^{j+})=F(x^{j+},z^{j+})+g(z^{j+})\leq F(x^{j+},z^{j})+g(z^{j})\\\ \leq
F(x^{j},z^{j})+g(z^{j})+\alpha\gamma_{j}\langle\nabla
p^{j},d^{j}\rangle=p(x^{j})+\alpha\gamma_{j}\langle\nabla p^{j},d^{j}\rangle.$

Then, finite termination of the inner loop follows from the bound
$\Phi_{j}\geq p(x^{j})$ valid for all $j\in\mathbb{N}$. The latter is obtained
by induction: if $j=0$, it is $\Phi_{j}=p(x^{j})$ by construction; if $j>0$,
the bound is due to Algorithm 4.1.

  * (ii)

The first inequality follows by combining (i) with the failure of the
condition at Algorithm 4.1. Together with Algorithm 4.1, this gives the lower
bound

$\Phi_{j+1}=(1-\nu)\Phi_{j}+\nu
p(x^{j+1})\geq(1-\nu)[p(x^{j+1})-\delta_{j}]+\nu p(x^{j+1})\\\
=p(x^{j+1})-(1-\nu)\delta_{j}\geq p(x^{j+1})$

as well as the upper bound

$\Phi_{j+1}=(1-\nu)\Phi_{j}+\nu
p(x^{j+1})\leq(1-\nu)\Phi_{j}+\nu(\Phi_{j}+\delta_{j})=\Phi_{j}+\nu\delta_{j}.$

  * (iii)

That $\Phi_{0}=p(x^{0})<\infty$ follows from the fact that
$z^{0}\in\mathbf{O}(x^{0})\in\operatorname{dom}g\neq\emptyset$. Then, the
assertion is a direct consequence of (ii), by telescoping
$p(x^{j})\leq\Phi_{j}\leq\ldots\leq\Phi_{0}$. ∎

We next consider an asymptotic analysis of the algorithm.

######  Lemma 4.3.

Let Assumptions 1 and 2 hold and consider a sequence of iterates
$\\{x^{j}\\}_{j\in\mathbb{N}}$ generated by Algorithm 4.1 with
$\varepsilon=0$. Then, the sequences $\\{\Phi_{j}\\}_{j\in\mathbb{N}}$ and
$\\{p(x^{j})\\}_{j\in\mathbb{N}}$ converge to a finite value $p_{\ast}\geq\inf
p$, the former from above.

###### Proof.

Regarding $\\{\Phi_{j}\\}_{j\in\mathbb{N}}$, the assertion follows from Lemma
4.2(ii) and Assumption 2. Then, rearranging the update rule at Algorithm 4.1,
one can observe that $p(x^{j+1})=\nu^{-1}(\Phi_{j+1}-\Phi_{j})+\Phi_{j}$
holds. Taking the limit for $j\to\infty$, by $\nu>0$ the convergence of
$\\{\Phi_{j}\\}_{j\in\mathbb{N}}$ implies that of
$\\{p(x^{j})\\}_{j\in\mathbb{N}}$ to the same finite value $p_{\ast}$. ∎

The following Theorem 4.4 provides a guarantee of global convergence to
$\Upsilon$-stationary points, namely that any accumulation point of the
iterates generated by Algorithm 4.1 is $\Upsilon$-stationary relative to the
minimization of $p$.

######  Theorem 4.4.

Let Assumptions 1 and 2 hold and consider a sequence of iterates
$\\{x^{j}\\}_{j\in\mathbb{N}}$ generated by Algorithm 4.1 with
$\varepsilon=0$. Let $x^{\ast}\in\mathbb{R}^{n}$ be an accumulation point of
$\\{x^{j}\\}_{j\in\mathbb{N}}$. Then, $x^{\ast}$ is $\Upsilon$-stationary for
$p$, namely there exists $z^{\ast}\in\mathbf{O}(x^{\ast})$ such that
$\nabla_{x}F(x^{\ast},z^{\ast})=0$.

###### Proof.

If Algorithm 4.1 terminates after finitely many iterations returning
$(x^{j},z^{j})$, then it must be that $\|\nabla p^{j}\|=0$, hence
$\nabla_{x}F(x^{j},z^{j})=0$. Since $z^{j}\in\mathbf{O}(x^{j})$ for all
$j\in\mathbb{N}$, then $x^{j}$ would be $\Upsilon$-stationary.

Now suppose instead that Algorithm 4.1 takes infinitely many iterations. Let
$\\{x^{j}\\}_{j\in J}$, $J:=\\{j_{0},j_{1},j_{2}\ldots\\}\subseteq\mathbb{N}$,
be a subsequence such that $x^{j}\to_{J}x^{\ast}$. Then, continuity of $p$ by
Proposition 4.1(v) gives $p(x^{j})\to_{J}p(x^{\ast})$. Owing to the linesearch
condition at Algorithm 4.1, since $j_{i+1}\geq j_{i}+1$, we have that

$p(x^{j_{i+1}})\leq
p(x^{j_{i}+1})\leq\Phi_{j_{i}}+\alpha\gamma_{j_{i}}\langle\nabla
p^{j_{i}},d^{j_{i}}\rangle<\Phi_{j_{i}}.$

Then, Lemma 4.3 implies that $\gamma_{j}\langle\nabla
p^{j},d^{j}\rangle\to_{J}0$, since $\alpha>0$, and, due to Algorithm 4.1, it
must be that $\gamma_{j}\|\nabla p^{j}\|\|d^{j}\|\to_{J}0$. As any term in the
limit could vanish, there are three (non exclusive) cases. Now we argue that,
up to extracting a subsequence, it is $\nabla p^{j}\to_{J}0$ in every case.
First, the claim is trivial if $\nabla p^{j}$ vanishes. Second, if
$d^{j}\to_{J}0$, conditions at Algorithm 4.1 imply that $\nabla p^{j}$
vanishes too, giving the claim. Third, if $\gamma_{j}\to_{J}0$, we proceed by
showing that there is a longer stepsize $\hat{\gamma}_{j}$ that also tends to
zero, along which the function did not decrease enough, and that this is
possible only if the gradient vanishes.

Let us start by noticing the failure of attempts at Algorithm 4.1. The
previous (hence rejected) stepsize
$\hat{\gamma}_{j}:=\beta^{-1}\gamma_{j}=\beta^{\ell_{j}-1}$ gives

$p(x^{j}+\hat{\gamma}_{j}d^{j}){}>{}\Phi_{j}+\alpha\hat{\gamma}_{i}\langle\nabla
p^{j},d^{j}\rangle{}\geq{}p(x^{j})+\alpha\hat{\gamma}_{j}\langle\nabla
p^{j},d^{j}\rangle,$ (20)

where the last inequality follows from the bound $\Phi_{j}\geq p(x^{j})$ of
Lemma 4.2(ii). In the case $\gamma_{j}\to_{J}0$, it must be
$\hat{\gamma}_{j}\to_{J}0$ too, since $\ell_{j}\to_{J}0$. Rearranging terms
and using the structure of $P$ in (17), the minimizing property of any
$z^{j}\in\mathbf{O}(x^{j})$ yields

$\alpha\langle\nabla
p^{j},d^{j}\rangle<\frac{p(x^{j}+\frac{\gamma_{j}}{\beta}d^{j})-p(x^{j})}{\frac{\gamma_{j}}{\beta}}\\\
\leq\frac{F(x^{j}+\frac{\gamma_{j}}{\beta}d^{j},z^{j})+g(z^{j})-F(x^{j},z^{j})-g(z^{j})}{\frac{\gamma_{j}}{\beta}}\\\
=\frac{F(x^{j}+\frac{\gamma_{j}}{\beta}d^{j},z^{j})-F(x^{j},z^{j})}{\frac{\gamma_{j}}{\beta}}.$

Defining $s^{j}:=\hat{\gamma}_{j}d^{j}$ for all $j\in\mathbb{N}$, we have that
$\|s^{j}\|\to_{J}0$. Possibly extracting a subsequence and redefining $J$, let
$s^{\ast}\in\mathbb{R}^{n}$ be such that $s^{j}/\|s^{j}\|\to_{J}s^{\ast}$.
Now, for all $j\in J$, by the mean value theorem there exists
$\zeta_{j}\in[0,1]$ such that

$\displaystyle\langle\nabla_{x}F(x^{j}+\zeta_{j}s^{j},z^{j}),s^{j}\rangle{}={}$
$\displaystyle F(x^{j}+s^{j},z^{j})-F(x^{j},z^{j})$ $\displaystyle{}={}$
$\displaystyle F(x^{j}+s^{j},z^{j})+g(z^{j})-F(x^{j},z^{j})-g(z^{j})$
$\displaystyle{}\geq{}$ $\displaystyle p(x^{j}+s^{j})-p(x^{j})$
$\displaystyle{}>{}$ $\displaystyle\alpha\langle\nabla
p^{j},s^{j}\rangle{}={}\alpha\langle\nabla_{x}F(x^{j},z^{j}),s^{j}\rangle,$

where the inequalities are due to the oracle’s optimality and (20). Dividing
by $\|s^{j}\|$ and rearranging give

$\frac{\langle\nabla_{x}F(x^{j}+\zeta_{j}s^{j},z^{j})-\alpha\nabla_{x}F(x^{j},z^{j}),s^{j}\rangle}{\|s^{j}\|}>0.$

Owing to $\zeta_{j}s^{j}\to_{J}0$ and continuity of $\nabla F$, taking the
limit yields

$\limsup\limits_{j\to_{J}\infty}\frac{\langle\nabla_{x}F(x^{j}+\zeta_{j}s^{j},z^{j})-\alpha\nabla_{x}F(x^{j},z^{j}),s^{j}\rangle}{\|s^{j}\|}\\\
{}={}\limsup\limits_{j\to_{J}\infty}(1-\alpha)\frac{\langle\nabla_{x}F(x^{j},z^{j}),s^{j}\rangle}{\|s^{j}\|}{}\leq{}0,$

where the inequality is due to
$\langle\nabla_{x}F(x^{j},z^{j}),s^{j}\rangle=\hat{\gamma}_{j}\langle\nabla
p^{j},d^{j}\rangle<0$ for all $j\in\mathbb{N}$, by Algorithm 4.1, and
$\alpha<1$. Recalling that $s^{j}/\|s^{j}\|\to_{J}s^{\ast}$, the inequalities
above demonstrate that it is necessarily
$\langle\nabla_{x}F(x^{j},z^{j}),s^{\ast}\rangle=\langle\nabla
p^{j},s^{\ast}\rangle\to_{J}0$. Then, by the conditions at Algorithm 4.1, it
must be $\nabla p^{j}\to_{J}0$.

Finally, owing to $\nabla p^{j}:=\nabla_{x}F(x^{j},z^{j})\to_{J}0$ in all
cases, $z^{j}\in\mathbf{O}(x^{j})$ for all $j\in J$, $x^{j}\to_{J}x^{\ast}$,
and continuity of $\nabla F$, we conclude that
$\nabla_{x}F(x^{\ast},z^{\ast})=0$ for some $z^{\ast}\in\mathbf{O}(x^{\ast})$,
proving the claim that $x^{\ast}$ is $\Upsilon$-stationary for $p$. ∎

###  4.3 Termination Criteria

Algorithm 4.1 would run indefinitely without Algorithm 4.1, generating an
infinite sequence of iterates. The merit and objective values tend to decrease
along the iterations and settle at some finite value, as established by Lemma
4.2(ii) and Lemma 4.3, respectively. Suitable termination criteria need to be
inserted for stopping and returning an iterate that, in some sense, is
satisfactory relatively to the minimization of $p$.

Since $z^{j}\in\mathbf{O}(x^{j})$ for all $j\in\mathbb{N}$, it suffices to
check the (approximate) stationarity of $x^{j}$. Thus, we equip Algorithm 4.1
with the termination condition
$\left\|\nabla_{x}F(x^{j},z^{j})\right\|\leq\varepsilon$, which encodes
$\varepsilon$-$\Upsilon$-stationarity of $x^{j}$ relative to $p$, by the
specific structure of $P$ in (17). Thanks to Theorem 4.4 we can now easily
infer termination of Algorithm 4.1 in finitely many steps for any tolerance
$\varepsilon>0$, which confirms that its output complies with the requirements
for the outer AL framework at Algorithm 3.1 of Algorithm 3.1.

######  Corollary 4.5.

Let Assumption 1 hold and consider a sequence of iterates
$\\{x^{j},z^{j}\\}_{j\in\mathbb{N}}$ generated by Algorithm 4.1 with
$\varepsilon>0$. Suppose that the sequence $\\{x^{j}\\}_{j\in\mathbb{N}}$
admits an accumulation point. Then, in finitely many steps an
$\varepsilon$-$\Upsilon$-stationary point $x^{\ast}\in\operatorname{dom}p$ is
returned with certificate $z^{\ast}\in\mathbf{O}(x^{\ast})$.

###### Proof.

That the algorithm terminates in finitely many iterates, say $j_{\ast}$ many,
follows from Theorem 4.4, since $\|\nabla_{x}F(x^{j},z^{j})\|=\|\nabla
p^{j}\|\to_{J}0$ along any convergent subsequence $\\{x^{j}\\}_{j\in J}$.
Since $z^{j}\in\mathbf{O}(x^{j})$ for all $j\in\mathbb{N}$, it follows that
the output $x^{\ast}:=x^{j_{\ast}}$ with certificate $z^{\ast}:=z^{j_{\ast}}$
satisfies $z^{\ast}\in\mathbf{O}(x^{\ast})$. Since the magnitude of
$\nabla_{x}F(x^{\ast},z^{\ast})$ is no more than $\varepsilon$ by the
termination criterion, the pair $(x^{\ast},z^{\ast})$ is
$\varepsilon$-$\Upsilon$-stationary for $p$. ∎

Note that the existence of limit points is assumed, and not guaranteed, in
Corollary 4.5. This can be satisfied by employing, in general, (level)
boundedness or coercivity arguments.

##  5 Numerical Example and Comparison

This section showcases the flexibility offered by the (nonconvex) generalized
programming framework and presents some numerical results that support the
implicit approach developed in Section 3. We consider as illustrative example
a well-known academic problem with vanishing constraints.

We have implemented the proposed approach in the Augmented Lagrangian Proximal
Solver (ALPS), as part of the open-source software package Bazinga.jl [12],
which includes the example discussed below. As our goal is to provide only a
motivation and numerical support for the implicit approach, we do not include
implementation details. Nonetheless, we shall point out that it closely
patterns that of ALS [15], an AL solver based on an explicit treatment of
slack variables. In particular, the AL subproblems at Algorithm 3.1 of
Algorithm 3.1 are solved by default using the accelerated proximal gradient
solver PANOC+ [17]. This choice is convenient in our numerical tests to obtain
comparable results across different approaches and formulations. As a proximal
gradient method, PANOC+ can handle the two-term structure of the explicit AL
function; it is in fact the default inner solver in ALS [15]. However, PANOC+
can act as a gradient-based descent method too, and thus can be adopted as
inner solver within Algorithm 3.1, in view of Section 4.1.

###  5.1 Problem with Vanishing Constraints

Let us consider the following two-dimensional problem, which arises from truss
topology optimization and is a well-known academic problem; it has been
thoroughly discussed in [25, Section 9.5.1].

$\displaystyle\operatorname*{minimize}_{x\in\mathbb{R}^{2}}\quad{}{}$
$\displaystyle 4x_{1}+2x_{2}$ (21)
$\displaystyle\operatorname{subject~{}to}\quad{}{}$ $\displaystyle x_{1}\geq
0,\quad x_{2}\geq 0$ $\displaystyle
x_{1}>0{}\quad\Rightarrow\quad{}x_{1}+x_{2}\geq 5\sqrt{2}$ $\displaystyle
x_{2}>0{}\quad\Rightarrow\quad{}x_{1}+x_{2}\geq 5$

Due to their structure, the logical implications give rise to the so called
vanishing constraints, usually reformulated by requiring
$x_{1}(x_{1}+x_{2}-5\sqrt{2}){}\geq{}0$ and $x_{2}(x_{1}+x_{2}-5){}\geq{}0$
along with $x_{1},x_{2}\geq 0$. The feasible set of (21) consists of the union
of an unbounded polyhedron
$\\{(x_{1},x_{2})\in\mathbb{R}_{+}^{2}\,|\,x_{1}+x_{2}\geq 5\sqrt{2}\\}$, of
an attached line segment $\\{(0,x_{2})\in\mathbb{R}^{2}\,|\,5\leq x_{2}\leq
5\sqrt{2}\\}$, and of the isolated point $(0,0)$. It admits a global minimum
at $x^{\star}=(0,0)$ and a local minimum at $x^{\diamond}=(0,5)$; see Figure
1.

There exist several ways to cast (21) into the form of (P). We are going to
consider three different formulations, with increasing degree of explicitness
for treating constraints: implicit, intermediate (or partially explicit,
avoiding slack variable components for trivial constraints), and (fully)
explicit. Let us define the cost function $f(x):=4x_{1}+2x_{2}$. Moreover, let
us conveniently denote by $D_{\text{VC}}\subset\mathbb{R}^{2}$ the set
associated to a vanishing constraint:
$D_{\text{VC}}:=\\{(a,b)\in\mathbb{R}^{2}\,|\,a\geq 0,ab\geq 0\\}$; it is
nonempty, closed, nonconvex, and its (set-valued) projection
$\operatorname{\Pi}_{D}$ can be easily evaluated.

  * •

An _implicit_ representation of (21) in the form (P) is given by the data
functions

$\displaystyle c_{\text{impl}}(x){}:={}$
$\displaystyle(x_{1},\,x_{1}+x_{2}-5\sqrt{2},\,x_{2},\,x_{1}+x_{2}-5)^{\top}$
$\displaystyle g_{\text{impl}}(c){}:={}$
$\displaystyle\operatorname{\iota}_{D_{\text{VC}}}(c_{1},c_{2})+\operatorname{\iota}_{D_{\text{VC}}}(c_{3},c_{4}).$

With size $(n,m)=(2,4)$, this seems the most compact formulation for (21).

  * •

An _explicit_ reformulation of (21) as (P) is obtained by introducing slack
variables as in (P${}_{\text{S}}$). This leads to a problem with size
$(n,m)=(6,8)$. Notice that this explicit reformulation is carried out
internally by ALS when fed with the implicit formulation, as a lifted
representation of the original problem [15, Section 3].

  * •

An _intermediate_ representation of (21) into the form of (P) is found by
avoiding slack variable components for trivial constraints. This intermediate
formulation has size $(n,m)=(4,6)$, with data functions

$\displaystyle c_{\text{inte}}(x){}:={}$
$\displaystyle(x_{1},\,x_{1}+x_{2}-5\sqrt{2}-x_{3},\,x_{2},\,x_{1}+x_{2}-5-x_{4},\,x_{3},\,x_{4})^{\top}$
$\displaystyle g_{\text{inte}}(c){}:={}$
$\displaystyle\operatorname{\iota}_{D_{\text{VC}}}(c_{1},c_{5})+\operatorname{\iota}_{D_{\text{VC}}}(c_{3},c_{6})+\operatorname{\iota}_{0}(c_{2})+\operatorname{\iota}_{0}(c_{4}).$

We execute ALPS on all problem representations. The outcome of ALS [15] on the
implicit formulation is also recorded, for comparison with the explicit one on
ALPS. We consider a grid of $51^{2}=2601$ starting points $x^{0}$ in
$[-5,20]^{2}$ and let the initial Lagrange multiplier estimate be $y^{0}=0$;
slack variables are initialized so as to satisfy the linear constraints
arising from the reformulation. Both solvers return a successful status for
all instances, always indicating either the global or local minimizer (within
$10^{-6}$ in Euclidean distance). In Figure 1 we depict the feasible region of
(21) and display these results in more detail. Each starting point $x^{0}$ is
given a mark, indicating the termination point which has been reached by
starting from $x^{0}$. Since the set $D_{\text{VC}}$ is nonconvex, the
iterates of Algorithm 3.1 depend on the arbitrary choices made when the
projection onto $D_{\text{VC}}$ is not a singleton. By selecting different
projections, we observed analogous numerical results, albeit not identical.

Table 1 collects statistics that summarize the numerical results. Considering
the _implicit_ formulation, ALPS returns the global minimizer
$x^{\star}=(0,0)$ in most cases ($93\%$). With the _intermediate_ formulation,
ALPS seems to consistently find $x^{\star}$ only when starting nearby
($28\%$), and even more so with the fully _explicit_ formulation ($19\%$) or
with ALS ($16\%$). Despite the fact that solutions and stationary points
remain unaffected by the lifting [15, Lemma 3.1], there is a significant gap
in terms of objective values attained in practice. But the advantages could go
beyond this aspect, as the implicit approach seems extremely effective in
reducing the (median) number of iterations and runtime too. We are drawn to
suspect that this is the oracle’s doing.

Figure 1: Comparison of solutions obtained using different methods and formulations. Feasible set (gray background) and grid of (2601 uniformly spaced) starting points and corresponding termination mark depending on the returned solution: global minimizer $x^{\star}=(0,0)$ (red circle) or local minimizer $x^{\diamond}=(0,5)$ (omitted mark). Table 1: Comparison of solutions obtained and performances using different methods and formulations. Statistics on the number of (outer and cumulative inner) iterations, elapsed runtime, and how often the global minimizer $x^{\star}$ is returned. |  | ALPS | ALPS | ALPS | ALS
---|---|---|---|---|---
|  | _implicit_ | _intermediate_ | _explicit_ |
outer | median | 10 | 10 | 10 | 10
iterations | max | 10 | 10 | 10 | 10
total inner | median | 10 | 77 | 117 | 117
iterations | max | 802 | 598 | 217 | 205
runtime | median | 0.17 | 2.1 | 1.9 | 3.0
[ms] | max | 155 | 130 | 46 | 36
global | abs | 2413 | 721 | 496 | 405
minimizer | rel | 92.77% | 27.72% | 19.07% | 15.57%

##  6 Concluding Remarks

The theoretical investigation of generalized optimization in the fully
nonconvex setting has led to a significant extension of the augmented
Lagrangian (AL) framework. With the help of a lifted reformulation, KKT-type
stationarity conditions have been derived in Lagrangian terms and complemented
with asymptotic and regularity notions. The implicit AL function and
subproblems have been defined, marginalizing their explicit counterparts with
respect to slack variables, seemingly sacrificing some regularity due to set-
valuedness. Then, we devised the concept of $\Upsilon$-stationarity, rooted on
a parametric interpretation, for characterizing the subproblems’ solutions.
Although weaker than M-stationarity for the implicit AL function,
$\Upsilon$-stationarity proved to be a strong enough tool for investigating
asymptotic properties. We considered a safeguarded AL algorithm and
established opportune convergence guarantees, on par with those in nonlinear
programming.

Embracing an implicit philosophy was not only a theoretical curiosity, though.
We developed and characterized a solver for the implicit AL subproblems, which
is indeed capable of generating approximate $\Upsilon$-stationary iterates.
Equipped with a proximal oracle, it handles fewer decision variables,
maintaining slack variables under the hood. This appears to produce a tighter
model and catalyze improved practical performance, compared to more explicit
approaches, as witness by numerical results on an illustrative problem.

Several issues and questions remain open. It is a common practice to assume
that the AL subproblem are well-posed, possibly based on coercivity or
boundedness arguments. In general, however, to overcome this concern one could
investigate the combination of AL schemes with trust-region strategies, as
hinted at in Remark 4.1. Another topic of interest is that of attentive
convergence, whose role should be closely inspected. It remains unclear
whether it can be discarded or not, in relation to asymptotic stationarity
conditions. Finally, future works should elaborate on a characterization of
stationary points as saddle points of the implicit AL function, a local
duality theory for generalized optimization, and the relationships with
proximal point algorithms.

### Acknowledgements

I wish to thank Andreas Themelis (Kyushu University) and Patrick Mehlitz (BTU
Cottbus-Senftenberg), for critical and constructive comments on an early draft
of this manuscript.

## References

  * [1] Roberto Andreani, Ernesto G. Birgin, José M. Martínez, and Maria Laura Schuverdt.  On augmented Lagrangian methods with general lower–level constraints.  SIAM Journal on Optimization, 18(4):1286–1309, 2008.
  * [2] Egon Balas.  Disjunctive Programming.  Springer, Cham, 2018.
  * [3] Matúš Benko and Patrick Mehlitz.  On implicit variables in optimization theory.  Journal of Nonsmooth Analysis and Optimization, 2:7215, 2021.
  * [4] Dimitri P. Bertsekas.  Constrained Optimization and Lagrange Multiplier Methods.  Athena Scientific, 1996.
  * [5] Ernesto G. Birgin and José M. Martínez.  Augmented Lagrangian method with nonmonotone penalty parameters for constrained optimization.  Computational Optimization and Applications, 51(3):941–965, 2012\.
  * [6] Ernesto G. Birgin and José M. Martínez.  Practical Augmented Lagrangian Methods for Constrained Optimization.  Society for Industrial and Applied Mathematics, Philadelphia, PA, 2014\.
  * [7] Stephen Boyd, Neal Parikh, Eric Chu, Borja Peleato, and Jonathan Eckstein.  Distributed optimization and statistical learning via the alternating direction method of multipliers.  Foundations and Trends in Machine Learning, 3(1):1–122, 2011.
  * [8] Xiaojun Chen, Lei Guo, Zhaosong Lu, and Jane J. Ye.  An augmented Lagrangian method for non-Lipschitz nonconvex programming.  SIAM Journal on Numerical Analysis, 55:168–193, 2017.
  * [9] Patrick L. Combettes and Jean-Christophe Pesquet.  Proximal splitting methods in signal processing, pages 185–212.  Springer, New York, 2011.
  * [10] Andrew R. Conn, Nicholas I. M. Gould, and Philippe L. Toint.  A globally convergent augmented Lagrangian algorithm for optimization with general constraints and simple bounds.  SIAM Journal on Numerical Analysis, 28(2):545–572, 1991.
  * [11] Andrew R. Conn, Nicholas I. M. Gould, and Philippe L. Toint.  Trust Region Methods.  Society for Industrial and Applied Mathematics, 2000.
  * [12] Alberto De Marchi.  Bazinga.jl: a toolbox for constrained structured optimization in Julia.  Online: https://github.com/aldma/Bazinga.jl, 2020\.  Accessed: 2023-01-21.
  * [13] Alberto De Marchi.  Augmented Lagrangian and Proximal Methods for Constrained Structured Optimization.  PhD thesis, University of the Bundeswehr Munich, 2021.
  * [14] Alberto De Marchi.  On a primal-dual Newton proximal method for convex quadratic programs.  Computational Optimization and Applications, 81(2):369–395, 3 2022\.
  * [15] Alberto De Marchi, Xiaoxi Jia, Christian Kanzow, and Patrick Mehlitz.  Constrained composite optimization and augmented Lagrangian methods.  Mathematical Programming, 201(1):863–896, 2023.
  * [16] Alberto De Marchi and Andreas Themelis.  An interior proximal gradient method for nonconvex optimization, 2022\.  arXiv:2208.00799v1.
  * [17] Alberto De Marchi and Andreas Themelis.  Proximal gradient algorithms under local Lipschitz gradient continuity.  Journal of Optimization Theory and Applications, 194(3):771–794, 2022.
  * [18] Neil K. Dhingra, Sei Zhen Khong, and Mihailo R. Jovanović.  The proximal augmented Lagrangian method for nonsmooth composite optimization.  IEEE Transactions on Automatic Control, 64(7):2861–2868, 2019.
  * [19] Brecht Evens, Puya Latafat, Andreas Themelis, Johan Suykens, and Panagiotis Patrinos.  Neural network training as an optimal control problem: An augmented Lagrangian approach.  In 2021 60th IEEE Conference on Decision and Control (CDC), pages 5136–5143. IEEE, 2021.
  * [20] Philip E. Gill and Daniel P. Robinson.  A primal-dual augmented Lagrangian.  Computational Optimization and Applications, 51(1):1–25, 2012.
  * [21] Geovani Nunes Grapiglia and Ya-xiang Yuan.  On the complexity of an augmented Lagrangian method for nonconvex optimization.  IMA Journal of Numerical Analysis, 41(2):1546–1568, 2020.
  * [22] Lei Guo and Jane J. Ye.  Necessary optimality conditions and exact penalization for non-Lipschitz nonlinear programs.  Mathematical Programming, 168(1):571–598, 2018.
  * [23] Ben Hermans, Andreas Themelis, and Panagiotis Patrinos.  QPALM: a proximal augmented Lagrangian method for nonconvex quadratic programs.  Mathematical Programming Computation, 14(3):497–541, 2022.
  * [24] Magnus R. Hestenes.  Multiplier and gradient methods.  Journal of Optimization Theory and Applications, 4(5):303–320, 1969\.
  * [25] Tim Hoheisel.  Mathematical Programs with Vanishing Constraints.  PhD thesis, University of Würzburg, Würzburg, Germany, 2009.
  * [26] Christian Kanzow, Daniel Steck, and Daniel Wachsmuth.  An augmented Lagrangian method for optimization problems in Banach spaces.  SIAM Journal on Control and Optimization, 56(1):272–291, 2018.
  * [27] Alexander Y. Kruger and Patrick Mehlitz.  Optimality conditions, approximate stationarity, and applications — a story beyond Lipschitzness.  ESAIM: Control, Optimisation and Calculus of Variations, 28:42, 2022\.
  * [28] Puya Latafat and Panagiotis Patrinos.  Asymmetric forward-backward-adjoint splitting for solving monotone inclusions involving three operators.  Computational Optimization and Applications, 68(1):57–93, 2017\.
  * [29] Patrick Mehlitz.  Asymptotic stationarity and regularity for nonsmooth optimization problems.  Journal of Nonsmooth Analysis and Optimization, 1:6575, 2020.
  * [30] Boris S. Mordukhovich.  Variational Analysis and Generalized Differentiation, Part I: Basic Theory, Part II: Applications.  Springer, Berlin, 2006.
  * [31] Jean Jacques Moreau.  Proximité et dualité dans un espace hilbertien.  Bulletin de la Société Mathématique de France, 93:273–299, 1965.
  * [32] Neal Parikh and Stephen Boyd.  Proximal algorithms.  Foundations and Trends in Optimization, 1(3):127–239, 2014.
  * [33] Michael J. D. Powell.  A method for nonlinear constraints in minimization problems, pages 283–298.  Academic Press, 1969.
  * [34] R. Tyrrell Rockafellar.  A dual approach to solving nonlinear programming problems by unconstrained optimization.  Mathematical Programming, 5(1):354–373, 1973.
  * [35] R. Tyrrell Rockafellar.  Augmented Lagrangians and applications of the proximal point algorithm in convex programming.  Mathematics of operations research, 1(2):97–116, 1976.
  * [36] R. Tyrrell Rockafellar.  Extended nonlinear programming.  In Gianni Di Pillo and Franco Giannessi, editors, Nonlinear Optimization and Related Topics, volume 36 of Applied Optimization, pages 381–399. Kluwer, 2000.
  * [37] R. Tyrrell Rockafellar.  Augmented Lagrangians and hidden convexity in sufficient conditions for local optimality.  Mathematical Programming, 198:159–194, 2023.
  * [38] R. Tyrrell Rockafellar.  Convergence of augmented Lagrangian methods in extensions beyond nonlinear programming.  Mathematical Programming, 199:375–420, 2023.
  * [39] R. Tyrrell Rockafellar and Roger J. B. Wets.  Variational Analysis, volume 317.  Springer, 1998.
  * [40] Andreas Themelis, Lorenzo Stella, and Panagiotis Patrinos.  Forward-backward envelope for the sum of two nonconvex functions: Further properties and nonmonotone linesearch algorithms.  SIAM Journal on Optimization, 28(3):2274–2303, 2018.