# An abstract convergence framework with application to inertial inexact
forward–backward methods

Silvia Bonettini  Dipartimento di Scienze Fisiche, Informatiche e Matematiche,
Università di Modena e Reggio Emilia, Via Campi 213/b, Modena, 41125, Italy
<EMAIL_ADDRESS><EMAIL_ADDRESS>Peter Ochs  Department of
Mathematics, University of Tübingen, Auf der Morgenstelle 10, Tübingen, 72076,
Germany
<EMAIL_ADDRESS>Marco Prato  Dipartimento di Scienze Fisiche,
Informatiche e Matematiche, Università di Modena e Reggio Emilia, Via Campi
213/b, Modena, 41125, Italy
<EMAIL_ADDRESS><EMAIL_ADDRESS>Simone Rebegoldi
Dipartimento di Ingegneria Industriale, Università di Firenze, Via di S. Marta
3, Firenze, 50139, Italy
<EMAIL_ADDRESS>
###### Abstract

In this paper we introduce a novel abstract descent scheme suited for the
minimization of proper and lower semicontinuous functions. The proposed
abstract scheme generalizes a set of properties that are crucial for the
convergence of several first-order methods designed for nonsmooth nonconvex
optimization problems. Such properties guarantee the convergence of the full
sequence of iterates to a stationary point, if the objective function
satisfies the Kurdyka–Łojasiewicz property. The abstract framework allows for
the design of new algorithms. We propose two inertial-type algorithms with
implementable inexactness criteria for the main iteration update step. The
first algorithm, i2Piano, exploits large steps by adjusting a local Lipschitz
constant. The second algorithm, iPila, overcomes the main drawback of line-
search based methods by enforcing a descent only on a merit function instead
of the objective function. Both algorithms have the potential to escape local
minimizers (or stationary points) by leveraging the inertial feature.
Moreover, they are proved to enjoy the full convergence guarantees of the
abstract descent scheme, which is the best we can expect in such a general
nonsmooth nonconvex optimization setup using first-order methods. The
efficiency of the proposed algorithms is demonstrated on two exemplary image
deblurring problems, where we can appreciate the benefits of performing a
linesearch along the descent direction inside an inertial scheme.

##  1 Introduction

The design of efficient first-order methods is vital for tackling composite
optimization problems of the form

$\min_{x\in\mathbb{R}^{n}}f(x),\quad f(x)=f_{0}(x)+f_{1}(x),$ (1)

where $f_{1}$ is convex and $f_{0}$ is continuously differentiable on an open
set containing the domain of $f_{1}$. Such problems are frequently encountered
in image processing and machine learning applications [10, 24, 26], where one
of the two terms is usually a data fidelity term and the other one encodes
some apriori information on the ground truth [10]. Popular and effective
first-order methods aiming at solving (1) include forward–backward (FB)
methods [31, 32, 40, 18], whose structure consists in the alternation of a
gradient step on $f_{0}$ followed by a proximal minimization step on $f_{1}$,
block coordinate methods [15, 20, 30, 34], Douglas-Rachford methods [31, 36]
and several others.

In recent years, the convergence of first-order descent methods in nonconvex
settings has been carefully addressed by relying on the so-called
Kurdyka–Łojasiewicz (KL) inequality [4, 14, 35]. This analytical property is
satisfied by a large number of objective functions arising in signal
processing and machine learning, such as real analytic or semialgebraic
functions (see e.g. [14, 12]), and, more generally, functions that are
definable in an $o$-minimal structure [33, 14, 13], thus making quite natural
to consider the KL inequality as a standard blanket assumption whenever the
objective function is nonconvex. The convergence of descent methods under the
KL assumption was first considered for gradient related methods in [1] and for
proximal methods in [2, 3], where the authors combine the KL inequality with
some crucial properties of descent methods to prove the convergence of the
iterates to a stationary point of the objective function, under some
boundedness assumption. In [4], the authors extend this seminal idea by
providing the first abstract descent scheme in the KL framework, namely a set
of abstract properties ensuring the convergence of a generic iterative scheme
to a stationary point if combined with the KL inequality. Specific algorithms
are then derived from the abstract scheme.
This approach, consisting in analyzing and devising new and existing
algorithms from an abstract scheme, has been successively adopted by several
authors [15, 34, 16, 40, 39, 21, 22]. In particular, in [40, 39], the authors
modify the abstract scheme proposed in [4] in order to include an inertial
term inside a classical FB scheme and devise the so-called iPiano (inertial
Proximal algorithm for nonconvex optimization) which can be considered as a
generalization of the Heavy-Ball method [42, 41]. Similarly, the authors in
[21, 22] adapt the abstract descent scheme to their proposed method VMILA
(Variable Metric Inexact Linesearch Algorithm); yet, unlike in [40], the
abstract scheme is designed to include an implementable inexactness criterion
for the computation of the proximal point.

The aim of this paper is to provide new theoretical and algorithmic tools able
for effectively solving problem (1) in very general settings, i.e., when
$f_{0}$ is not necessarily convex and/or the proximity operator of $f_{1}$ is
not necessarily available in closed form. These difficulties are very common
in the framework of inverse problems, for example when the acquisition model
is nonlinear or needs to be blindly estimated [5], leading to a nonconvex data
fidelity term, and/or the regularization term involves Total-Variation-like
functionals or sparsity in nontrivial spaces, so that its proximal operator
cannot be computed in closed form or is too costly to compute [6, 18, 23, 50].

More specifically, we give two different types of contributions: (i) From the
theoretical point of view, we propose a novel abstract descent scheme for
proving the convergence of iterative methods under the KL assumption, which
includes as special cases the abstract schemes in [39] and in [22], which are
currently considered as state-of-the-art. Notably, our proposed abstract
scheme provides the theoretical foundation for FB algorithms that include both
an inertial term in the iteration rule and an inexact computation of the
proximity operator. (ii) In the second part of the paper, we materialize this
new theoretical framework into two novel inexact inertial FB algorithms, which
improve the modelling flexibility and the practical performance, while
providing the same convergence guarantees, compared to other FB-type
algorithms. The first algorithm, denominated i2Piano (inertial inexact
Proximal algorithm for nonconvex optimization), can be considered as an
inexact version of the proximal Heavy-Ball algorithm [40, 39] equipped with a
backtracking procedure based on a local version of the Descent Lemma. The
second one is denoted iPila (inertial Proximal inexact line–search algorithm),
which features an inertial–like step followed by a linesearch procedure along
the descent direction of a suitable merit function. Such a linesearch strategy
allows to compute the inexact proximal point only once per iteration, unlike
the backtracking procedure of i2Piano.

The paper is organized as follows. In Section 2 some basic notions on
variational analysis and the definition of the KL property are reported. In
Section 3 the proposed abstract scheme is presented and its convergence
properties analysed. Section 4 is devoted to the design of the two algorithms
i2Piano and iPila. Their analysis within the abstract framework presented in
Section 3 is performed in Section A. Finally, two numerical tests on image
deblurring problems are reported in Section 5.

##  2 Preliminaries

In the remainder of the paper, we denote with
$\overline{\mathbb{R}}=\mathbb{R}\cup\\{-\infty,+\infty\\}$ the extended real
numbers set and $\mathbb{R}^{n\times n}$ the set of $n\times n$ real-valued
matrices, while $\|\cdot\|$ denotes the Euclidean norm. Given a function
${\mathcal{F}}:\mathbb{R}^{n}\to\overline{\mathbb{R}}$ and denoting with
$\operatorname{dom}({\mathcal{F}})=\\{x\in\mathbb{R}^{n}:\
{\mathcal{F}}(x)<+\infty\\}$ the domain of ${\mathcal{F}}$, we say that
${\mathcal{F}}$ is proper if $\operatorname{dom}({\mathcal{F}})\neq\emptyset$
and ${\mathcal{F}}$ is finite on $\operatorname{dom}({\mathcal{F}})$. The
distance operator of a point $x\in\mathbb{R}^{n}$ to a set
$\Omega\subset\mathbb{R}^{n}$ is defined as

$\mathrm{dist}(x,\Omega)=\inf_{y\in\Omega}\|x-y\|.$

Observe that, if $\Omega=\Omega_{1}\times\Omega_{2}$, where
$\Omega_{1}\subset\mathbb{R}^{n_{1}}$, $\Omega_{2}\subset\mathbb{R}^{n_{2}}$,
with $n_{1}+n_{2}=n$, then for all $x=(x_{1},x_{2})\in\mathbb{R}^{n}$, with
$x_{1}\in\mathbb{R}^{n_{1}}$, $x_{2}\in\mathbb{R}^{n_{2}}$, we have

$\mathrm{dist}(x,\Omega)=\sqrt{\mathrm{dist}(x_{1},\Omega_{1})^{2}+\mathrm{dist}(x_{2},\Omega_{2})^{2}}.$
(2)

This follows by observing that
$\|x-z\|^{2}=\|x_{1}-z_{1}\|^{2}+\|x_{2}-z_{2}\|^{2}$ for all
$z=(z_{1},z_{2})\in\Omega$.

######  Definition 1.

Given a proper, lower semicontinuous function
${\mathcal{F}}:\mathbb{R}^{n}\rightarrow\overline{\mathbb{R}}$, the _Fréchet
subdifferential_ of ${\mathcal{F}}$ at
$\overline{z}\in\mathrm{dom}({\mathcal{F}})$ is defined as the set [45,
Definition 8.3(a)]

$\hat{\partial}{\mathcal{F}}(\overline{z})=\left\\{w\in\mathbb{R}^{n}:\liminf_{u\to\overline{z},u\neq\overline{z}}\frac{{\mathcal{F}}(u)-{\mathcal{F}}(\overline{z})-(u-\overline{z})^{T}w}{\|u-\overline{z}\|}\geq
0\right\\}.$

Furthermore, the _limiting subdifferential_ of $\mathcal{F}$ at $\overline{z}$
is given by [45, Definition 8.3(b)]

$\partial\mathcal{F}(\overline{z})=\\{w\in\mathbb{R}^{n}:\exists\
z^{(k)}\to\overline{z},\ \mathcal{F}(z^{(k)})\to{\mathcal{F}}(\overline{z}),\
w^{(k)}\in\hat{\partial}{\mathcal{F}}(z^{(k)})\to w\ {\text{as}}\
k\to\infty\\}.$

######  Definition 2.

Let ${\mathcal{F}}:\mathbb{R}^{n}\rightarrow\overline{\mathbb{R}}$ be a
proper, lower semicontinuous function. A point $\overline{z}\in\mathbb{R}^{n}$
is stationary for ${\mathcal{F}}$ if
$0\in\partial{\mathcal{F}}(\overline{z})$.

Let us introduce the lazy slope of $\mathcal{F}$ at $\overline{z}$, which is
given by [34, p. 877]

$\|\partial{\mathcal{F}}(\overline{z})\|_{-}=\inf_{v\in\partial{\mathcal{F}}(\overline{z})}\|v\|.$

It is then easy to prove the following sufficient criterion for establishing
if a point $\overline{z}\in\mathbb{R}^{n}$ is stationary for the function
$\mathcal{F}$.

######  Lemma 1.

[34, Lemma 2.1] Let
${\mathcal{F}}:\mathbb{R}^{n}\rightarrow\overline{\mathbb{R}}$ be a proper,
lower semicontinuous function and $\overline{z}\in\mathbb{R}^{n}$. If there
exists $\\{z^{(k)}\\}_{k\in\mathbb{N}}\subset\mathbb{R}^{n}$ such that
$z^{(k)}\rightarrow\overline{z}$,
${\mathcal{F}}(z^{(k)})\rightarrow{\mathcal{F}}(\overline{z})$ and
$\displaystyle\liminf_{k\rightarrow\infty}\|\partial\mathcal{F}(z^{(k)})\|_{-}=0$,
then $0\in\partial{\mathcal{F}}(\overline{z})$.

######  Definition 3.

Let ${\mathcal{F}}:\mathbb{R}^{n}\longrightarrow\overline{\mathbb{R}}$ be a
proper, lower semicontinuous function. The function ${\mathcal{F}}$ is said to
have the KL property at $\overline{z}\in\mathrm{dom}(\partial{\mathcal{F}})$
if there exist $\upsilon\in(0,+\infty]$, a neighborhood $U$ of $\overline{z}$,
a continuous concave function $\phi:[0,\upsilon)\longrightarrow[0,+\infty)$
with $\phi(0)=0$, $\phi\in C^{1}(0,\upsilon)$, $\phi^{\prime}(s)>0$ for all
$s\in(0,\upsilon)$, such that the following inequality is satisfied

$\phi^{\prime}({\mathcal{F}}(z)-{\mathcal{F}}(\overline{z}))\|\partial{\mathcal{F}}(z)\|_{-}\geq
1$

for all $z\in
U\cap\\{z\in\mathbb{R}^{n}:{\mathcal{F}}(\overline{z})<{\mathcal{F}}(z)<{\mathcal{F}}(\overline{z})+\upsilon\\}$.
If ${\mathcal{F}}$ satisfies the KL property at each point of
$\mathrm{dom}({\mathcal{F}})$, then ${\mathcal{F}}$ is called a KL function.

The KL property is always satisfied around non stationary points [3, Remark
4(b)], whereas it might fail to hold around stationary points. If we assume
$\mathcal{F}$ is continuously differentiable with
$\mathcal{F}(\overline{z})=0$, then the KL property can be rewritten as

$\|\nabla(\phi\circ\mathcal{F})(z)\|\geq 1$ (3)

around the stationary point $\overline{z}$, which means that the values of
$\mathcal{F}$ can be reparametrized by a so-called _desingularization
function_ $\phi$ such that a “singular region”, i.e., a region in which the
gradients are arbitrarily small, turns into a “regular region”, i.e., a
neighbourhood of $\overline{z}$ where the gradients are bounded away from zero
[4].

In [15, Lemma 6], the following _uniformized_ version of the KL property is
introduced, where the KL inequality holds with the same desingularization
function for all points in a suitable neighborhood of a compact set where the
function is constant.

######  Lemma 2.

Let ${\mathcal{F}}:\mathbb{R}^{n}\to\overline{\mathbb{R}}$ be a proper, lower
semicontinuous function and $X\subset\mathbb{R}^{n}$ a compact set. Suppose
that ${\mathcal{F}}$ satisfies the KL property at each point belonging to $X$
and that ${\mathcal{F}}$ is constant over $X$, i.e.,
${\mathcal{F}}(\bar{x})=\bar{\mathcal{F}}\in\mathbb{R}$ for all $\bar{x}\in
X$. Then, there exists $\mu,\upsilon>0$ and a function $\phi$ as in Definition
3 such that

$\phi^{\prime}({\mathcal{F}}(z)-\bar{\mathcal{F}})\|\partial{\mathcal{F}}(z)\|_{-}\geq
1,\quad\forall z\in\bar{B}$ (4)

where the set $\bar{B}$ is defined as

$\bar{B}=\\{z\in\mathbb{R}^{n}:\mathrm{dist}(z,X)<\mu\mbox{ and
}\bar{\mathcal{F}}<{\mathcal{F}}(z)<\bar{\mathcal{F}}+\upsilon\\}.$ (5)

##  3 Abstract algorithm scheme

In the following, we are interested in proving the convergence of an abstract
descent algorithm to a stationary point of a proper, lower semicontinuous
function ${\mathcal{F}}$. Such abstract algorithm is defined through a
specific set of properties that are shared by several first-order methods
designed for nonsmooth nonconvex optimization, including gradient descent
methods [1], forward–backward methods [18, 29, 40] and block coordinate
methods [4, 20, 34]. Similarly to other abstract descent algorithms in the KL
framework, the two main ingredients guaranteeing the convergence of our scheme
are the sufficient decrease condition and the relative error condition, the
latter being related to the minimization subproblem that one has to
(inexactly) solve at each iteration of a first-order method. However, unlike
in previous works in the literature, we require that the relative error
condition is satisfied at a point that might be different from the actual
iterate generated by the method. As we will see in Section 4, this simple
modification allows to circumvent the issue of the actual implementation of
the relative error condition, which seems rather difficult to impose in
practice [21, 38], allowing to include inexact forward-backward methods
equipped with an implementable inexactness criterion for the solution of the
minimization subproblem.

######  Condition 3 (Abstract algorithm scheme).

Let ${\mathcal{F}}:\mathbb{R}^{n}\times\mathbb{R}^{m}\to\overline{\mathbb{R}}$
be a proper, lower semicontinuous function and
$\Phi:\mathbb{R}^{n}\times\mathbb{R}^{q}\to\overline{\mathbb{R}}$ a proper,
lower semicontinuous, bounded from below function. Consider two sequences
$\\{{x}^{(k)}\\}_{k\in\mathbb{N}}$, $\\{u^{(k)}\\}_{{k\in\mathbb{N}}}$ in
$\mathbb{R}^{n}$, a sequence $\\{\rho^{(k)}\\}_{k\in\mathbb{N}}$ in
$\mathbb{R}^{m}$, a sequence $\\{s^{(k)}\\}_{{k\in\mathbb{N}}}$ in
$\mathbb{R}^{q}$ and a sequence of nonnegative real numbers
$\\{d_{k}\\}_{k\in\mathbb{N}}$ such that the following relations are
satisfied.

  1. [H1]

There exists a sequence of positive real numbers
$\\{a_{k}\\}_{k\in\mathbb{N}}$ such that

$\Phi({x}^{(k+1)},s^{(k+1)})+a_{k}d_{k}^{2}\leq\Phi({x}^{(k)},s^{(k)}),\quad\forall\
k\geq 0.$

  2. [H2]

There exists a sequence of nonnegative real numbers
$\\{r_{k}\\}_{k\in\mathbb{N}}$ with $\lim\limits_{k\to\infty}r_{k}=0$ such
that

$\Phi({x}^{(k+1)},s^{(k+1)})\leq{\mathcal{F}}(u^{(k)},\rho^{(k)})\leq\Phi({x}^{(k)},s^{(k)})+r_{k},\quad\forall\
k\geq 0.$

  3. [H3]

There exist $b>0$, a sequence of positive real numbers
$\\{b_{k}\\}_{k\in\mathbb{N}}$, a summable sequence of nonnegative real
numbers $\\{\zeta_{k}\\}_{k\in\mathbb{N}}$, a non-empty finite index set
$I\subset\mathbb{Z}$ and $\theta_{i}\geq 0$, $i\in I$ with $\sum_{i\in
I}\theta_{i}=1$ such that, setting $d_{j}=0$ for $j\leq 0$, we have

$b_{k+1}\|\partial{\mathcal{F}}(u^{(k)},\rho^{(k)})\|_{-}\leq b\sum_{i\in
I}\theta_{i}d_{k+1-i}+\zeta_{k+1},\quad\forall\ k\geq 0.$

  4. [H4]

