∎

11institutetext: Lahcen El Bourkhissi 22institutetext: Department of Automatic
Control and Systems Engineering, National University for Science and
Technology Politehnica Bucharest, Bucharest, 060042, Romania.
22email<EMAIL_ADDRESS>33institutetext: Ion Necoara 44institutetext:
Department of Automatic Control and Systems Engineering, National University
for Science and Technology Politehnica Bucharest, 060042 Bucharest. Gheorghe
Mihoc-Caius Iacob Institute of Mathematical Statistics and Applied Mathematics
of the Romanian Academy, 050711 Bucharest, Romania.
44email<EMAIL_ADDRESS>
# Complexity of linearized quadratic penalty for optimization with nonlinear
equality constraints ††thanks: The research leading to these results has
received funding from: the European Union’s Horizon 2020 research and
innovation programme under the Marie Sklodowska-Curie grant agreement no.
953348; UEFISCDI, Romania, PN-III-P4-PCE-2021-0720, under project L2O-MOC, no.
70/2022.

Lahcen El Bourkhissi  Ion Necoara

(Received: date / Accepted: date)

###### Abstract

In this paper we consider a nonconvex optimization problem with nonlinear
equality constraints. We assume that both, the objective function and the
functional constraints, are locally smooth. For solving this problem, we
propose a linearized quadratic penalty method, i.e., we linearize the
objective function and the functional constraints in the penalty formulation
at the current iterate and add a quadratic regularization, thus yielding a
subproblem that is easy to solve, and whose solution is the next iterate.
Under a dynamic regularization parameter choice, we derive convergence
guarantees for the iterates of our method to an $\epsilon$ first-order optimal
solution in $\mathcal{O}(1/{\epsilon^{3}})$ outer iterations. Finally, we show
that when the problem data satisfy Kurdyka-Lojasiewicz property, e.g., are
semialgebraic, the whole sequence generated by our algorithm converges and we
derive convergence rates. We validate the theory and the performance of the
proposed algorithm by numerically comparing it with the existing methods from
the literature.

###### Keywords:

Nonconvex optimization nonlinear functional constraints linearized quadratic
penalty convergence analysis.

###### MSC:

68Q25 90C06 90C30.

##  1 Introduction

In many fields, including machine learning, matrix optimization, statistics,
control and signal processing, a variety of applications can be reformulated
as nonconvex optimization problems with nonlinear functional equality
constraints, see e.g., Fes:20 ; LukSab:19 ; HonHaj:17 ; Roy:19 . In this
paper, we propose a solution to these optimization problems using a penalty
approach.

Related work: The penalty framework holds a key role in the realm of
theoretical and numerical optimization, with its history tracing back at least
to Cou:43 . Numerous works have examined penalty methods for a wide range of
problems, see e.g., PolTre:73 ; Ber:76 ; KonMel:18 ; IzmSol:23 ; Pol:87 ;
NocWri:06 ; Fle:87 ; FiaMcC:68 ; CarGou:11 ; LinMa:22 . For example, in
IzmSol:23 a general class of penalty functions of the form
$\|\cdot\|_{q}^{q}$, with $q>0$, is considered and bounds on the distance of
the solution of the penalty subproblem to the solution of the original problem
are derived in terms of the penalty parameter $\rho$. In particular, under
strict Mangasarian–Fromovitz constraint qualification and second-order
sufficiency, a bound of the form $\mathcal{O}({\rho^{1-q}})$ is derived and it
becomes zero for $q\in(0,1]$, provided that $\rho$ is sufficiently large. In
CarGou:11 , using a Lipschitz penalty function (e.g., Euclidean norm
$\|\cdot\|$), a new algorithm is proposed based on the linearization of both,
the objective function and the functional constraints, in the penalty
subproblem and, additionally, incorporating a dynamic quadratic
regularization. It is shown that this algorithm requires at most
$\mathcal{O}({\epsilon^{-2}})$ functions and first derivatives evaluations to
obtain a reduction in the size of a first-order criticality measure below some
accuracy $\epsilon$. In contrast to CarGou:11 , in KonMel:18 a quadratic
penalty method is used when dealing with nonconvex composite programs subject
to linear constraints. The study demonstrates that the algorithm, which
involves sequentially minimizing the quadratic penalty function, converges to
an $\epsilon$ stationary point when using the accelerated composite gradient
(ACG) method for solving the subproblem. Notably, this convergence is achieved
within $\mathcal{O}(\epsilon^{-3})$ ACG iterations.

Furthermore, in BirGar:16 a two-phase scaled penalty type algorithm is
presented for general nonlinear programs. The first phase finds an approximate
feasible point, while the second phase involves generating a $p$ Taylor
approximation of a quadratic penalty function and the next iterate is the
solution that satisfies a decrease in this penalty function. This algorithm
requires $\mathcal{O}(\epsilon^{1-2(p+1)/p})$ number of evaluations of the
problem’s functions and derivatives up to order $p$ to achieve an $\epsilon$
approximate first-order critical point for the original problem. It is
noteworthy that when $p=1$, no hessian information is used and the rate is
$\mathcal{O}(\epsilon^{-3})$. Finally, in LinMa:22 an inexact proximal-point
penalty method is proposed, where at each iteration one solves inexactly the
quadratic penalty function plus a proximal term. It is proved that when the
objective function and constraints are nonconvex, this method converges to an
$\epsilon$ first-order solution in $\mathcal{O}({\epsilon^{-3}})$ iterations,
while for nonconvex objective and convex functional constraints this rate can
be improved to $\mathcal{O}({\epsilon^{-2.5}})$.

Another promising approach when dealing with nonlinear constraints is the
augmented Lagrangian (AL), which can be seen as a generalization of the
quadratic penalty method. This approach has been extensively explored in the
literature for nonlinear programs, as demonstrated in various works BerTsi:03
; BirMar:14 ; BirMar:20 ; BoyPar:11 ; XieWri:21 ; CohHal:21 ; Yas:22 ;
SahEft:19 . In particular, in BirMar:14 an in-depth analysis of practical AL
methods is given and connections between quadratic penalty and AL methods are
established. Moreover, in the recent study XieWri:21 , a proximal augmented
Lagrangian method is examined for nonconvex, yet smooth, optimization
problems. It is proved that when an approximate first- (second-) order
solution of the subproblem is found inexactly, with an error approaching
asymptotically zero, then an $\epsilon$ first- (second-) order solution of the
original problem is obtained within $\mathcal{O}(\epsilon^{\eta-2})$ outer
iterations, where $\eta\in[0,2]$ is an user-defined parameter.

Finally, Sequential Convex Programming (SCP) or inexact restoration
techniques, offer an alternative to the penalty/augmented Lagrangian-based
approaches MesBau:21 ; TraDie:10 ; BueMar:20 . More, specifically, SCP solves
a sequence of convex approximations of the original problem by linearizing the
nonconvex parts of the objective and functional constraints while preserving
structures that can be exploited by convex optimization techniques, see e.g.,
GraBoy:14 ; NecKva:15 . However, to the best of our knowledge, SCP methods
converge under mild assumptions only locally MesBau:21 ; TraDie:10 . On the
other hand, inexact restoration methods are based on a two phase procedure,
first phase focuses on seeking a feasible point, while the second phase aims
to generate a point ensuring a decrease in the objective function. In
BueMar:20 it has been proven that such strategy yields an approximate KKT
point in $\mathcal{O}\left({\epsilon^{-2}}\right)$ iterations. Nevertheless,
the two-phase nature of the algorithm leads to increased computational effort
in each iteration.

Contributions: Our approach, the linearized quadratic penalty method (LQP),
addresses several key limitations present in the previous works. Notably, in
CarGou:11 , the subproblem is nondifferentiable due to the adoption of a
nonsmooth penalty function, and in KonMel:18 , only linear constraints can be
handled. Furthermore, the Proximal AL method in XieWri:21 incurs high
computational costs for solving nonconvex subproblems, and the SCP schemes
MesBau:21 ; TraDie:10 provide only local convergence guarantees. Specifically,
our main contributions are:

  * (i)

We linearize at each iteration the cost function and the nonlinear functional
equality constraints in the quadratic penalty formulation and we consider a
dynamic regularization term. This leads to a new algorithm, called linearized
quadratic penalty (LQP) method, which at each iteration requires solving a
simple unconstrained quadratic subproblem, that is easy to solve.  Note that
algorithms based on removing some parts of the second derivatives have been
proposed previously in the optimization literature, e.g., Gauss-Newton,
quadratic penalty and augmented Lagrangian methods have employed similar
strategies NocWri:06 ; KreMar:00 ; CarGou:11 ; BirGar:16 . However, it is
essential to highlight that these methods typically address specific problem
classes, such as those with a null objective function, or generate the next
iterate as an approximate stationary point of the quadratic penalty/augmented
Lagrangian function, which is computationally expensive. One contribution of
this paper consists of application of such a technique to a broader class of
problems, encompassing scenarios beyond those covered by previous papers,
while generating the next iterate with the only requirement that the penalty
function decreases, which is computationally cheap. Additionally, we derive a
complete convergence analysis. To the best of our knowledge, our paper stands
as the first to present a complete complexity analysis for a quadratic
penalty-type approach, specifically when second derivatives of the objective
function and functional constraints are neglected, while still retaining
minimal information about the Hessian of the penalty term through the
computation of only first derivatives.

  * (ii)

We provide rigorous proofs of global asymptotic convergence, proving that the
iterates converge to a critical point of the penalty function. Additionally,
our method guarantees convergence to an $\epsilon$ first-order solution of the
original problem in at most $\mathcal{O}(1/{\epsilon^{3}})$ iterations.

  * (iii)

Under the Kurdyka-Lojasiewicz (KL) property, we demonstrate the convergence of
the entire sequence generated by our algorithm and derive improved convergence
rates that depend on the KL parameter.

  * (iv)

In comparison with CarGou:11 , where a Lipschitz penalty function is adopted
and their rate may seem better than ours, our approach considerably simplifies
the subproblem as our update reduces to solving a linear system. In addition,
our method uses an adaptive and consequently potentially faster strategy to
update the dynamic regularization term. Still comparing the complexity of the
subproblems, the algorithms in XieWri:21 ; LinMa:22 are difficult to implement
in practice due to their highly nonconvexity caused by the presence of
nonlinear constraints in the subproblem from each iteration. Moreover, unlike
KonMel:18 , our LQP method can handle nonlinear equality constraints. Finally,
our LQP algorithm benefits from global convergence guarantees, which is an
advantage over the SCP methods, for which only local convergence has been
established.

These contributions collectively make our LQP method a powerful and versatile
tool for efficiently solving smooth nonconvex problems with nonlinear equality
constraints, surpassing some limitations of existing techniques. The
performance of LQP is also validated through numerical simulations.

The paper is structured as follows. In Section 2, we introduce our problem of
interest and some notions necessary for our analysis. In Section 3, we present
our algorithm, followed in Section 4 by its convergence analysis. Finally, in
Section 5, we compare numerically our method with existing algorithms.

##  2 Problem formulation and preliminaries

In this paper, we consider the following nonlinear optimization problem:

$\displaystyle\underset{x\in\mathbb{R}^{n}}{\min}$ $\displaystyle f(x)$ (1)
s.t. $\displaystyle\hskip 1.9919ptF(x)=0,$

where $f:\mathbb{R}^{n}\to{\mathbb{R}}$ and
$F(x)\triangleq{(f_{1}(x),...,f_{m}(x))}^{T}$, with
$f_{i}:\mathbb{R}^{n}\to{\mathbb{R}}$ for all $i=1:m$. We assume the functions
$f,f_{i}$ are continuously differentiable for all $i=1:m$, $f$ is nonconvex
and $F$ is nonlinear. Moreover, we assume that the problem is well-posed,
i.e., the feasible set is nonempty and the optimal value is finite.  Problem
(LABEL:eq1) encompasses a broad spectrum of optimization problems, including
those with inequality constraints. In fact, any nonlinear program with
inequality constraints of the form:

$\displaystyle\underset{y\in\mathbb{R}^{p}}{\min}$ $\displaystyle g(y)$ (2)
s.t. $\displaystyle G(y)=0,\;\;H(y)\leq 0,$

can be transformed into problem (LABEL:eq1) by introducing slack variables as
follows:

$H(y)+\text{diag}(s)s=0.$

By defining $n=m+p$ and $x=(y,s)^{T}\in\mathbb{R}^{n}$, and putting
$f(x)=g(y)$, and

$F(x)=\begin{bmatrix}G(y)\\\ H(y)+\text{diag}(s)s\end{bmatrix},$

problem (LABEL:eq_gen) is equivalent to problem (LABEL:eq1). It is important
to note that transforming inequality constraints into equality constraints
through this method is not always the best approach, especially when the
inequality constraints are convex. However, we do not consider convex
constraints in this paper, deferring that to future work.  Before introducing
the main assumptions for our analysis, we would like to clarify some
notations. We use $\|.\|$ to denote the 2-norm of a vector or of a matrix,
respectively. For a differentiable function $f:\mathbb{R}^{n}\to\mathbb{R}$,
we denote by $\nabla f(x)\in\mathbb{R}^{n}$ its gradient at a point $x$. For a
differentiable vector function $F:\mathbb{R}^{n}\to\mathbb{R}^{m}$, we denote
its Jacobian at a given point $x$ by ${J}_{F}(x)\in\mathbb{R}^{m\times n}$.
Let us now present the main assumptions considered for problem (LABEL:eq1):

###### Assumption 1

Assume that there exists $\rho_{0}\geq 0$ such that
$f(x)+\frac{\rho_{0}}{2}\|F(x)\|^{2}$ has a compact level set, i.e., for an
appropriate value of $\alpha\in\mathbb{R}$, the following set is nonempty and
compact:

$\mathcal{S}_{\alpha}^{0}\triangleq{\\{x:\,f(x)+\frac{\rho_{0}}{2}\|F(x)\|^{2}\leq\alpha\\}}.$

###### Assumption 2

Given a compact set $\mathcal{S}\subseteq\mathbb{R}^{n}$, there exist positive
constants $M_{F},\sigma,L_{f},L_{F}$ such that $f$ and $F$ satisfy the
following conditions:

  1. (i)

$\|\nabla f(x)-\nabla f(y)\|\leq L_{f}\|x-y\|\;\text{ for all
}x,y\in\mathcal{S}$.

  2. (ii)

$\|J_{F}(x)\|\leq M_{F},\hskip 14.22636pt\|J_{F}(x)-J_{F}(y)\|\leq
L_{F}\|x-y\|\;\text{ for all }x,y\in\mathcal{S}$.

  3. (iii)

Problem (LABEL:eq1) satisfies Linear Independance Constraint Qualification
(LICQ) condition  for all  $x\in\mathcal{S}$.

###### Assumption 3

There exists finite $\bar{\alpha}$ such that $f(x)\leq\bar{\alpha}$ for all
$x\in\\{x:\,\|F(x)\|\leq 1\\}$.