If $\\{(x^{(k_{j})},\rho^{(k_{j})})\\}_{j\in\mathbb{N}}$ is a subsequence of
$\\{({x}^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$ converging to some
$(x^{*},\rho^{*})\in\mathbb{R}^{n}\times\mathbb{R}^{m}$, then we have for
$\\{u^{(k_{j})}\\}_{j\in\mathbb{N}}$:

$\lim_{j\to\infty}\|u^{(k_{j})}-x^{(k_{j})}\|=0,\quad\lim_{j\to\infty}{\mathcal{F}}(u^{(k_{j})},\rho^{(k_{j})})={\mathcal{F}}(x^{*},\rho^{*}).$

  5. [H5]

There exists a positive real number $p>0$ and $k^{\prime}\in\mathbb{Z}$ such
that

$\|{x}^{(k+1)}-{x}^{(k)}\|\leq pd_{k+k^{\prime}},\quad\forall\ k\geq 0.$

  6. [H6]

The sequences $\\{a_{k}\\}_{k\in\mathbb{N}}$, $\\{b_{k}\\}_{k\in\mathbb{N}}$
satisfy the following conditions

$\sum_{k=0}^{+\infty}b_{k}=+\infty,\quad\sup_{k\in\mathbb{N}}\frac{1}{b_{k}a_{k}}<+\infty,\quad\inf_{k\in\mathbb{N}}a_{k}>0.$

The abstract scheme given in Conditions 3 can be seen as a further extension
of the one proposed in [39], which is indeed recovered by setting
$\Phi={\mathcal{F}}$, $u^{(k)}={x}^{(k+1)}$ and $\rho^{(k)}\equiv s^{(k)}$. In
the following, we discuss in detail conditions [H1]-[H6] and their relation
with the abstract scheme in [39].

  * •

Condition [H1] is the sufficient decrease condition, which is imposed on the
proper, lower semicontinuous function $\Phi$. The quantity $a_{k}d_{k}^{2}$
measures the amount of the decrease. Note that, in earlier works based on the
KL property [4, 15, 18, 34], condition [H1] is usually presented by setting
$d_{k}=\|x^{(k+1)}-{x}^{(k)}\|_{2}$, $s^{(k)}\equiv 0$ and $\Phi(x,s)=f(x)$,
being $f$ the function to minimize. The generalized condition reported here is
the same one introduced in the recent work [39], although here it is imposed
on a function $\Phi$ that is different from the function ${\mathcal{F}}$
appearing in condition [H3], whereas in [39] the two functions were identical.

  * •

Condition [H2] is crucial in the convergence proof of Theorem 5. Indeed it
ensures that the sequence
$\\{{\mathcal{F}}(u^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$ converges to a
limit value ${\mathcal{F}}^{*}\in\mathbb{R}$, thus allowing to apply the
uniformized KL property at the point $(u^{(k)},\rho^{(k)})$ for all
sufficiently large $k$ and, furthermore, it enables the combination of [H3]
and [H1] with the KL inequality. Imposing condition [H2] is required only when
$\Phi\neq{\mathcal{F}}$, this is why it does not appear in [39].

  * •

Condition [H3] is the so-called relative error condition, which is related to
the (possibly) inexact solution of the minimization subproblem performed at
each iteration of a first-order method. In the previous literature [4, 15, 18,
34], such condition is usually employed by setting $u^{(k)}=x^{(k+1)}$,
$d_{k}=\|x^{(k+1)}-{x}^{(k)}\|_{2}$, $\rho^{(k)}\equiv 0$, $I=\\{1\\}$,
$\theta_{1}=1$ and ${\mathcal{F}}(u,\rho)=f(u)$, being $f$ the function to
minimize. In [39], a general positive term $d_{k}$, a finite index set $I$, a
variable parameter $s^{(k)}$ and a generic merit function ${\mathcal{F}}$ are
employed, while keeping $u^{(k)}=x^{(k+1)}$ and $\rho^{(k)}\equiv 0$. Here we
also allow the sequence $\\{u^{(k)}\\}_{k\in\mathbb{N}}$ to be distinct from
$\\{x^{(k)}\\}_{k\in\mathbb{N}}$ and the parameters
$\\{\rho^{(k)}\\}_{k\in\mathbb{N}}$ to vary at each iteration. The reason to
do so comes from the fact that condition [H3] is hard to enforce
algorithmically on $x^{(k+1)}$ when the minimization subproblem is solved
inexactly, as noted in [21, 22, 38]. However, if a specific, implementable
inexactness criterion is adopted for the solution of the subproblem, then the
same condition holds for a merit function ${\mathcal{F}}$ evaluated at a
different iterate $(u^{(k)},\rho^{(k)})$. For instance, this is observed in
the convergence analysis of VMILA algorithm [18, 21, 22]. In [18], VMILA is
included in the KL framework by setting
$\Phi(x,s)={\mathcal{F}}(x,\rho)=f(x)$, $u^{(k)}=x^{(k+1)}$,
$d_{k}=\|x^{(k+1)}-{x}^{(k)}\|_{2}$, $\rho^{(k)}\equiv s^{(k)}\equiv 0$ and
noting that, in so doing, the relative error condition holds only by exactly
computing the proximal operator. In [21], the authors include VMILA in the
abstract scheme in a different way, by using $\Phi(x,s)=f(x)$, a merit
function ${\mathcal{F}}$ defined upon the concept of forward–backward envelope
of $f$ [49], the iterate $u^{(k)}$ as the inexact proximal-gradient point
$\tilde{y}^{(k)}$, $d_{k}=\|x^{(k+1)}-{x}^{(k)}\|_{2}$ and $\rho^{(k)}$ as the
error parameter due to the computation of $\tilde{y}^{(k)}$. Finally, in [22],
VMILA is framed by setting $\Phi(x,s)=f(x)$,
${\mathcal{F}}(u,\rho)=f(u)+\rho^{2}/2$, the iterate $u^{(k)}$ as the exact
proximal–gradient point $y^{(k)}$ and $d_{k}=-h^{(k)}(\tilde{y}^{(k)})$, where
$h^{(k)}$ is the function to minimize when computing the approximation
$\tilde{y}^{(k)}$ of the exact point $y^{(k)}$. In other words, from the
analysis in [22], it turns out that we are able to enforce the relative error
condition at the exact point $y^{(k)}$, which we do not need to compute
explicitly, provided that the approximation $\tilde{y}^{(k)}$ is computed
using a specific criterion.

  * •

Condition [H4] is the analogue of the so-called continuity condition in [39].
Here we impose the property for all converging subsequences
$(x^{(k_{j})},\rho^{(k_{j})})$, whereas in [39] it is only required the
existence of one such subsequence. This is because, unlike in [39], we need to
ensure that the distance between $\\{(u^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$
and the limit set of $\\{({x}^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$ converges
to $0$ (see Lemma 4(iii)).

  * •

Condition [H5] is also called the distance condition. It states the connection
between the general term $d_{k}$ and the Euclidean norm, which is fundamental
in order to prove the finite length of the sequence
$\\{{x}^{(k)}\\}_{k\in\mathbb{N}}$. Note that the distance condition given in
[39] is slightly more general than [H5]; however, that condition alone allows
to prove only the finite length of the sequence
$\\{d_{k}\\}_{k\in\mathbb{N}}$, which in general does not imply the
convergence of the sequence $\\{{x}^{(k)}\\}_{k\in\mathbb{N}}$. In order to
obtain the strongest result, condition [H5] is then imposed in [39, Theorem
10].

  * •

Condition [H6] is the same as the parameter condition in [39]. These
requirements on the sequences $\\{a_{k}\\}_{k\in\mathbb{N}}$,
$\\{b_{k}\\}_{k\in\mathbb{N}}$ were first introduced in [34] in order to
generalize the abstract descent scheme in [4].

In the remainder of this section, we will denote with
$\\{{x}^{(k)}\\}_{k\in\mathbb{N}}$, $\\{u^{(k)}\\}_{{k\in\mathbb{N}}}$,
$\\{\rho^{(k)}\\}_{k\in\mathbb{N}}$, $\\{s^{(k)}\\}_{k\in\mathbb{N}}$ the
sequences complying with Conditions 3. Furthermore, let us define the set of
all limit points of the sequence
$\\{({x}^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$:

$\Omega^{*}(x^{(0)},\rho^{(0)})=\\{(x^{*},\rho^{*})\in\mathbb{R}^{n}\times\mathbb{R}^{m}:\exists\
\\{k_{j}\\}_{j\in\mathbb{N}}\subset\mathbb{N}\mbox{ such that
}({x}^{(k_{j})},\rho^{(k_{j})})\rightarrow(x^{*},\rho^{*})\\}.$

Note that the set $\Omega^{*}(x^{(0)},\rho^{(0)})$ can be written as

$\Omega^{*}(x^{(0)},\rho^{(0)})=X^{*}(x^{(0)})\times R^{*}(\rho^{(0)})$

where $X^{*}(x^{(0)})=\\{x^{*}\in\mathbb{R}^{n}:\exists\
\\{k_{j}\\}_{j\in\mathbb{N}}\subset\mathbb{N}\mbox{ such that
}{x}^{(k_{j})}\rightarrow x^{*}\\}\subset\mathbb{R}^{n}$,
$R^{*}(\rho^{(0)})=\\{\rho^{*}\in\mathbb{R}^{m}:\exists\
\\{k_{j}\\}_{j\in\mathbb{N}}\subset\mathbb{N}\mbox{ such that
}\rho^{(k_{j})}\rightarrow\rho^{*}\\}\subset\mathbb{R}^{m}$.

######  Lemma 4.

Let Conditions 3 be satisfied. Suppose that
$\\{({x}^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$ is a bounded sequence. Then
the following facts hold true.

  1. (i)

$\Omega^{*}(x^{(0)},\rho^{(0)})$ is nonempty and compact.

  2. (ii)

There exists ${\mathcal{F}}^{*}\in\mathbb{R}$ such that
$\displaystyle\lim_{k\to\infty}\Phi({x}^{(k)},s^{(k)})=\lim_{k\to\infty}{\mathcal{F}}(u^{(k)},\rho^{(k)})={\mathcal{F}}^{*}$.

  3. (iii)

We have

$\lim_{k\to\infty}\mathrm{dist}(({x}^{(k)},\rho^{(k)}),\Omega^{*}({x}^{(0)},\rho^{(0)}))=\lim_{k\to\infty}\mathrm{dist}((u^{(k)},\rho^{(k)}),\Omega^{*}({x}^{(0)},\rho^{(0)}))=0.$

  4. (iv)

We have ${\mathcal{F}}(x^{*},\rho^{*})={\mathcal{F}}^{*}$, $\forall\
(x^{*},\rho^{*})\in\Omega^{*}(x^{(0)},\rho^{(0)})$.

###### Proof.

(i) Since the sequence $\\{({x}^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$ is
bounded, it admits at least a limit point and, hence,
$\Omega^{*}(x^{(0)},\rho^{(0)})$ is nonempty. Compactness is proved in [15,
Lemma 5].
(ii) From [H1] we have that the sequence
$\\{\Phi({x}^{(k)},s^{(k)})\\}_{k\in\mathbb{N}}$ is nonincreasing and, since
$\Phi$ is bounded from below, there exists ${\mathcal{F}}^{*}\in\mathbb{R}$
such that

$\lim_{k\to\infty}\Phi({x}^{(k)},s^{(k)})={\mathcal{F}}^{*}.$

The previous relation combined with [H2] proves Part (ii).
(iii) Since $\\{({x}^{(k)},\rho^{(k)})\\}_{{k\in\mathbb{N}}}$ is bounded and
by definition of $\Omega^{*}(x^{(0)},\rho^{(0)})$, we have

$\lim_{k\to\infty}\mathrm{dist}(({x}^{(k)},\rho^{(k)}),\Omega^{*}({x}^{(0)},\rho^{(0)}))=0.$

Observing that $\Omega^{*}(x^{(0)},\rho^{(0)})=X^{*}({x}^{(0)})\times
R^{*}(\rho^{(0)})$ and recalling (2), the previous limit implies

$\lim_{k\to\infty}\mathrm{dist}({x}^{(k)},X^{*}({x}^{(0)}))=0,\ \ \ \
\lim_{k\to\infty}\mathrm{dist}(\rho^{(k)},R^{*}(\rho^{(0)}))=0.$

Combining the boundedness of $\\{{x}^{(k)}\\}_{{k\in\mathbb{N}}}$ with the
definition of $X^{*}({x}^{(0)})$ and property [H4], we obtain
$\displaystyle\lim_{k\to\infty}\mathrm{dist}(u^{(k)},X^{*}({x}^{(0)}))=0$,
which together with the second limit above yields

$\lim_{k\to\infty}\mathrm{dist}((u^{(k)},\rho^{(k)}),\Omega^{*}({x}^{(0)},\rho^{(0)}))=0.$

(iv) This point follows directly from part (ii) and [H4]. $\square$

######  Theorem 5.

Let Conditions 3 be satisfied, and suppose that
$\\{({x}^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$ is a bounded sequence and that
${\mathcal{F}}$ is a KL function. Then the following statements are true.

  * (i)

The sequence $\\{d_{k}\\}_{k\in\mathbb{N}}$ is summable, i.e., it satisfies

$\sum_{k=0}^{+\infty}d_{k}<+\infty.$

  * (ii)

The sequence $\\{{x}^{(k)}\\}_{k\in\mathbb{N}}$ has finite length, i.e., it
satisfies

$\sum_{k=0}^{+\infty}\|{x}^{(k+1)}-{x}^{(k)}\|<+\infty$

and thus $\\{{x}^{(k)}\\}_{k\in\mathbb{N}}$ is a convergent sequence.

  * (iii)

If also $\\{\rho^{(k)}\\}_{k\in\mathbb{N}}$ converges, then the sequence
$\\{({x}^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$ converges to a stationary
point for ${\mathcal{F}}$.

###### Proof.

(i) By Lemma 4(i)-(iv), the function ${\mathcal{F}}$ is constant over the
compact set $\Omega^{*}(x^{(0)},\rho^{(0)})$, therefore we can apply Lemma 2.
Let $\upsilon,\mu,\phi,\bar{B}$ be as in Lemma 2. Thanks to Lemma 4(ii)-(iii)
and [H2], there exists a positive integer $k_{0}$ such that

$\Phi({x}^{(k)},s^{(k)})+r_{k}<{\mathcal{F}}^{*}+\upsilon,\quad\mathrm{dist}((u^{(k)},\rho^{(k)}),\Omega^{*}(x^{(0)},\rho^{(0)}))<\mu$
(6)

for all $k\geq k_{0}$. Without loss of generality, up to a translation of the
iteration index, we can assume $k_{0}=0$.
Let us now set $c=\sup_{k}1/(a_{k}b_{k})$, where $c<+\infty$ due to [H6],
$\zeta_{k}^{\prime}=\zeta_{k}/b$ and

$\phi_{k}=\frac{b}{c}(\phi(\Phi({x}^{(k)},s^{(k)})-{\mathcal{F}}^{*})-\phi(\Phi({x}^{(k+1)},s^{(k+1)})-{\mathcal{F}}^{*}))$

and prove that

$2d_{k}\leq\phi_{k}+\sum_{i\in
I}\theta_{i}d_{k-i}+\zeta_{k}^{\prime},\quad\forall\ k\geq 1.$ (7)

We first observe that the definition of $\phi_{k}$ is well posed, since
${\mathcal{F}}^{*}$ is the limit of the nonincreasing sequence
$\\{\Phi({x}^{(k)},s^{(k)})\\}_{k\in\mathbb{N}}$ (see Lemma 4) and, hence, we
have $\Phi({x}^{(k)},s^{(k)})\geq{\mathcal{F}}^{*}$, $\forall\
k\in\mathbb{N}$. Moreover, the monotonicity of both function $\phi$ and
sequence $\\{\Phi({x}^{(k)},s^{(k)})\\}_{k\in\mathbb{N}}$ implies that
$\phi_{k}\geq 0$, $\forall\ k\in\mathbb{N}$.
Let us now consider the two cases $d_{k}=0$ and $d_{k}>0$ separately.
If $d_{k}=0$, inequality (7) holds trivially. Otherwise, for any iteration
index $k\geq 1$ such that $d_{k}>0$, taking into account [H1] and [H2], we can
write

${\mathcal{F}}^{*}\leq\Phi({x}^{(k+1)},s^{(k+1)})<\Phi({x}^{(k)},s^{(k)})\leq{\mathcal{F}}(u^{(k-1)},\rho^{(k-1)}).$

On the other hand, the rightmost inequality in [H2] gives

${\mathcal{F}}(u^{(k-1)},\rho^{(k-1)})\leq\Phi({x}^{(k-1)},s^{(k-1)})+r_{k-1},$

which, in view of (6), implies that $(u^{(k-1)},\rho^{(k-1)})\in\bar{B}$.
Then, we can write the KL inequality related to the function ${\mathcal{F}}$
at $(u^{(k-1)},\rho^{(k-1)})$:

$\phi^{\prime}({\mathcal{F}}(u^{(k-1)},\rho^{(k-1)})-{\mathcal{F}}^{*})\geq\frac{1}{\|\partial{\mathcal{F}}(u^{(k-1)},\rho^{(k-1)})\|_{-}}.$

Furthermore, combining the previous inequality with [H3] yields

$\phi^{\prime}({\mathcal{F}}(u^{(k-1)},\rho^{(k-1)})-{\mathcal{F}}^{*})\geq\frac{1}{\frac{b}{b_{k}}\sum\limits_{i\in
I}\theta_{i}d_{k-i}+\frac{1}{b_{k}}\zeta_{k}}.$

Since $\phi$ is concave, $\phi^{\prime}$ is nonincreasing. Therefore, [H2]
implies

$\phi^{\prime}(\Phi({x}^{(k)},s^{(k)})-{\mathcal{F}}^{*})\geq\phi^{\prime}({\mathcal{F}}(u^{(k-1)},\rho^{(k-1)})-{\mathcal{F}}^{*}).$

Exploiting again the concavity of $\phi$, we have

$\displaystyle\phi(\Phi({x}^{(k)},s^{(k)})-{\mathcal{F}}^{*})$
$\displaystyle-\phi(\Phi({x}^{(k+1)},s^{(k+1)})-{\mathcal{F}}^{*})\geq$
$\displaystyle\phi^{\prime}(\Phi({x}^{(k)},s^{(k)})-{\mathcal{F}}^{*})(\Phi({x}^{(k)},s^{(k)})-\Phi({x}^{(k+1)},s^{(k+1)})).$

Combining the last three relations with [H1] leads to

$\displaystyle\phi(\Phi({x}^{(k)},s^{(k)})-{\mathcal{F}}^{*})-\phi(\Phi({x}^{(k+1)},s^{(k+1)})-{\mathcal{F}}^{*})$
$\displaystyle\geq$
$\displaystyle\frac{\Phi({x}^{(k)},s^{(k)})-\Phi({x}^{(k+1)},s^{(k+1)})}{\frac{b}{b_{k}}\sum\limits_{i\in
I}\theta_{i}d_{k-i}+\frac{1}{b_{k}}\zeta_{k}}$ $\displaystyle\geq$
$\displaystyle\frac{a_{k}d_{k}^{2}}{\frac{b}{b_{k}}\sum\limits_{i\in
I}\theta_{i}d_{k-i}+\frac{1}{b_{k}}\zeta_{k}}.$

Recalling the definition of $\phi_{k}$ and $\zeta_{k}^{\prime}$, the above
inequality implies the following one

$d_{k}^{2}\leq\phi_{k}\left(\sum_{i\in
I}\theta_{i}d_{k-i}+\zeta_{k}^{\prime}\right).$

Taking the square root of both sides and using the inequality $2\sqrt{uv}\leq
u+v$ on the right-hand-side, we obtain (7).
Summing (7) from 1 to $k$ leads to

$2\sum_{j=1}^{k}d_{j}\leq\sum_{j=1}^{k}\phi_{j}+\sum_{j=1}^{k}\sum_{i\in
I}\theta_{i}d_{j-i}+\sum_{j=1}^{k}\zeta^{\prime}_{j}.$ (8)

We now observe that

$\displaystyle\sum_{j=1}^{k}\phi_{j}$ $\displaystyle=$
$\displaystyle\frac{b}{c}(\phi(\Phi(x^{(1)},s^{(1)})-{\mathcal{F}}^{*})-\phi(\Phi(x^{(k+1)},s^{(k+1)})-{\mathcal{F}}^{*}))$
$\displaystyle\leq$
$\displaystyle\frac{b}{c}\phi(\Phi(x^{(1)},s^{(1)})-{\mathcal{F}}^{*}),$

where the rightmost inequality follows from the positive sign of $\phi$.
Furthermore, the second sum in the right-hand side of (8) can be rewritten as
below

$\displaystyle\sum_{j=1}^{k}\sum_{i\in I}\theta_{i}d_{j-i}$
$\displaystyle=\sum_{i\in I}\sum_{j=1}^{k}\theta_{i}d_{j-i}=\sum_{i\in
I}\sum_{r=1-i}^{k-i}\theta_{i}d_{r}$ $\displaystyle\leq\sum_{i\in
I}\sum_{r=1-i}^{0}\theta_{i}d_{r}+\left(\sum_{i\in
I}\theta_{i}\right)\sum_{r=1}^{k}d_{r}+\sum_{i\in
I}\sum_{r=k+1}^{k-i}\theta_{i}d_{r}$ $\displaystyle=\sum_{i\in
I}\sum_{r=1-i}^{0}\theta_{i}d_{r}+\sum_{r=1}^{k}d_{r}+\sum_{i\in
I}\sum_{r=k+1}^{k-i}\theta_{i}d_{r}$

where we have used the change of variable $r=j-1$ and the property $\sum_{i\in
I}\theta_{i}=1$. Note that the sums appearing in the previous relation are
assumed to be zero whenever the start index of the summation is larger than
the termination index. Therefore we can write

$2\sum_{j=1}^{k}d_{j}\leq\sum_{i\in
I}\sum_{r=1-i}^{0}\theta_{i}d_{r}+\sum_{r=1}^{k}d_{r}+\sum_{i\in
I}\sum_{r=k+1}^{k-i}\theta_{i}d_{r}+\frac{b}{c}\phi(\Phi(x^{(1)},s^{(1)})-{\mathcal{F}}^{*})+\sum_{j=1}^{k}\zeta^{\prime}_{j}$

which clearly implies

$\sum_{j=1}^{k}d_{j}\leq\sum_{i\in
I}\sum_{r=1-i}^{0}\theta_{i}d_{r}+\sum_{i\in
I}\sum_{r=k+1}^{k-i}\theta_{i}d_{r}+\frac{b}{c}\phi(\Phi(x^{(1)},s^{(1)})-{\mathcal{F}}^{*})+\sum_{j=1}^{k}\zeta^{\prime}_{j}.$
(9)

At this point, observe that the first two sums in the right-hand side of (9)
are finite linear combinations of the terms $\\{d_{r}\\}_{r\in\mathbb{N}}$.
Conditions [H1] and [H6] ensure that $d_{k}\rightarrow 0$, hence those sums
are converging to $0$ for $k\to\infty$. Noting also that
$\\{\zeta^{\prime}_{k}\\}_{k\in\mathbb{N}}$ is summable and taking the limit
of (9) for $k\to\infty$, we obtain

$\sum_{k=0}^{\infty}d_{k}<+\infty.$ (10)

(ii) Combining (10) with [H5], we also obtain

$\sum_{k=0}^{\infty}\|{x}^{(k+1)}-{x}^{(k)}\|\leq
p\sum_{k=0}^{\infty}d_{k+k^{\prime}}<+\infty.$

which implies that the sequence $\\{{x}^{(k)}\\}_{k\in\mathbb{N}}$ converges
to a point $x^{*}\in\mathbb{R}^{n}$.
(iii) Let $(x^{*},\rho^{*})\in\mathbb{R}^{n}\times\mathbb{R}^{m}$ be the
unique limit point of the sequence
$\\{({x}^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$, namely
$(x^{(k)},\rho^{(k)})\rightarrow(x^{*},\rho^{*})$. By using [H4], it follows
that $(u^{(k)},\rho^{(k)})\rightarrow(x^{*},\rho^{*})$ and
${\mathcal{F}}(u^{(k)},\rho^{(k)})\rightarrow{\mathcal{F}}(x^{*},\rho^{*})$.
Furthermore, summing [H3] for $k=0,\ldots,K$ yields

$\sum_{k=0}^{K}b_{k+1}\|\partial{\mathcal{F}}(u^{(k)},\rho^{(k)})\|_{-}\leq
b\sum_{k=0}^{K}\sum_{i\in I}\theta_{i}d_{k+1-i}+\sum_{k=0}^{K}\zeta_{k}.$

Taking the limit for $K\rightarrow\infty$, using (10) and recalling that
$\\{\zeta_{k}\\}_{k\in\mathbb{N}}$ is summable, we obtain

$\sum_{k=0}^{\infty}b_{k+1}\|\partial{\mathcal{F}}(u^{(k)},\rho^{(k)})\|_{-}<+\infty.$

Since [H6] requires $\sum_{k}b_{k}=+\infty$, the previous relation implies
that

$\liminf_{k\rightarrow\infty}\|\partial{\mathcal{F}}(u^{(k)},\rho^{(k)})\|_{-}=0.$

In conclusion, the sequence $\\{(u^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$
satisfies all the hypotheses of Lemma 1, which means that
$0\in\partial\mathcal{F}(x^{*},\rho^{*})$.

$\square$

##  4 Applications of the abstract scheme

In this section, we show how we can devise some brand new
forward–backward–type algorithms satisfying Conditions 3 and, hence, guarantee
their convergence to a stationary point in virtue of Theorem 5. In particular,
from now on, we address the problem

$\min_{x\in\mathbb{R}^{n}}f(x),\ \ \ f(x)=f_{0}(x)+f_{1}(x),$ (11)

where we assume that $f_{0},f_{1}$ are as follows:

  1. [A1]

$f_{1}:\mathbb{R}^{n}\to\overline{\mathbb{R}}$ is a proper, lower
semicontinuous, convex function;

  2. [A2]

$f_{0}:\mathbb{R}^{n}\to\mathbb{R}$ is continuously differentiable on an open
set $\Omega_{0}\supset\overline{\mathrm{dom}(f_{1})}$;

  3. [A3]

$f_{0}$ has $L-$Lipschitz continuous gradient on $\mathrm{dom}(f_{1})$, i.e.,

$\|\nabla f_{0}(x)-\nabla f_{0}(y)\|\leq L\|x-y\|,\ \ \forall
x,y\in\mathrm{dom}(f_{1}),$

for some $L>0$.

  4. [A4]

$f$ is bounded from below.

Under the above assumptions, for any $z\in\mathrm{dom}(f_{1})$, the following
subdifferential calculus rules hold [45, Proposition 8.12, Exercise 8.8(c)]

$\displaystyle\partial f_{1}(z)$
$\displaystyle=\\{w\in\mathbb{R}^{n}:f_{1}(y)\geq f_{1}(z)+\langle
w,y-z\rangle,\ \forall y\in\mathbb{R}^{n}\\}$ $\displaystyle\partial f(z)$
$\displaystyle=\\{\nabla f_{0}(z)\\}+\partial f_{1}(z).$ (12)

One of the most popular strategies for accelerating first order methods
consists in including an inertial, or heavy ball, term to the iteration rule.
This idea was originally proposed in [42] for the gradient descent method, and
consists in combining the gradient direction and the direction obtained from
the last two iterates. This approach has been further developed in the seminal
papers [37, 8] and successively realized in a variety of algorithms. In this
paper, in the general nonsmooth nonconvex setup of (11), we consider the
therefrom inspired inertial proximal-gradient method proposed in [40, 39],
which is defined by the iteration

${x}^{(k+1)}={\mathrm{prox}}_{\alpha_{k}f_{1}}({x}^{(k)}-\alpha_{k}\nabla
f_{0}({x}^{(k)})+\beta_{k}({x}^{(k)}-{x}^{(k-1)})),$ (13)

where $\alpha_{k},\beta_{k}$ are suitably chosen parameters. By definition of
the proximity operator, the above updating rule consists in defining the new
point ${x}^{(k+1)}$ as the unique solution of the minimization problem

$\min_{y\in\mathbb{R}^{n}}f_{1}(y)-f_{1}({x}^{(k)})+\langle\nabla
f_{0}({x}^{(k)})-\frac{\beta_{k}}{\alpha_{k}}({x}^{(k)}-{x}^{(k-1)}),y-{x}^{(k)}\rangle+\frac{1}{2\alpha_{k}}\|y-{x}^{(k)}\|^{2}.$
(14)

We will refer to the minimizer of this problem as the _inertial proximal
gradient point_. If $\beta_{k}=0$, we recover the standard proximal gradient
point, otherwise the inertial step ${x}^{(k)}-{x}^{(k-1)}$ is included in the
argument of the proximal gradient operator, with the aim of improving the
convergence behaviour of the overall method.

However, in several practical situations the exact minimization of (14) is
infeasible, while high precision approximations can be computed efficiently.
Therefore, in the following we address the key challenge of designing
algorithms that inexactly compute the inertial proximal gradient point with
implementable conditions that still preserve the convergence guarantees of
Theorem 5. More precisely, we propose two new inexact inertial–type
algorithms, where the second one also features a linesearch procedure along a
descent direction of a suitable merit function. The convergence analysis of
both algorithms can be performed in the abstract framework provided by
Conditions 3. We stress that, due to the inexactness in the computation of the
inertial proximal gradient point, our algorithms cannot be cast in the
abstract frameworks proposed in previous works.

###  4.1 Inexactness criterion of the inertial proximal gradient point

We start our presentation by defining the inexactness criterion for the
inertial proximal gradient point, which generalizes the one proposed in [17,
18, 21].

Given two positive parameters $\alpha_{k},\beta_{k}$, consider the function
$h^{(k)}:\mathbb{R}^{n}\times\mathbb{R}^{n}\times\mathbb{R}^{n}\to\mathbb{R}\cup\\{+\infty\\}$
defined as follows:

$h^{(k)}(y;x,s)=f_{1}(y)-f_{1}(x)+\langle\nabla
f_{0}(x)-\frac{\beta_{k}}{\alpha_{k}}(x-s),y-x\rangle+\frac{1}{2\alpha_{k}}\|y-x\|^{2}.$
(15)

Clearly, the inertial proximal gradient point (13) is the minimizer of the
above function with respect to the first argument, with $x={x}^{(k)}$,
$s={x}^{(k-1)}$.
Given $({x}^{(k)},s^{(k)})$, we denote by $\hat{y}^{(k)}$ the (exact)
minimizer of the function in (15)

$\hat{y}^{(k)}=\operatornamewithlimits{argmin}\limits_{y\in\mathbb{R}^{n}}h^{(k)}(y;x^{(k)},s^{(k)})={\mathrm{prox}}_{\alpha_{k}f_{1}}(x^{(k)}-\alpha_{k}\nabla
f_{0}(x^{(k)})+\beta_{k}(x^{(k)}-s^{(k)})).$ (16)

The point $\hat{y}^{(k)}$ is the unique point satisfying the optimality
condition

$\displaystyle 0\in\partial h^{(k)}(\hat{y}^{(k)};{x}^{(k)},s^{(k)})$
$\displaystyle\Updownarrow$ (17)
$\displaystyle-\frac{1}{\alpha}_{k}(\hat{y}^{(k)}-x^{(k)}+\alpha_{k}\nabla
f_{0}(x^{(k)})-\beta_{k}(x^{(k)}-s^{(k)}))\in\partial f_{1}(\hat{y}^{(k)}).$

Borrowing the ideas in [17, 22], we define an approximation of $\hat{y}^{(k)}$
as any point $\tilde{y}^{(k)}\in\mathrm{dom}(f_{1})$ such that

$h^{(k)}(\tilde{y}^{(k)};x^{(k)},s^{(k)})-h^{(k)}(\hat{y}^{(k)};x^{(k)},s^{(k)})\leq-\frac{\tau}{2}h^{(k)}(\tilde{y}^{(k)};x^{(k)},s^{(k)}),$
(18)

for a given constant $\tau\geq 0$ not depending on $k$. Since $\hat{y}^{(k)}$
is the unique minimizer of $h^{(k)}(\,\cdot\,;x^{(k)},s^{(k)})$, for $\tau=0$
the previous inequality implies $\tilde{y}^{(k)}=\hat{y}^{(k)}$, therefore the
inexactness occurs for $\tau>0$. The above condition is equivalent to the
following one:

$h^{(k)}(\tilde{y}^{(k)};x^{(k)},s^{(k)})\leq\left(\frac{2}{2+\tau}\right)h^{(k)}(\hat{y}^{(k)};x^{(k)},s^{(k)})\leq
0,$ (19)

where the rightmost inequality is a consequence of the fact that
$\hat{y}^{(k)}$ is a minimizer of $h^{(k)}(\cdot;x^{(k)},s^{(k)})$ and
$h^{(k)}(x^{(k)};x^{(k)},s^{(k)})=0$. Therefore, we have
$h^{(k)}(\tilde{y}^{(k)};x^{(k)},s^{(k)})\leq 0$ and condition (18) can be
rewritten in equivalent way as

$0\in\partial_{\epsilon_{k}}h^{(k)}(\tilde{y}^{(k)};x^{(k)},s^{(k)}),\ \
\mbox{with }\epsilon_{k}=-\frac{\tau}{2}h(\tilde{y}^{(k)};x^{(k)},s^{(k)}),$
(20)

where $\partial_{\epsilon}h(\tilde{y};x,s)=\\{w\in\mathbb{R}^{n}:h(y;x,s)\geq
h(\tilde{y};x,s)+\langle w,y-\tilde{y}\rangle-\epsilon,\ \forall
y\in\mathbb{R}^{n}\\}$ is the $\epsilon$-subdifferential in the first argument
of the convex function $h(\cdot;x,s)$ at point $\tilde{y}$ [51, p. 82].
Therefore, the point $\tilde{y}^{(k)}$ is defined upon a relaxation of the
optimality condition (4.1), where the subdifferential of
$h^{(k)}(\cdot;x^{(k)},s^{(k)})$ is replaced by the
$\epsilon_{k}$-subdifferential and the accuracy parameter $\epsilon_{k}$ is
chosen in a specific way, which is crucial for preserving the theoretical
convergence properties.

Even if the inclusion (20) is implicit, a point
$\tilde{y}^{(k)}\in\mathrm{dom}(f_{1})$ satisfying (20) can be actually
computed in practice in some special cases with a well defined, explicit
primal–dual procedure, as explained in the following section.

####  4.1.1 Computation of the inexact inertial proximal gradient point

In this section we describe a procedure for computing a point
$\tilde{y}^{(k)}$ satisfying criterion (20). Our discussion is mainly based on
the approach presented in [23], which is a development of [18, 17, 47, 50].

We consider the case when the convex term in (11) has the structure

$f_{1}(x)=\sum_{i=1}^{p}{g_{i}}(M_{i}x)+\xi(x),$ (21)

where, for all $i\in\\{1,\ldots,p\\}$, $M_{i}\in\mathbb{R}^{m_{i}\times n}$
and ${g_{i}}:\mathbb{R}^{m_{i}}\rightarrow\overline{\mathbb{R}}$,
$\xi:\mathbb{R}^{n}\rightarrow\overline{\mathbb{R}}$ are proper convex
functions that are continuous on their domain. Moreover, we assume that
$\mathrm{dom}({g_{i}}\circ M_{i})\supseteq\mathrm{dom}(\xi)$, i.e.,
$\mathrm{dom}(f_{1})=\mathrm{dom}(\xi)$, and $g_{i}$, $\xi$ admit a closed
form formula for computing their associated proximity operators.

With the aim of describing a procedure for computing a point satisfying (20),
we rewrite the convex minimization problem in (14), omitting the iteration
index for simplicity, as

$\min_{z\in\mathbb{R}^{n}}\
h(z;x,s)\equiv\sum_{i=1}^{p}{g_{i}}(M_{i}z)+\xi(z)+\frac{1}{2\alpha}\|z-\bar{x}\|^{2}+c,$

where $\bar{x}=x-\alpha\nabla f_{0}(x)+\beta(x-s)$ and
$c=-\frac{\alpha}{2}\|\nabla f_{0}(x)-\beta/\alpha(x-s)\|^{2}-f_{1}(x)$ depend
on $x$ and $s$ but not on the optimization variable $z$. Recalling the
relation ${g_{i}}(M_{i}z)=\max_{w_{i}\in\mathbb{R}^{m_{i}}}\langle
w_{i},M_{i}z\rangle-{g_{i}}^{*}(w_{i})$, where
${g_{i}}^{*}:\mathbb{R}^{m_{i}}\rightarrow\overline{\mathbb{R}}$ is the convex
conjugate of $g_{i}$, and plugging it in the previous equation, results in the
following primal–dual formulation

$\min_{z\in\mathbb{R}^{n}}\max_{w\in\mathbb{R}^{m}}\langle
w,Mz\rangle-\sum_{i=1}^{p}{g_{i}}^{*}(w_{i})+\xi(z)+\frac{1}{2\alpha}\|z-\bar{x}\|^{2}+c,$

where $M=(M_{1}^{T}\ M_{2}^{T}\ \cdots M_{p}^{T})^{T}\in\mathbb{R}^{m\times
n}$, $m=m_{1}+\ldots+m_{p}$, and
$w=(w_{1}^{T},w_{2}^{T},\cdots,w_{p}^{T})^{T}\in\mathbb{R}^{m}$. Some further
simple manipulations lead to the following equivalent reformulation

$\displaystyle\max_{w\in\mathbb{R}^{m}}\min_{z\in\mathbb{R}^{n}}\ \xi(z)$
$\displaystyle+\frac{1}{2\alpha}\|z-(\bar{x}-\alpha
M^{T}w)\|^{2}-\frac{1}{2\alpha}\|\bar{x}-\alpha M^{T}w\|^{2}+$
$\displaystyle+\frac{1}{2\alpha}\|\bar{x}\|^{2}-\sum_{i=1}^{p}{g_{i}}^{*}(w_{i})+c.$

Minimizing with respect to $z$ and using the definition of proximity operator
yields the dual problem

$\max_{w\in\mathbb{R}^{m}}\psi(w;x,s)\equiv
G(w;x,s)-\sum_{i=1}^{p}{g_{i}}^{*}(w_{i}),$ (22)

where

$\displaystyle G(w;x,s)$
$\displaystyle=\xi({\mathrm{prox}}_{\alpha\xi}(\bar{x}-\alpha M^{T}w))$
$\displaystyle+\frac{1}{2\alpha}\|{\mathrm{prox}}_{\alpha\xi}(\bar{x}-\alpha
M^{T}w)-(\bar{x}-\alpha M^{T}w)\|^{2}$
$\displaystyle-\frac{1}{2\alpha}\|\bar{x}-\alpha
M^{T}w\|^{2}+\frac{1}{2\alpha}\|\bar{x}\|^{2}+c.$

The function $G(\,\cdot\,;x,s)$ is concave, continuously differentiable with
$\nabla G(w;x,s)=M{\mathrm{prox}}_{\alpha\xi}(\bar{x}-\alpha M^{T}w)$ [7,
Proposition 12.30], and its domain is the whole dual space $\mathbb{R}^{m}$.
Note also that, by definition,

$h(z;x,s)\geq\psi(w;x,s)\ \ \forall z\in\mathbb{R}^{n},w\in\mathbb{R}^{m}$

and the equality holds if and only if $z=\hat{y}$ and $w=\hat{w}$, where
$\hat{w}$ denotes a solution of the dual problem (22). In particular, the
following relations hold:

$\displaystyle\hat{y}$ $\displaystyle=$
$\displaystyle{\mathrm{prox}}_{\alpha\xi}(\bar{x}-\alpha M^{T}\hat{w})$
$\displaystyle\psi(\hat{w};x,s)$ $\displaystyle=$ $\displaystyle
h(\hat{y};x,s)$ $\displaystyle\psi(w;x,s)$ $\displaystyle\leq$ $\displaystyle
h(\hat{y};x,s)\ \ \forall w\in\mathbb{R}^{m}.$

Then, if a primal-dual pair
$(\tilde{y},\tilde{w})\in\mathbb{R}^{n}\times\mathbb{R}^{m}$ satisfies the
inequality

$h(\tilde{y};x,s)\leq\frac{2}{2+\tau}\psi(\tilde{w};x,s),$ (23)

the point $\tilde{y}$ complies with (20). A pair $(\tilde{y},\tilde{w})$
satisfying (23) can be computed by proceeding as follows:

  1. 1.

apply an iterative method to the dual problem (22) generating a sequence
$\\{w^{(\ell)}\\}_{\ell\in\mathbb{N}}\subset\mathrm{dom}(\psi)$ such that
$\\{w^{(\ell)}\\}_{\ell\in\mathbb{N}}$ converges to $\hat{w}$ and
$\\{\psi(w^{(\ell)};x,s)\\}_{\ell\in\mathbb{N}}$ converges to
$\psi(\hat{w};x,s)$;

  2. 2.

define the corresponding primal sequence
$\tilde{y}^{(\ell)}={\mathrm{prox}}_{\alpha\xi}(\bar{x}-\alpha
M^{T}w^{(\ell)})$;

  3. 3.

stop the dual iterations when
$h(\tilde{y}^{(\ell)};x,s)\leq\eta\psi(\tilde{w}^{(\ell)};x,s)$, where
$\eta=2/(2+\tau)$;

  4. 4.

set $\tilde{y}=\tilde{y}^{(\ell)}$.

Note that the described procedure is well defined: indeed, by continuity,
$\tilde{y}^{(\ell)}$ converges to $\hat{y}$, and since ${g_{i}},\xi$ are
assumed to be continuous, $h(\tilde{y}^{(\ell)};x,s)$ converges to
$h(\hat{y};x,s)$. If $\mathrm{dom}({g_{i}}\circ
M_{i})\supseteq\mathrm{dom}(\xi)$, as we assumed at the beginning of this
section, then all the points $\tilde{y}^{(\ell)}$ generated by this procedure
(including its output $\tilde{y}$) belong to the domain of $f_{1}$, which
coincides with the domain of $h(\cdot;x,s)$. Feasibility is required during
the dual iterations for checking the stopping criterion, as well as for the
inexactness criterion (19) to be fulfilled.

We remark that the objective function $-\psi(w;x,z)$ of the dual problem (22)
is the sum of the smooth convex term $-G(w;x,z)$ and the convex function
$\sum_{i=1}^{p}{g_{i}}^{*}(w_{i})$, whose proximity operator is easy to
compute, as it can be computed by separately evaluating the proximity
operators associated to ${g_{i}}^{*}$ ($i=1,...,p$), for which Moreau’s
identity and, consequently, the proximity operator of $g_{i}$, can be
exploited. This means that, for example, any (exact) forward–backward method
can be applied to the dual problem (22) for generating the sequence
$\\{w^{(\ell)}\\}_{\ell\in\mathbb{N}}$.

###  4.2 i2Piano: inertial inexact proximal algorithm for nonconvex
optimization

In this section we propose a generalization of the inertial method iPiano,
first proposed in [40] and further developed in [39], introducing the
possibility of an inexact computation of the inertial proximal gradient point.

Our proposed method is reported in Algorithm 1 and denoted as i2Piano
(inertial inexact Proximal algorithm for nonconvex optimization). Let us
describe the i2Piano iteration in detail. STEP 1–4 determine the stepsize
$\alpha_{k}$ and inertial parameter $\beta_{k}$ at iteration $k$. Given the
parameters $\alpha_{k},\beta_{k}$, STEP 5 seeks to find a possibly inexact
inertial proximal point, i.e., an inexact minimizer of the function (15).
According to (20), the i2Piano iterate is any point
${x}^{(k+1)}=\tilde{y}^{(k)}$ such that

$0\in\partial_{\epsilon_{k}}h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)}),\ \ \
\mbox{ with
}\epsilon_{k}=-\frac{\tau}{2}h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})$

for some fixed constant $\tau\geq 0$. When $\tau=0$, we recover the exact
inertial proximal gradient point provided by the iPiano method [39, 40]. If
$\tau>0$, as explained in Section 4.1.1, Step 5 of i2Piano can be practically
implemented with an inner loop consisting of an iterative optimization method
applied to the dual of problem
$\min_{y\in\mathbb{R}^{n}}h^{(k)}(y;{x}^{(k)},{x}^{(k-1)})$, until the
stopping condition (23) is met. In the implementation of i2Piano, besides the
stepsize $\alpha_{k}$ and the inertial parameter $\beta_{k}$, a further
parameter, $L_{k}$, is introduced (cf. STEP 6) with the aim of estimating a
local Lipschitz constant of $\nabla f_{0}$ that allows us to take larger
steps. In particular, $L_{k}$ is successively increased by a factor $\eta>1$
until the following descent condition holds

$f_{0}({x}^{(k+1)})\leq f_{0}(x^{(k)})+\langle\nabla
f_{0}(x^{(k)}),{x}^{(k+1)}-x^{(k)}\rangle+\frac{L_{k}}{2}\|{x}^{(k+1)}-x^{(k)}\|^{2}.$
(24)

Clearly, if the Lipschitz constant $L$ is known, the previous inequality is
satisfied with $L_{k}\geq L$ for all $k$, and the loop between STEP 2 and STEP
6 can be skipped. In this case, also the choice of the inertial parameter
$\beta_{k}$ and the corresponding steplength $\alpha_{k}$ could be simplified.
For example, $\\{\beta_{k}\\}_{k\in\mathbb{N}}$ can be set as a constant
sequence $\beta_{k}\equiv\beta<1$, or as a prefixed sequence such that
$\beta_{k}<1$ and $\lim_{k}\beta_{k}=1$, whereas $\alpha_{k}$ can be selected
as in Step 4 with $L_{k}=L$.

The choice of the parameters in Algorithm i2Piano is very similar to the one
proposed in [39], and aims at adaptively estimating the Lipschitz constant,
which is quite common in forward–backward methods [18, 25, 40, 39, 48]. We
stress that the main novelty in Algorithm i2Piano is not on the parameters
choice, but in the possibility of inexactly computing the proximity operator
according to an implementable inexactness criterion, which, as far as we know,
is new in the framework of inertial methods for nonconvex optimization.

Algorithm 1 i2Piano: inertial inexact proximal algorithm for nonconvex
optimization.

Choose $x^{(-1)},x^{(0)}\in\mathrm{dom}(f_{1})$, $\delta\geq\gamma>0$,
$\eta>1$, $0<L_{min}\leq L_{max}$, $\tau\geq 0$. Set
$\theta=2/(\sqrt{2+\tau}+\sqrt{\tau})^{2}$ and choose $\omega\in[0,1)$ if
$\tau>0$, $\omega\in[0,1]$ if $\tau=0$.
FOR $k=0,1,\ldots$

  * Step 1.

Choose $L_{k}\in[L_{min},L_{max}]$.

Step 2.

Set $\displaystyle b_{k}=\frac{L_{k}+2\delta}{L_{k}+2\gamma}$.

Step 3.

Set
$\beta_{k}=\displaystyle\frac{1+\theta\omega}{2}\cdot\frac{b_{k}-1}{b_{k}-\frac{1}{2}}$.

Step 4.

Set $\alpha_{k}=\displaystyle\frac{1+\theta\omega-2\beta_{k}}{L_{k}+2\gamma}$.

Step 5.

Compute $\tilde{y}^{(k)}$ such that

$0\in\partial_{\epsilon_{k}}h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},{x}^{(k-1)}),\ \
\ \mbox{ with
}\epsilon_{k}=-\frac{\tau}{2}h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},{x}^{(k-1)}).$

Step 6.

Check the local descent:

If

$\displaystyle f_{0}(\tilde{y}^{(k)})\leq f_{0}(x^{(k)})+\langle\nabla
f_{0}(x^{(k)}),\tilde{y}^{(k)}-x^{(k)}\rangle+\frac{L_{k}}{2}\|\tilde{y}^{(k)}-x^{(k)}\|^{2}$

-

Set ${x}^{(k+1)}=\tilde{y}^{(k)}$.

Else

-

Set $L_{k}=\eta L_{k}$.

-

Go to Step 2.

END

Under some standard boundedness assumptions and within the Kurdyka–Lojasiewicz
framework employed within the abstract framework of Section 3, we can prove
that the sequence generated by i2Piano converges to a stationary point, as
stated below.

######  Theorem 6.

Let ${\mathcal{F}}:\mathbb{R}^{n}\times\mathbb{R}\to\bar{\mathbb{R}}$ be the
merit function defined as

${\mathcal{F}}(u,\rho)=f(u)+\frac{1}{2}\rho^{2},\quad\forall
u\in\mathbb{R}^{n},\ \rho\in\mathbb{R}.$ (25)

Suppose that $\mathcal{F}$ is a KL function and assume that the sequence
$\\{x^{(k)}\\}_{k\in\mathbb{N}}$ generated by i2Piano is bounded. Then,
$\\{x^{(k)}\\}_{k\in\mathbb{N}}$ converges to a stationary point of $f$.

In order not to excessively slowing down the reading of this section, we
postpone the proof of Theorem 6 to Section A.2.

######  Remark 1.

We underline that $\mathcal{F}$ is a KL function if, for instance, $f$ and
$\frac{1}{2}\|\cdot\|^{2}$ are definable in the same $o-$minimal structure
[13, Definition 7]. Indeed functions definable in an $o-$minimal structure
satisfy the KL property on their domain [13, Theorem 11] and $o-$minimal
structures are closed with respect to the sum, see [13, Remark 5] and
references therein. Examples of functions definable in an $o-$minimal
structure are semialgebraic, real analytic functions and subanalytic functions
which are continuous in their closed domain.

######  Remark 2.

Theorem 6 requires the boundedness of the iterates as hypothesis. A standard
way to assert such a condition is when the Lyapunov function $\Phi$ defined in
(55) is coercive, since this assumption combined with the descent property
(56) guarantees that the sequence $\\{(x^{(k)},x^{(k-1)})\\}_{k\in\mathbb{N}}$
is included in a (bounded) level set of the coercive function $\Phi$.

###  4.3 iPila: inertial proximal inexact line–search algorithm

In the following we introduce a novel algorithm combining a line–search along
the descent direction and an inertial proximal-gradient step as a special case
of our abstract scheme.

A line–search procedure for the objective function $f$ requires a descent
direction $d\in\mathbb{R}^{n}$, i.e., a vector such that the directional
derivative $f^{\prime}(x;d)=\lim_{\lambda\downarrow
0}(f(x+\lambda)-f(x))/\lambda$ is negative. If $\beta_{k}=0$ and
$\tilde{y}^{(k)}$ satisfies (18), then the vector $\tilde{y}^{(k)}-{x}^{(k)}$
is a descent direction for $f$ at $x^{(k)}$, as explained in [17, 18] In fact,
if $\beta_{k}=0$, the inequality
$h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})\leq 0$ would be enough to
guarantee the descent property of the vector $\tilde{y}^{(k)}-{x}^{(k)}$, even
though it would not be sufficient to guarantee the convergence of the
iterates, as shown in [17]. Unfortunately, none of the two previous conditions
can guarantee that $\tilde{y}^{(k)}-{x}^{(k)}$ is a descent direction for $f$
at $x^{(k)}$ when $\beta_{k}>0$, i.e., when inertia is incorporated in the
iterative scheme.

In this section we show that in the general case $\beta_{k}\geq 0$, the point
$\tilde{y}^{(k)}$ can still be used to define a descent direction for a
suitable merit function. Then, we propose a line–search procedure along this
direction that enable us to define a descent algorithm such that the merit
function monotonically decreases along the iterates. Differently from the
backtracking procedure in i2Piano, the proposed line–search requires to solve
the minimization subproblem (14) only once per iteration. Finally, we show
that the new algorithm can be analyzed in the framework of Section 3, in order
to prove the convergence of the iterates to a stationary point of $f$. In this
case, unlike i2Piano, one of the two merit functions involved in the abstract
scheme will play an active role in the algorithm, as it will be explicitly
computed at each iteration to determine the new point.

We define the merit function $\Phi$ appearing in [H1]–[H2] as follows:

$\Phi:\mathbb{R}^{n}\times\mathbb{R}^{n}\to\mathbb{R}\cup\\{+\infty\\},\ \ \
\Phi(x,s)=f(x)+\frac{1}{2}\|x-s\|^{2},$ (26)

where the variable $s$ will be considered as an actual optimization variable,
independent on $x$. The function $\Phi$ can be decomposed as
$\Phi(x,s)=\Phi_{0}(x,s)+\Phi_{1}(x,s)$, where

$\Phi_{0}(x,s)=f_{0}(x)+\frac{1}{2}\|x-s\|^{2},\ \ \Phi_{1}(x,s)=f_{1}(x).$

The function $\Phi_{0}$ is differentiable with gradient

$\nabla\Phi_{0}(x,s)=\begin{pmatrix}\nabla_{x}\Phi_{0}(x,s)\\\
\nabla_{s}\Phi_{0}(x,s)\end{pmatrix}=\begin{pmatrix}\nabla f_{0}(x)+x-s\\\
s-x\end{pmatrix}.$

It is easy to see that $\nabla\Phi_{0}(x,s)$ is Lipschitz continuous and in
particular it holds that

$\|\nabla\Phi_{0}(x,s)-\nabla\Phi_{0}(\bar{x},\bar{s})\|\leq{L_{\Phi_{0}}}\left\|\begin{pmatrix}x-\bar{x}\\\
s-\bar{s}\end{pmatrix}\right\|,\ \ \forall x,\bar{x},s,\bar{s},$ (27)

with ${L_{\Phi_{0}}}=L+2$, where $L$ is the Lipschitz constant of $\nabla
f_{0}$.
Given a vector $d\in\mathbb{R}^{n}\times\mathbb{R}^{n}$

$d=\begin{pmatrix}d_{x}\\\ d_{s}\end{pmatrix},$ (28)

the directional derivative of $\Phi$ at the point $(x,s)$ with respect to the
direction $d$ can be written as

$\displaystyle\Phi^{\prime}(x,s;d_{x},d_{s})$ $\displaystyle=$
$\displaystyle\Phi_{0}^{\prime}(x,s;d_{x},d_{s})+\Phi_{1}^{\prime}(x,s;d_{x},d_{s})=f_{1}^{\prime}(x;d_{x})+\langle\nabla\Phi_{0}(x,s),d\rangle$
$\displaystyle=$ $\displaystyle f_{1}^{\prime}(x;d_{x})+\langle\nabla
f_{0}(x)+x-s,d_{x}\rangle+\langle s-x,d_{s}\rangle,$

which always exists thanks to the convexity of $f_{1}$.
A vector $d\in\mathbb{R}^{2n}$ is called a descent direction for $\Phi$ at
$(x,s)$ when

$\Phi^{\prime}(x,s;d_{x},d_{s})<0.$

Assume now that the vector $d_{x}$ has the form $d_{x}=y-x$, where $y$ is a
point belonging to the domain of $f_{1}$; then, from [44, Theorem 23.1] we
have

$\Phi^{\prime}(x,s;y-x,d_{s})\leq f_{1}(y)-f_{1}(x)+\langle\nabla
f_{0}(x)+x-s,d_{x}\rangle+\langle s-x,d_{s}\rangle.$ (29)

The above inequality holds independently on the form of $d_{s}$.

####  4.3.1 Defining a descent direction for the merit function

Assume that $({x}^{(k)},s^{(k)})$ is a given point in
$\mathrm{dom}(f_{1})\times\mathbb{R}^{n}$, while $\alpha_{k}$, $\beta_{k}$ are
two given parameters. Let the function $h^{(k)}(y;x,s)$ be defined as in (15).
Given a fixed tolerance parameter $\tau\geq 0$, we denote by $\tilde{y}^{(k)}$
any point in $\mathrm{dom}(f_{1})$ satisfying (20). Given the parameter
$\gamma_{k}\geq 0$, consider

$d^{(k)}=\begin{pmatrix}d_{x}^{(k)}\\\ d_{s}^{(k)}\end{pmatrix},$ (30)

where

$\displaystyle d_{x}^{(k)}$ $\displaystyle=$
$\displaystyle\tilde{y}^{(k)}-{x}^{(k)}$ (31) $\displaystyle d_{s}^{(k)}$
$\displaystyle=$
$\displaystyle\left(1+\frac{\beta_{k}}{\alpha_{k}}\right)(\tilde{y}^{(k)}-{x}^{(k)})+\gamma_{k}({x}^{(k)}-s^{(k)}).$
(32)

It can be shown (see Lemma 16 in Section A.3) that $d^{(k)}$ is a descent
direction for $\Phi$ at $({x}^{(k)},s^{(k)})$ and, therefore, we can introduce
a backtracking procedure along it to seek for a sufficient decrease of the
merit function.

In particular, given a point $({x}^{(k)},s^{(k)})$ and the direction $d^{(k)}$
(31)–(32), our proposed line–search algorithm computes a positive parameter
$\lambda_{k}^{+}$ satisfying the following generalized Armijo inequality.

$\Phi({x}^{(k)}+\lambda^{+}_{k}d_{x}^{(k)},s^{(k)}+\lambda^{+}_{k}d_{s}^{(k)})\leq\Phi({x}^{(k)},s^{(k)})+\sigma\lambda^{+}_{k}\Delta_{k},\quad\sigma\in(0,1),$
(33)

where

$\Delta_{k}=h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})-\gamma_{k}\|{x}^{(k)}-s^{(k)}\|^{2}.$

The implementation of this rule via a backtracking procedure is given in
Algorithm 2.

Algorithm 2 Armijo line–search

INPUT: $({x}^{(k)},s^{(k)})\in\mathbb{R}^{n}\times\mathbb{R}^{n}$, $d^{(k)}$,
$\Delta_{k}$ as in (68), $\sigma,\delta\in(0,1)$

  * Set $\lambda^{+}=1$.

  * WHILE $\Phi({x}^{(k)}+\lambda^{+}d_{x}^{(k)},s^{(k)}+\lambda^{+}d_{s}^{(k)})>\Phi({x}^{(k)},s^{(k)})+\sigma\lambda^{+}\Delta_{k}$

    * Set $\lambda^{+}=\delta\lambda^{+}$

  * Set $\lambda_{k}^{+}=\lambda^{+}$

END
OUTPUT: $\lambda_{k}^{+}$.

In Section A.3 we show that Algorithm 2 is well posed, i.e., it terminates in
a finite number of steps.

Algorithm 3 iPila: inertial proximal inexact line–search algorithm

INPUT: $(x^{(0)},s^{(0)})\in\mathrm{dom}(f_{1})\times\mathbb{R}^{n}$,
$\sigma\in(0,1)$, $0<\alpha_{min}\leq\alpha_{max}$, $\beta_{max}>0$,
$0<\gamma_{min}\leq\gamma_{max}$, $\tau\geq 0$.
FOR $k=0,1,\ldots$

  * STEP 1.

Choose $\alpha_{k}\in[\alpha_{min},\alpha_{max}]$,
$\beta_{k}\in[0,\beta_{max}]$

STEP 2.

Compute $\tilde{y}^{(k)}$ such that

$0\in\partial_{\epsilon_{k}}h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)}),\ \ \
\mbox{ with
}\epsilon_{k}=-\frac{\tau}{2}h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)}).$

STEP 3.

Choose $\gamma_{k}\in[\gamma_{min},\gamma_{max}]$.

STEP 4.

Compute
$\Delta_{k}=h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})-\gamma_{k}\|{x}^{(k)}-s^{(k)}\|^{2}$.

STEP 5.

Compute the search direction

$\displaystyle d_{x}^{(k)}$ $\displaystyle=$
$\displaystyle\tilde{y}^{(k)}-{x}^{(k)}$ $\displaystyle d_{s}^{(k)}$
$\displaystyle=$
$\displaystyle\left(1+\frac{\beta_{k}}{\alpha_{k}}\right)(\tilde{y}^{(k)}-{x}^{(k)})+\gamma_{k}({x}^{(k)}-s^{(k)})$

STEP 6.

Compute $\lambda_{k}\in(0,1]$ such that

$\Phi({x}^{(k)}+\lambda_{k}d_{x}^{(k)},s^{(k)}+\lambda_{k}d_{s}^{(k)})\leq\Phi({x}^{(k)},s^{(k)})+\sigma\lambda_{k}\Delta_{k}$

with the line–search backtracking algorithm.

STEP 7.

Define the new point as

$({x}^{(k+1)},s^{(k+1)})=\left\\{\begin{array}[]{l}(\tilde{y}^{(k)},{x}^{(k)})\
\ \mbox{ if
}\Phi(\tilde{y}^{(k)},{x}^{(k)})\leq\Phi({x}^{(k)},s^{(k)})+\sigma\lambda_{k}\Delta_{k}\\\
({x}^{(k)}+\lambda_{k}d_{x}^{(k)},s^{(k)}+\lambda_{k}d_{s}^{(k)})\ \ \mbox{
otherwise}\end{array}\right.$

END

The descent direction and the backtracking procedure described above are at
the basis of the new algorithm, named iPila (inertial Proximal inexact
line–search algorithm), which formally consists in a descent method for the
merit function $\Phi(x,s)$ and exploits the inertial inexact proximal gradient
point for defining the search direction. In particular, it generates a
sequence of iterates $\\{({x}^{(k)},s^{(k)})\\}_{k\in\mathbb{N}}$ and a
sequence of steplength parameters $\\{\lambda_{k}\\}_{{k\in\mathbb{N}}}$
fulfilling the following decrease condition

$\Phi({x}^{(k+1)},s^{(k+1)})\leq\Phi({x}^{(k)},s^{(k)})+\sigma\lambda_{k}\Delta_{k}\
\ \ \mbox{ and }\ \ \
\Phi({x}^{(k+1)},s^{(k+1)})\leq\Phi(\tilde{y}^{(k)},x^{(k)}).$

The connection with the inertial methods is in fact that, when
$({x}^{(k+1)},s^{(k+1)})=(\tilde{y}^{(k)},{x}^{(k)})$ is selected at STEP 7,
the following iteration will consist of an actual inertial step. In practice,
the condition at STEP 7 can be considered as an alternative acceptance rule
for the inexact inertial proximal gradient point, having a similar role than
the condition at STEP 6 of i2Piano (see also (24)). The main difference is
that here the acceptance condition is based on the Armijo inequality, while
the one in i2Piano is based on the Descent Lemma.

Notice that the inexact evaluation of the proximity operator in iPila is
required only once per iteration, unlike in i2Piano, where it is needed at
each step of the loop for selecting the parameter $L_{k}$, until inequality
(24) is satisfied. The Armijo condition also results in a larger freedom of
choosing the parameters $\alpha_{k},\beta_{k}$, which here satisfy very
minimal assumptions. A possible strategy to choose these parameters preserving
both the theoretical prescriptions and the benefits deriving from the presence
of an inertial step is described in Section 5.

The Armijo line–search strategy is very well established in optimization,
however, to the best of our knowledge, its use in combination with an
inertial/heavy-ball step is completely new. Moreover, even if the use of a
merit function is quite common in the theoretical analysis of optimization
methods, here we propose to explicitly compute it in the algorithm
implementation, enforcing its decrease by means of the line–search procedure.

Under the same assumptions stated for i2Piano, we can show the convergence of
iPila to a stationary point (the proof is postponed in Section A.3).

######  Theorem 7.

Suppose that the function $\mathcal{F}$ defined in (25) is a KL function.
Moreover, assume that the sequence $\\{x^{(k)}\\}_{k\in\mathbb{N}}$ generated
by iPila is bounded. Then, $\\{x^{(k)}\\}_{k\in\mathbb{N}}$ converges to a
stationary point of $f$.

We refer the reader to Remark 1-2 for conditions on $f$ guaranteeing that
$\mathcal{F}$ is a KL function and the sequence
$\\{x^{(k)}\\}_{k\in\mathbb{N}}$ is bounded.

##  5 Numerical illustration

The aim of this section is to apply our proposed methods to a couple of
difficult problems arising in image restoration, by implementing them in
accordance with the theoretical guarantees stated in Section A.

The goal of image restoration is to recover a good quality image from a noisy
blurred one. Following the variational approach, the clean image is obtained
by solving an optimization problem with the structure (11), where the
objective function includes a measure of the data fidelity and a
regularization/penalization term, incorporating all the a priori information
on the desired solution. In addition, a nonnegativity constraint is often
imposed for physical reasons. Then, the variational model has the following
form:

$\min_{x\in\mathbb{R}^{n}_{\geq 0}}{\mathcal{D}}(Hx,g)+{\mathcal{R}}(x),$ (34)

where $g\in\mathbb{R}^{n}$ is the noisy blurred data, $H\in\mathbb{R}^{n\times
n}$ represents the blurring operator,
$\mathcal{R}:\mathbb{R}^{n}\to\mathbb{R}$ is the regularization term, and
${\mathcal{D}}(\cdot,\cdot)$ is the data discrepancy function. Moreover,
$\mathbb{R}^{n}_{\geq 0}=\\{x\in\mathbb{R}^{n}:x_{i}\geq 0,i=1,\ldots,n\\}$ is
the nonnegative orthant. In the following, we will consider two instances of
the image restoration model (34), simulating two different kinds of noise,
impulse noise and signal dependent Gaussian noise.

###  5.1 Image denoising and deblurring in presence of impulse noise

When data suffers of impulse noise, the more suitable function to measure the
data discrepancy is the $\ell_{1}-$norm:

${\mathcal{D}}(Hx,g)=\|Hx-g\|_{1}.$

On the other side, with the aim to preserve some sharpness in the restored
image, an edge preserving regularization term has to be included in the
variational model. In this section we consider as regularization function the
one proposed in [28, 27], i.e.,

${\mathcal{R}}(x)=\rho\sum_{\ell=1}^{q}{\theta_{\ell}}\sum_{i=1}^{n}\log(1+(K_{\ell}x)_{i}^{2}),$

where the matrices $K_{\ell}\in\mathbb{R}^{n\times n}$ correspond to a
convolution with a given filter $k_{\ell}$, while $\rho,\theta_{\ell}$ are
positive parameters. In particular, the set of $48$ filters $k_{\ell}$ of size
$7\times 7$ and corresponding coefficients $\theta_{\ell}$ have been computed
with the approach proposed in [27]. Finally, it is possible to prove that
${\mathcal{R}}(x)$ has Lipschitz-continuous gradient, using the same arguments
as in [22]. Then, setting $f_{0}(x)={\mathcal{R}}(x)$ complies with
assumptions [A2]-[A3]. The nonnegativity constraint, expressed by means of the
indicator function of the nonnegative orthant $\iota_{\geq 0}(x)$, can be
included in the convex, nonsmooth term of the objective function, i.e.,
$f_{1}(x)=\|Hx-g\|_{1}+\iota_{\geq 0}(x)$. The regularization parameter $\rho$
has been manually tuned in order to have a good quality restoration. Its value
has been set equal to 0.08 for all the runs. Note that $f_{1}$ does not have a
closed-form proximal operator, therefore it is necessary to employ
implementable inexactness criteria as the one proposed in Section 4.1. Due to
the nonconvexity and the lack of a closed form formula for the proximity
operator, this problem is really challenging.

The proposed algorithms have been implemented in Matlab R2019a on a laptop
equipped with a 2.60 GHz Intel Core i7-4510U processor and 8 GB of RAM; the
Matlab code is available online at [19]. The parameters of Algorithm i2Piano
have been set as $\delta=0.5$, $\gamma=10^{-5}$, $\eta=1.5$, $\omega=0.95$.
The estimate of the Lipschitz constant $L_{k}$ is updated in a nondecreasing
way. In particular, the initial value $L_{0}$ is set as an input parameter;
then, at Step 1 of each iteration, the first tentative value is set as
$L_{k}=L_{k-1}$. This value is possibly increased until inequality (24) is
met. Actually, more sophisticated updating rules for this parameter could be
adopted; however the objective function of the considered image restoration
problem is very costly to evaluate, therefore a more conservative parameters
selection rule has shown to be more convenient. As for Algorithm iPila, the
parameters settings aim to mimic that of the inertial method i2Piano. Indeed,
introducing the additional parameters $\delta,\gamma,L_{k},b_{k}>0$, with
$b_{k}=\frac{L_{k}+2\delta}{L_{k}+2\gamma}$, we set
$\alpha_{k},\beta_{k},\gamma_{k}$ as follows

$\beta_{k}=\frac{b_{k}-1}{b_{k}-\frac{1}{2}},\qquad\qquad\alpha_{k}=2\frac{1-\beta_{k}}{L_{k}+2\gamma},\qquad\qquad\gamma_{k}=\gamma.$

This choice is motivated by the following arguments. If $L_{k}$ is a good
local approximation of the Lipschitz constant satisfying condition

$f_{0}(\tilde{y}^{(k)})\leq f_{0}(x^{(k)})+\langle\nabla
f_{0}(x^{(k)}),\tilde{y}^{(k)}-x^{(k)}\rangle+\frac{L_{k}}{2}\|\tilde{y}^{(k)}-x^{(k)}\|^{2},$
(35)

then reasoning as in the proof of Proposition 12, and choosing $\delta=0.5$
and $\gamma_{k}=\gamma$, we obtain

$\Phi(\tilde{y}^{(k)},{x}^{(k)})\leq\Phi({x}^{(k)},s^{(k)})+\Delta_{k}+\frac{1}{2\alpha_{k}}\|\tilde{y}^{(k)}-{x}^{(k)}\|^{2}.$

Hence, if $L_{k}$ satisfies (35), the point $(\tilde{y}^{(k)},{x}^{(k)})$ will
be likely accepted at Step 7, as also confirmed by the numerical experience.
This reasoning suggests to implement algorithm iPila as follows. We check the
condition
$\Phi(\tilde{y}^{(k)},{x}^{(k)})\leq\Phi({x}^{(k)},s^{(k)})+\sigma\Delta_{k}$
right after Step 4: if the condition holds, then the steps from 4 to 7 are
skipped in order to avoid unnecessary computations, and the next point is
directly defined as $({x}^{(k+1)},s^{(k+1)})=(\tilde{y}^{(k)},{x}^{(k)})$;
otherwise, the value $L_{k}$ is increased by a factor $\eta=1.5$, and the
line–search in steps 5–7 is performed in order to compute the next point. By
possibly increasing $L_{k}$, we aim at improving the chances that
$(\tilde{y}^{(k+1)},x^{(k+1)})$ is accepted at the next iteration, thus
reducing the computational time due to the line–search reductions steps. We
point out that this procedure for computing the parameters in iPila complies
with all the theoretical prescriptions in Section A. For iPila, the parameter
in the Armijo condition, is set to $\sigma=10^{-4}$, while
$\gamma_{k}=\gamma=10^{-5}$ for all $k$. The initial estimate $L_{0}$ of the
Lipschitz constant has been set equal to one for both algorithms.

For both i2Piano and iPila, the inexact proximal point $\tilde{y}^{(k)}$ is
computed by approximately solving the dual of problem
$\min_{y\in\mathbb{R}^{n}}h^{(k)}(y;{x}^{(k)},s^{(k)})$, which is a quadratic
problem with simple constraints, with FISTA (more details can be found in [22]
and references therein). The accuracy of the approximation is controlled by
the parameter $\tau$: in our experiments we set $\tau=10^{6}$, which
corresponds to a good balancing of the computational complexity among inner
and outer iterations. An extensive performance assessment of the algorithms
with respect to this and other parameters is out of the scope of this paper,
and it will be subject of future research.

The deblurring test problem has been obtained by first artificially blurring a
good quality image, then simulating impulse noise on the 15% of the pixels
with `imnoise`. The clean and the noisy image are reported in Figure 1 (a) and
(b). Assuming reflective boundary conditions, matrix-vector multiplications
involving $H$ and $H^{T}$ can be implemented efficiently with the DCT
transform. In particular, each inner (dual) iteration requires the computation
of two matrix-vector products of this kind. In fact, in our experiments, only
one or two inner iterations per outer iteration are, in general, needed to
satisfy the inner stopping criterion.

We compare our proposed algorithms to an inexact version of the standard
forward–backward method with backtracking in [8] (denoted as iISTA in the
following), and with the variable metric line–search based method denominated
VMILAn [18, 22] equipped with its standard parameters settings. The inexact
computation of the proximal gradient point in iISTA and VMILAn is implemented
exactly in the same way as for i2Piano and iPila. Moreover, the initial
estimate of the Lipschitz constant in iISTA is set to one, i.e., it is the
same choice made for i2Piano and iPila.

For comparing different algorithms on a nonconvex problem, we adopt the same
approach employed in [29, 39]: we first numerically estimated the “optimal”
value $f^{*}$ as the smallest function value among the ones obtained by
running each algorithm for a huge number of iterations. Then, we run again the
three algorithms, computing at each iteration the relative difference of the
objective function value with respect the reference value $f^{*}$. The result
of this comparison is depicted in Figure 2.

|  |
---|---|---
(a) | (b) | (c)
|  |
(d) | (e) | (f)
Figure 1: Image deblurring test problem. First row, impulse noise: (a)
Original image ($512\times 512$ pixels); (b) Noisy blurred image, PSNR =
13.19; (c) Restored image, PSNR = 23.72. Second row, signal dependent Gaussian
noise: (d) Original image ($256\times 256$ pixels); (e) Noisy blurred image,
PSNR = 23.38; (f) Restored image, PSNR = 29.09

Panel (a) reports the values $(f({x}^{(k)})-f^{*})/f^{*}$ provided by the
algorithms i2Piano, iPila, VMILAn, and iISTA with respect to the computational
time.

###  5.2 Image denoising and deblurring in presence of signal dependent
Gaussian noise

We now consider a second image deblurring problem, where we assume that the
data contains signal dependent Gaussian noise. In this case the discrepancy
functional is given by [29]

${\mathcal{D}}(Hx,g)=\frac{1}{2}\sum_{i=1}^{n}\frac{((Hx)_{i}-g_{i})^{2}}{a_{i}(Hx)_{i}+c_{i}}+\log(a_{i}(Hx)_{i}+c_{i}),$

where $a_{i},c_{i}$ are positive parameters. The functional above is
nonconvex, smooth, and assuming that the entries of the blurring matrix $H$
are nonnegative, its domain contains the nonnegative orthant. Moreover, its
gradient is Lipschitz continuous, even if an estimation of the Lipschitz
constant is difficult to compute. As regularization term, we adopt the Total
Variation function [46]

${\mathcal{R}}(x)=\rho\sum_{i=1}^{n}\|\nabla_{i}x\|,$

where the two components of $\nabla_{i}x\in\mathbb{R}^{2}$ contains the
differences of the pixel $i$ with its vertical and horizontal neighbours,
respectively. In this case, we can split the smooth and nonsmooth part of the
objective function as $f=f_{0}+f_{1}$ with $f_{0}(x)={\mathcal{D}}(Hx,g)$ and
$f_{1}(x)={\mathcal{R}}(x)+\iota_{\geq 0}(x)$. For the minimization of the
functional described above, we compare i2Piano and iPila with the variable
metric forward backward algorithm VMFB proposed in [29], and again with iISTA.
In particular, as a benchmark test, we adopt the same test problem provided in
[43], with the same setting for the parameters $a_{i},b_{i},c_{i},\rho$. The
true image and the blurred noisy one are reported in Figure 1 (d) and (e).

As done in the previous section, we numerically compute an approximation of
the optimal value $f^{*}$ by running all algorithms for a huge number of
iterations, keeping the last function value for each algorithm, and then
retaining the smallest value among them. Then, we evaluate the optimization
capability of the four algorithms in terms of the quantity
$(f({x}^{(k)})-f^{*})/f^{*}$ (see Figure 2 (b)). The parameters in i2Piano,
iPila and iISTA are set in the same way as for the impulse noise test problem,
while we run VMFB using the implementation released by the same authors [43].

The numerical results show that our proposed methods i2Piano and iPila are
able to effectively solve challenging problems and that the presence of an
inertial step can improve the effectiveness of standard forward–backward
methods. On the one hand, iPila outperforms i2Piano on both test problems,
which suggests that combining an inertial step with a linesearch along the
descent direction (rather than along the arc) may be extremely beneficial in
terms of computational times. On the other hand, iPila performs well also in
comparison with VMILAn and VMFB, showing that inertial forward–backward
algorithms can be competitive with variable metric approaches, provided that a
sensible linesearch strategy is adopted.

|
---|---
(a) | (b)
Figure 2: Relative decrease of the objective function with respect to the
computational time for the two test problems. Panel (a): impulse noise. Panel
(b): signal dependent Gaussian noise.

##  6 Conclusions

We proposed two novel inertial-type forward–backward algorithms for solving
nonsmooth nonconvex optimization problems. Both algorithms are equipped with
implementable inexactness criteria for computing the proximal operator of the
convex part. The first algorithm, i2Piano, performs a classical backtracking
procedure based on a local version of the Descent Lemma, whereas the second
algorithm, iPila, is based on an innovative linesearch procedure along the
descent direction of a suitable merit function. We showed that both algorithms
converge to a stationary point of the problem, under some standard assumptions
on the iterates sequence and the objective function. The convergence analysis
is cast into an abstract framework that generalizes related work and thereby
unfolds the formulation and convergence guarantees of our novel algorithms. We
demonstrated the efficiency of the proposed algorithms on a couple of
challenging image restoration problems, showing that the innovative approach
employed in iPila may be beneficial in terms of computational times. Future
work could be devoted to the design and analysis of novel effective rules for
selecting the parameters in the iPila algorithm, in order to further improve
its performance with respect to more traditional forward–backward algorithms.
Another possible development could be the adoption of a Bregman distance in
the computation of the inexact proximal–gradient point of i2Piano and iPila,
as done for other first-order methods in the KL framework [9, 16, 17].

##  Appendix A Appendix: convergence Analysis of i2Piano and iPila

In this section we analyze the convergence properties of algorithms i2Piano
and iPila, showing that they can be both considered as special cases of the
abstract scheme presented in Section 3.

###  A.1 Preliminary results

We collect here new basic results concerning the inexactness criterion
introduced in Section 4.1, which is incorporated in our proposed algorithms.
The following lemma is a consequence of the strong convexity of the function
$h(\cdot;x,s)$ defined in (15), and will be often employed in the following.
In order to simplify the notation, here we omit the iteration index $k$.

######  Lemma 8.

Suppose that Assumptions [A1]–[A2] hold true. For a given pair
$(x,s)\in\mathrm{dom}(f_{1})\times\mathbb{R}^{n}$, let $\hat{y},\tilde{y}$ be
defined as in (16),(18). Then, the following inequalities hold.

$\displaystyle\frac{1}{2\alpha}\|\hat{y}-x\|^{2}$ $\displaystyle\leq$
$\displaystyle\left(1+\frac{\tau}{2}\right)(-h(\tilde{y};x,s))$ (36)
$\displaystyle\frac{1}{2\alpha}\|\tilde{y}-\hat{y}\|^{2}$ $\displaystyle\leq$
$\displaystyle\frac{\tau}{2}(-h(\tilde{y};x,s))$ (37)
$\displaystyle\frac{\theta}{2\alpha}\|\tilde{y}-x\|^{2}$ $\displaystyle\leq$
$\displaystyle(-h(\tilde{y};x,s)),\ \ \mbox{ with
}\theta=1/\left(\sqrt{1+\frac{\tau}{2}}+\sqrt{\frac{\tau}{2}}\right)^{2}\leq
1.$ (38)

###### Proof.

Inequalities (36)–(37) follow by combining the strong convexity of the
function $h(\cdot;x,s)$ and condition (18) as in [22, Lemma 2]. As for (38),
we have

$\displaystyle\frac{1}{2\alpha}\|\tilde{y}-x\|^{2}$ $\displaystyle=$
$\displaystyle\frac{1}{2\alpha}\|\tilde{y}-\hat{y}+\hat{y}-x\|^{2}=\frac{1}{2\alpha}\|\tilde{y}-\hat{y}\|^{2}+\frac{1}{2\alpha}\|\hat{y}-x\|^{2}+\frac{1}{\alpha}\langle\tilde{y}-\hat{y},\hat{y}-x\rangle$
$\displaystyle\leq$
$\displaystyle\frac{1}{2\alpha}\|\tilde{y}-\hat{y}\|^{2}+\frac{1}{2\alpha}\|\hat{y}-x\|^{2}+\frac{1}{\alpha}\|\tilde{y}-\hat{y}\|\cdot\|\hat{y}-x\|$
$\displaystyle\leq$
$\displaystyle\left(1+\frac{\tau}{2}\right)(-h(\tilde{y};x,s))+\frac{\tau}{2}(-h(\tilde{y};x,s))+2\sqrt{1+\frac{\tau}{2}}\sqrt{\frac{\tau}{2}}(-h(\tilde{y};x,s))$
$\displaystyle=$
$\displaystyle\left(\sqrt{1+\frac{\tau}{2}}+\sqrt{\frac{\tau}{2}}\right)^{2}(-h(\tilde{y};x,s)),$

where the last inequality follows from the application of (36)-(37). $\square$

The next lemma provides a subgradient $\hat{v}\in\partial f(\hat{y})$ whose
norm is bounded from above by a quantity containing
$\sqrt{-h(\tilde{y};x,s)}$. Its proof is omitted since it is almost identical
to the one of Lemma 3 in [22].

######  Lemma 9.

Suppose Assumptions [A1]–[A3] hold true. Let $x$ be a point in
$\mathrm{dom}(f_{1})$ and let $\hat{y},\tilde{y}$ be defined as in (16)–(18).
Moreover, assume that $\alpha\in[\alpha_{min},\alpha_{max}]$, with
$0<\alpha_{min}\leq\alpha_{max}$ and $\beta\in[0,\beta_{max}]$, with
$\beta_{max}\geq 0$. Then, there exists a subgradient $\hat{v}\in\partial
f(\hat{y})$ such that

$\displaystyle\|\hat{v}\|$ $\displaystyle\leq$ $\displaystyle
p(\|\hat{y}-x\|+\|x-s\|)$ (39) $\displaystyle\leq$ $\displaystyle
q(\sqrt{-h(\tilde{y};x,s)}+\|x-s\|),$ (40)

where the two constants $p,q$ depend only on
$\alpha_{min},\alpha_{max},\beta_{max}$ and on the Lipschitz constant $L$.

The following lemma is the equivalent of Lemma 4-5 in [22].

######  Lemma 10.

Suppose Assumptions [A1]–[A3] hold true and assume
$0<\alpha_{min}\leq\alpha\leq\alpha_{max}$, $\beta\in[0,\beta_{max}]$. Let
$(x,s)$ be a point in $\mathrm{dom}(f_{1})\times\mathbb{R}^{n}$ and let
$\hat{y},\tilde{y}$ be defined as in (16)–(18), for some $\tau\geq 0$. Then,
there exists $c,d,\bar{c},\bar{d}\in\mathbb{R}$ depending only on
$\alpha_{min},\alpha_{max},\beta_{max},\tau$ such that

$\displaystyle f(\hat{y})$ $\displaystyle\geq$ $\displaystyle
f(\tilde{y})+ch(\tilde{y};x,s)-d\|x-s\|^{2}$ (41) $\displaystyle f(\hat{y})$
$\displaystyle\leq$ $\displaystyle
f(x)-\bar{c}h(\tilde{y};x,s)+\bar{d}\|x-s\|^{2}.$ (42)

###### Proof.

From the Descent Lemma [11, Proposition A.24] we have

$f_{0}(\hat{y})\geq f_{0}(\tilde{y})-\langle\nabla
f_{0}(\hat{y}),\tilde{y}-\hat{y}\rangle-\frac{L}{2}\|\tilde{y}-\hat{y}\|^{2}.$
(43)

The inclusion $0\in\partial_{\epsilon}h(\tilde{y};x,s)$ in (20) implies that
there exists a vector $e\in\mathbb{R}^{n}$ with

$\frac{1}{2\alpha}\|e\|^{2}\leq\epsilon$ (44)

such that

$-\frac{1}{\alpha}(\tilde{y}-x+\alpha\nabla
f_{0}(x)-\beta(x-s)+e)\in\partial_{\epsilon}f_{1}(\tilde{y})$

(see [18] and references therein). The definition of
$\epsilon$-subdifferential implies

$f_{1}(\hat{y})\geq
f_{1}(\tilde{y})-\frac{1}{\alpha}\langle\tilde{y}-x,\hat{y}-\tilde{y}\rangle+\frac{\beta}{\alpha}\langle\hat{y}-\tilde{y},x-s\rangle-\langle\nabla
f_{0}(x),\hat{y}-\tilde{y}\rangle-\frac{1}{\alpha}\langle
e,\hat{y}-\tilde{y}\rangle-\epsilon.$ (45)

Summing inequalities (43) and (45) yields

$\displaystyle f(\hat{y})$ $\displaystyle\geq$ $\displaystyle
f(\tilde{y})-\langle\nabla f_{0}(x)-\nabla
f_{0}(\hat{y}),\hat{y}-\tilde{y}\rangle+\frac{\beta}{\alpha}\langle\hat{y}-\tilde{y},x-s\rangle$
(46)
$\displaystyle-\frac{1}{\alpha}\langle\tilde{y}-x,\hat{y}-\tilde{y}\rangle-\frac{1}{\alpha}\langle
e,\hat{y}-\tilde{y}\rangle-\frac{L}{2}\|\tilde{y}-\hat{y}\|^{2}-\epsilon.$

Now we consider each term at the right-hand-side in the above inequality so as
to obtain a lower bound. Using the Cauchy-Schwarz inequality, Assumption [A3],
(36) and (37) we obtain

$\displaystyle\langle\nabla f_{0}(x)-\nabla
f_{0}(\hat{y}),\hat{y}-\tilde{y}\rangle$ $\displaystyle\leq$
$\displaystyle\|\nabla f_{0}(x)-\nabla f_{0}(\hat{y})\|\|\hat{y}-\tilde{y}\|$
(47) $\displaystyle\leq$ $\displaystyle
L\alpha_{max}\sqrt{2\tau\left(1+\frac{\tau}{2}\right)}(-h(\tilde{y};x,s)).$

Similarly, using again the Cauchy-Schwarz inequality, (37) and (38), we can
write

$\displaystyle\frac{1}{\alpha}\langle\tilde{y}-x,\hat{y}-\tilde{y}\rangle$
$\displaystyle\leq$
$\displaystyle\frac{1}{\alpha_{min}}\|\tilde{y}-x\|\|\hat{y}-\tilde{y}\|\leq\frac{\alpha_{max}}{\alpha_{min}}\sqrt{\frac{2\tau}{\theta}}(-h(\tilde{y};x,s)).$
(48)

Moreover, from (44) and (20) we obtain
$\|e\|\leq\sqrt{2\alpha_{max}\epsilon}\leq\sqrt{\alpha_{max}\tau(-h(\tilde{y};x,s))}$
which, using also (37), yields

$\displaystyle\frac{1}{\alpha}\langle e,\hat{y}-\tilde{y}\rangle$
$\displaystyle\leq$
$\displaystyle\frac{1}{\alpha}\|e\|\|\hat{y}-\tilde{y}\|\leq\frac{\alpha_{max}\tau}{\alpha_{min}}(-h(\tilde{y};x,s)).$
(49)

Finally, using (37), we can also write

$\frac{\beta}{\alpha}\langle\hat{y}-\tilde{y},x-s\rangle\geq-\frac{\beta}{2\alpha}(\|\hat{y}-\tilde{y}\|^{2}+\|x-s\|^{2})\geq-\frac{\beta_{max}}{2\alpha_{min}}\left(-\alpha_{max}\tau
h(\tilde{y};x,s)+\|x-s\|^{2}\right).$ (50)

Combining (46) with (38), (47), (48), (49), (50) and (20), gives (41) with

$c=L\alpha_{max}\sqrt{2\tau\left(1+\frac{\tau}{2}\right)}+\frac{\alpha_{max}}{\alpha_{min}}\sqrt{\frac{2\tau}{\theta}}+\frac{\alpha_{max}\tau}{\alpha_{min}}+\frac{L\alpha_{max}\tau}{2}+\frac{\tau}{2}+\frac{\beta_{max}\alpha_{max}\tau}{2\alpha_{min}},\
\ d=\frac{\beta_{max}}{2\alpha_{min}}.$

As for (42), using the Descent Lemma we obtain

$f_{0}(y)\leq f_{0}(x)+\langle\nabla
f_{0}(x),y-x\rangle+\frac{L}{2}\|y-x\|^{2},$

for all $x,y\in\mathrm{dom}(f_{1})$. Summing $f_{1}(y)$ on both sides yields

$\displaystyle f(y)$ $\displaystyle\leq$ $\displaystyle
f(x)+f_{1}(y)-f_{1}(x)+\langle\nabla
f_{0}(x),y-x\rangle+\frac{L}{2}\|y-x\|^{2}$ $\displaystyle\leq$ $\displaystyle
f(x)+h(y;x,s)+\frac{L}{2}\|y-x\|^{2}+\frac{\beta}{\alpha}\langle
x-s,y-x\rangle$ $\displaystyle\leq$ $\displaystyle
f(x)+h(y;x,s)+\frac{L}{2}\|y-x\|^{2}+\frac{\beta}{2\alpha}(\|x-s\|^{2}+\|y-x\|^{2}).$

From the previous inequality with $y=\hat{y}$, recalling that
$h(\hat{y};x,s)\leq 0$ and combining with (36) yields (42), where the
constants are set as $\bar{c}=(L\alpha_{max}+\beta_{max})(1+{\tau}/{2})$,
$\bar{d}={\beta_{max}}/{(2\alpha_{min})}$. $\square$

###  A.2 Convergence analysis of i2Piano

Our aim now is to frame i2Piano in the abstract scheme defined by Conditions 3
to enjoy the favourable convergence guarantees that are provided by Theorem 5.
The line of the proof developed in this section is based on an extension of
the arguments in [39]. We start the convergence analysis by showing that
i2Piano is well-posed and that its parameters satisfy some useful relations.

######  Lemma 11.

The loop between STEP 2 and STEP 6 terminates in a finite number of steps. In
particular, there exists $\overline{L}>0$ such that $L_{k}\leq\overline{L}$,
$\forall\ k\geq 0$. Moreover, we have

$0\leq\beta_{k}\leq\frac{1+\theta\omega}{2},\ \ \forall\ k\geq 0$ (51)

and there exist two positive constants $\alpha_{min},\alpha_{max}$ with
$0<\alpha_{min}\leq\alpha_{max}$ such that
$\alpha_{k}\in[\alpha_{min},\alpha_{max}]$, $\forall k\geq 0$. We also have

$\displaystyle\frac{1+\theta\omega}{2\alpha_{k}}-\frac{L_{k}}{2}-\frac{\beta_{k}}{2\alpha_{k}}$
$\displaystyle=$ $\displaystyle\delta$ (52)
$\displaystyle\delta-\frac{\beta_{k}}{2\alpha_{k}}$ $\displaystyle=$
$\displaystyle\gamma.$ (53)

###### Proof.

Since $\eta>1$, after a finite number of steps the tentative value of $L_{k}$
satisfies $L_{k}\geq L$, where $L$ is the Lipschitz constant of $\nabla
f_{0}$. Then, from the Descent Lemma, the inequality at Step 6 is satisfied.
From $\delta\geq\gamma$ we have $b_{k}\geq 1$, which implies (51). A simple
inspection shows that the following equalities hold:

$b_{k}=\frac{1+\theta\omega-{\beta_{k}}}{1+\theta\omega-2\beta_{k}}\Rightarrow\frac{L_{k}+2\delta}{L_{k}+2\gamma}=\frac{1+\theta\omega-{\beta_{k}}}{1+\theta\omega-2\beta_{k}},$

which leads to rewriting the parameter $\alpha_{k}$ as

$\alpha_{k}=\frac{1+\theta\omega-\beta_{k}}{L_{k}+2\delta}.$ (54)

Then there holds $\alpha_{k}\geq\alpha_{min}$ with
$\alpha_{min}=(1+\theta\omega)/(2(\overline{L}+2\delta))$ and, since
$\theta\omega\leq 1$, we also have $\alpha_{k}\leq\alpha_{max}$ with
$\alpha_{max}=2/L_{min}$. Moreover, we have

$\frac{1+\theta\omega}{\alpha_{k}}-\frac{L_{k}}{2}-\frac{\beta_{k}}{2\alpha_{k}}=\frac{1+\theta\omega-\beta_{k}}{2\alpha_{k}}-\frac{L_{k}}{2}=\frac{L_{k}+2\delta}{2}-\frac{L_{k}}{2}=\delta$

and

$\delta-\frac{\beta_{k}}{2\alpha_{k}}=\frac{1+\theta\omega}{2\alpha_{k}}-\frac{L_{k}}{2}-\frac{\beta_{k}}{\alpha_{k}}=\frac{1+\theta\omega-2\beta_{k}}{2\alpha_{k}}-\frac{L_{k}}{2}=\frac{L_{k}+2\gamma}{2}-\frac{L_{k}}{2}=\gamma.$

$\square$
Notice that the case $\tau=0$, which corresponds to the exact computation of
the inertial proximal gradient point at Step 5, implies $\theta=1$ and,
choosing $\omega=1$, the parameters settings in i2Piano are exactly the same
as in [40]. The need of introducing the parameter $\omega$ is mainly technical
and will be explained in the following.
We now prove that condition [H1] holds for i2Piano when the corresponding
merit function $\Phi$ is defined as follows:

$\displaystyle\Phi:\mathbb{R}^{n}\times\mathbb{R}^{n}\rightarrow\overline{\mathbb{R}},\
\ \ \ \Phi(x,s)=f(x)+\delta\|x-s\|^{2}.$ (55)

######  Proposition 12.

Let $\\{x^{(k)}\\}_{k\in\mathbb{N}}$ be the sequence generated by i2Piano.
Then there holds

$\Phi(x^{(k+1)},x^{(k)})\leq\Phi(x^{(k)},x^{(k-1)})-\gamma\|{x}^{(k)}-{x}^{(k-1)}\|^{2}+(1-\omega)h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)}).$
(56)

Consequently, there exist $\\{s^{(k)}\\}_{k\in\mathbb{N}}$,
$\\{d_{k}\\}_{k\in\mathbb{N}}$, $\\{a_{k}\\}_{k\in\mathbb{N}}$ such that [H1]
holds with $\Phi$ as in (55).

###### Proof.

By summing the quantity $f_{1}({x}^{(k+1)})$ to both sides of inequality (24)
we obtain

$\displaystyle f({x}^{(k+1)})$ $\displaystyle\leq$ $\displaystyle
f({x}^{(k)})\\!+\\!f_{1}({x}^{(k+1)})-f_{1}({x}^{(k)})\\!+\\!\langle\nabla
f_{0}({x}^{(k)}),{x}^{(k+1)}-{x}^{(k)}\rangle\\!+\\!\frac{L_{k}}{2}\|{x}^{(k+1)}-{x}^{(k)}\|^{2}$
$\displaystyle=$ $\displaystyle
f({x}^{(k)})+h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})-\left(\frac{1}{2{\alpha_{k}}}-\frac{L_{k}}{2}\right)\|{x}^{(k+1)}-{x}^{(k)}\|^{2}$
$\displaystyle+\frac{\beta_{k}}{\alpha_{k}}\langle{x}^{(k+1)}-{x}^{(k)},{x}^{(k)}-{x}^{(k-1)}\rangle$
$\displaystyle\leq$ $\displaystyle
f({x}^{(k)})+h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})-\left(\frac{1}{2{\alpha_{k}}}-\frac{L_{k}}{2}\right)\|{x}^{(k+1)}-{x}^{(k)}\|^{2}$
$\displaystyle+\frac{\beta_{k}}{2\alpha_{k}}(\|{x}^{(k+1)}-{x}^{(k)}\|^{2}+\|{x}^{(k)}-{x}^{(k-1)}\|^{2})$
$\displaystyle\leq$ $\displaystyle
f({x}^{(k)})+(1-\omega)h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})-\frac{\theta\omega}{2\alpha_{k}}\|{x}^{(k+1)}-{x}^{(k)}\|^{2}$
$\displaystyle-\left(\frac{1}{2{\alpha_{k}}}-\frac{L_{k}}{2}\right)\|{x}^{(k+1)}-{x}^{(k)}\|^{2}+\frac{\beta_{k}}{2\alpha_{k}}(\|{x}^{(k+1)}-{x}^{(k)}\|^{2}+\|{x}^{(k)}-{x}^{(k-1)}\|^{2})$
$\displaystyle=$ $\displaystyle
f({x}^{(k)})-\left(\frac{1+\theta\omega}{2{\alpha_{k}}}-\frac{L_{k}}{2}-\frac{\beta_{k}}{2\alpha_{k}}\right)\|{x}^{(k+1)}-{x}^{(k)}\|^{2}$
$\displaystyle+(1-\omega)h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})+\frac{\beta_{k}}{2\alpha_{k}}\|{x}^{(k)}-{x}^{(k-1)}\|^{2}$