Note that these assumptions are standard in the nonconvex optimization
literature, in particular in penalty-type methods, see e.g., CarGou:11 ;
XieWri:21 ; CohHal:21 . In fact, these assumptions, except LICQ (although it
is also standard in nonconvex optimization), are not too restrictive as they
need to hold only locally. Any twice-differentiable function is Lipschitz on a
compact set, and its gradient is also Lipschitz on any compact set. Indeed,
large classes of functions satisfy these assumptions, as discussed below.

###### Remark 1

Assumption 1 holds e.g., when $f(x)+\rho\|F(x)\|^{2}$ is coercive for some
$\rho\geq 0$, or when $f(x)$ is strongly convex or $f(x)$ is bounded bellow
and $F(x)=x^{T}x-1$, as is the case in dictionary learning applications. It
also holds when $f(x)=\frac{1}{2}x^{T}Qx-p^{T}x,F(x)=Ax-b$ and $Q$ is a
positive definite matrix on $\text{null}(A):=\\{x:\,Ax=0\\}$.

###### Remark 2

Assumption 2 allows general classes of problems. In particular, condition (i)
holds if $f(x)$ is differentiable and $\nabla f(x)$ is locally Lipschitz
continuous on a neighborhood of $\mathcal{S}$. Conditions (ii) hold when
$F(x)$ is differentiable on a neighborhood of $\mathcal{S}$ and $J_{F}(x)$ is
locally Lipschitz continuous on $\mathcal{S}$. Finally, the LICQ assumption
guarantees the existence of dual multipliers and is commonly used in nonconvex
optimization, see e.g., NocWri:06 ; XieWri:21 .

###### Remark 3

For Assumption 3 to hold, it is sufficient the set $\\{x:\,\|F(x)\|\leq 1\\}$
to be compact. In fact, we do not need this assumption if we can choose the
starting point $x_{0}$ such that $F(x_{0})=0$, i.e., the initial point is
feasible.

The following lemma is an immediate consequence of Assumption 1.

###### Lemma 1

If Assumption 1 holds, then $f(x)+\frac{\rho_{0}}{2}\|F(x)\|^{2}$ is lower
bounded:

$\bar{L}\triangleq{\inf_{x\in\mathbb{R}^{n}}\\{f(x)+\frac{\rho_{0}}{2}\|F(x)\|^{2}\\}}>-\infty.$
(3)

We are interested in (approximate) first-order (KKT) solution of (LABEL:eq1).
Hence, let us introduce the following definitions:

###### Definition 1

[First-order solution and $\epsilon$ first-order solution of (LABEL:eq1)]. The
vector $x^{*}$ is a first-order (KKT) solution of (LABEL:eq1) if
$\exists\lambda^{*}\in\mathbb{R}^{m}$ such that:

$\nabla f(x^{*})+J_{F}(x^{*})^{T}\lambda^{*}=0\hskip 8.5359pt\text{and}\hskip
8.5359ptF(x^{*})=0.$

Moreover, $x^{*}$ is an $\epsilon$ first-order solution of (LABEL:eq1) if
$\exists\lambda^{*}\in\mathbb{R}^{m}$ such that:

$\|\nabla f(x^{*})+J_{F}(x^{*})^{T}\lambda^{*}\|\leq\epsilon\hskip
8.5359pt\text{and}\hskip 8.5359pt\|F(x^{*})\|\leq\epsilon.$

Finally, let us introduce the Kurdyka-Lojasiewicz (KL) property, which leads
to improvements in the convergence rates of our algorithm. In particular, its
incorporation will prove useful in improving the global performance of our
algorithm. Let $\Phi:\mathbb{R}^{d}\to\bar{\mathbb{R}}$ be a proper lower
semicontinuous function. For $-\infty<\tau_{1}<\tau_{2}\leq+\infty$, we define
$[\tau_{1}<\Phi<\tau_{2}]=\\{x\in\mathbb{R}^{d}:\tau_{1}<\Phi(x)<\tau_{2}\\}$.
Let $\tau\in(0,+\infty]$. We denote by $\Psi_{\tau}$ the set of all continuous
concave functions $\varphi:[0,\tau]\to[0,+\infty)$ such that $\varphi(0)=0$
and $\varphi$ is continuously differentiable on $(0,\tau)$, with
$\varphi^{\prime}(s)>0$ over $(0,\tau)$.

###### Definition 2

Let $\Phi:\mathbb{R}^{d}\to\bar{\mathbb{R}}$ be a proper lower semicontinuous
function that takes constant value on a set $\Omega$. We say that $\Phi$
satisfies the KL property on $\Omega$ if there exists $\epsilon>0,\tau>0$, and
$\varphi\in\Psi_{\tau}$ such that for every $x^{*}\in\Omega$ and every element
$x$ in the intersection $\\{x\in\mathbb{R}^{d}:\text{
dist}(x,\Omega)<\epsilon\\}\cap[\Psi(x^{*})<\Psi(x)<\Psi(x^{*})+\tau]$, we
have:

$\varphi^{\prime}\big{(}\Phi(x)-\Phi(x^{*})\big{)}\text{dist}\big{(}0,\partial\Phi(x)\big{)}\geq
1.$

This definition covers many classes of functions arising in practical
optimization problems. For example, if $f$ is a proper closed semialgebraic
function, then $f$ is a KL function with exponent $\nu\in[0,1)$, see AttBol:13
. The function $g(Ax)$, where $g$ is strongly convex on a compact set and
twice differentiable, and $A\in\mathbb{R}^{m\times n}$, is a KL function.
Convex piecewise linear/quadratic functions such as
$\|x\|_{1},\|x\|_{0},\gamma\sum_{i=1}^{k}|x_{[i]}|$, where $|x_{[i]}|$ is the
$i$th largest entry in $x,\;k\leq n$ and $\gamma\in(0,1]$; the indicator
function $\delta_{\Delta}(x)$, where
$\Delta=\\{x\in\mathbb{R}^{n}:e^{T}x=1,x\geq 0\\}$; least-squares problems
with the Smoothly Clipped Absolute Deviation (SCAD) Fan:97 ; Minimax Concave
Penalty (MCP) regularized functions Zha:10 are all KL functions.

##  3 A linearized quadratic penalty method

In this section, we propose a new algorithm for solving nonconvex problem
(LABEL:eq1) using the quadratic penalty framework. Let us first introduce few
notations. The penalty function associated with the problem (LABEL:eq1) is

$\mathcal{P}_{\rho}(x)=f(x)+\frac{\rho}{2}{\|F(x)\|^{2}}.$ (4)

The gradient of $\mathcal{P}_{\rho}$ is: $\nabla\mathcal{P}_{\rho}(x)={\nabla
f(x)+J_{F}(x)}^{T}\left(\rho F(x)\right)$. In the rest of this paper, for the
sake of clarity, we provide the proofs of all the lemmas in Appendix. The
following lemma shows that $\mathcal{P}_{\rho}$ is locally smooth, (i.e., it
has Lipschitz continuous gradient locally).

###### Lemma 2

[Smoothness of $\mathcal{P}_{\rho}$] If Assumption 2 holds on a given compact
set $\mathcal{S}$, then there exists $L_{\mathcal{P}_{\rho}}>0$ such that:

$\|\nabla\mathcal{P}_{\rho}(x)-\nabla\mathcal{P}_{\rho}(x^{\prime})\|\leq{L_{\mathcal{P}_{\rho}}}\left\|x-x^{\prime}\right\|\hskip
14.22636pt\forall x,x^{\prime}\in\mathcal{S},$

where
$L_{\mathcal{P}_{\rho}}=L_{f}+\rho(L_{F}\sup_{x\in\mathcal{S}}\|F(x)\|+M_{F}^{2})$.

###### Proof

See Appendix.

Further, let us denote the following function derived from linearization of
the objective function and the functional constraints, at a given point
$\bar{x}$, in the penalty function:

$\displaystyle\bar{\mathcal{P}}_{\rho}(x;\bar{x})=f(\bar{x})+\langle\nabla
f(\bar{x}),x-\bar{x}\rangle+\frac{\rho}{2}{\|F(\bar{x})+J_{F}(\bar{x})(x-\bar{x})\|^{2}}.$

To solve the optimization problem (LABEL:eq1) we propose the following
Linearized quadratic penalty (LQP) algorithm, i.e., we linearize the objective
function and the functional constraints in the penalty function at the current
iterate and add an adaptive quadratic regularization.

Algorithm 1 Linearized quadratic penalty (LQP) method

1:$\textbf{Initialization: }x_{0},\underaccent{\bar}{\beta}>0\text{ and
}\rho>0$.

2:$k\leftarrow 0$

3:while stopping criterion is not satisfied  do

4: $\text{generate a proximal parameter
}\beta_{k+1}\geq\underaccent{\bar}{\beta}$ such that

5:
$x_{k+1}\leftarrow\operatorname*{arg\,min}_{x\in\mathbb{R}^{n}}{\bar{\mathcal{P}}_{\rho}(x;x_{k})+\frac{\beta_{k+1}}{2}{\|x-x_{k}\|}^{2}}$
satisfies the following descent:
$\mathcal{P}_{\rho}(x_{k+1})\leq\mathcal{P}_{\rho}(x_{k})-\frac{\beta_{k+1}}{2}\|x_{k+1}-x_{k}\|^{2}.$
(5)

6: $k\leftarrow k+1$

7:end while

To the best of our knowledge LQP algorithm is new and its convergence
behaviour has not been analyzed before in the literature. Note that the
objective function in the subproblem of step 5 of Algorithm 1 is always
unconstrained, strongly convex and quadratic. Therefore, finding a solution of
the subproblem in step 5 is equivalent to solving a linear system of
equalities for which efficient solution methods exist, i.e.:

$x_{k+1}=x_{k}-(\rho
J_{F}(x_{k})^{T}J_{F}(x_{k})+\beta_{k+1}I_{n})^{-1}(\nabla f(x_{k})+\rho
J_{F}(x_{k})^{T}F(x_{k})).$

###### Remark 4

Note that, from the previous relation, when the objective function is null and
$\rho=1$, the update of $x_{k+1}$ takes the form:

$x_{k+1}=x_{k}-(J_{F}(x_{k})^{T}J_{F}(x_{k})+\beta_{k+1}I_{n})^{-1}(J_{F}(x_{k})^{T}F(x_{k})),$

which corresponds to the Levenberg-Marquardt method NocWri:06 . Moreover, when
$\beta_{k+1}=0$, we obtain the following update:

$x_{k+1}=x_{k}-(J_{F}(x_{k})^{T}J_{F}(x_{k}))^{-1}J_{F}(x_{k})^{T}F(x_{k}),$

which yields the Gauss-Newton method NocWri:06 .

Let us show that we can always choose $\beta_{k+1}$ guaranteeing the descent
property (5). Indeed, according to Lemma 2 the penalty function
$\mathcal{P}_{\rho}$ is smooth, and thus one can always choose
$\underaccent{\bar}{\beta}\leq\beta_{k+1}\leq L_{\mathcal{P}_{\rho}}$ such
that:

$\displaystyle\mathcal{P}_{\rho}(x_{k+1})-\mathcal{P}_{\rho}(x_{k})-\langle\nabla\mathcal{P}_{\rho}(x_{k}),x_{k+1}-x_{k}\rangle\leq\frac{\beta_{k+1}}{2}\|x_{k+1}-x_{k}\|^{2}.$
(6)

Next lemma shows that requiring the newly generated iterate $x_{k+1}$ to
satisfy the descent (5) is achievable.

###### Lemma 3

[Existence of $\beta_{k+1}$] If the sequence, $\\{x_{k}\\}_{k\geq 0}$
generated by Algorithm 1 is in some compact set $\mathcal{S}$ on which
Assumption 2 holds and we generate $\beta_{k+1}$ as in (6), then the descent
property (5) is valid.

###### Proof

See Appendix.

Lemma 3 proves that our algorithm is implementable, since choosing
$\beta_{k+1}$ as in (6) already guarantees the descent in (5). In fact, a
value of $\beta_{k+1}$ generated by (5) might be much smaller than the one
generated by (6). In practice, the regularization parameter $\beta_{k}$ can be
determined using a backtracking scheme see Algorithm 2.

Algorithm 2 Backtracking procedure

1:$\textbf{Initialization: }\mu>1\text{ and
}\beta^{0}_{k+1}=\max\\{\beta_{k}/\mu,\underline{\beta}\\}$.

2:$i\leftarrow 0$

3:while (5) is not satisfied for $\beta_{k+1}=\beta^{i}_{k+1}$ do

4: $\beta^{i+1}_{k+1}\leftarrow\mu\beta^{i}_{k+1}$

5:
$x^{i+1}_{k+1}\leftarrow\operatorname*{arg\,min}_{x\in\mathbb{R}^{n}}{\bar{\mathcal{P}}_{\rho}(x;x_{k})+\frac{\beta^{i+1}_{k+1}}{2}{\|x-x_{k}\|}^{2}}$

6: $i\leftarrow i+1$

7:end while

8:$\beta_{k+1}\leftarrow\beta^{i+1}_{k+1},\quad x_{k+1}\leftarrow
x^{i+1}_{k+1}$

Specifically, at iteration $k$ we start with an initial value
$\beta^{0}_{k+1}=\max\\{\beta_{k}/\mu,\underline{\beta}\\}$ for some $\mu>1$
and generate the new iterate $x_{k+1}$ by solving the subproblem in Step 5 of
Algorithm 1. Next, we test whether (5) is satisfied or not. If it is
satisfied, we proceed to Step 6 of Algorithm 1. However, if (5) is not
satisfied, we geometrically increase $\beta^{0}_{k+1}$ by a geometric factor
(with inner values of the form $\beta^{i}_{k+1}:=\mu^{i}\beta^{0}_{k+1}$, with
$i\geq 0$), until (5) is satisfied. This process finishes in a finite number
of steps, as proved in Lemma 3. The value of $\beta^{i}_{k+1}$ for which (5)
is satisfied, denoted $\beta^{i_{0}}_{k+1}$, is then chosen as the value of
$\beta_{k+1}$.

###### Remark 5

It is clear that when the sequence $\\{\beta_{k}\\}_{k\geq 1}$ is generated
such that it satisfies (5) or (6), then $\\{\beta_{k}\\}_{k\geq 1}$ is upper
bounded by the Lipschitz constant $L_{\mathcal{P}_{\rho}}$, i.e.,

$\beta_{k}\leq L_{\mathcal{P}_{\rho}}\hskip 14.22636pt\forall k\geq 1.$

There are good reasons for adopting a dynamic regularization $\beta_{k}$. One
of them is that, in general, $L_{\mathcal{P}_{\rho}}$ is not known a priori
since the parameters of the problem are usually unknown. Another reason is
that, since $\beta_{k+1}$ is generated based on the inequalities (5) or (6),
which must be tested only for the point $x_{k+1}$, this might results in
$\beta_{k}$’s that are much smaller than the Lipschitz constant
$L_{\mathcal{P}_{\rho}}$. In the sequel we also denote:

$\Delta x_{k}=x_{k}-x_{k-1}\hskip 8.5359pt\forall k\geq 1.$

##  4 Convergence analysis

This section is dedicated to the convergence analysis of Algorithm 1. This
analysis holds significance because, to the best of our knowledge, the
proposed method differs from the conventional quadratic penalty-based methods.
Specifically, the update step for regularized penalty methods with quadratic
regularization is computed, given the current point $\bar{x}_{k}$, as follows
(see e.g., KonMel:18 ; LinMa:22 ):

$\bar{x}_{k+1}\;\text{approximate stationary point
of}\;\left\\{f(x)+\frac{\rho}{2}\|F(x)\|^{2}+\frac{\beta_{k+1}}{2}\|x-\bar{x}_{k}\|^{2}\right\\}.$

In contrast, our update step, given the current point $x_{k}$, takes the form:

$x_{k+1}\leftarrow\arg\min_{x}\left\\{l_{f}(x;x_{k})+\frac{\rho}{2}\|l_{F}(x;x_{k})\|^{2}+\frac{\beta_{k+1}}{2}\|x-x_{k}\|^{2}\right\\},$

where $l_{f}(x;x_{k})=f(x_{k})+\langle\nabla f(x_{k}),x-x_{k}\rangle$ and
$l_{F}(x;x_{k})=F(x_{k})+J_{F}(x_{k})(x-x_{k})$. Note that it is very
difficult to compute $\bar{x}_{k+1}$, as the corresponding subproblem is
highly nonconvex, while computing $x_{k+1}$ is equivalent to solving a linear
system. Moreover, our algorithm utilizes only first-order derivatives to
construct a surrogate model at the current point $x_{k}$, retaining
information about the Hessian of the penalty term through $J_{F}^{T}J_{F}$.
Clearly, the algorithms are distinct; in fact, for the first update, an
approximate stationary point of the penalty function must be computed
KonMel:18 ; LinMa:22 , which requires an additional inner iterative process
for solving the subproblem, while for the second one there is no need for
inner iterations since $x_{k+1}$ can be computed explicitly and $\beta_{k+1}$
can be chosen such that the decrease (5) holds. Consequently, the convergence
results of the proposed method cannot be inferred from those of the standard
quadratic penalty methods. Therefore, a new convergence analysis is necessary
for the proposed method. In this section, we derive the asymptotic convergence
of LQP algorithm (Algorithm 1) and prove its efficiency to obtain an
$\epsilon$ first-order solution for problem (LABEL:eq1). Moreover, we provide
better convergence rates under the KL condition. A notable distinction in our
analysis, as compared to other penalty methods from literature, e.g.,
CarGou:11 ; KonMel:18 , lies in the way we prove convergence. While penalty
methods often rely on criticality measures directly linked to the choice of
regularization parameter, which can necessitate numerous inner iterations, we
take inspiration from augmented Lagrangian methods to prove convergence, see
e.g., XieWri:21 . More specifically, we use a Lyapunov analysis to prove the
convergence of our algorithm. Bellow, we are using the quadratic penalty
function, $\mathcal{P}_{\rho}$, as a Lyapunov function. The evaluation of the
Lyapunov function along the iterates of LQP algorithm is denoted by:

$P_{k}=\mathcal{P}_{\rho}(x_{k})\hskip 8.5359pt\forall k\geq 0.$ (7)

It is clear from Lemma 3 that $\\{P_{k}\\}_{k\geq 0}$ is decreasing and later
we will prove that it is bounded from bellow. In the sequel, we assume that
$x_{0}$ is chosen such that:

$\|F(x_{0})\|^{2}\leq\min\left\\{1,\frac{2c_{0}}{\rho}\right\\}\hskip
19.91684pt\text{ for some }c_{0}>0,$ (8)

and that $f(x_{0})\leq\bar{\alpha}$. Let us define:

$\alpha\triangleq\bar{\alpha}+c_{0}.$ (9)

Moreover, let us choose:

$\displaystyle\rho\geq\max\left\\{\\!{1},3\rho_{0}\right\\}.$ (10)

Using the definition of $\mathcal{P}_{\rho}$, we have:

$\displaystyle\mathcal{P}_{\rho}(x_{0})$
$\displaystyle=f(x_{0})+\frac{\rho}{2}\|F(x_{0})\|^{2}{\overset{{\eqref{eq4}}}{\leq}}\bar{\alpha}+c_{0}.$
(11)

It then follows, after some re-arrangements, that:

$\displaystyle\bar{\alpha}+c_{0}-\bar{L}$ $\displaystyle\geq
f(x_{0})+\frac{\rho}{2}\|F(x_{0})\|^{2}-\bar{L}$ (12)
$\displaystyle{\overset{{(\rho\geq
3\rho_{0})}}{\geq}}f(x_{0})+\frac{\rho_{0}}{2}\|F(x_{0})\|^{2}-\bar{L}+\frac{\rho}{3}\|F(x_{0})\|^{2}{\overset{{\eqref{lem1}}}{\geq}}\frac{\rho}{3}\|F(x_{0})\|^{2}\geq
0.$

The following lemma shows that if the sequence $\\{x_{k}\\}_{k\geq 0}$
generated by Algorithm 1 is bounded, then the Lyapunov sequence
$\\{P_{k}\\}_{k\geq 0}$ is also bounded.

###### Lemma 4

Consider Algorithm 1 and let $\\{P_{k}\\}_{k\geq 0}$ as defined in (7). If the
sequence $\\{x_{k}\\}_{k\geq 0}$ generated by Algorithm 1 is in some compact
set $\mathcal{S}$ on which Assumptions 1, 2 and 3 hold and, moreover, $\rho$
is chosen as in (10) and $x_{0}$ is chosen as in (8), then we have the
following:

$\displaystyle\bar{L}\leq P_{k}\leq\alpha\hskip 14.22636pt\forall k\geq 0,$
(13a)

where $\bar{L}$ is defined in (3) and $\alpha$ defined in (9) for any fixed
constant $c_{0}$.

###### Proof

See Appendix.

Let us now bound the gradient of the penalty function.

###### Lemma 5

[Boundedness of $\nabla\mathcal{P}_{\rho}$] If the sequence
$\\{x_{k}\\}_{k\geq 0}$ generated by Algorithm 1 is in some compact set
$\mathcal{S}$ on which Assumption 2 holds, then there exists
$\Gamma_{k+1}:=L_{f}+\rho(L_{F}\Delta_{F}+4M_{F}^{2})+\beta_{k+1}$ for all
$k\geq 0$, where
$0\leq\Delta_{F}:=\sup_{x\in\mathcal{S}}\\{\|F(x)\|\\}<\infty$ such that:

$\|\nabla\mathcal{P}_{\rho}(x_{k+1})\|\leq\Gamma_{k+1}\|x_{k+1}-x_{k}\|.$

###### Proof

See Appendix.

Since $\\{\beta_{k}\\}_{k\geq 1}$ is bounded (see Remark 5), it follows that
$\\{\Gamma_{k}\\}_{k\geq 1}$ is finite. Hence, in the sequel we define this
bound as follows:

$\Gamma_{\text{max}}:=\sup_{k\geq 1}\\{\Gamma_{k}\\}.$ (14)

###  4.1 Global asymptotic convergence

In this section we prove global convergence for the iterates generated by
Algorithm 1 and also convergence rates to an $\epsilon$ first-order (KKT)
solution. Based on the previous lemmas, we are now ready to present the global
asymptotic convergence of the iterates of Algorithm 1.

###### Theorem 4.1

[Limit points] If the sequence, $\\{x_{k}\\}_{k\geq 0}$, generated by
Algorithm 1 is in some compact set $\mathcal{S}$ on which Assumptions 1, 2 and
3 hold and, moreover, $\rho$ is chosen as in (10) and $x_{0}$ is chosen as in
(8), then any limit point $x_{\rho}^{*}$ of the sequence $\\{x_{k}\\}_{k\geq
0}$ is a critical point of the penalty function $\mathcal{P}_{\rho}(\cdot)$
defined in (4), i.e., $\nabla\mathcal{P}_{\rho}(x_{\rho}^{*})=0$.
Equivalently:

$\nabla f(x_{\rho}^{*})+{J_{F}(x_{\rho}^{*})}^{T}\lambda^{*}=0,\hskip
8.5359pt\text{where}\hskip 8.5359pt\lambda^{*}=\rho F(x_{\rho}^{*}).$

Additionally, if the penalty parameter $\rho$ is chosen to be large than
$\epsilon^{-1}$ (i.e., $\rho\geq\epsilon^{-1}$), then we are also $\epsilon$
feasible, meaning that $\|F(x^{*}_{\rho})\|\leq\mathcal{O}(\epsilon)$.

###### Proof

Using Lemma 3 and the fact that $\beta_{k}\geq\underaccent{\bar}{\beta}>0$, we
have:

$\displaystyle\frac{\beta_{k+1}}{2}\|$ $\displaystyle\Delta x_{k+1}\|^{2}\leq
P_{k}-P_{k+1}\hskip 8.5359pt\forall k\geq 0.$

Let $k\geq 0$, by summing up the above inequality from $i=0$ to $i=k$, we
obtain:

$\displaystyle\sum_{i=0}^{k}{\frac{\beta_{i+1}}{2}\|\Delta x_{i+1}\|^{2}}$
$\displaystyle\leq P_{0}-P_{k+1}{\overset{{\text{Lemma
}\ref{bbound}}}{\leq}}\alpha-\bar{L}.$ (15)

Since (15) holds for any $k\geq 0$, we have:

$\sum_{i=1}^{\infty}{\frac{\beta_{i+1}}{2}\|\Delta x_{i+1}\|^{2}}<\infty.$

This, together with the fact that $\beta_{k}\geq\underaccent{\bar}{\beta}>0$,
yields that:

$\lim_{k\to\infty}{\|\Delta x_{k}\|}=0.$ (16)

Since the sequence $\\{x_{k}\\}_{k\geq 0}$ is bounded, then there exists a
convergent subsequence, let us say $\\{x_{k}\\}_{k\in\mathcal{K}}$, with the
limit $x_{\rho}^{*}$. From Lemma 5 and (14), we have:

$\|\nabla\mathcal{P}_{\rho}(x_{\rho}^{*})\|=\lim_{k\in\mathcal{K}}{\|\nabla\mathcal{P}_{\rho}(x_{k})\|}\leq\Gamma_{\text{max}}\lim_{k\in\mathcal{K}}\|\Delta
x_{k}\|=0.$

Therefore, $\nabla\mathcal{P}_{\rho}(x_{\rho}^{*})=0$. Hence,

$\nabla f(x_{\rho}^{*})+{J_{F}(x_{\rho}^{*})}^{T}\lambda^{*}=0,\hskip
8.5359pt\text{where}\hskip 8.5359pt{\lambda^{*}}=\rho F(x_{\rho}^{*}).$

Additionally, provided that $\rho\geq\epsilon^{-1}$, we get:

$\|F(x_{\rho}^{*})\|=\frac{\|\lambda^{*}\|}{\rho}\leq\|\lambda^{*}\|\epsilon.$

Bellow, we show that $\lambda^{*}$ is bounded. Indeed, let $x^{*}$ be a KKT
point close to $x_{\rho}^{*}$. From LICQ, it follows that there exists a
$y^{*}$ such that the pair $(x^{*},y^{*})$ satisfies the following KKT
conditions for problem (LABEL:eq1):

$\nabla f(x^{*})+{J_{F}(x^{*})}^{T}y^{*}=0\hskip 8.5359pt\text{and}\hskip
8.5359ptF(x^{*})=0.$

It then follows that:

$y^{*}=-{J_{F}(x^{*})}^{+}\nabla f(x^{*}).$

Moreover, considering:

$\nabla f(x_{\rho}^{*})+{J_{F}(x_{\rho}^{*})}^{T}\lambda^{*}=0,\hskip
8.5359pt\text{where}\hskip 8.5359pt\lambda^{*}=\rho F(x_{\rho}^{*}),$

we arrive at the following relation:

$\lambda^{*}=-{J_{F}(x_{\rho}^{*})}^{+}\nabla f(x_{\rho}^{*}).$

Subsequently, it follows that:

$\|\lambda^{*}-y^{*}\|=\|{J_{F}(x^{*})}^{+}\nabla
f(x^{*})-{J_{F}(x_{\rho}^{*})}^{+}\nabla f(x_{\rho}^{*})\|.$

Given the continuity of $J_{F}^{+}$ and $\nabla f$, along with the fact that
$x^{*}$ and $x_{\rho}^{*}$ belong to the compact set $\mathcal{S}$, we
conclude that there exists a constant $M\geq 0$ such that:

$\|\lambda^{*}-y^{*}\|\leq M.$

Hence:

$\|\lambda^{*}\|\leq M+\|y^{*}\|,$

and recall that from LICQ, $\|y^{*}\|$ is bounded. This proves the statement
of our theorem. ∎

For the remainder of this paper, we define:

$\bar{\beta}:=\sup_{k\geq 1}\\{\beta_{k}\\}\leq L_{\mathcal{P}_{\rho}}\hskip
8.5359pt\text{ and }\hskip 8.5359pt\underaccent{\bar}{\beta}:=\inf_{k\geq
1}\\{\beta_{k}\\}>0.$ (17)

Let us now present one of the main results of this paper, which derives the
complexity of Algorithm 1 to find an $\epsilon$ first-order solution of
problem (LABEL:eq1).

###### Theorem 4.2

[First-order complexity] Consider Algorithm 1 and let $\\{P_{k}\\}_{k\geq 0}$
be defined as in (7). If the sequence $\\{x_{k}\\}_{k\geq 0}$ generated by
Algorithm 1 is in some compact set $\mathcal{S}$ on which Assumptions 1, 2 and
3 hold and, moreover, $\rho$ is chosen as in (10) and $x_{0}$ is chosen as in
(8), then for any $\epsilon>0$, Algorithm 1 yields an $\epsilon$ stationary
point of (4) after $K=\mathcal{O}(\epsilon^{-2})$ iterations, i.e., there
exists $k^{*}\in\\{1,...,K\\}$ such that
$\|\nabla\mathcal{P}_{\rho}(x_{k^{*}})\|\leq\epsilon$. Additionally, if
$k^{*}$ is suffiently large such that
$\|x_{k^{*}+1}-x_{k^{*}}\|\leq\mathcal{O}(\epsilon)$ and
$\rho\geq\max\\{1,3\rho_{0},\mathcal{O}(\epsilon^{-1})\\}$, then $x_{k^{*}}$
is an $\epsilon$ first-order solution of (LABEL:eq1) after
$K=\mathcal{O}(\epsilon^{-3})$ iterations.

###### Proof

From Lemma 3 and Lemma 5, we can deduce the following inequality:

$P_{k+1}-P_{k}\leq-\frac{\beta_{k+1}}{2\Gamma^{2}_{k+1}}\|\nabla\mathcal{P}_{\rho}(x_{k+1})\|^{2}\hskip
8.5359pt\forall k\geq 0.$

By considering Lemma 2 and (6), it becomes clear that there exist constants
$a_{k+1}$ and $b_{k+1}$ such that $\beta_{k+1}=a_{k+1}\rho+b_{k+1}$.
Additionally, since $\underaccent{\bar}{\beta}\leq\beta_{k}\leq\bar{\beta}$,
we can express $\beta_{k+1}$ and $\Gamma_{k+1}$ as follows:

$\beta_{k+1}=\mathcal{O}(\rho),\hskip
14.22636pt\Gamma_{k+1}=\mathcal{O}(\rho).$