where the first equality is obtained by adding and subtracting to the right-
hand-side the quantity
$\|{x}^{(k+1)}-{x}^{(k)}\|^{2}/(2{\alpha_{k}})+\beta_{k}/\alpha_{k}\langle{x}^{(k+1)}-{x}^{(k)},{x}^{(k)}-{x}^{(k-1)}\rangle$,
the subsequent inequality follows from the basic relation $2\langle
a,b\rangle\leq\|a\|^{2}+\|b\|^{2}$ and the next one from (38).
Recalling (52), the above inequality can be conveniently rewritten as

$\displaystyle f(x^{(k+1)})+\delta\|x^{(k+1)}$
$\displaystyle-x^{(k)}\|^{2}\leq f(x^{(k)})+\delta\|x^{(k)}-x^{(k-1)}\|^{2}$
$\displaystyle\ \
+\left(\frac{\beta_{k}}{2\alpha_{k}}-\delta\right)\|x^{(k)}-x^{(k-1)}\|^{2}+(1-\omega)h({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)}).$

Finally, exploiting (53), we obtain condition [H1] with $\Phi$ given in (55),
$a_{k}=1$ and

$\displaystyle s^{(k)}$ $\displaystyle=$ $\displaystyle{x}^{(k-1)}$ (57)
$\displaystyle d_{k}^{2}$ $\displaystyle=$
$\displaystyle\gamma\|{x}^{(k)}-{x}^{(k-1)}\|^{2}-(1-\omega)h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)}).$
(58)

$\square$
Under Assumption [A4], $\Phi$ is bounded from below, hence, condition [H1]
implies

$\lim_{k\to\infty}\|x^{(k)}-x^{(k-1)}\|=0$ (59)

and, if $\omega<1$, also

$\lim_{k\to\infty}h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})=0.$ (60)

The choice $\omega<1$ is enforced when $\tau>0$ in order to obtain (60). If we
take $\omega=1$, with the same arguments as above we still obtain (59). It is
also worth noticing that, for large values of $\tau$, that is when a coarser
accuracy is allowed in the computation of $\tilde{y}^{(k)}$, the parameter
$\theta$ can be very small and this also influences the choice of $\alpha_{k}$
and $\beta_{k}$.

In order to prove condition [H2], let us now introduce the second merit
function as follows

${\mathcal{F}}:\mathbb{R}^{n}\times\mathbb{R}\to\overline{\mathbb{R}},\ \ \ \
{\mathcal{F}}(u,\rho)=f(u)+\frac{1}{2}\rho^{2}.$ (61)

The following result is proved using similar arguments as the ones used in
Lemma 4-5 in [22].

######  Proposition 13.

Let $\\{{x}^{(k)}\\}_{k\in\mathbb{N}}$ be the sequence generated by i2Piano
and, for each $k$, let the point $\hat{x}^{(k+1)}$ be defined as