Consequently, the above inequality can be further simplified to:

$P_{k+1}-P_{k}\leq-\mathcal{O}(\rho^{-1})\|\nabla\mathcal{P}_{\rho}(x_{k+1})\|^{2}\hskip
8.5359pt\forall k\geq 0.$

Let $k\geq 0$. By summing up the above inequality from $i=0$ to $i=K$, we get:

$\mathcal{O}(\rho^{-1})\sum_{i=0}^{K}{\|\nabla\mathcal{P}_{\rho}(x_{i+1})\|^{2}}\leq\alpha-\bar{L}.$

This implies that there exists $k^{*}\in{0,...,K}$ such that:

$\|\nabla\mathcal{P}_{\rho}(x_{k^{*}+1})\|^{2}\leq\frac{\alpha-\bar{L}}{K+1}\mathcal{O}(\rho).$

Given that $\alpha-\bar{L}$ is a constant, we can express this as:

$\|\nabla\mathcal{P}_{\rho}(x_{k^{*}+1})\|\leq{\frac{1}{\sqrt{K+1}}}\mathcal{O}(\rho^{\frac{1}{2}}).$

It follows that for any $\epsilon>0$,
$\|\nabla\mathcal{P}_{\rho}(x_{k^{*}+1})\|\leq\epsilon$ provided that
$K=\mathcal{O}(\frac{\rho}{\epsilon^{2}})$ i.e.:

$\displaystyle\|\nabla
f(x_{k^{*}})+J_{F}(x_{k^{*}})^{T}\lambda_{k^{*}}\|\leq\epsilon,\text{ where
}\lambda_{k^{*}}=\rho F(x_{k^{*}}).$ (18)

Moreover, if $\|x_{k^{*}+1}-x_{k^{*}}\|\leq\mathcal{O}(\epsilon)$, then it
follows from the optimality condition for $x_{k^{*}+1}$:

$\nabla f(x_{k^{*}})+\rho
J_{F}(x_{k^{*}})^{T}l_{F}(x_{k^{*}+1};x_{k^{*}})+\beta_{k^{*}+1}\Delta
x_{k^{*}+1}=0,$

that:

$\rho J_{F}(x_{k^{*}})^{T}F(x_{k^{*}})=-\nabla f(x_{k^{*}})-\rho
J_{F}(x_{k^{*}})^{T}J_{F}(x_{k^{*}})\Delta x_{k^{*}+1}-\beta_{k^{*}+1}\Delta
x_{k^{*}+1}.$

This implies that:

$\displaystyle\|J_{F}(x_{k^{*}})^{T}F(x_{k^{*}})\|\leq\frac{\|\nabla
f(x_{k^{*}})\|}{\rho}+\left\|J_{F}(x_{k^{*}})^{T}J_{F}(x_{k^{*}})+\frac{\beta_{k^{*}+1}}{\rho}I\right\|\|\Delta
x_{k^{*}+1}\|.$

Since $\|x_{k+1}-x_{k}\|\to 0$ as $k\to\infty$ (see equation (16)) and
$\|\nabla f(x_{k^{*}})\|$ and
$\left\|J_{F}(x_{k^{*}})^{T}J_{F}(x_{k^{*}})+\frac{\beta_{k^{*}+1}}{\rho}I\right\|$
are bounded (due to continuity of $\nabla f$ and $J_{F}$ and boundedness of
the sequence $\\{x_{k}\\}_{k\geq 0}$), then considering
$\rho\geq\mathcal{O}(\epsilon^{-1})$ and
$\|x_{k^{*}+1}-x_{k^{*}}\|\leq\mathcal{O}(\epsilon)$ and that
$J_{F}(x_{k^{*}})$ has full row rank (see Assumption 2), we further obtain:

$\sigma_{\text{min}}(J_{F}(x_{k^{*}}))\|F(x_{k^{*}})\|\leq\|J_{F}(x_{k^{*}})^{T}F(x_{k^{*}})\|\leq\mathcal{O}(\epsilon).$

This implies that:

$\displaystyle\|F(x_{k^{*}+1})$
$\displaystyle\|\leq\|F(x_{k^{*}+1})-F(x_{k^{*}})\|+\|F(x_{k^{*}})\|$
$\displaystyle\leq\|F(x_{k^{*}})\|+M_{F}\|x_{k^{*}+1}-x_{k^{*}}\|\leq\mathcal{O}(\epsilon).$

Therefore, combining previous relation with (18), it follows that $x_{k^{*}}$
is an $\epsilon$ first-order solution of problem (LABEL:eq1) and this is
achieved after
$K=\mathcal{O}(\frac{\rho}{\epsilon^{2}})=\mathcal{O}(\frac{1}{\epsilon^{3}})$
iterations. This concludes our proof. ∎

From previous theorems, one can see that in addition to its straightforward
implementation, our algorithm also enjoys convincing global theoretical
convergence results, which gives it an advantage over approaches where only
local convergence can be guaranteed, such as SCP schemes MesBau:21 . One of
its key advantages lies in its avoidance of calling complicated subroutines,
as the subproblem is smooth convex quadratic and unconstrained, making it
remarkably efficient compared to CarGou:11 , where the subproblem is not
smooth, and compared to XieWri:21 ; LinMa:22 , where the subproblem is
nonconvex. In fact, solving the subproblem in our approach simplifies to
solving a linear system of equalities. With global theoretical convergence to
an $\epsilon$ first-order solution in at most
$\mathcal{O}(\frac{1}{\epsilon^{3}})$ iterations, our method ensures the
reliable finding of optimal solutions for a broad spectrum of nonconvex
optimization problems, even those with nonlinear constraints, unlike KonMel:18
which can only handle linear equality constraints. Hence, its simplicity and
effectiveness make it an attractive algorithm for a wide range of practical
applications.

###### Remark 6

It is essential to emphasize the significance of conducting a complexity
analysis study for a given scheme. This analysis not only provides worst-case
bounds, but also yields valuable insights into the tuning of algorithmic
parameters. For instance, the above study of the convergence of the proposed
method plays a pivotal role in determining the choice of $\beta_{k}$. This
choice might have been made differently without the insights derived from the
complexity analysis. Similarly, Theorem 2 addresses the proper selection of
$\rho$.

###  4.2 Improved rates under KL

In this section, under the KL condition, we provide better convergence rates
for the iterates of Algorithm 1. In particular, we prove that the whole
sequence $\\{x_{k}\\}_{k\geq 0}$ converges. Let’s first derive several
inequalities that will be useful for future proofs. From lemma 5, we get the
following:

$\|\nabla\mathcal{P}_{\rho}(x_{k+1})\|^{2}\leq\Gamma_{\text{max}}^{2}\|\Delta
x_{k+1}\|^{2}.$ (19)

Then, it follows from (19) and Lemma 3, that:

$P_{k+1}-P_{k}\leq-\frac{\underaccent{\bar}{\beta}}{2\Gamma_{\text{max}}^{2}}\left\|\nabla\mathcal{P}_{\rho}(x_{k+1})\right\|^{2}.$
(20)

We denote by $\text{crit }\mathcal{P}_{\rho}$ the set of critical points of
the function $\mathcal{P}_{\rho}(\cdot)$ defined in (7). Furthermore, we
denote $\mathcal{E}_{k}=P_{k}-P^{*}$, where $P^{*}=\lim_{k\to\infty}{P_{k}}$
(recall that the sequence $\\{P_{k}\\}_{k\geq 0}$ is decreasing and bounded
from bellow according to Lemma 3 and Lemma 4, respectively, hence it is
convergent). Let us denote the set of limit points of $\\{x_{k}\\}_{k\geq 0}$
by:

$\Omega:=\\{x_{\rho}^{*}\;:\;\exists\text{ a convergent
subsequence}\;\\{x_{k}\\}_{k\in\mathcal{K}}\;\text{such
that}\lim_{k\in\mathcal{K}}{x_{k}}=x_{\rho}^{*}\\}.$

Let’s now introduce the following lemma, which highlights the relationship
between the set of limit points of the generated sequence, denoted $\Omega$,
and the critical points of the Lyapunov function, denoted by
$\text{crit}\,\mathcal{P}_{\rho}$.

###### Lemma 6

Consider Algorithm 1 and let $\\{P_{k}\\}_{k\geq 0}$ be defined as in (7). If
the sequence, $\\{x_{k}\\}_{k\geq 0}$, generated by Algorithm 1 is in some
compact set $\mathcal{S}$ on which Assumptions 1, 2 and 3 hold and, moreover,
$\rho$ is chosen as in (10) and $x_{0}$ is chosen as in (8), then the
following statements hold:

  1. (i)

$\Omega$ is a compact subset of crit $\mathcal{P}_{\rho}$ and
$\lim_{k\to\infty}{\text{dist}(x_{k},\Omega)}=0$.

  2. (ii)

For any $x\in\Omega,$ we have $\mathcal{P}_{\rho}(x)=P^{*}$.

###### Proof

See Appendix.

Let us now show that the sequence $\|\Delta x_{k}\|$ has finite length,
provided that a KL condition holds.

###### Lemma 7

If the sequence, $\\{x_{k}\\}_{k\geq 0}$, generated by Algorithm 1 is in some
compact set $\mathcal{S}$ on which Assumptions 1, 2 and 3 hold and assume that
$\mathcal{P}_{\rho}(\cdot)$ satisfies the KL property on $\Omega$, then
$\\{x_{k}\\}_{k\geq 0}$ satisfies the finite length property, i.e.,

$\sum_{k=1}^{\infty}{\|\Delta x_{k}\|}<\infty,$

and consequently converges to a critical point of $\mathcal{P}_{\rho}(\cdot)$.

###### Proof

See Appendix.

Lemma 7 shows that under the KL property the set of limit points of the
sequence $\\{x_{k}\\}_{k\geq 0}$ generated by Algorithm 1 is a singleton,
denoted by $x_{\rho}^{*}$. In the next theorem, we will prove how fast the
iterates converge to this limit point (see also Yas:22 for a similar
reasoning).

###### Theorem 4.3

[Convergence of the whole sequence $\\{x_{k}\\}_{k\geq 0}$] Let Assumptions 1,
2 and 3 hold and $x_{\rho}^{*}$ is the limit point of the sequence
$\\{x_{k}\\}_{k\geq 0}$ generated by Algorithm 1. If
$\mathcal{P}_{\rho}(\cdot)$ satisfies KL property at $x_{\rho}^{*}$, then
there exists $k_{1}\geq 1$ such that for all $k\geq k_{1}$ we have:

$\|x_{k}-x_{\rho}^{*}\|\leq
C\max\\{\varphi(\mathcal{E}_{k}),\sqrt{\mathcal{E}_{k-1}}\\},$

where $C>0$ and $\varphi\in\Psi_{\tau}$, with $\tau>0$, denotes a
desingularizing function.

###### Proof

From Lemma 3, we get:

$\displaystyle
P_{k+1}-P_{k}{\overset{{}}{\leq}}-\frac{\underaccent{\bar}{\beta}}{2}\|x_{k+1}-x_{k}\|^{2}.$
(21)

Based on our choice of $\rho$ and $\beta_{k}$, the sequence
$\\{P_{k}\\}_{k\geq 0}$ is monotonically decreasing, see Lemma 3, and
consequently $\\{\mathcal{E}_{k}\\}_{k\geq 0}$ is monotonically decreasing.
Using (21) and the fact that $\\{\mathcal{E}_{k}\\}_{k\geq 0}$ is nonnegative,
we have for all $k\geq 0$:

$\|\Delta
x_{k+1}\|\leq\sqrt{\frac{2}{\underaccent{\bar}{\beta}}}\sqrt{\mathcal{E}_{k}}.$
(22)

Since $P_{k}\to P^{*}$, ${x_{k}}\to x_{\rho}^{*}$ and
$\mathcal{P}_{\rho}(\cdot)$ satisfies the KL property at $x_{\rho}^{*}$, then
there exists $k_{1}=k_{1}(\epsilon,\tau)\geq 1$ such that $\forall k>k_{1}$,
we have $\|x_{k}-x_{\rho}^{*}\|\leq\epsilon$ and $P^{*}<P_{k}<P^{*}+\tau$, and
the following KL property holds:

$\varphi^{\prime}(\mathcal{E}_{k})\|\nabla\mathcal{P}_{\rho}(x_{k})\|\geq 1.$
(23)

Since $\varphi$ is concave function, we have
$\varphi(\mathcal{E}_{k})-\varphi(\mathcal{E}_{k+1})\geq\varphi^{\prime}(\mathcal{E}_{k})(\mathcal{E}_{k}-\mathcal{E}_{k+1})$.
Therefore, from (21) and (23), we get:

$\displaystyle\|x_{k+1}-x_{k}\|^{2}$
$\displaystyle\leq\varphi^{\prime}(\mathcal{E}_{k})\|x_{k+1}-x_{k}\|^{2}\|\nabla\mathcal{P}_{\rho}(x_{k})|$
$\displaystyle\leq\frac{2}{\underaccent{\bar}{\beta}}\varphi^{\prime}(\mathcal{E}_{k})(\mathcal{E}_{k}-\mathcal{E}_{k+1})\|\nabla\mathcal{P}_{\rho}(x_{k})\|$
$\displaystyle\leq\frac{2}{\underaccent{\bar}{\beta}}\Big{(}\varphi(\mathcal{E}_{k})-\varphi(\mathcal{E}_{k+1})\Big{)}\|\nabla\mathcal{P}_{\rho}(x_{k})\|.$

Note that for given $a,b,c\geq 0$, if we have ${a^{2}}\leq 2b\times c$, and
recognizing that we always have $2b\times c\leq(b+c)^{2}$, then
${a^{2}}\leq(b+c)^{2}$, which in turn implies $a\leq b+c$. It follows that for
any $\theta>0$, we have (by taking $a=\|\Delta x_{k+1}\|$,
$b=\frac{\theta}{\underaccent{\bar}{\beta}}\Big{(}\varphi(\mathcal{E}k)-\varphi(\mathcal{E}{k+1})\Big{)}$,
and $c=\frac{1}{\theta}\|\nabla\mathcal{P}_{\rho}(x_{k})\|$):

$\displaystyle{\|\Delta x_{k+1}\|}\leq$
$\displaystyle\frac{\theta}{\underaccent{\bar}{\beta}}\Big{(}\varphi(\mathcal{E}_{k})-\varphi(\mathcal{E}_{k+1})\Big{)}+\frac{1}{\theta}\|\nabla\mathcal{P}_{\rho}(x_{k})\|.$
(24)

Furthermore, we have:

$\displaystyle\|\nabla\mathcal{P}_{\rho}(x_{k})\|{\overset{{}}{\leq}}\Gamma_{\text{max}}\|\Delta
x_{k}\|,$

where the last inequality follows from the definition of $\mathcal{P}_{\rho}$
and the properties of the derivative together with Lemma 5 and (14). Then,
(24) becomes:

$\displaystyle\|\Delta x_{k+1}\|\leq$
$\displaystyle\frac{\theta}{\underaccent{\bar}{\beta}}\Big{(}\varphi(\mathcal{E}_{k})-\varphi(\mathcal{E}_{k+1})\Big{)}+\frac{\Gamma_{\text{max}}}{\theta}\|\Delta
x_{k}\|.$