$\hat{x}^{(k+1)}=\operatornamewithlimits{argmin}\limits_{y\in\mathbb{R}^{n}}h^{(k)}(y;{x}^{(k)},{x}^{(k-1)}).$
(62)

Then, there exist $\\{\rho^{(k)}\\}_{k\in\mathbb{N}}$,
$\\{r_{k}\\}_{k\in\mathbb{N}}$ such that condition [H2] holds with $\Phi$
defined in (55), $\mathcal{F}$ defined in (61), $s^{(k)}$ defined in (57) and
$u^{(k)}=\hat{x}^{(k+1)}$.

###### Proof.

When $\omega=1$ we necessarily have $\tau=0$, which means
${x}^{(k+1)}=\hat{x}^{(k+1)}$. Therefore, [H2] direcly follows from [H1] with
$u^{(k)}={x}^{(k+1)}$, $\rho^{(k)}=\sqrt{2\delta}\|{x}^{(k)}-{x}^{(k-1)}\|$,
$r_{k}=0$. Consider now the case $\omega<1$. From (41) and (42), we directly
obtain

$\displaystyle f(\hat{x}^{(k+1)})$ $\displaystyle\geq$ $\displaystyle
f({x}^{(k+1)})+ch^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})-d\|{x}^{(k)}-{x}^{(k-1)}\|^{2}$
$\displaystyle f(\hat{x}^{(k+1)})$ $\displaystyle\leq$ $\displaystyle
f({x}^{(k)})-\bar{c}h({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})+\bar{d}\|{x}^{(k)}-{x}^{(k-1)}\|^{2},$

where $c,d,\bar{c},\bar{d}$ are defined as in Lemma 10 and do not depend on
$k$. Combining the two inequalities above we obtain

$\displaystyle f($
$\displaystyle{x}^{(k+1)})+\delta\|{x}^{(k+1)}-{x}^{(k)}\|^{2}\leq$
$\displaystyle\leq
f(\hat{x}^{(k+1)})+\delta\|{x}^{(k+1)}-{x}^{(k)}\|^{2}{-ch^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})}+d\|{x}^{(k)}-{x}^{(k-1)}\|^{2}$
$\displaystyle\leq
f({x}^{(k)})-(c+\bar{c})h({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})+(d+\bar{d})\|{x}^{(k)}-{x}^{(k-1)}\|^{2}+\delta\|{x}^{(k+1)}-{x}^{(k)}\|^{2}.$

Recalling the definition of $\mathcal{F}$ in (61), the above inequalities can
be rewritten as

$f({x}^{(k+1)})+\delta\|{x}^{(k+1)}-{x}^{(k)}\|^{2}\leq{\mathcal{F}}(\hat{x}^{(k+1)},\rho^{(k)})\leq
f({x}^{(k)})+\delta\|{x}^{(k)}-{x}^{(k-1)}\|^{2}+r_{k}$ (63)

where $\rho^{(k)},r_{k}$ are given by

$\displaystyle\rho^{(k)}$ $\displaystyle=$
$\displaystyle\sqrt{2}(\delta\|{x}^{(k+1)}-{x}^{(k)}\|^{2}{-ch^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})}+d\|{x}^{(k)}-{x}^{(k-1)}\|^{2})^{\frac{1}{2}}$
$\displaystyle r_{k}$ $\displaystyle=$
$\displaystyle-(c+\bar{c})h({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})+(d+\bar{d})\|{x}^{(k)}-{x}^{(k-1)}\|^{2}+\delta\|{x}^{(k+1)}-{x}^{(k)}\|^{2}.$

From (59)–(60) we obtain $\displaystyle\lim_{k\to\infty}r_{k}=0$, hence,
assumption [H2] is satisfied with the above settings and with
$u^{(k)}=\hat{x}^{(k+1)}$. $\square$
Next we show that condition [H3] holds for i2Piano. The following result
combines elements of Lemma 3 in [22] and Lemma 17 in [39].

######  Proposition 14.

There exist $b>0$, $I\subset\mathbb{Z}$, $\\{\theta_{i}\\}_{i\in I}$ such that
condition [H3] holds with
$\\{u^{(k)}\\}_{k\in\mathbb{N}},\\{\rho^{(k)}\\}_{k\in\mathbb{N}},\\{d_{k}\\}_{k\in\mathbb{N}}$
defined as in Proposition 13 and in (58), $\zeta_{k+1}=0$ and $b_{k+1}=1$.