It follows that:

$\displaystyle\left(1-\frac{\Gamma_{\text{max}}}{\theta}\right)\|\Delta
x_{k+1}\|\leq$
$\displaystyle\frac{\theta}{\underaccent{\bar}{\beta}}\Big{(}\varphi(\mathcal{E}_{k})-\varphi(\mathcal{E}_{k+1})\Big{)}+\frac{\Gamma_{\text{max}}}{\theta}\Big{(}\|\Delta
x_{k}\|-\|\Delta x_{k+1}\|\Big{)}.$

By summing up the above inequality over $k\geq k_{1}$, we get:

$\displaystyle\left(1-\frac{\Gamma_{\text{max}}}{\theta}\right)\sum_{k\geq
k_{1}}{\|\Delta x_{k+1}\|}\leq$
$\displaystyle\frac{\theta}{\underaccent{\bar}{\beta}}\varphi(\mathcal{E}_{{k_{1}}})+\frac{\Gamma_{\text{max}}}{\theta}\|\Delta
x_{k_{1}}\|.$

Let us now choose $\theta>0$ such that
$0<\frac{\Gamma_{\text{max}}}{\theta}<1$ and define $\delta_{0}$ as
$\delta_{0}=1-\frac{\Gamma_{\text{max}}}{\theta}>0$. Then, we have:

$\displaystyle\sum_{k\geq k_{1}}{\|\Delta x_{k+1}\|}\leq$
$\displaystyle\frac{\theta}{\underaccent{\bar}{\beta}\delta_{0}}\varphi(\mathcal{E}_{{k_{1}}})+\frac{\Gamma_{\text{max}}}{\theta\delta_{0}}\|\Delta
x_{k_{1}}\|.$

Hence, using the triangle inequality, we get for any $k\geq k_{1}$:

$\displaystyle\|x_{k}-x_{\rho}^{*}\|\leq\sum_{l\geq
k}{\|x_{l}-x_{l+1}\|}\leq\frac{\theta}{\underaccent{\bar}{\beta}\delta_{0}}\varphi(\mathcal{E}_{{k}})+\frac{\Gamma_{\text{max}}}{\theta\delta_{0}}\|\Delta
x_{k}\|.$

Further, using (22), it follows that:

$\displaystyle\|x_{k}-x_{\rho}^{*}\|\leq\frac{\theta}{\underaccent{\bar}{\beta}\delta_{0}}\varphi(\mathcal{E}_{{k}})+\frac{\Gamma_{\text{max}}}{\theta\delta_{0}}\sqrt{\frac{2}{\underaccent{\bar}{\beta}}}\sqrt{\mathcal{E}_{k-1}}\leq
C\max\\{\varphi(\mathcal{E}_{{k}}),\sqrt{\mathcal{E}_{k-1}}\\},$

where
$C=\frac{\theta}{\underaccent{\bar}{\beta}\delta_{0}}+\frac{\Gamma_{\text{max}}}{\theta\delta_{0}}\sqrt{\frac{2}{\underaccent{\bar}{\beta}}}.$
This concludes our proof.∎

The next theorem derives the convergence rates of the sequence generated by
Algorithm 1 when the Lyapunov function satisfies the KL property with the
desingularizing function $\varphi$:

$\varphi:[0,\tau)\to[0,+\infty),\;\varphi(s)=s^{1-\nu},\text{ where
}\nu\in[0,1).$

###### Theorem 4.4

[Convergence rates of $\\{x_{k}\\}_{k\geq 0}$] Let the assumptions of Theorem
4.3 hold and $x_{\rho}^{*}$ be the limit point of the sequence
$\\{x_{k}\\}_{k\geq 0}$ generated by Algorithm 1. If
$\mathcal{P}_{\rho}(\cdot)$ satisfies KL property at $x_{\rho}^{*}$ with
desingularizing function

$\varphi:[0,\tau)\to[0,+\infty),\;\varphi(s)=s^{1-\nu},\text{ where
}\nu\in[0,1),$

then the following rates hold:

  1. 1.

If $\nu=0$, then $x_{k}$ converges to $x_{\rho}^{*}$ in a finite number of
iterations.

  2. 2.

If
$\nu\in(0,\frac{1}{2}{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}]}$,
then there exists $k_{1}$ such that for all $k\geq k_{1}$, we have the
following linear convergence:

$\|x_{k}-x_{\rho}^{*}\|\leq
C\frac{\sqrt{\mathcal{E}_{k_{1}}}}{\sqrt{(1+\bar{c}\mathcal{E}_{k_{1}}^{2\nu-1}})^{k-k_{1}}},\quad\text{where}\;\;\bar{c}=\frac{\underaccent{\bar}{\beta}}{(1-\nu)^{2}\Gamma_{\text{max}}^{2}}.$

  3. 3.

If $\nu\in(\frac{1}{2},1)$, then there exists $k_{1}$ such that for all
$k>k_{1}$, we have the following sublinear convergence:

$\|x_{k}-x_{\rho}^{*}\|\leq
C\left(\frac{1}{\mu(k-k_{1})+\mathcal{E}_{k_{1}}^{1-2\nu}}\right)^{\frac{1-\nu}{2\nu-1}}.$

###### Proof

Let $\nu\in[0,1)$ and for all $s\in[0,\tau),\varphi(s)=s^{1-\nu}$ and
$\varphi^{\prime}(s)=(1-\nu)s^{-\nu}$. From Theorem 4.3, it follows that there
exists $k_{1}\geq 1$ such that $\forall k\geq k_{1}$, we have:

$\|x_{k}-x_{\rho}^{*}\|\leq
C\max\\{\mathcal{E}_{k}^{1-\nu},\sqrt{\mathcal{E}_{k-1}}\\}.$ (25)

Furthermore, (23) yields:

${\mathcal{E}_{k}}^{\nu}\leq(1-\nu)\|\nabla\mathcal{P}_{\rho}(x_{k})\|\hskip
14.22636pt\forall k\geq k_{1}.$

Moreover, from (20), we have for any $k\geq 1$:

$\|\nabla\mathcal{P}_{\rho}(x_{k})\|^{2}\leq\frac{2\Gamma_{\text{max}}^{2}}{\underaccent{\bar}{\beta}}(\mathcal{E}_{k-1}-\mathcal{E}_{k}).$

Hence,
${\mathcal{E}_{k}}^{2\nu}\leq\frac{2(1-\nu)^{2}\Gamma_{\text{max}}^{2}}{\underaccent{\bar}{\beta}}(\mathcal{E}_{k-1}-\mathcal{E}_{k})\hskip
14.22636pt\forall k>k_{1}.$ Setting
$\bar{c}=\frac{\underaccent{\bar}{\beta}}{2(1-\nu)^{2}\Gamma_{\text{max}}^{2}}>0,$
we get:

$\bar{c}{\mathcal{E}_{k}}^{2\nu}\leq\mathcal{E}_{k-1}-\mathcal{E}_{k}\hskip
14.22636pt\forall k>k_{1}.$ (26)

  1. 1.

Let $\nu=0$. If $\mathcal{E}_{k}>0$ for any $k>k_{1}$, we have
$\bar{c}\leq\mathcal{E}_{k-1}-\mathcal{E}_{k}$. As $k$ goes to infinity, the
right hand side approaches zero. Then, $0<\bar{c}\leq 0$ which is a
contradiction. Hence, there exists $k>k_{1}$ such that $\mathcal{E}_{k}=0.$
Then, $\mathcal{E}_{k}\to 0$ in a finite number of steps and from (25),
$x_{k}\to x_{\rho}^{*}$ in a finite number of steps.

  2. 2.

Let
$\nu\in(0,\frac{1}{2}{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}]}$.
Then,
$2\nu-1{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color<EMAIL_ADDRESS>Let $k>k_{1}$. Since $\\{\mathcal{E}_{i}\\}_{i\geq k_{1}}$ is monotonically
decreasing, then $\mathcal{E}_{i}\leq\mathcal{E}_{k_{1}}$ for any
$i\in\\{k_{1}+1,k_{1}+2,...,k\\}$ and

$\bar{c}{\mathcal{E}_{k_{1}}}^{2\nu-1}\mathcal{E}_{k}\leq\mathcal{E}_{k-1}-\mathcal{E}_{k}\hskip
14.22636pt\forall k>k_{1}.$

Rearranging this, we get for all $k>k_{1}$:

$\mathcal{E}_{k}\leq\frac{\mathcal{E}_{k-1}}{1+\bar{c}{\mathcal{E}_{k_{1}}}^{2\nu-1}}\leq\frac{\mathcal{E}_{k-2}}{(1+\bar{c}{\mathcal{E}_{k_{1}}}^{2\nu-1})^{2}}\leq...\leq\frac{\mathcal{E}_{k_{1}}}{(1+\bar{c}{\mathcal{E}_{k_{1}}}^{2\nu-1})^{k-k_{1}}}.$

Since $\\{\mathcal{E}_{k}\\}_{k\geq k_{1}}$ is monotonically decreasing, then
$\mathcal{E}_{k}\leq\mathcal{E}_{k-1}$, and since $\nu\in(0,\frac{1}{2}]$,
then $1-\nu\geq\frac{1}{2}$. Hence,

$\mathcal{E}_{k}^{1-\nu}\leq\mathcal{E}_{k}^{\frac{1}{2}}\leq\mathcal{E}_{k-1}^{\frac{1}{2}}.$

Then, it follows that:

$\max\\{\mathcal{E}_{k}^{1-\nu},\sqrt{\mathcal{E}_{k-1}}\\}=\sqrt{\mathcal{E}_{k-1}}.$

Then, we have
$\max\\{{\mathcal{E}_{k}}^{1-\nu},\sqrt{\mathcal{E}_{k-1}}\\}=\sqrt{\mathcal{E}_{k-1}}.$
It then follows that:

$\|x_{k}-x_{\rho}^{*}\|\leq
C\frac{\sqrt{\mathcal{E}_{k_{1}}}}{\sqrt{(1+\bar{c}\mathcal{E}_{k_{1}}^{2\nu-1}})^{k-k_{1}}},$

  3. 3.

Let $\nu\in(\frac{1}{2},1)$. From (26), we have:

$\bar{c}\leq(\mathcal{E}_{k-1}-\mathcal{E}_{k}){\mathcal{E}^{-2\nu}_{k}}\hskip
14.22636pt\forall k>k_{1}.$ (27)

Let $h:\mathbb{R}^{*}_{+}\to\mathbb{R}$ be defined as $h(s)=s^{-2\nu}$ for any
$s\in\mathbb{R}^{*}_{+}$. It is clear that $h$ is monotonically decreasing and
$\forall s\in\mathbb{R}_{+},h^{\prime}(s)=-2\nu s^{-(1+2\nu)}<0$. Since
$\mathcal{E}_{k}\leq\mathcal{E}_{k-1}$ for all $k>k_{1}$, then
$h(\mathcal{E}_{k-1})\leq h(\mathcal{E}_{k})$ for all $k>k_{1}$. We consider
two scenarios: in the first case, for any $k>k_{1}$, there exists an
$r_{0}\in(1,+\infty)$ such that: $h(\mathcal{E}_{k})\leq
r_{0}h(\mathcal{E}_{k-1})$; in the second case, for some $k>k_{1}$, regardless
of the choice of $r_{0}\in(1,+\infty)$, we have
$h(\mathcal{E}_{k})>r_{0}h(\mathcal{E}_{k-1})$.

Case 1: There exists $r_{0}\in(1,+\infty)$ such that $h(\mathcal{E}_{k})\leq
r_{0}h(\mathcal{E}_{k-1}),\;\forall k>k_{1}.$ Then, from (27) we get:

$\displaystyle\bar{c}\leq{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}(\mathcal{E}_{k-1}-\mathcal{E}_{k})h(\mathcal{E}_{k})}\leq(\mathcal{E}_{k-1}-\mathcal{E}_{k})r_{0}h(\mathcal{E}_{k-1})=r_{0}h(\mathcal{E}_{k-1})\int_{\mathcal{E}_{k}}^{\mathcal{E}_{k-1}}{1\,ds}$
$\displaystyle{\overset{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}h(\cdot)\text{
decreases}}}{\leq}}r_{0}\int_{\mathcal{E}_{k}}^{\mathcal{E}_{k-1}}{h(s)\,ds}=r_{0}\int_{\mathcal{E}_{k}}^{\mathcal{E}_{k-1}}{s^{-2\nu}\,ds}=\frac{r_{0}}{1-2\nu}({\mathcal{E}^{1-2\nu}_{k-1}}-{\mathcal{E}^{1-2\nu}_{k}}).$

Since $\nu>\frac{1}{2}$, then $2\nu-1>0$ and multiplying both sides of the
above inequality by $\frac{2\nu-1}{r_{0}}$, we get that:

$0<\frac{\bar{c}(2\nu-1)}{r_{0}}\leq{\mathcal{E}^{1-2\nu}_{k}}-{\mathcal{E}^{1-2\nu}_{k-1}}.$

Let us define $\hat{c}=\frac{\bar{c}(2\nu-1)}{r_{0}}$ and
$\hat{\nu}=1-2\nu<0$. We then get:

$0<\hat{c}\leq{\mathcal{E}^{\hat{\nu}}_{k}}-{\mathcal{E}^{\hat{\nu}}_{k-1}}\hskip
14.22636pt\forall k>k_{1}.$ (28)

Case 2: For any $r_{0}\in(1,+\infty)$, there exists some $k>k_{1}$ such that
$h(\mathcal{E}_{k})>r_{0}h(\mathcal{E}_{k-1})$, which can be expressed as
${\mathcal{E}_{k}^{-2\nu}}\geq r_{0}{\mathcal{E}_{k-1}^{-2\nu}}$. This
implies: $r_{0}{\mathcal{E}_{k-1}^{2\nu}}\geq{\mathcal{E}_{k}^{2\nu}}.$ By
setting $q=\frac{1}{r_{0}^{\frac{1}{2\nu}}}\in(0,1)$, we obtain,

$q\mathcal{E}_{k-1}\geq\mathcal{E}_{k}.$

Since $\hat{\nu}=1-2\nu<0$ we have
$q^{\hat{\nu}}{\mathcal{E}^{\hat{\nu}}_{k-1}}\leq{\mathcal{E}^{\hat{\nu}}_{k}}$
and then, it follows that:

$(q^{\hat{\nu}}-1){\mathcal{E}^{\hat{\nu}}_{k-1}}\leq{\mathcal{E}^{\hat{\nu}}_{k}}-{\mathcal{E}^{\hat{\nu}}_{k-1}}.$

Since $q\in(0,1)$ and $\hat{\nu}<0$, we have $q^{\hat{\nu}}-1>0$ and since,
moreover, $\mathcal{E}_{k}\to 0^{+}$ as $k\to\infty$ which is equivalent to
$\mathcal{E}^{\hat{\nu}}_{k}\to+\infty$ as $k\to\infty$ , then there exists
$\tilde{c}$ such that
$(q^{\hat{\nu}}-1){\mathcal{E}^{\hat{\nu}}_{k-1}}\geq\tilde{c}$ for all
$k>k_{1}$ (it is sufficient to take
$\tilde{c}=(q^{\hat{\nu}}-1){\mathcal{E}^{\hat{\nu}}_{k_{1}}}$). Therefore, it
follows that:

$0<\tilde{c}\leq(q^{\hat{\nu}}-1){\mathcal{E}^{\hat{\nu}}_{k-1}}\leq{\mathcal{E}^{\hat{\nu}}_{k}}-{\mathcal{E}^{\hat{\nu}}_{k-1}}.$
(29)

Hence, regardless of the cases involved (case 1 or case 2), one can always
combine (28) and (29) by choosing ${\mu}=\min\\{\hat{c},\tilde{c}\\}>0$, to
obtain

$0<{\mu}\leq{\mathcal{E}^{\hat{\nu}}_{k}}-{\mathcal{E}^{\hat{\nu}}_{k-1}}\hskip
14.22636pt\forall k>k_{1}.$

Summing the above inequality from $k_{1}+1$ to some $k>k_{1}+1$ gives

$\mu\sum_{i=k_{1}+1}^{k}{1}\leq\sum_{i=k_{1}+1}^{k}{\left({\mathcal{E}^{\hat{\nu}}_{i}}-{\mathcal{E}^{\hat{\nu}}_{i-1}}\right)}.$

By simplifying the summation, we get:

$\mu(k-k_{1})+{\mathcal{E}^{\hat{\nu}}_{k_{1}}}\leq{\mathcal{E}^{\hat{\nu}}_{k}}.$

Moreover, since $\hat{\nu}=1-2\nu<0$, it follows that

$\mathcal{E}_{k}\leq(\mu(k-k_{1})+{\mathcal{E}^{\hat{\nu}}_{k_{1}}})^{\frac{1}{\hat{\nu}}}=(\mu(k-k_{1})+{\mathcal{E}^{1-2{\nu}}_{k_{1}}})^{\frac{1}{1-2{\nu}}}.$

Furthermore, since $\\{\mathcal{E}_{k}\\}_{k\geq k_{1}}$ is monotonically
decreasing, then $\mathcal{E}_{k}\leq\mathcal{E}_{k-1}$, and since
$\nu\in(\frac{1}{2},1)$, then $0<1-\nu<\frac{1}{2}$. Hence

$\mathcal{E}_{k}^{1-\nu}\leq\mathcal{E}_{k-1}^{1-\nu}\leq\quad\text{and}\quad\mathcal{E}_{k-1}^{\frac{1}{2}}\leq\mathcal{E}_{k-1}^{1-\nu}.$

Then, it follows that:

$\max\\{\mathcal{E}_{k}^{1-\nu},\sqrt{\mathcal{E}_{k-1}}\\}\leq\max\\{\mathcal{E}_{k-1}^{1-\nu},\sqrt{\mathcal{E}_{k-1}}\\}={\mathcal{E}^{1-\nu}_{k-1}}.$

Thus, (25) becomes:

$\|x_{k}-x_{\rho}^{*}\|\leq
C\left(\frac{1}{\mu(k-k_{1}-1)+\mathcal{E}_{k_{1}}^{1-2\nu}}\right)^{\frac{1-\nu}{2\nu-1}},\hskip
14.22636pt\forall k>k_{1}.$

This concludes our proof.∎

Note that in Theorems 4.3 and 4.4, $x_{\rho}^{*}$ represents a critical point
of the Lyapunov function $\mathcal{P}_{\rho}(\cdot)$, which is obtained for
any choice of $\rho$ that satisfies (10). However, as indicated in Theorem
4.2, it is essential to choose the penalty parameter $\rho$ large enough for
$x_{\rho}^{*}$ to qualify as an $\epsilon$ first order solution to the problem
$\eqref{eq1}$. Consequently, by following the prescribed choice of $\rho$ in
Theorem 4.2 and considering the KL condition, our algorithm guarantees
convergence to an $\epsilon$ first order solution at the (sub)linear rates
determined by Theorem 4.4.

###  4.3 Selection of the penalty parameter $\rho$

The above results, showcasing the total number of iterations required to find
an $\epsilon$ first-order solution to the problem, rely on the assumption that
the penalty parameter $\rho$ exceeds a certain threshold. However, determining
this threshold beforehand poses challenges as it depends on unknown parameters
of the functions within the problem and the algorithm’s parameters. To
overcome this challenge, we propose a scheme that allows us to determine a
sufficiently large value of $\rho$ without needing specific parameter
information. Inspired by Algorithm 3 in XieWri:21 , our approach involves
repeatedly calling Algorithm 1 as an inner loop. If Algorithm 1 fails to
converge within a predetermined number of iterations, we gradually increase
the penalty parameter $\rho$ by a constant multiple in the outer loop. The
detailed implementation of this scheme is outlined in Algorithm 3. Remarkably,
in certain applications, such as Model Predictive Control, we need to solve a
series of problems with a desired level of accuracy $\epsilon$. Interestingly,
the optimal penalty parameter for the initial optimization problem often works
well for subsequent problems as well. This implies that investing more
computational resources into the first optimization problem can lead to
considerable savings for the remaining optimization problems. Thus, this
adaptive approach ensures the effectiveness and efficiency of our algorithm,
even in scenarios where precise parameter information is unavailable. This
adaptability makes our method a practical and robust solution for a wide range
of optimization problems, with potential benefits in various real-world
applications.

Algorithm 3 LQP method with trial value of $\rho$

1:$\textbf{Initialization: }\tau>1,\epsilon>0\text{ and }\rho_{0}>0;$

2:$k\leftarrow 0$

3:while Infeasible  do

4: $\rho_{k+1}\leftarrow\tau\rho_{k}$

5: Call Algorithm 1 using warm start

6: $k\leftarrow k+1$

7:end while

Algorithm 3 is well-defined and terminates in a finite number of iterations,
provided that the parameter $\tau>1$ and that $\lambda^{*}$ is bounded by some
positive constant $M$ (guaranted by LICQ). In fact, during the $k$-th
iteration of Algorithm 3, we have $\rho_{k+1}=\tau^{k+1}\rho_{0}$.
Consequently, from Theorem 4.1, $\lambda^{*}=\rho F(x^{*}_{\rho})$, and we
have:

$\|F(x_{\rho}^{*})\|=\frac{\|\lambda^{*}\|}{\tau^{k+1}\rho_{0}}\leq\frac{M}{\rho_{0}\tau^{k+1}}.$

Hence, $\|F(x_{\rho}^{*})\|\leq\epsilon$ provided that
$\frac{M}{\rho_{0}\tau^{k+1}}\leq\epsilon$, which is equivalent to
$\tau^{k+1}\geq\frac{M}{\rho_{0}\epsilon}$. Thus, we require:

$(k+1)\log\tau\geq\log\frac{M}{\rho_{0}}+\log\frac{1}{\epsilon}.$

This leads to the following condition, considering that
$\tau=\frac{1}{\epsilon}$:

$k\geq\frac{\log\frac{M}{\rho_{0}}}{\log\frac{1}{\epsilon}}.$

Therefore, for $\tau=\frac{1}{\epsilon}$, one has to increase $\rho$ at most
$\frac{\log\frac{M}{\rho_{0}}}{\log\frac{1}{\epsilon}}$ times to reach a value
that guarantees the desired precision.

##  5 Numerical results

In this section we numerically compare Algorithm 1 (LQP) with SCP algorithm
MesBau:21 and IPOPT WacBie:06 and with Inexact Restoration algorithm, which we
refer to as the 2-phases method BueMar:20 , on nonconvex optimization problems
with nonlinear equality constraints. The simulations are implemented in Python
and executed on a PC with (CPU 2.90GHz, 16GB RAM). Since one cannot guarantee
that the SCP iterates converge to a first-order (KKT) point, we choose the
following stopping criteria: we stop the algorithms when the difference
between two consecutive values of the objective function is less than a
tolerance $\epsilon_{1}=10^{-3}$ and the norm of constraints is less than a
tolerance $\epsilon_{2}=10^{-5}$. For the implementation of our method, we
choose $\beta_{k}$ constant and equal to $10^{4}$, and we select the penalty
parameter to be constant for a given problem, specifically $\rho$ is chosen
from the interval $[10^{7},10^{9}]$. A problem is considered to be solved by a
method if the stopping criteria are achieved in less than 30 minutes and if
the method in question finds a viable point whose objective value is
sufficiently close to the value found by the best algorithm. The numerical
results are illustrated in Table 1 and Figure 1.

(n,m) Alg LQP SCP IPOPT 2-phasesBueMar:20 # iter cpu # iter cpu # iter cpu #
iter cpu $f^{*}$ $\|F\|$ $f^{*}$ $\|F\|$ $f^{*}$ $\|F\|$ $f^{*}$ $\|F\|$
OPTCTRL3 5 0.11 5 0.16 7 7.40 - - (119,80) 2047.99 1.43e-06 2048.01 4.52e-10
2048.01 1.84e-08 - - OPTCTRL3 7 0.35 7 1.55 10 11.99 - - (1199,800) 18459.99
4.89e-06 18460.22 1.84e-09 18460.22 6.33e-09 - - OPTCTRL3 24 10.63 24 105.08
11 26.95 - - (4499,3000) 74464.82 1.43e-06 74465.03 6.87e-09 74465.03 1.09e-08
- - DTOC4 3 0.11 4 0.14 3 7.61 4 0.22 (299,198) 2.94 2.33e-07 2.94 3.48e-09
2.94 1.02e-08 2.94 1.29e-06 DTOC4 3 0.29 3 0.87 3 12.45 6 1.39 (1497,998) 2.87
1.16e-06 2.88 1.16e-06 2.88 1.20e-08 2.88 6.28e-06 DTOC4 4 1.07 4 5.81 3 23.51
3 2.33 (2997,1998) 2.87 2.832-07 2.87 2.83e-07 2.87 9.33e-09 2.87 2.64e-07
DTOC4 4 2.87 3 16.74 3 29.02 4 6.68 (4497,2998) 2.87 3.02e-07 2.87 4.87e-10
2.87 3.66e-08 2.87 3.00e-06 DTOC4 4 32.30 4 566.80 3 146.73 - - (14997,9998)
2.87 1.05e-07 2.87 1.05e-07 2.86 4.49e-09 - - DTOC5 5 0.12 5 0.21 4 10.63 13
0.36 (198,99) 1.53 6.88e-08 1.53 6.64e-08 1.53 4.99e-11 1.53 7.38e-06 DTOC5 7
0.29 7 0.98 3 12.06 5 0.73 (998,499) 1.53 3.42e-06 1.53 3.45e-06 1.53 7.76e-07
1.53 1.47e-06 DTOC5 10 1.33 10 4.52 3 18.74 6 1.95 (1998,999) 1.53 1.56e-06
1.53 1.56e-06 1.53 6.88e-08 1.53 4.02e-06 DTOC5 24 47.03 24 799.07 3 75.25 16
137.53 (9998,4999) 1.54 1.97e-07 1.54 1.96e-07 1.53 2.49e-07 1.53 3.84e-07
DTOC6 30 0.30 30 0.72 8 8.22 81 2.10 (200,100) 729.00 2.20e-06 729.00 1.01e-09
727.98 3.52e-10 728.29 1.14e-06 DTOC6 32 1.17 32 4.56 9 13.10 106 14.96
(1000,500) 6848.01 2.97-06 6848.06 2.7e-10 6846.61 7.56e-10 6847.45 1.92e-06
DTOC6 33 3.88 33 15.25 9 18.85 118 33.22 (2000,1000) 17177.66 8.95e-06
17178.07 1.46e-09 17176.45 5.24e-09 17177.56 1.25e-06 ORTHREGA 54 0.52 18 0.47
66 7.37 - - (133,64) 350.30 2.42e-06 414.53 1.87e-06 350.30 3.02e-10 - -
ORTHREGA 55 1.22 39 1.73 76 10.14 - - (517,256) 1414.05 1.08e-06 1664.80
1.24e-06 1414.05 6.19e-10 - - ORTHREGA 94 14.24 67 31.78 14 23.99 - -
(2053,1024) 5661.43 3.19e-06 6654.78 2.07e-06 5661.43 9.25e-07 - - MSS1 17
0.22 12 0.15 53 13.52 42 0.58 (90, 73) -16.00 6.04e-06 -8.71e-08 1.76e-06
-16.00 4.17e-08 -15.48 3.93e-06 MSS2 19 2.64 21 8.05 7 14.65 15 7.75 (756,
703)) -121.03 6.20e-06 -2.53e-10 6.12e-06 -26.97 5.96e-08 -0.37 3.55e-06 MSS3
68 78.73 22 135.15 - - 23 155.51 (2070, 1981) -339.69 8.17e-06 -5.29e-09
7.76e-06 - - -1e-03 7.98e-06 OPTCTRL6 - - 24 13.46 13 27.42 - - (4499, 3000) -
- 74465.03 3.27e-09 74465.03 2.32e-09 - - OPTCDEG2 4 1.31 3 2.71 4 7.92 3 4.87
(4499, 3000) 59.15 3.07e-07 59.08 5.85e-08 227.72 6.55e-08 59.08 9.39e-10
ORTHREGC 101 85.75 29 20.93 16 31.14 6 19.62 (5005, 2500) 1681.78 8.65e-06
94.81 8.42e-06 94.81 7.52e-07 133.86 8.68e-08 EIGENB2 11 18.76 5 5.28 27 56.61
9 53.24 (2550, 1275) 110.50 3.72e-13 110.50 1.61e-14 0.00 5.45e-09 110.49
1.12e-08 EIGENC2 6 3.12 6 4.68 13 24.93 10 83.99 (2652, 1326) 11162.34
3.25e-14 11162.75 4.64e-16 0.00 8.43e-10 11150.50 2.12e-10 EIGENACO 19 9.98 8
1.75 - - 14 10.52 (2550, 1275) 10731.01 1.00e-15 22425.04 2.37e-18 - -
10731.25 1.54e-09 EIGENBCO 6 11.80 5 1.23 9 19.58 11 50.84 (2550, 1275) 49.50
4.19e-14 49.50 5.79e-16 0.00 3.18e-17 59.50 8.11e-07 EIGENCCO 7 4.06 7 1.55 13
42.88 7 125.49 (2652, 1326) 11100.50 7.56e-11 11100.51 1.99e-10 0.00 2.75e-12
11100.50 2.01e-10 DTOC1NA 50 50.67 4 3.86 5 12.08 14 40.81 (5994, 3996) 47.64
5.95e-11 47.66 5.03e-13 4.15 7.44e-11 47.64 2.96e-09 DTOC1NB 46 46.65 4 3.89 5
11.71 13 38.32 (5994, 3996) 48.45 6.94e-11 48.47 1.68e-14 7.13 6.31e-12 48.44
2.46e-08 DTOC1NC 27 28.36 6 5.76 3 8.12 6 17.78 (5994, 3996) 58.63 1.43e-09
58.64 3.20e-11 35.21 5.80e-10 58.64 1.67e-07 DTOC1ND 22 24.04 7 6.70 3 12.79 6
17.75 (5994, 3996) 66.65 1.81e-11 66.66 3.38e-11 47.63 7.85e-09 66.67 6.36e-07
SPINOP 11 4.25 - - - - 7 6.23 (1327, 1325) 2998.78 4.10e-07 - - - - 3000.33
1.06e-10 DIXCHLNV 114 54.32 - - 28 47.31 8 3.29 (1000, 500) 4.87 8.94e-07 - -
0.00 3.58e-08 12.53 8.28e-08 GLIDER 69 136.18 - - 53 107.53 - - (5207, 4808)
-543.46 8.27e-07 - - -1247.97 9.01e-08 - - DTOC2 11 12.16 21 19.73 11 22.64 15
44.98 (5994, 3996) 0.52 3.39e-06 0.91 9.66e-06 0.50 6.21e-09 0.52 9.51e-06
ROBOTARM 12 9.18 - - 7 10.92 - - (4400, 3202) 18.34 2.88e-07 - - 9.14 2.05e-08
- - ROCKET 251 142.92 - - 5 7.49 11 6.90 (2403, 2002) -1.00 2.24e-06 - - -1.00
2.73e-07 -0.99 3.23e-06 CATMIX 12 2.65 9 2.84 3 4.36 4 2.30 (2401, 1600) -0.03
7.32e-08 -0.03 5.71e-09 -0.04 4.88e-09 -0.34 3.19e-07