###### Proof.

From the separable structure of ${\mathcal{F}}$, if $v\in\partial f(u)$ and
$\rho\in\mathbb{R}$, we have that $(v,\rho)\in\partial{\mathcal{F}}(u,\rho)$.
In particular,

$\|\partial{\mathcal{F}}(u,\rho)\|_{-}\leq\|v\|+\lvert\rho\rvert,\mbox{ for
all }v\in\partial f(u),\rho\in\mathbb{R}.$ (65)

If $\omega=1$, we are in the case $\tau=0$. This means that
${x}^{(k+1)}=\hat{x}^{(k+1)}$,
$\rho^{(k)}=\sqrt{2\delta}\|{x}^{(k)}-{x}^{(k-1)}\|$. Moreover, (58) reduces
to
$d_{k}=\sqrt{\gamma}\|{x}^{(k)}-{x}^{(k-1)}\|=\sqrt{\gamma/(2\delta)}\rho^{(k)}$.
From (39), we have that there exists a subgradient $\hat{v}^{(k+1)}\in\partial
f(\hat{x}^{(k+1)})$ and a positive constant $p$ such that

$\displaystyle\|\hat{v}^{(k+1)}\|$ $\displaystyle\leq$ $\displaystyle
p(\|{x}^{(k+1)}-{x}^{(k)}\|+\|{x}^{(k)}-{x}^{(k-1)}\|)=\frac{p}{\sqrt{\gamma}}(d_{k+1}+d_{k}).$

Hence, $\|\hat{v}^{(k+1)}\|+\rho^{(k)}\leq
pd_{k+1}/{\sqrt{\gamma}}+d_{k}(p+\sqrt{2\delta})/{\sqrt{\gamma}}$ and,
recalling (65), [H3] follows with $b=(p+\sqrt{2\delta})/{\sqrt{\gamma}}$,
$I=\\{0,1\\}$, $\theta_{0}=\theta_{1}=\frac{1}{2}$.
Consider now the case $\omega<1$. From (40), there exists a subgradient
$\hat{v}^{(k+1)}\in\partial f(\hat{x}^{(k+1)})$ and a positive constant $q$
such that

$\displaystyle\|\hat{v}^{(k+1)}\|$ $\displaystyle\leq$ $\displaystyle
q\sqrt{-h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})}+q\|{x}^{(k)}-{x}^{(k-1)}\|$
(66) $\displaystyle\leq$ $\displaystyle
q\sqrt{-h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})}+q\sqrt{\frac{2\alpha_{max}}{\theta}}\sqrt{-h^{(k-1)}({x}^{(k)};{x}^{(k-1)},x^{(k-2)})}$
$\displaystyle\leq$
$\displaystyle\frac{q}{\sqrt{1-\omega}}\sqrt{-(1-\omega)h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})+\gamma\|{x}^{(k)}-{x}^{(k-1)}\|^{2}}$
$\displaystyle+q\sqrt{\frac{2\alpha_{max}}{\theta(1-\omega)}}\sqrt{-(1-\omega)h^{(k-1)}({x}^{(k)};{x}^{(k-1)},x^{(k-2)})+\gamma\|{x}^{(k-1)}-x^{(k-2)}\|^{2}}$
$\displaystyle\leq$
$\displaystyle\frac{q}{\sqrt{1-\omega}}d_{k}+q\sqrt{\frac{2\alpha_{max}}{\theta(1-\omega)}}d_{k-1}.$

From (LABEL:sigmadef) and (38) we have

$\displaystyle\rho^{(k)}$ $\displaystyle\leq$
$\displaystyle\sqrt{2}\left(-(2\delta\alpha_{max}/\theta+c)h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})+d\|{x}^{(k)}-{x}^{(k-1)}\|^{2}\right)^{\frac{1}{2}}$
(67) $\displaystyle\leq$ $\displaystyle rd_{k}$

where
$r=\sqrt{2}\max\left\\{(2\delta\alpha_{max}/\theta+c)/(1-\omega),d/\gamma\right\\}^{\frac{1}{2}}$.
Then, combining (66) with (67), in view of (65) we obtain

$\|\partial{\mathcal{F}}(\hat{x}^{(k+1)},\rho^{(k)})\|_{-}\leq\frac{b}{2}(d_{k}+d_{k-1}),\
\ \mbox{ with
}b=2\max\left\\{\frac{q}{\sqrt{1-\omega}}+r,q\sqrt{\frac{2\alpha_{max}}{\theta(1-\omega)}}\right\\}$

and the thesis follows with $I=\\{1,2\\}$,
$\theta_{1}=\theta_{2}=\frac{1}{2}$. $\square$
The following lemma holds for all methods whose iterates satisfy [H1],[H2],
[H3] with ${\mathcal{F}}$ defined as in (61), and it is crucial to ensure
condition [H4] for i2Piano.

######  Proposition 15.

Let Assumptions [A1]–[A4] be satisfied and assume that
$\\{{x}^{(k)}\\}_{k\in\mathbb{N}}$ satisfies [H1],[H2], [H3] with
${\mathcal{F}}$ defined as in (61). If
$\\{(x^{(k_{j})},\rho^{(k_{j})})\\}_{j\in\mathbb{N}}$ is a subsequence of
$\\{({x}^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$ converging to some
$(x^{*},\rho^{*})\in\mathbb{R}^{n}\times\mathbb{R}^{m}$ such that
$\lim_{j\to\infty}\|u^{(k_{j})}-x^{(k_{j})}\|=0$, then $\rho^{*}=0$ and
$\lim_{j\to\infty}{\mathcal{F}}(u^{(k_{j})},\rho^{(k_{j})})={\mathcal{F}}(x^{*},\rho^{*})$.

###### Proof.

Recalling that [H1] implies $\displaystyle\lim_{k\to\infty}d_{k}=0$, from [H3]
and thanks to the separable structure of $\mathcal{F}$, we obtain that there
exists $\hat{v}^{(k)}\in\partial f(u^{(k)})$ such that
$\displaystyle\lim_{k\to\infty}\|\hat{v}^{(k)}\|=\displaystyle\lim_{k\to\infty}\rho^{(k)}=0$.
In particular, in view of (12), we can write $\hat{v}^{(k)}=\nabla
f_{0}(u^{(k)})+w^{(k)}$, where $w^{(k)}\in\partial f_{1}(u^{(k)})$. Therefore,
by continuity of $\nabla f_{0}$, the following implication holds

$\lim_{k\to\infty}\hat{v}^{(k)}=0\Rightarrow\lim_{j\to\infty}w^{(k_{j})}=-\nabla
f_{0}(x^{*}).$

Adding the quantity $\frac{1}{2}(\rho^{(k_{j})})^{2}$ to both sides of the
subgradient inequality yields

$\displaystyle f_{1}(x^{*})+\frac{1}{2}(\rho^{(k_{j})})^{2}$
$\displaystyle\geq$ $\displaystyle f_{1}(u^{(k_{j})})+\langle
w^{(k_{j})},x^{*}-u^{(k_{j})}\rangle+\frac{1}{2}(\rho^{(k_{j})})^{2}$
$\displaystyle=$
$\displaystyle{\mathcal{F}}(u^{(k_{j})},\rho^{(k_{j})})-f_{0}(u^{(k_{j})})+\langle
w^{(k_{j})},x^{*}-u^{(k_{j})}\rangle$

where the last equality is obtained by adding and subtracting
$f_{0}(u^{(k_{j})})$ to the right-hand-side. Taking limits on both sides we
obtain

$f_{1}(x^{*})+\frac{1}{2}(\rho^{*})^{2}\geq\lim_{j\to\infty}{\mathcal{F}}(u^{(k_{j})},\rho^{(k_{j})})-f_{0}(x^{*}),$

which, rearranging terms, gives
$\displaystyle\lim_{j\to\infty}{\mathcal{F}}(u^{(k_{j})},\rho^{(k_{j})})\leq{\mathcal{F}}(x^{*},\rho^{*})$.
On the other side, by assumption, $\mathcal{F}$ is lower semicontinuous,
therefore
$\displaystyle\lim_{j\to\infty}{\mathcal{F}}(u^{(k_{j})},\rho^{(k_{j})})\geq{\mathcal{F}}(x^{*},\rho^{*})$,
which completes the proof. $\square$
We are now ready to prove Theorem 6, which states the convergence of the
i2Piano iterates to a stationary point.

###### Proof of Theorem 6.

By Propositions 12-13-14, we know that conditions [H1]-[H2]-[H3] hold for
i2Piano. Furthermore, if $\omega=1$, then
$u^{(k)}=\hat{x}^{(k+1)}={x}^{(k+1)}$, and (59) directly implies that
$\displaystyle\lim_{k\to\infty}\|u^{(k)}-{x}^{(k)}\|=0$. If $\omega<1$, using
(36) and (60) we have

$\displaystyle\lim_{k\to\infty}\|u^{(k)}-{x}^{(k)}\|$
$\displaystyle=\lim_{k\to\infty}\|\hat{x}^{(k+1)}-{x}^{(k)}\|$
$\displaystyle\leq\lim_{k\to\infty}\sqrt{-2\alpha_{max}\left(1+\frac{\tau}{2}\right)h^{(k)}({x}^{(k+1)};{x}^{(k)},{x}^{(k-1)})}=0.$

Then, in both cases, the assumptions of Proposition 15 are satisfied and,
therefore, condition [H4] holds. Finally, condition [H5] follows from (58),
while condition [H6] is trivially satisfied, since both sequences
$\\{a_{k}\\}_{k\in\mathbb{N}}$ and $\\{b_{k}\\}_{k\in\mathbb{N}}$ are
constant. Then, Theorem 5 applies and guarantees that the sequence
$\\{(x^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$ converges to a stationary point
$(x^{*},\rho^{*})$ of $\mathcal{F}$. Note that, since $\mathcal{F}$ is the sum
of separable functions, its subdifferential can be written as
$\partial\mathcal{F}(x,\rho)=\partial f(x)\times\\{\rho\\}$. Then,
$(x^{*},\rho^{*})$ is stationary for $\mathcal{F}$ if and only if $\rho^{*}=0$
and $0\in\partial f(x^{*})$. Hence, $x^{*}$ is a stationary point for $f$ and
$\\{x^{(k)}\\}_{k\in\mathbb{N}}$ converges to it. $\square$

###  A.3 Convergence analysis of iPila

This section aims to develop the convergence framework for iPila. The first
issue to be addressed is the well posedness of the line–search algorithm. To
this end, we prove first the following Lemma.

######  Lemma 16.

Let $d^{(k)}\in\mathbb{R}^{2n}$ be defined according to (31)-(32), where
$0<{\alpha_{k}}\leq\alpha_{max}$, $\beta_{k}\geq 0$,
$\gamma_{k}\geq\gamma_{min}$, for some given positive real constants
$\alpha_{max}>0$, $\gamma_{min}\geq 0$. Define $\Delta_{k}\in\mathbb{R}_{\leq
0}$ as

$\Delta_{k}=h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})-\gamma_{k}\|{x}^{(k)}-s^{(k)}\|^{2}.$
(68)

Then, we have

$\displaystyle\Phi^{\prime}(x^{(k)},s^{(k)};d_{x}^{(k)},d_{s}^{(k)})$
$\displaystyle\leq\Delta_{k}$ (69)
$\displaystyle\leq-a\|\tilde{y}^{(k)}-x^{(k)}\|^{2}-\gamma_{min}\|x^{(k)}-s^{(k)}\|^{2}$
(70)

where $a>0$ is a positive constant. Therefore,
$\Phi^{\prime}({x}^{(k)},s^{(k)};d_{x}^{(k)},d_{s}^{(k)})<0$ whenever
$\tilde{y}^{(k)}\neq{x}^{(k)}$ or ${x}^{(k)}\neq s^{(k)}$.

###### Proof.

We first observe that (29) with $x={x}^{(k)}$, $s=s^{(k)}$,
$y=\tilde{y}^{(k)}$, $d_{x}=d_{x}^{(k)}$, $d_{s}=d_{s}^{(k)}$, gives:

$\displaystyle\Phi^{\prime}$
$\displaystyle(x^{(k)},s^{(k)};d_{x}^{(k)},d_{s}^{(k)})$ $\displaystyle\leq
f_{1}(\tilde{y}^{(k)})-f_{1}(x^{(k)})+\langle\nabla
f_{0}(x^{(k)})+x^{(k)}-s^{(k)},d_{x}^{(k)}\rangle+\langle
s^{(k)}-x^{(k)},d_{s}^{(k)}\rangle$ (71)
$\displaystyle=f_{1}(\tilde{y}^{(k)})-f_{1}({x}^{(k)})+\langle\nabla
f_{0}({x}^{(k)}),\tilde{y}^{(k)}-{x}^{(k)}\rangle+\langle{x}^{(k)}-s^{(k)},\tilde{y}^{(k)}-{x}^{(k)}\rangle+$
$\displaystyle\phantom{=}+\left(1+\frac{\beta_{k}}{\alpha_{k}}\right)\langle
s^{(k)}-{x}^{(k)},\tilde{y}^{(k)}-{x}^{(k)}\rangle-\gamma_{k}\|{x}^{(k)}-s^{(k)}\|^{2}$
$\displaystyle=f_{1}(\tilde{y}^{(k)})-f_{1}({x}^{(k)})+\langle\nabla
f_{0}({x}^{(k)})-\frac{\beta_{k}}{\alpha_{k}}({x}^{(k)}-s^{(k)}),\tilde{y}^{(k)}-{x}^{(k)}\rangle-\gamma_{k}\|{x}^{(k)}-s^{(k)}\|^{2}$
$\displaystyle\leq
h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})-\gamma_{k}\|{x}^{(k)}-s^{(k)}\|^{2}$
(72)
$\displaystyle\leq-\frac{\theta}{2\alpha_{k}}\|\tilde{y}^{(k)}-{x}^{(k)}\|^{2}-\gamma_{k}\|{x}^{(k)}-s^{(k)}\|^{2}$

where the last inequality follows from (38). Then, the thesis follows from
$\alpha_{k}\leq\alpha_{max}$ with $a=\theta/{2\alpha_{max}}$. $\square$

The well posedness of the line–search procedure and its main properties are
summarized in the following lemma.

######  Lemma 17.

Let the assumptions of Lemma 16 be satisfied with, in addition,
$\beta_{k}\leq\beta_{max}$, $\gamma\leq\gamma_{max}$, ${\beta_{max}\geq 0}$,
$\gamma_{max}\geq 0$. Then, the Armijo backtracking line–search algorithm
terminates in a finite number of steps, and there exists $\lambda_{min}>0$
such that the parameter $\lambda_{k}^{+}$ computed with the line–search
algorithm satisfies

$\lambda_{k}^{+}\geq\lambda_{min}.$ (73)

###### Proof.

Let us first prove that

$\Delta_{k}\leq-C\|d^{(k)}\|^{2},\ \ \forall k\in\mathbb{N}.$ (74)

for a positive constant $C>0$. Setting
$\delta_{k}=1+\frac{\beta_{k}}{\alpha_{k}}$, the bounds on the parameters
imply that $1\leq\delta_{k}\leq\bar{\delta}$, with
$\bar{\delta}=1+\frac{\beta_{max}}{\alpha_{min}}$. By definition of $d^{(k)}$
in (30) we have

$\displaystyle\|d^{(k)}\|^{2}$
$\displaystyle=\|d_{x}^{(k)}\|^{2}+\|d_{s}^{(k)}\|^{2}$
$\displaystyle=\|\tilde{y}^{(k)}-{x}^{(k)}\|^{2}+\|\delta_{k}(\tilde{y}^{(k)}-{x}^{(k)})+\gamma_{k}({x}^{(k)}-s^{(k)})\|^{2}$
$\displaystyle=(1+\delta_{k}^{2})\|\tilde{y}^{(k)}-{x}^{(k)}\|^{2}+\gamma_{k}^{2}\|{x}^{(k)}-s^{(k)}\|^{2}+2\delta_{k}\gamma_{k}\langle\tilde{y}^{(k)}-{x}^{(k)},{x}^{(k)}-s^{(k)}\rangle$
$\displaystyle\leq(1+\delta_{k}^{2}+\delta_{k}\gamma_{k})\|\tilde{y}^{(k)}-{x}^{(k)}\|^{2}+(\gamma_{k}^{2}+\delta_{k}\gamma_{k})\|{x}^{(k)}-s^{(k)}\|^{2}$
$\displaystyle\leq(1+\bar{\delta}^{2}+\bar{\delta}\gamma_{max})\|\tilde{y}^{(k)}-{x}^{(k)}\|^{2}+\gamma_{k}(\gamma_{max}+\bar{\delta})\|{x}^{(k)}-s^{(k)}\|^{2}$
$\displaystyle\leq\frac{1}{C}(a\|\tilde{y}^{(k)}-{x}^{(k)}\|^{2}+\gamma_{k}\|{x}^{(k)}-s^{(k)}\|^{2})$

where
$C=1/\max\\{(1+\bar{\delta}^{2}+\bar{\delta}\gamma_{max})/a,\gamma_{max}+\bar{\delta}\\}$.
Multiplying both sides of the last inequality above by $C$ and combining with
(70) gives (74).
Since $\Phi_{0}$ has $M$-Lipschitz continuous gradient, with $M=L+2$ (see
(27)), we can apply the Descent Lemma obtaining

$\displaystyle\Phi_{0}({x}^{(k)}+\lambda d_{x}^{(k)},s^{(k)}+\lambda
d_{s}^{(k)})$
$\displaystyle\leq\Phi_{0}({x}^{(k)},s^{(k)})+\lambda\langle\nabla\Phi_{0}({x}^{(k)},s^{(k)}),d^{(k)}\rangle+\frac{L_{\Phi_{0}}}{2}\lambda^{2}\|d^{(k)}\|^{2}$
$\displaystyle\leq\Phi_{0}({x}^{(k)},s^{(k)})+\lambda\langle\nabla\Phi_{0}({x}^{(k)},s^{(k)}),d^{(k)}\rangle-\frac{L_{\Phi_{0}}}{2C}\lambda^{2}\Delta_{k},$
(75)

where the second inequality follows from (74). From the Jensen’s inequality
applied to the convex function $f_{1}$ we also obtain

$\displaystyle\Phi_{1}({x}^{(k)}+\lambda d_{x}^{(k)},s^{(k)}+\lambda
d_{s}^{(k)})$ $\displaystyle=$ $\displaystyle f_{1}({x}^{(k)}+\lambda
d_{x}^{(k)})=f_{1}(\lambda\tilde{y}^{(k)}+(1-\lambda){x}^{(k)})$ (76)
$\displaystyle\leq$ $\displaystyle(1-\lambda)f_{1}({x}^{(k)})+\lambda
f_{1}(\tilde{y}^{(k)}).$

Summing (75) with (76) gives

$\displaystyle\Phi({x}^{(k)}$ $\displaystyle+\lambda
d_{x}^{(k)},s^{(k)}+\lambda d_{s}^{(k)})\leq$
$\displaystyle\leq\Phi({x}^{(k)},s^{(k)})+\lambda\left(f_{1}(\tilde{y}^{(k)})-f_{1}({x}^{(k)})+\langle\nabla\Phi_{0}({x}^{(k)},s^{(k)}),d^{(k)}\rangle\right)-\frac{L_{\Phi_{0}}}{2C}\lambda^{2}\Delta_{k}$
$\displaystyle\leq\Phi({x}^{(k)},s^{(k)})+\lambda\Delta_{k}-\frac{L_{\Phi_{0}}}{2C}\lambda^{2}\Delta_{k},$

where the last inequality follows from (71)–(72). The above relation implies

$\Phi({x}^{(k)}+\lambda d_{x}^{(k)},s^{(k)}+\lambda
d_{s}^{(k)})\leq\Phi({x}^{(k)},s^{(k)})+\lambda(1-\rho\lambda)\Delta_{k},\ \
\forall\lambda\in[0,1]$

with $\rho=\frac{L_{\Phi_{0}}}{2C}$. Moreover, comparing the above inequality
with the Armijo condition (33) shows that the last one is surely fulfilled
when $\lambda_{k}^{+}$ satisfies $1-\rho\lambda_{k}^{+}\geq\sigma$, that is
when $\lambda_{k}^{+}\leq(1-\sigma)/\rho$.
Since $\lambda_{k}^{+}$ in the backtracking procedure is obtained starting
from 1 and by successive reductions of a factor $\delta<1$, we have
$\lambda_{k}^{+}\geq\delta^{L_{\Phi_{0}}}$, where $L_{\Phi_{0}}$ is the
smallest nonnegative integer such that
$\delta^{L_{\Phi_{0}}}\leq(1-\sigma)/\rho$. Therefore, (73) is satisfied with
$\lambda_{min}=\delta^{L_{\Phi_{0}}}$. $\square$