Table 1: Comparison between LQP, SCP and IPOPT on selected problems from the
CUTEst collection.

In Table 1, we compare the number of iterations, CPU time (sec), objective
value, and feasibility violation (the Euclidean norm of functional
constraints) for LQP, SCP, IPOPT, and Inexact Restoration on real problems
selected from the CUTEst collection GouOrb:15 . We have chosen all problems
from CUTEst with nonlinear equality constraints, excluding obvious cases
(small sizes and problems with only linear equality constraints). Note that
for the test cases from CUTEst, our algorithm LQP is capable of yielding
solutions much faster than the other methods. However, for few problems, the
optimal values yielded by LQP, SCP and the 2-phases methods are not the best,
as those given by IPOPT. We attribute this phenomenon to the problem
structure, as observed in the test group ”EIGEN” and certain instances in
”DTOC” test group. Another possible reason, specifically for LQP, could be the
optimal objective value being much smaller than the constraint penalty, as we
suspect is the case for the problem instance ”ORTHREGC”. Moreover, from Table
1, it can be seen that our method fails only in one case, ’OPTCTRL6’, where it
cannot provide some information about the solution. In contrast, the other
methods fail in more than one case, thus proving that LQP is robust compared
to SCP, IPOPT, and the 2-phases methods.

Figure 1: Performance profiles for the computation time (left) and the number
of iterations (right).

Figure 1 shows the performance profiling for computation time (left) and
number of iterations (right) for the 3 algorithms. In the performance
profiling, the vertical axis $P(r_{p,s}\leq\tau)$ ($P(r_{p,s}\leq k)$)
represent the proportion of problems in the numerical experiments for which
$r_{p,s}$ does not exceed $\tau$ ($k$), respectively, where $r_{p,s}$ is the
ratio of the computational time (number of iterations) that the $s$ solver
takes to solve problem $p$ to the shortest computational time (number of
iterations) among the three algorithms to solve problem $p$, respectively. It
is clear from the computational time profile, Figure 1 (left), that the
proposed algorithm, LQP, approaches the limit line $1$ faster than SCP, IPOPT,
2-phases. However, for the number of iterations this is not always the case.

Based on these preliminary experiments, one can see that our new algorithm,
LQP, stands out as an efficient and viable method. Its ease of implementation,
theoretical guarantees for convergence and practical performance make it
suitable for addressing hard optimization problems with nonlinear equality
constraints.

##  6 Conclusions

In this paper, we have proposed a linearized quadratic penalty method for
solving (locally) smooth optimization problems with nonlinear equality
constraints. In this method we have linearized the cost function and
functional constraints within the quadratic penalty function and added a
regularization term. By dynamically generating the regularization (proximal)
parameter, we have derived convergence rate to an $\epsilon$ first-order
optimal solution and improved convergence rates under the KL condition.
Moreover, we have numerically shown that the proposed algorithm is efficient,
comparing it with three known algorithms, SCP, IPOPT and Inexact Restoration.

Furthermore, a potential avenue for extending this work involves the
incorporation of full second-order derivatives of both the objective function
and constraints. We believe that the analysis of such a scheme is viable,
following a similar rationale as presented in this paper, but the subproblem
may become nonconvex in this case. Finally, it would be interesting to extend
this method to the augmented Lagrangian framework. This extension would
preserve the simplicity of the subproblems and would allow for controlling
feasibility through dual variables rather than the penalty parameter,
potentially leading to smaller values for this parameter.

## Conflict of interest

The authors declare that they have no conflict of interest.

## Data availability

It is not applicable.

## Appendix

Proof of Lemma 2. Using the definition of $\mathcal{P}_{\rho}$ and
$\nabla{\mathcal{P}_{\rho}}$, we have:

$\displaystyle\|\nabla\mathcal{P}_{\rho}(x)-\nabla\mathcal{P}_{\rho}(x^{\prime})\|$
$\displaystyle=$ $\displaystyle\|\nabla f(x)+{J_{F}(x)}^{T}(\rho F(x))-\nabla
f(x^{\prime})-{J_{F}(x^{\prime})}^{T}(\rho F(x^{\prime}))\|$ $\displaystyle=$
$\displaystyle\left\|\nabla f(x)-\nabla
f(x^{\prime})+\Big{(}{J_{F}(x)}-{J_{F}(x^{\prime})}\Big{)}^{T}\Big{(}\rho
F(x)\Big{)}+{J_{F}(x^{\prime})}^{T}\Big{(}\rho\big{(}F(x)-F(x^{\prime})\big{)}\Big{)}\right\|$
$\displaystyle\leq$ $\displaystyle\|\nabla f(x)-\nabla
f(x^{\prime})\|+\|{J_{F}(x)}-{J_{F}(x^{\prime})}\|\|\rho
F(x)\|+\|{J_{F}(x^{\prime})}\|\Big{(}\rho\|F(x)-F(x^{\prime})\|\Big{)}$
$\displaystyle{\overset{\text{Ass. }\ref{assump2}}{\leq}}$
$\displaystyle\Big{(}L_{f}+\rho(L_{F}\|F(x)\|+M_{F}^{2})\Big{)}\|x-x^{\prime}\|\leq
L_{\mathcal{P}_{\rho}}\left\|x-x^{\prime}\right\|,$

where
$L_{\mathcal{P}_{\rho}}=L_{f}+\rho(L_{F}\sup_{x\in\mathcal{S}}\|F(x)\|+M_{F}^{2}).$
This concludes our proof. ∎

Proof of Lemma 3. Using the optimality of $x_{k+1}$, we have:

$\nabla f(x_{k})+\rho
J_{F}^{T}(x_{k})(F(x_{k})+J_{F}(x_{k})(x_{k+1}-x_{k}))+\beta_{k+1}(x_{k+1}-x_{k})=0.$

On the other hand, we have:

$\displaystyle\langle\nabla\mathcal{P}_{\rho}(x_{k}),x_{k+1}-x_{k}\rangle$
$\displaystyle=\langle\nabla f(x_{k})+\rho J_{F}^{T}(x_{k})F(x_{k}),\Delta
x_{k+1}\rangle$ $\displaystyle=-\rho\langle J_{F}^{T}(x_{k})J_{F}(x_{k})\Delta
x_{k+1},\Delta x_{k+1}\rangle-\beta_{k+1}\|\Delta x_{k+1}\|^{2}$
$\displaystyle=-\rho\|J_{F}(x_{k})\Delta x_{k+1}\|^{2}-\beta_{k+1}\|\Delta
x_{k+1}\|^{2}\leq-\beta_{k+1}\|\Delta x_{k+1}\|^{2},$ (30)

where in the second equality we used the previous relation related to
optimality of $x_{k+1}$. By combining the above inequality with (6), we get:

$\mathcal{P}_{\rho}(x_{k+1})-\mathcal{P}_{\rho}(x_{k}){\overset{{\eqref{eq_assu}}}{\leq}}\frac{\beta_{k+1}}{2}\|\Delta
x_{k+1}\|^{2}+\langle\nabla\mathcal{P}_{\rho}(x_{k}),\Delta
x_{k+1}\rangle{\overset{{\eqref{used_just_after}}}{\leq}}-\frac{\beta_{k+1}}{2}\|\Delta
x_{k+1}\|^{2}.$

This proves our statement.∎

Proof of Lemma 4. First, we prove, $P_{k}\leq\alpha\;\forall k\geq 0$, using
induction arguments. Since $x_{0},x_{1}\in\mathcal{S}$, then from Lemma 3 for
$k=0$, it follows that:

$\displaystyle
f(x_{1})+\frac{\rho}{2}\|F(x_{1})\|^{2}+\frac{\beta_{1}}{2}\|x_{1}-x_{0}\|^{2}{\overset{{}}{\leq}}f(x_{0})+\frac{\rho}{2}\|F(x_{0})\|^{2}{\overset{{\eqref{ine10}}}{\leq}}\bar{\alpha}+c_{0}.$
(31)

Furthermore, exploiting the definition of $P_{k}$ for $k=0$ and $k=1$, we
have:

$\displaystyle P_{1}$
$\displaystyle=\mathcal{P}_{\rho}(x_{1})=f(x_{1})+\frac{\rho}{2}\|F(x_{1})\|^{2}\leq\mathcal{P}_{\rho}(x_{0})=P_{0}{\overset{{\eqref{ine111}}}{\leq}}\bar{\alpha}+c_{0}=\alpha.$
(32)

It then follows that for $k=0$ and $k=1$, $(P_{k}\leq\alpha)$ is verified.
Now, assume that $(P_{k}\leq\alpha)$ holds for some $k\geq 1$ (induction
hypothesis) and we will prove that it continues to hold for $k+1$. Since
$x_{k},x_{k+1}\in\mathcal{S}$ and from Lemma 3, we have:

$P_{k+1}-P_{k}\leq-\frac{\beta_{k+1}}{4}\|\Delta x_{k+1}\|^{2}\leq 0.$

Together with the induction hypothesis, we obtain:

$P_{k+1}\leq P_{k}{\overset{{}}{\leq}}\alpha.$

It remains to prove that $P_{k}\geq\bar{L},\forall k\geq 0$. Using (7), we
have:

$\displaystyle P_{k}=f(x_{k})+\frac{\rho}{2}\|F(x_{k})\|^{2}\geq
f(x_{k})+\frac{\rho_{0}}{2}\|F(x_{k})\|^{2}{\overset{{(\text{Lemma
}\ref{lem1})}}{\geq}}\bar{L}.$

It follows that the sequence $\\{P_{k}\\}_{k\geq 0}$ is bounded from below.
Finally, (13) is proved, which completes our proof.∎

Proof of Lemma 5. Using the optimality condition, we have:

$\nabla f(x_{k})+{J_{F}(x_{k})}^{T}(\rho
l_{F}(x_{k+1};x_{k}))+\beta_{k+1}(x_{k+1}-x_{k})=0.$

It then follows, by exploiting the definition of $\mathcal{P}_{\rho}$ and the
properties of the derivative, that:

$\displaystyle\nabla\mathcal{P}_{\rho}(x_{k+1})=$ $\displaystyle\nabla
f(x_{k+1})+{J_{F}(x_{k+1})}^{T}\big{(}\rho F(x_{k+1})\big{)}$ $\displaystyle=$
$\displaystyle\nabla f(x_{k+1})-\nabla
f(x_{k})+\rho\big{(}{J_{F}(x_{k+1})}-{J_{F}(x_{k})}\big{)}^{T}l_{F}(x_{k+1};x_{k})-\beta_{k+1}\Delta
x_{k+1}$
$\displaystyle+\rho{J_{F}(x_{k+1})}^{T}\big{(}F(x_{k+1})-l_{F}(x_{k+1;x_{k}})\big{)}.$

Using basic properties of the Euclidean norm, we further get:

$\displaystyle\|\nabla\mathcal{P}_{\rho}(x_{k+1})\|\leq$
$\displaystyle\|\nabla f(x_{k+1})-\nabla
f(x_{k})\|+\rho\|F(x_{k})\|\|{J_{F}(x_{k+1})}-{J_{F}(x_{k})}\|+\beta_{k+1}\|\Delta
x_{k+1}\|$
$\displaystyle+\rho\|J_{F}(x_{k})\|\big{(}\|{J_{F}(x_{k+1})}\|+\|{J_{F}(x_{k})}\|\big{)}\|\Delta
x_{k+1}\|$
$\displaystyle+\rho\|{J_{F}(x_{k+1})}\|\|F(x_{k+1})-F(x_{k})\|+\rho\|{J_{F}(x_{k+1})}\|^{2}\|\Delta
x_{k+1}\|$ $\displaystyle{\overset{{\text{Ass. }\ref{assump2}}}{\leq}}$
$\displaystyle\big{(}L_{f}+\rho(L_{F}\|F(x_{k})\|+4M_{F}^{2})+\beta_{k+1}\big{)}\|\Delta
x_{k+1}\|$ $\displaystyle\leq$
$\displaystyle(L_{f}+\rho(L_{F}\Delta_{F}+4M_{F}^{2})+\beta_{k+1})\|\Delta
x_{k+1}\|.$

This proves our claim.∎

Proof of Lemma 6. (i) Since $\\{x_{k}\\}_{k\geq 0}$ is bounded, then there
exists a convergent subsequence $\\{x_{k}\\}_{k\in\mathcal{K}}$ such that
$\lim_{k\in\mathcal{K}}{x_{k}}=x_{\rho}^{*}$. Hence $\Omega$ is nonempty.
Moreover, $\Omega$ is compact since it is bounded and closed. On the other
hand, for any $x_{\rho}^{*}\in\Omega$, there exists a sequence of increasing
integers $\mathcal{K}$ such that $\lim_{k\in\mathcal{K}}{x_{k}}=x_{\rho}^{*}$
and using Lemma 5 and (16), it follows that:

$\|\nabla\mathcal{P}_{\rho}(x_{\rho}^{*})\|=\lim_{k\in\mathcal{K}}{\|\nabla\mathcal{P}_{\rho}(x_{k})\|}=0.$

Hence, $x_{\rho}^{*}\in\text{crit }\mathcal{P}_{\rho}$ and
$0\leq\lim_{k\to\infty}{\text{dist}(x_{k},\Omega)}\leq\lim_{k\in\mathcal{K}}{\text{dist}(x_{k},\Omega)}=\text{dist}(x_{\rho}^{*},\Omega)=0$.
(ii) Since $\mathcal{P}_{\rho}(\cdot)$ is continuous and
$\\{\mathcal{P}_{\rho}(x_{k})=P_{k}\\}_{k\geq 0}$ converges to $P^{*}$, then
any subsequence $\\{\mathcal{P}_{\rho}(x_{k})=P_{k}\\}_{k\in\mathcal{K}}$ that
converges, converge to the same limit $P^{*}$. This concludes our proof.∎

Proof of Lemma 7. From Lemma 3 and (17), we have:

$\displaystyle P_{k+1}-P_{k}$ $\displaystyle{\overset{{\text{Lemma
}\ref{lemma3}}}{\leq}}-\frac{\underaccent{\bar}{\beta}}{2}\|\Delta
x_{k+1}\|^{2}.$ (33)

Since $P_{k}\to P^{*}$ and $\\{P_{k}\\}_{k\geq 0}$ is monotonically decreasing
to $P^{*}$, then it follows that the error sequence
$\\{\mathcal{E}_{k}\\}_{k\geq 0}$, is nonnegative, monotonically decreasing
and converges to $0$. We distinguish two cases.

Case 1: There exists $k_{1}\geq 1$ such that $\mathcal{E}_{k_{1}}=0$. Then,
$\mathcal{E}_{k}=0\;\forall k\geq k_{1}$ and it follows that:

$\|x_{k+1}-x_{k}\|^{2}\leq\frac{2}{\underaccent{\bar}{\beta}}(\mathcal{E}_{k}-\mathcal{E}_{k+1})=0\hskip
8.5359pt\forall k\geq k_{1}.$

From Lemma 4 the sequence $\\{x_{k}\\}_{k\geq 0}$ is bounded, and thus:

$\sum_{k=1}^{\infty}{\|\Delta x_{k}\|}=\sum_{k=1}^{k_{1}}{\|\Delta
x_{k}\|}{\overset{{}}{<}}\infty.$

Case 2: The error $\mathcal{E}_{k}>0\;\forall k\geq 0$. Then, there exists
$k_{1}=k_{1}(\epsilon,\tau)\geq 1$ such that $\forall k\geq k_{1}$ we have
$\text{dist}(x_{k},\Omega)\leq\epsilon$,
$P^{*}<\mathcal{P}_{\rho}(x_{k})<P^{*}+\tau$ and

$\varphi^{\prime}(\mathcal{E}_{k})\|\nabla\mathcal{P}_{\rho}(x_{k})\|\geq 1,$
(34)

where $\epsilon>0,\tau>0$ and $\varphi\in\Psi_{\tau}$ are well defined and
correspond to those in Definition 2, recall that $\mathcal{P}_{\rho}(\cdot)$
satisfies the KL property on $\Omega$. Since $\varphi$ is concave, we have
$\varphi(\mathcal{E}_{k})-\varphi(\mathcal{E}_{k+1})\geq\varphi^{\prime}(\mathcal{E}_{k})(\mathcal{E}_{k}-\mathcal{E}_{k+1})$.
Then, from (33) and (34) we get:

$\displaystyle\|x_{k+1}-x_{k}\|^{2}$
$\displaystyle\leq\varphi^{\prime}(\mathcal{E}_{k})\|x_{k+1}-x_{k}\|^{2}\|\nabla\mathcal{P}_{\rho}(x_{k})\|$
$\displaystyle\leq\frac{2}{\underaccent{\bar}{\beta}}\varphi^{\prime}(\mathcal{E}_{k})(\mathcal{E}_{k}-\mathcal{E}_{k+1})\|\nabla\mathcal{P}_{\rho}(x_{k})\|\leq\frac{2}{\underaccent{\bar}{\beta}}\Big{(}\varphi(\mathcal{E}_{k})-\varphi(\mathcal{E}_{k+1})\Big{)}\|\nabla\mathcal{P}_{\rho}(x_{k})\|.$

Note that for given $a,b,c\geq 0$, if we have ${a^{2}}\leq 2b\times c$, and
recognizing that we always have $2b\times c\leq(b+c)^{2}$, then
${a^{2}}\leq(b+c)^{2}$, which in turn implies $a\leq b+c$. It follows that for
any $\theta>0$, we have (by taking $a=\|\Delta x_{k+1}\|$,
$b=\frac{\theta}{\underaccent{\bar}{\beta}}\Big{(}\varphi(\mathcal{E}k)-\varphi(\mathcal{E}{k+1})\Big{)}$,
and $c=\frac{1}{\theta}\|\nabla\mathcal{P}_{\rho}(x_{k})\|$):

$\displaystyle{\|\Delta x_{k+1}\|}\leq$
$\displaystyle\frac{\theta}{\underaccent{\bar}{\beta}}\Big{(}\varphi(\mathcal{E}_{k})-\varphi(\mathcal{E}_{k+1})\Big{)}+\frac{1}{\theta}\|\nabla\mathcal{P}_{\rho}(x_{k})\|.$
(35)

Furthermore, we have:
$\|\nabla{P}(x_{k})\|{\overset{{}}{\leq}}\Gamma_{\text{max}}\|\Delta x_{k}\|.$
Then, (35) becomes:

$\displaystyle{\|\Delta x_{k+1}\|}\leq$
$\displaystyle\frac{\theta}{\underaccent{\bar}{\beta}}\Big{(}\varphi(\mathcal{E}_{k})-\varphi(\mathcal{E}_{k+1})\Big{)}+\frac{\Gamma_{\text{max}}}{\theta}\|\Delta
x_{k}\|.$

Let us now choose $\theta>0$ so that $0<\frac{\Gamma_{\text{max}}}{\theta}<1$
and define the parameter $\delta_{0}$ as:
$\delta_{0}=1-\frac{\Gamma_{\text{max}}}{\theta}>0$. Then, by summing up the
above inequality from $k=\underaccent{\bar}{k}$ to $k=K$ and using the
property: $\sum_{k=\underaccent{\bar}{k}}^{K}{\|\Delta
x_{k}\|}=\sum_{k=\underaccent{\bar}{k}}^{K}{\|\Delta x_{k+1}\|}+\|\Delta
x_{\underaccent{\bar}{k}}\|-\|\Delta x_{{K+1}}\|$, we get:

$\displaystyle\sum_{k=\underaccent{\bar}{k}}^{K}$ $\displaystyle\|\Delta
x_{k+1}\|\leq\frac{\theta}{\underaccent{\bar}{\beta}\delta_{0}}\Big{(}\varphi(\mathcal{E}_{\underaccent{\bar}{k}})-\varphi(\mathcal{E}_{K+1})\Big{)}+\frac{\Gamma_{\text{max}}}{\theta\delta_{0}}\|\Delta
x_{\underaccent{\bar}{k}}\|-\frac{\Gamma_{\text{max}}}{\theta\delta_{0}}\|\Delta
x_{{K+1}}\|.$

Using the fact that $\\{\mathcal{E}_{k}\\}_{k\geq k_{1}}$ is monotonically
decreasing and that the function $\varphi$ is positive and increasing, which
yields $\varphi(\mathcal{E}_{k})\geq\varphi(\mathcal{E}_{k+1})>0$, then:

$\displaystyle\sum_{k=\underaccent{\bar}{k}}^{K}{\|\Delta x_{k+1}\|}\leq$
$\displaystyle\frac{\theta}{\underaccent{\bar}{\beta}\delta_{0}}\varphi(\mathcal{E}_{\underaccent{\bar}{k}})+\frac{\Gamma_{\text{max}}}{\theta\delta_{0}}\|\Delta
x_{\underaccent{\bar}{k}}\|.$

It is clear that the right-hand side of the above inequality is bounded for
any $K\geq\underaccent{\bar}{k}$. Letting $K\to\infty$, we get that:
$\sum_{k=\underaccent{\bar}{k}}^{\infty}{\|\Delta x_{k+1}\|}<\infty.$ From
Lemma 4, the sequence $\\{x_{k}\\}_{k\geq 0}$ is bounded. Then, it follows
that: $\sum_{k=1}^{\underaccent{\bar}{k}}{\|\Delta x_{k}\|}<\infty.$ Hence:
$\sum_{k=1}^{\infty}{\|\Delta x_{k}\|}<\infty$. Let $m,n\in\mathcal{Z}_{+}$
such that $n\geq m$, we have:

$\displaystyle\|x_{n}-x_{m}\|=\|\sum_{k=m}^{n-1}{\Delta
x_{k+1}}\|\leq\sum_{k=m}^{n-1}{\|\Delta x_{k+1}\|}.$

Since $\sum_{k=0}^{\infty}{\|\Delta x_{k+1}\|}<\infty$, it follows that
$\forall\varepsilon>0,\exists N\in\mathcal{Z}_{+}$ such that $\forall m,n\geq
N$ where $n\geq m$, we have: $\|x_{n}-x_{m}\|\leq\varepsilon$. This implies
that $\\{x_{k}\\}_{k\geq 0}$ is a Cauchy sequence and converges. Moreover, by
Theorem 4.1, $\\{x_{k}\\}_{k\geq 0}$ converges to a critical point of
$\mathcal{P}_{\rho}(\cdot)$. This concludes our proof.∎

## References

  * (1) H. Attouch, J. Bolte and B. Svaiter, _Convergence of descent methods for semi-algebraic and tame problems: proximal algorithms, forward-backward splitting, and regularized Gauss–Seidel methods_ , Mathematical Programming, 137, 91–129, 2013.
  * (2) D.P. Bertsekas _On penalty and multiplier methods for constrained minimization_ , SIAM Journal on Control and Optimization, 14, 216–235, 1976.
  * (3) D.P. Bertsekas and J.N. Tsitsiklis, _Parallel and distributed computation: numerical methods_ , Athena Scientific, 2003.
  * (4) E.G. Birgin, J.L. Gardenghi, J.M. Martínez, S.A. Santos and Ph.L. Toint, _Evaluation complexity for nonlinear constrained optimization using unscaled KKT conditions and high-order models_ , SIAM Journal on Optimization 26, 951-967, 2016.
  * (5) E.G. Birgin and J.M. Martínez, _Practical augmented Lagrangian methods for constrained optimization_ , SIAM, Philadelphia, 2014.
  * (6) E.G. Birgin and J.M. Martínez, _Complexity and performance of an augmented Lagrangian algorithm_ , Optimization Methods and Software, 35(5), 885-920, 2020.
  * (7) S. Boyd, N. Parikh, E. Chu, B. Peleato and J. 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) L.F. Bueno and J.M. Martínez, _On the complexity of an inexact restoration method for constrained optimization_ , SIAM Journal on Optimization, 30, 80-101, 2020.
  * (9) C. Cartis, N. Gould and P. Toint, _On the evaluation complexity of composite function minimization with applications to nonconvex nonlinear programming_ , SIAM Journal on Optimization, 21, 1721-1739, 2011.
  * (10) E. Cohen, N. Hallak and M. Teboulle, _A Dynamic alternating direction of multipliers for nonconvex minimization with nonlinear functional equality constraints_ , Journal of Optimization Theory and Applications, 193, 324–353, 2022.
  * (11) R. Courant, _Variational methods for the solution of problems of equilibrium and vibration_ , Bulletin of the American Mathematical Society, 49, 1–23, 1943.
  * (12) J. Fan, _Comments on wavelets in statistics: a review by A. Antoniadis_ , Journal of the Italian Statistical Society, 6, 131–138, 1997.
  * (13) J.A. Fessler, _Optimization methods for magnetic resonance image reconstruction, key models and optimization algorithms_ , IEEE Signal Processing Magazine, 37(1), 33–40, 2020.
  * (14) A.V. Fiacco and G. McCormick, _Nonlinear programming: sequential unconstrained minimization techniques_ , Wiley, 1968.
  * (15) R. Fletcher, _Practical methods of optimization_ , 2nd edition, Wiley, 1987.
  * (16) N. Gould, D. Orban and P. Toint, _CUTEst: a constrained and unconstrained testing environment with safe threads for mathematical optimization_ , Computational Optimization and Applications, 60, 545–557, 2015.
  * (17) M. Grant and S. Boyd, _CVX: Matlab software for disciplined convex programming_ , version 2.1, 2014.
  * (18) M. Hong, D. Hajinezhad and M.M. Zhao, _Prox-PDA: The proximal primal-dual algorithm for fast distributed nonconvex optimization and learning over networks_ , International Conference on Machine Learning, 70, 1529–1538, 2017.
  * (19) A.F. Izmailov and M.V. Solodov, _Convergence rate estimates for penalty methods revisited_ , Computational Optimization and Applications, 85, 973–992, 2023.
  * (20) W. Kong, J. Melo and R. Monteiro, _Complexity of a quadratic penalty accelerated inexact proximal point method for solving linearly constrained nonconvex composite programs_ , SIAM Journal on Optimization, 29(4), 2566–2593, 2018.
  * (21) N. Krejić, J.M. Martínez, M. Mello and E. Pilotta, _Validation of an augmented Lagrangian algorithm with a Gauss-Newton hessian approximation using a set of hard-spheres problems_ , Computational Optimization and Applications, 16. 247-263, 2000
  * (22) Q. Lin, R. Ma and Y. Xu, _Complexity of an inexact proximal-point penalty method for constrained smooth non-convex optimization_ , Computational Optimization and Applications, 82, 175–224, 2022.
  * (23) D.R. Luke, S. Sabach and M. Teboulle, _Optimization on spheres: models and proximal algorithms with computational performance comparisons_ , SIAM Journal on Mathematics of Data Science, 1(3), 408–445, 2019.
  * (24) M.F. Sahin, A. Eftekhari, A. Alacaoglu, F. Latorre and V. Cevher, _An inexact augmented Lagrangian framework for nonconvex optimization with nonlinear constraints_ , Neural Information Processing Systems, 2019.
  * (25) F. Messerer, K. Baumgärtner and M. Diehl, _Survey of sequential convex programming and generalized Gauss-Newton methods_ , ESAIM: Proceedings and Surveys, 2021.
  * (26) I. Necoara and S. Kvamme, _DuQuad: A toolbox for solving convex quadratic programs using dual (augmented) first order algorithms_ , Conference on Decision and Control, 2043-2048, 2015.
  * (27) J. Nocedal and S. J. Wright, _Numerical Optimization_ , Springer, 2006.
  * (28) B.T. Polyak, _Introduction to optimization_ , Optimization Software Inc., Publications Division, 1987.
  * (29) B.T. Polyak and N. V. Tret’yakov, _The method of penalty estimates for conditional extremum problems_ , USSR Computational Mathematics and Mathematical Physics, 13, 42–58, 1973.
  * (30) J.O. Royset, _Variational analysis in modern statistics_ , Mathematical Programming, 174, 2019.
  * (31) Q. Tran-Dinh and M. Diehl, _Local convergence of sequential convex programming for nonconvex optimization_ , In M. Diehl, F. Glineur, E. Jarlebring, and W. Michiels, editors, _Recent advances in optimization and its application in engineering_ , 93–103, Springer, 2010.
  * (32) A. Wächter and L.T. Biegler, _On the implementation of a primal-dual interior point filter line search algorithm for large-scale nonlinear programming_ , Mathematical Programming 106(1), 25-57, 2006.
  * (33) Y. Xie and S.J. Wright, _Complexity of proximal augmented Lagrangian for nonconvex optimization with nonlinear equality constraints_ , Journal of Scientific Computing, 86, 2021.
  * (34) M. Yashtini, _Convergence and rate analysis of a proximal linearized ADMM for nonconvex nonsmooth optimization_ , Journal of Global Optimization, 84, 913–939, 2022.
  * (35) C. Zhang, _Nearly unbiased variable selection under minimax concave penalty_ , Annals of Statistics, 38, 894–942, 2010.