In the remaining of this section we show that iPila can be cast in the
framework of the abstract scheme. We first show that conditions [H1]–[H3] are
satisfied.

######  Proposition 18.

Let $\\{({x}^{(k)},s^{(k)})\\}_{k\in\mathbb{N}}$ be the sequence generated by
iPila. Then, condition [H1] holds with $d_{k}=\sqrt{-\Delta_{k}}$,
$a_{k}=\sigma\lambda_{min}$. Moreover, under Assumption [A4], we also have

$0=\lim_{k\to\infty}\|{x}^{(k)}-s^{(k)}\|=\lim_{k\to\infty}h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})=\lim_{k\to\infty}\|\tilde{y}^{(k)}-{x}^{(k)}\|.$
(77)

###### Proof.

From the updating rule at STEP 7 and from Lemma 17, we have

$\displaystyle\Phi({x}^{(k+1)},s^{(k+1)})$ $\displaystyle\leq$
$\displaystyle\Phi({x}^{(k)},s^{(k)})+\sigma\lambda_{k}\Delta_{k}\leq\Phi({x}^{(k)},s^{(k)})+\sigma\lambda_{min}\Delta_{k}.$

Then, condition [H1] is satisfied with $d_{k}^{2}=-\Delta_{k}$,
$a_{k}=\sigma\lambda_{min}$. Since from Assumption [A4] $f$ is bounded from
below, $\Phi$ is bounded from below as well. Therefore, [H1] implies
$\displaystyle-\sum_{k=0}^{\infty}\Delta_{k}<\infty$ which, in turn, yields
$\displaystyle\lim_{k\to\infty}\Delta_{k}=0$. Recalling (70), this implies
(77). $\square$
In the following we describe the setup for proving [H2], with the second
auxiliary function defined as in (61).

######  Proposition 19.

Let $\\{({x}^{(k)},s^{(k)})\\}_{{k\in\mathbb{N}}}$ be the sequence generated
by Algorithm iPila with $\gamma_{min}>0$ and let ${\mathcal{F}}$ be defined as
in (61). Then, there exist
$\\{\rho^{(k)}\\}_{k\in\mathbb{N}},\\{r_{k}\\}_{k\in\mathbb{N}}\subset\mathbb{R}$,
with $\lim_{k\to\infty}r_{k}=0$, such that [H2] holds with
$u^{(k)}=\hat{y}^{(k)}$, where $\hat{y}^{(k)}$ is the exact minimizer of
$h^{(k)}(y;{x}^{(k)},s^{(k)})$, i.e.,

$\hat{y}^{(k)}=\operatornamewithlimits{argmin}\limits_{y\in\mathbb{R}^{n}}h^{(k)}(y;{x}^{(k)},s^{(k)}).$

###### Proof.

From STEP 7, we have

$\displaystyle\Phi({x}^{(k+1)},s^{(k+1)})$ $\displaystyle\leq$
$\displaystyle\Phi(\tilde{y}^{(k)},{x}^{(k)})=f(\tilde{y}^{(k)})+\frac{1}{2}\|\tilde{y}^{(k)}-{x}^{(k)}\|^{2}$
$\displaystyle\leq$ $\displaystyle
f(\hat{y}^{(k)})-\left(c+\frac{\alpha_{max}}{\theta}\right)h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})+d\|{x}^{(k)}-s^{(k)}\|^{2},$

where the last inequality follows from (41) and (38). Setting

$\rho^{(k)}=\sqrt{2}\left(-\left(c+\frac{\alpha_{max}}{\theta}\right)h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})+d\|{x}^{(k)}-s^{(k)}\|^{2}\right)^{\frac{1}{2}},$
(78)

we obtain
$\Phi({x}^{(k+1)},s^{(k+1)})\leq{\mathcal{F}}(\hat{y}^{(k)},\rho^{(k)})$,
which represents the left-most inequality in [H2], with
$u^{(k)}=\hat{y}^{(k)}$. On the other hand, from inequality (42) we obtain

${\mathcal{F}}(\hat{y}^{(k)},\rho^{(k)})=f(\hat{y}^{(k)})+\frac{1}{2}(\rho^{(k)})^{2}\leq
f({x}^{(k)})-\bar{c}h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})+\bar{d}\|{x}^{(k)}-s^{(k)}\|^{2}+\frac{1}{2}(\rho^{(k)})^{2}.$

Setting
$r_{k}=(\rho^{(k)})^{2}/2-\bar{c}h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})+(\bar{d}-\frac{1}{2})\|{x}^{(k)}-s^{(k)}\|^{2}$,
we can write

${\mathcal{F}}(\hat{y}^{(k)},\rho^{(k)})\leq
f({x}^{(k)})+\frac{1}{2}\|{x}^{(k)}-s^{(k)}\|^{2}+r_{k}=\Phi({x}^{(k)},s^{(k)})+r_{k}.$

From (77) we have that $\lim_{k\to\infty}r_{k}=0$ and this proves [H2].
$\square$

######  Proposition 20.

Let $\\{({x}^{(k)},s^{(k)})\\}_{{k\in\mathbb{N}}}$ be the sequence generated
by iPila with $\gamma_{min}>0$. Then, there exists a positive constant $b$
such that [H3] is satisfied with $I=\\{1\\}$, $\theta_{1}=1$, $\zeta_{k}=0$.

###### Proof.

From (40) we know that there exists a subgradient $\hat{v}^{(k)}\in\partial
f(\hat{y}^{(k)})$ such that

$\|\hat{v}^{(k)}\|\leq
q\sqrt{-h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})}+q\|{x}^{(k)}-s^{(k)}\|$
(79)

and, reasoning as in the proof of Proposition 15, it follows that

$\|\partial{\mathcal{F}}(\hat{y}^{(k)},\rho^{(k)})\|_{-}\leq\left\|\begin{pmatrix}\hat{v}^{(k)}\\\
\rho^{(k)}\end{pmatrix}\right\|\leq\|\hat{v}^{(k)}\|+\lvert\rho^{(k)}\rvert.$
(80)

Let us analyze the two terms at the right-hand side of the inequality above,
showing that both can be bounded from above with a multiple of
$\sqrt{-\Delta_{k}}$. From (79) we obtain

$\displaystyle\|\hat{v}^{(k)}\|$ $\displaystyle\leq$ $\displaystyle
q\sqrt{-h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})+\gamma_{min}\|{x}^{(k)}-s^{(k)}\|^{2}}$
$\displaystyle+\frac{q}{\sqrt{\gamma_{min}}}\sqrt{\gamma_{min}\|{x}^{(k)}-s^{(k)}\|^{2}-h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})}.$

which, setting $A=q\left(1+1/\sqrt{\gamma_{min}}\right)$ and using (70),
yields

$\|\hat{v}^{(k)}\|\leq
A\sqrt{-h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})+\gamma_{min}\|{x}^{(k)}-s^{(k)}\|^{2}}\leq
A\sqrt{-\Delta_{k}}.$ (81)

On the other hand, from definition (78)

$\displaystyle\rho^{(k)}$ $\displaystyle\leq$ $\displaystyle
B\sqrt{-h^{(k)}(\tilde{y}^{(k)};{x}^{(k)},s^{(k)})+\gamma_{k}\|{x}^{(k)}-s^{(k)}\|^{2}}=B\sqrt{-\Delta_{k}}$

where
$B=\sqrt{2}\max\left\\{c+\alpha_{max}/{\theta},d/{\gamma_{min}}\right\\}^{\frac{1}{2}}$.
Therefore, combining the last inequality above with (81) and (80), yields
$\|\partial{\mathcal{F}}(\hat{y}^{(k)},\rho^{(k)})\|_{-}\leq(A+B)\sqrt{-\Delta_{k}}$
which proves that [H3] is satisfied with $I=\\{1\\}$, $\theta_{1}=1$,
$\zeta_{k}=0$, $b=A+B$, $b_{k}=1$. $\square$
We are now ready for proving Theorem 7, which states the main convergence
result for Algorithm iPila.

###### Proof of Theorem 7.

By Proposition 18, 19, and 20, we know that conditions [H1]-[H2]-[H3] hold for
iPila. From (36) we have
$\displaystyle\lim_{k\to\infty}\|u^{(k)}-{x}^{(k)}\|=\displaystyle\lim_{k\to\infty}\|\hat{y}^{(k)}-{x}^{(k)}\|=0$.
Hence we can apply Proposition 15 and conclude that [H4] holds. Moreover,
condition [H5] holds as a consequence of Lemma 16, since
$d_{k}=\sqrt{-\Delta_{k}}$ and
$\|\tilde{y}^{(k)}-x^{(k)}\|\geq\|x^{(k+1)}-x^{(k)}\|/\lambda_{k}\geq\|x^{(k+1)}-x^{(k)}\|$
(see STEP 7). Finally, condition [H6] is trivially satisfied, since both
sequences $\\{a_{k}\\}_{k\in\mathbb{N}}$ and $\\{b_{k}\\}_{k\in\mathbb{N}}$
are constant. Then Theorem 5 applies and guarantees that the sequence
$\\{(x^{(k)},\rho^{(k)})\\}_{k\in\mathbb{N}}$ converges to a stationary point
$(x^{*},\rho^{*})$ of $\mathcal{F}$. Note that, since $\mathcal{F}$ is the sum
of separable functions, its subdifferential can be written as
$\partial\mathcal{F}(x,\rho)=\partial f(x)\times\\{\rho\\}$. Then,
$(x^{*},\rho^{*})$ is stationary for $\mathcal{F}$ if and only if $\rho^{*}=0$
and $0\in\partial f(x^{*})$. Hence $x^{*}$ is a stationary point for $f$ and
$\\{x^{(k)}\\}_{k\in\mathbb{N}}$ converges to it. $\square$

##### Funding and data availability statement

Silvia Bonettini, Marco Prato and Simone Rebegoldi are members of the INdAM
research group GNCS. Peter Ochs acknowledges funding by the German Research
Foundation (DFG Grant OC 150/3-1).
The test image barbara in the experiments described in Section 5.1 is included
in the software available at [19], while the jetplane one used in Section 5.2
can be downloaded from [43].

## References

  * [1] P. A. Absil, R. Mahony, and B. Andrews.  Convergence of the iterates of descent methods for analytic cost functions.  SIAM J. Optim., 16(2):531–547, 2005.
  * [2] H. Attouch and J. Bolte.  On the convergence of the proximal algorithm for nonsmooth functions involving analytic features.  Mathematical Programming, 116(1-2):5–16, 2009.
  * [3] H. Attouch, J. Bolte, P. Redont, and A. Soubeyran.  Proximal alternating minimization and projection methods for nonconvex problems: an approach based on the Kurdyka-łojasiewicz inequality.  Math. Oper. Res., 35(2):438–457, May 2010.
  * [4] H. Attouch, J. Bolte, and B. F. Svaiter.  Convergence of descent methods for semi-algebraic and tame problems: proximal algorithms, forward–backward splitting, and regularized Gauss–Seidel methods.  Math. Program., 137(1–2):91–129, February 2013.
  * [5] G. R. Ayers and J. C. Dainty.  Iterative blind deconvolution method and its applications.  Opt. Lett., 13(7):547–549, July 1988.
  * [6] F. Bach, R. Jenatton, J. Mairal, and G. Obozinski.  Structured sparsity through convex optimization.  Stat. Sci., 27(4):450–468, 2012.
  * [7] H. H. Bauschke and P. L. Combettes.  Convex analysis and monotone operator theory in Hilbert spaces.  CMS Books on Mathematics. Springer, New York, 2011.
  * [8] A. Beck and M. Teboulle.  A fast iterative shrinkage-thresholding algorithm for linear inverse problems.  SIAM J. Imaging Sci., 2(1):183–202, 2009.
  * [9] M. Benning, M.M. Betcke, M.J. Ehrhardt, and C.-B. Schönlieb.  Choose Your Path Wisely: Gradient Descent in a Bregman Distance Framework.  SIAM J. Imaging Sci., 14(2):814–843, 2021.
  * [10] M. Bertero, P. Boccacci, and V. Ruggiero.  Inverse Imaging with Poisson Data.  IOP Publishing, Bristol, 2018.
  * [11] D. Bertsekas.  Nonlinear programming.  Athena Scientific, Belmont, 1999.
  * [12] J. Bolte, A. Daniilidis, O. Ley, and L. Mazet.  Characterizations of Łojasiewicz inequalities: Subgradient flows, talweg, convexity.  Trans. Am. Math. Soc., 362(6):3319–3363, 2010.
  * [13] J. Bolte, A. Daniilidis, and Shiota M.  Clarke subgradients of stratifiable functions.  SIAM J. Optim., 10:556–572, 2007.
  * [14] J. Bolte, A. Danilidis, and A. Lewis.  The Łojasiewicz inequality for nonsmooth subanalytic functions with applications to subgradient dynamical systems.  SIAM J. Optim., 17(4):1205–1223, 2007.
  * [15] J. Bolte, S. Sabach, and M. Teboulle.  Proximal alternating linearized minimization for nonconvex and nonsmooth problems.  Math. Program., 146(1–2):459–494, August 2014.
  * [16] J. Bolte, S. Sabach, M. Teboulle, and Y. Vaisbourd.  First order methods beyond convexity and Lipschitz gradient continuity with applications to quadratic inverse problems.  SIAM J. Optim., 28(3):2131–2151, 2018.
  * [17] S. Bonettini, I. Loris, F. Porta, and M. Prato.  Variable metric inexact line–search based methods for nonsmooth optimization.  SIAM J. Optim., 26(2):891–921, 2016.
  * [18] S. Bonettini, I. Loris, F. Porta, M. Prato, and S. Rebegoldi.  On the convergence of a linesearch based proximal-gradient method for nonconvex optimization.  Inverse Probl., 33(5):055005, 2017.
  * [19] S. Bonettini, P. Ochs, M. Prato, and S. Rebegoldi.  inertial inexact Proximal algorithm for nonconvex optimization (i2Piano) and inertial Proximal inexact linesearch algorithm (iPila) software.  http://www.oasis.unimore.it/site/home/software.html, 2021.
  * [20] S. Bonettini, M. Prato, and S. Rebegoldi.  A block coordinate variable metric linesearch based proximal gradient method.  Comput. Optim. Appl., 71(1):5–52, 2018.
  * [21] S. Bonettini, M. Prato, and S. Rebegoldi.  Convergence of inexact forward-backward algorithms using the forward-backward envelope.  SIAM J. Optim., 30(4):3069–3097, 2020.
  * [22] S. Bonettini, M. Prato, and S. Rebegoldi.  New convergence results for the inexact variable metric forward-backward method.  Appl. Math. Comput., 392:125719, 2021.
  * [23] S. Bonettini, S. Rebegoldi, and V. Ruggiero.  Inertial variable metric techniques for the inexact forward-backward algorithm.  SIAM J. Sci. Comput., 40(5):A3180–A3210, 2018.
  * [24] L. Bottou, F. E. Curtis, and J. Nocedal.  Optimization methods for large-scale machine learning.  SIAM Rev., 60(2):223–311, 2018.
  * [25] L. Calatroni and A. Chambolle.  Backtracking strategies for accelerated descent methods with smooth composite objectives.  SIAM J. Optim., 29(3):1772–1798, 2019.
  * [26] A. Chambolle and T. Pock.  An introduction to continuous optimization for imaging.  Acta Numer., 25:161–319, 2016.
  * [27] Y. Chen, R. Ranftl, and T. Pock.  Insights into analysis operator learning: From patch-based sparse models to higher order mrfs.  IEEE Trans. Image Process., 23(3):1060–1072, 2014.
  * [28] Yunjin Chen, Thomas Pock, René Ranftl, and Horst Bischof.  Revisiting loss-specific training of filter-based MRFs for image restoration.  In Joachim Weickert, Matthias Hein, and Bernt Schiele, editors, Pattern Recognition, pages 271–281, Berlin, Heidelberg, 2013. Springer Berlin Heidelberg.
  * [29] E. Chouzenoux, J.-C. Pesquet, and A. Repetti.  Variable metric forward-backward algorithm for minimizing the sum of a differentiable function and a convex function.  J. Optim. Theory Appl., 162(1):107–132, July 2014.
  * [30] E. Chouzenoux, J.-C. Pesquet, and A. Repetti.  A block coordinate variable metric forward-backward algorithm.  J. Global Optim., 66(3):457–485, 2016.
  * [31] P. L. Combettes and J.-C. Pesquet.  Proximal splitting methods in signal processing.  In H. H. Bauschke, R. S. Burachik, P. L. Combettes, V. Elser, D. R. Luke, and H. Wolkowicz, editors, Fixed-point algorithms for inverse problems in science and engineering, Springer Optimization and Its Applications, pages 185–212. Springer, New York, NY, 2011.
  * [32] P.L. Combettes and V. R. Wajs.  Signal recovery by proximal forward-backward splitting.  Multiscale Model. Simul., 4(4):1168–1200, 2005.
  * [33] L. Van den Dries.  Tame topology and $\o$-minimal structures.  150 184. Cambridge University Press, 1998.
  * [34] P. Frankel, G. Garrigos, and J. Peypouquet.  Splitting methods with variable metric for Kurdyka–Łojasiewicz functions and general convergence rates.  J. Opt. Theory Appl., 165:874–900, 2015.
  * [35] K. Kurdyka.  On gradients of functions definable in o-minimal structures.  Ann. Inst. Fourier, 48(3):769–783, 1998.
  * [36] G. Li and T. K. Pong.  Douglas–Rachford splitting for nonconvex optimization with application to nonconvex feasibility problems.  Math. Program., 159:371–401, 2016.
  * [37] Y. Nesterov.  Smooth minimization of non-smooth functions.  Math. Program., 103(1):127–152, May 2005.
  * [38] D. Noll.  Convergence of non-smooth descent methods using the Kurdyka–Łojasiewicz inequality.  J. Opt. Theory Appl., 160(2):553–572, 2014.
  * [39] P. Ochs.  Unifying abstract inexact convergence theorems and block coordinate variable metric iPiano.  SIAM J. Optim., 29(1):541–570, 2019.
  * [40] P. Ochs, Y. Chen, T. Brox, and T. Pock.  iPiano: Inertial proximal algorithm for non-convex optimization.  SIAM J. Imaging Sci., 7(2):1388–1419, 2014.
  * [41] B. Polyak.  Introduction to optimization.  Optimization Software - Inc., Publication Division, New York, 1987.
  * [42] B. T. Polyak.  Some methods of speeding up the convergence of iteration methods.  USSR Comput. Math. Math. Phys., 4:1–17, 1964.
  * [43] A. Repetti and E. Chouzenoux.  RestoVMFB Lab: Matlab Toolbox for image restoration with the variable metric forward-backward algorithm.  http://www-syscom.univ-mlv.fr/~chouzeno/Logiciel.html, 2013.
  * [44] R. T. Rockafellar.  Convex Analysis.  Princeton University Press, Princeton, NJ, 1970.
  * [45] R. T. Rockafellar, R. J.-B. Wets, and M. Wets.  Variational Analysis, volume 317 of Grundlehren der Mathematischen Wissenschaften.  Springer, Berlin, 1998.
  * [46] L.I. Rudin, S. Osher, and E. Fatemi.  Nonlinear total variation based noise removal algorithms.  J. Phys. D., 60(1–4):259–268, November 1992.
  * [47] S. Salzo and S. Villa.  Inexact and accelerated proximal point algorithms.  J. Convex Anal., 19(4):1167–1192, December 2012.
  * [48] K. Scheinberg, D. Goldfarb, and X. Bai.  Fast first–order methods for composite convex optimization with backtracking.  Found. Comput. Math., 14:389–417, 2014.
  * [49] L. Stella, A. Themelis, and P. Patrinos.  Forward-backward quasi-Newton methods for nonsmooth optimization problems.  Comput. Optim. Appl., 67:443–487, 2017.
  * [50] S. Villa, S. Salzo, L. Baldassarre, and A. Verri.  Accelerated and inexact forward-backward algorithms.  SIAM J. Optim., 23(3):1607–1633, 2013.
  * [51] A. Zalinescu.  Convex analysis in general vector spaces.  World Scientific Publishing, Singapore, 2002.