# Stochastic Monotone Inclusion with Closed Loop Distributions

Hamza Ennaji, Jalal M. Fadili, Hedy Attouch

<figure><img src="image_1.png" /></figure>

## To cite this version:

Hamza Ennaji, Jalal M. Fadili, Hedy Attouch. Stochastic Monotone Inclusion with Closed Loop Distributions. Evolution Equations and Control Theory, 2025, 10.3934/eect.2025022. hal-04974921

HAL Id: hal-04974921

https://hal.science/hal-04974921v1

Submitted on 4 Mar 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
1 STOCHASTIC MONOTONE INCLUSION WITH CLOSED LOOP DISTRIBUTIONS

2 HAMZA ENNAJI<sup>c,*</sup>, JALAL FADILI<sup>‡</sup>, AND HEDY ATTOUCH<sup>◦</sup>

*Dedicated to the memory of Hedy Attouch, outstanding mathematician and beloved collaborator.*

**ABSTRACT.** In this paper, we study in a Hilbertian setting, first and second-order monotone inclusions related to stochastic optimization problems with decision-dependent distributions. The studied dynamics are formulated as monotone inclusions governed by Lipschitz perturbations of maximally monotone operators where the concept of equilibrium plays a central role. We discuss the relationship between the $W_1$-Wasserstein Lipschitz behavior of the distribution and the so-called coarse Ricci curvature. As an application, we consider the monotone inclusions associated with stochastic optimisation problems involving the sum of a smooth function with Lipschitz gradient, a proximable function and a composite term.

3 1. INTRODUCTION

4 Recently, many problems in machine learning and risk management come in form of stochastic optimisation problems. Such problems aim to learn a decision rule from a data sample. This can be formulated in terms of optimization problems of the form

$$ \min_{x \in \mathbb{R}^n} \mathbb{E}_{\xi \sim m}(f(x, \xi)) + g(x), \qquad (1) $$

7 where $m$ is a probability measure, $\mathbb{E}_{\xi \sim m}$ is the expectation operator with respect to the measure $m$, $f(x, \xi)$ is a loss function of the decision $x$ at data point $\xi$ and $g$ is a regularizer.

8 In this work, we are interested in (1) in the case where the distribution $m$ depends itself on the decision $x$, i.e., problems of the form

$$ \min_{x \in \mathbb{R}^n} \mathbb{E}_{\xi \sim m_x}(f(x, \xi)) + g(x), \qquad (2) $$

11 In this case, one tries to learn a decision rule from a decision-dependent data distribution $m_x$. Problems of the form (2) were addressed in the framework of performative prediction proposed in [42, 46] and discussed with further algorithmic aspects in [34]. A typical example concerns prediction of loan default risks, that is the chance that a borrower won't be able to repay their loan. More precisely, banks take into account several parameters, including the default risk, to decide whether to accept a consumer's loan application and, if so, what interest rate will apply. It is clear that a high default risk implies a high interest rate, but a high interest rate increases the consumer's default risk. Thus, the predictive performance of the bank's model is not calibrated with respect to future results obtained by acting on

<footnote><sup>c</sup> CORRESPONDING AUTHOR.</footnote>

<footnote><sup>*</sup> UNIV. GRENOBLE ALPES, CNRS, GRENOBLE INP*, LJK, 38000 GRENOBLE, FRANCE.</footnote>

<footnote><sup>‡</sup> ENSICAEN, NORMANDIE UNIVERSITÉ, CNRS, GREYC, FRANCE.</footnote>

<footnote><sup>◦</sup> IMAG, UNIVERSITÉ MONTPELLIER, CNRS, FRANCE.</footnote>

<footnote>E-mail addresses: hamza.ennaji@univ-grenoble-alpes.fr, jalal.fadili@ensicaen.fr, hedy.attouch@umontpellier.fr.</footnote>

<footnote n="2020">2020 Mathematics Subject Classification. Primary 34G25, 37N40, 46N10, 47H05, 49M30, 60J20.</footnote>

<footnote>Key words and phrases. Monotone inclusions; Dynamical systems; Gradient flows; Inertial dynamics; Ollivier-Ricci curvature.</footnote>
1 the model. Another example concerns navigation apps, such as Google Maps (see, e.g., [32, 40]), which suggest routes with low travel time to users. This influences users' decisions to pick such routes and consequently, increases traffic on these routes, impacting travel times. 2
4 Further applications and illustrations can be found in [46, Appendix A]. 5

6 In general, problems of the form (2) are difficult to solve. However, a natural approach consists in performing a repeated minimization procedure, i.e., throughout iterations, one solves 7

$$x_{t+1} \in \underset{x \in \mathbb{R}^n}{\operatorname{argmin}} \mathbb{E}_{\xi \sim m_{x_t}} (f(x, \xi)) + g(x), \quad (3)$$

8 and then updates the distribution $m_{x_{t+1}}$. Under suitable assumptions that will be specified later, the sequence $(x_t)_t$ generated by the repeated minimization procedure (3) admits a fixed point $\bar{x}$. Such a point turns out to be an equilibrium with respect to the distribution $m_{(.)}$ in the following sense: 9
10

11

$$\bar{x} \in \underset{x \in \mathbb{R}^n}{\operatorname{argmin}} \mathbb{E}_{\xi \sim m_{\bar{x}}} (f(x, \xi)) + g(x), \quad (4)$$

12 that is, $\bar{x}$ solves (2) for the induced distribution $m_{\bar{x}}$. So instead of solving (2), we look at an 13 equilibrium point in the sense of (4). Notice that in terms of operators, (4) can be written 14 (formally, for instance) as the monotone inclusion:

14

$$A(\bar{x}) + B(\bar{x}) \ni 0 \quad (5)$$

15 with $A(x) = \partial g(x)$ and $B(x) = \mathbb{E}_{\xi \sim m_{\bar{x}}} (\nabla f(x, \xi))$ ($\nabla$ being always the gradient w.r.t. to the variable $x$), where $A+B$ is monotone. That is $\bar{x}$ is a zero of the sum of two monotone operators. One strategy to solve problems of the form (5) is to consider some continuous and discrete dynamical systems whose trajectories may converge, under suitable assumptions, to an element in $(A+B)^{-1}(0)$, the zero set of the sum $A+B$. For instance, when $B \equiv 0$ and $A = \partial g$ where $g$ is a proper convex lower semicontinuous function on $\mathbb{R}^n$, it is well known since the works of Brézis, Baillon and Bucker [21, 28], that each trajectory of the subgradient flow 16
17
18
19
20
21
22

16

17

18

19

20

21

22

$$\dot{x}(t) + \partial g(x) \ni 0 \text{ with } x(0) = x_0 \in \mathbb{R}^n, \quad (6)$$

23

converges to a minimizer of $g$, and thus a zero of $\partial g$, provided $\operatorname{argmin} g \neq \emptyset$.

24

Designing algorithms and dynamical systems with rapid convergence properties to solve monotone inclusions is at the core of many fields in modern optimization, partial differential equation, game theory, etc. The literature is extensive, and to name only a few, the reader is referred to [1, 2, 7, 9, 17, 25, 37, 39, 50] and the references therein.

25

26

27

28

In this work, we address closed-loop differential inclusions of the form

29

$$\dot{x}(t) + A(x(t)) + \mathbb{E}_{\xi \sim m_{x(t)}} (B(x(t), \xi)) \ni 0, \quad (7)$$

30

where $A$ is a maximally monotone operator, and $B(\cdot, \xi)$ is a single-valued mapping for all $\xi \in \Xi$. The particularity of such a dynamic is, of course, the presence of the random operator $\mathbb{E}_{\xi \sim m_{(.)}} (B(\cdot, \xi))$ where the random variable $\xi$ has a trajectory-dependent distribution $m_{x(t)}$. Thus, it is not straightforward how to address (7) within the classical framework (see, e.g., [22, 26]). Yet, a clever reformulation of (7), based on the notion of equilibrium, as a monotone inclusion governed by a Lipschitz perturbation of a maximally monotone operator will allow us to tackle this issue. Then we investigate inertial dynamics related to (7). Indeed, since the work of Polyak [47], who considered a system of the form

31

32

33

34

35

36

37

$$\ddot{x}(t) + \gamma \dot{x}(t) + \nabla f(x(t)) = 0, \quad (\text{HBF})$$
1 where $\gamma > 0$ is called the viscous damping coefficient, the introduction of inertial dynamics to accelerate optimization methods has gained a lot of attention and led to many developments (see, e.g., [3, 10, 15, 18, 19, 50] and the references therein).

4 In this paper, we then consider second-order dynamics of the form

$$ \ddot{x}(t) + \gamma(t)\dot{x}(t) + \nabla f_{m_{\bar{x}}}(x(t)) + \omega\nabla^2 f_{m_{\bar{x}}}(x(t))\dot{x}(t) + e_{\bar{x}}(x(t)) + \omega \frac{d}{dt} e_{\bar{x}}(x(t)) = 0, \quad (8) $$

5 where $f_{m_{\bar{x}}}(x) = E_{\xi \sim m_{\bar{x}}}(f(x, \xi))$, $e_{\bar{x}}: \mathcal{H} \to \mathcal{H}$ is a perturbation operator and $\omega$ is the so-called Hessian-driven damping coefficient. When $\gamma(t) \equiv \gamma$, $f_{m_{\bar{x}}} = f$ (i.e., without a stochastic structure) and $e_{\bar{x}}(x) = 0$, systems of the form (8) were first studied in [3]. Later, this system was combined with an asymptotic vanishing damping $\gamma(t) = \frac{\alpha}{t}$, for $\alpha > 0$ in [19]. Several recent studies have been devoted to this topic (see, e.g., [9, 25, 37, 39, 49]).

10 1.1. Statement of the problem. Throughout, $\mathcal{H}$ is a real Hilbert space endowed with the scalar product $\langle \cdot, \cdot \rangle$ and induced norm $\| \cdot \|$, and $\Xi$ is a Polish space, i.e., separable and completely metrizable.

13 We consider the closed-loop differential inclusion

$$ \left\{ \begin{array}{l} \dot{x}(t) + A(x(t)) + E_{\xi \sim m_{x(t)}}(B(x(t), \xi)) \ni 0, \text{ a.e } t > t_0 > 0, \\ x(t_0) = x_0 \in \overline{\text{dom}(A)}, \end{array} \right. \qquad (9) $$

14 where $m_x$ is a family of probability distributions on $\Xi$ indexed by $x \in \mathcal{H}$. We will work under the standing assumption:

16 Assumption 1.

17 $A: \mathcal{H} \Rightarrow \mathcal{H}$ is a set-valued maximal monotone operator such that $\text{int}(\text{dom}(A)) \neq \emptyset$;

18 $B: \mathcal{H} \times \Xi \to \mathcal{H}$ is single-valued with $\xi \mapsto B(x, \xi)$ $m_x$-measurable, and $\exists \beta > 0$ such that $\xi \mapsto B(x, \xi)$ is $\beta$-Lipschitz continuous for every $x \in \mathcal{H}$;

20 Observe that when $x \mapsto m_x(C)$ is a measurable function on $\mathcal{H}$ for each fixed $C \in \mathcal{B}$, where $\mathcal{B}$ is a countably generated $\sigma$-field on $\Xi$, $(m_x)_{x \in \mathcal{H}}$ can be viewed as a random walk on $(\mathcal{H} \times \Xi, \mathcal{F} \otimes \mathcal{B})$ where $\mathcal{F}$ is a $\sigma$-field on $\mathcal{H}$; see Section 5.

23 ▶ Example 1.1. 1. Typically, (9) covers the case of stochastic optimization problems with a state-dependent distribution studied in [34, 46] by taking $A = \partial g$ and $B = \nabla f$ for $g \in \Gamma_0(\mathcal{H})$ and $f \in C^1(\mathcal{H} \times \Xi)$ whose gradient $\nabla f(x, \cdot)$ is $\beta$-Lipschitz continuous for every $x \in \mathcal{H}$. The last assumption ensures that $E_{\xi \sim m}(f(\cdot, \xi))$ is $C^1(\mathcal{H})$ whose gradient is $E_{\xi \sim m}(\nabla f(\cdot, \xi))$ which is $\beta$-Lipschitz continuous.

28 2. Taking $A = N_K$ the normal cone of a nonempty closed convex set $K$ of admissible decisions, we recover the framework of variational inequalities addressed recently in [31].

31 To simplify the presentation, we set, for any measure $m \in \mathcal{P}(\Xi)$

$$ B_m(x) = E_{\xi \sim m}(B(x, \xi)) \quad \text{and} \quad F_m(x) = A(x) + B_m(x). \qquad (10) $$

32 Using this notation, the system (9) can be simply rewritten as

$$ \left\{ \begin{array}{l} \dot{x}(t) + F_{m_{x(t)}}(x(t)) \ni 0, \text{ a.e } t > t_0, \\ x(t_0) = x_0 \in \overline{\text{dom}(A)}. \end{array} \right. \qquad (\text{SMI}) $$

33 The acronym (SMI) for the above dynamic stands for Stochastic Monotone Inclusion. Though "monotone" may seem as an abuse of terminology because the measure $m_\cdot$ depends on the trajectory, and $B_m$ is not even monotone. We will show later that this terminology is
1. still justified as (SMI) can be reformulated as a Lipschitzian perturbation of a monotone inclusion (see Section 3.2).

2. **Contributions and organization of the paper.** The paper is organized as follows. In Section 3 we address first-order monotone inclusions with closed-loop distributions. We prove the existence of equilibria as well as the well-posedness of the dynamics and convergence properties of the trajectories. In Section 4 we study asymptotic convergence properties of the trajectories of second-order dynamics with closed-loop distributions via viscous and Hessian damping. This allows us in particular to cover problems of the form (2). Section 5 contains a discussion concerning the Lipschitz behavior of the family $(m_x)_{x \in \mathcal{H}}$ with respect to the Wasserstein distance and some consequences in the framework of Markov chains on metric random walk spaces. In Section 6 we discuss the inertial primal-dual algorithm as an application of our results. Finally, Section 7 contains some conclusions and discusses some future works.

## 2. NOTATION AND PRELIMINARIES

In this section we fix some notation and present some notions and results that will be used.

2.1. **Convex analysis.** The domain of a function $g$ on $\mathcal{H}$ is defined by $\mathrm{dom}(g) = \{x \in \mathcal{H} : g(x) < \infty\}$. We denote by $\Gamma_0(\mathcal{H})$ the class of proper (bounded from below and $\mathrm{dom}(g) \neq \emptyset$), lower semicontinuous (l.s.c) and convex functions on $\mathcal{H}$ with values in $\mathbb{R} \cup \{+\infty\}$. We say that $g$ is $\alpha$-strongly convex, for $\alpha > 0$, if $g - \frac{\alpha}{2} \|\cdot\|^2$ is convex.

The subdifferential of $g$ is defined as

$$ \partial g : x \in \mathcal{H} \mapsto \{v \in \mathcal{H} : g(y) \geq g(x) + \langle v, y - x \rangle\}. $$

We recall the following Fermat's optimality condition for $g \in \Gamma_0(\mathcal{H})$,

$$ 0 \in \partial g(x^*) \Leftrightarrow x^* \in \mathrm{argmin}_{\mathcal{H}} g(\mathcal{H}). $$

**Definition 1** (Differentiability). Let $g: \mathcal{H} \to \mathbb{R} \cup \{\infty\}$ and $x \in \mathrm{int}(\mathrm{dom}(g))$. We say that $g$ is (Fréchet) differentiable at $x$ if there exists $v \in \mathcal{H}$ such that

$$ \lim_{h \to 0} \frac{g(x+h) - g(x) - \langle v, h \rangle}{\|h\|} = 0. $$

The unique vector $v$ satisfying this condition is the gradient of $g$ at $x$ denoted by $\nabla g(x)$.

If $g \in \Gamma_0(\mathcal{H})$ and differentiable at $x$, then $\partial g(x) = \{\nabla g(x)\}$.

We also recall the following.

**Definition 2** (L-smoothness). Let $L \ge 0$ and $g: \mathcal{H} \to \mathbb{R} \cup \{\infty\}$. We say that $g$ is L-smooth over $D \subset \mathcal{H}$ if it is differentiable over $D$ and

$$ \|\nabla g(x) - \nabla g(y)\| \le L \|x - y\| \quad \text{for any } x, y \in D. $$

We denote by $C_L^{1,1}(D)$ the class of L-smooth functions over $D$.

Given a function $g: \mathcal{H} \to \mathbb{R} \cup \{+\infty\}$, it proximal mapping is defined through

$$ \mathrm{prox}_f(x) = \underset{y \in \mathcal{H}}{\operatorname{argmin}} \left\{ g(y) + \frac{1}{2} \|y - x\|^2 \right\} \quad \text{for any } x \in \mathcal{H}. $$

When $g = \iota_K$ the indicator function of a nonempty closed convex set $K \subset \mathcal{H}$, then $\mathrm{prox}_f = \mathrm{proj}_K$, the projector onto $K$. For further details and notion, we refer the reader to [23].
2.2. Operator theory. The domain of the set-valued operator $A: \mathcal{H} \Rightarrow \mathcal{H}$ is dom($A$) = $\{x \in \mathcal{H}: A(x) \neq \emptyset\}$, its graph is gra($A$) = $\{[x, u] \in \mathcal{H} \times \mathcal{H}: u \in Ax\}$ and its zeros set is zer($A$) = $\{x \in \mathcal{H}: 0 \in A(x)\} := A^{-1}(0)$. A selection of $A$ is an operator $T: \text{dom } A \to \mathcal{H}$ such that, $Tx \in Ax$ for any $x \in \text{dom } A$. We write $A: \mathcal{H} \to \mathcal{H}$ to indicate that $A$ is single-valued. In the following, we gather some main properties that are essential for the rest of the paper.

**Definition 3.** • We say that $A$ is $\beta$-Lipschitz continuous if it is single-valued over its domain and

$$ \|Ax - Ay\| \le \beta\|x - y\|, \quad \forall x, y \in \text{dom } A. \qquad (11) $$

• We say that $A: \mathcal{H} \Rightarrow \mathcal{H}$ is monotone if

$$ \langle x - y, u - v \rangle \ge 0, \quad \forall [x, u], [y, v] \in \text{gra } A. \qquad (12) $$

• We say that $A$ is maximal monotone if there exists no monotone operator $B$, i.e., satisfying (12), such that $\text{gra } A \subset \text{gra } B$.

• We say that $A$ is uniformly monotone with modulus $\phi: [0, \infty) \to [0, \infty)$ if $\phi$ is increasing, $\phi(0) = 0$, $\lim_{t \to \infty} \phi(t) = \infty$ and

$$ \langle x - y, u - v \rangle \ge \|x - y\|\phi(\|x - y\|), \quad \forall [x, u], [y, v] \in \text{gra } A. \qquad (13) $$

• We say that $A$ is $\mu$-strongly monotone, with $\mu > 0$, if

$$ \langle x - y, u - v \rangle \ge \mu \|x - y\|^2, \quad \forall [x, u], [y, v] \in \text{gra } A. \qquad (14) $$

**Remark 2.1.** • Note that if $A$ is $\mu$-strongly monotone is equivalent to saying that $A - \mu\text{Id}$ is monotone.

• The definition of uniform monotonicity given by (13) is slightly different from the one in [23, Definition 22.1].

• If $A$ is $\mu$-strongly monotone, then it is uniformly monotone with modulus $\phi(t) = \mu t$.

**Example 2.1.** • The typical example of a maximal monotone operator is the subdifferential $\partial g$ of a function $g \in \Gamma_0(\mathcal{H})$. We usually refer to such an operator as a subpotential maximal monotone operator.

2.3. Monotone inclusions. Let $A$ be a maximal monotone operator on $\mathcal{H}$ and a single-valued mapping $D: [t_0, +\infty[ \times \text{dom}(A) \to \mathcal{H}$ and consider the following differential inclusion

$$ \begin{cases} \dot{x}(t) + A(x(t)) + D(t, x(t)) \ni 0, & t \in [t_0, T] \\ x(t_0) = x_0 \in \text{dom}(A). \end{cases} \qquad (15) $$

**Definition 4.** We will say that $x: [t_0, T] \to \mathcal{H}$ is a strong solution trajectory on $[t_0, T]$ of (15) if the following properties are satisfied:

(a) $x$ is continuous on $[t_0, T]$ and absolutely continuous on any compact subset of $]t_0, T[$ (hence almost everywhere differentiable);

(b) $x(t) \in \text{dom}(A)$ for almost every $t \in ]t_0, T[$, and (15) is verified for almost every $t \in ]t_0, T[$.

A trajectory $x: [t_0, +\infty[ \to \mathcal{H}$ is a strong global solution of (15) if it is a strong solution on $[t_0, T]$ for any $T > t_0$.

For further details, we refer the reader to the classical monographs [26] or [22].
2.4. **Transportation distance.** Denote $\mathcal{P}(\Xi)$ the space of probability measures on $\Xi$. For $\mathbf{m}_1, \mathbf{m}_2 \in \mathcal{P}(\Xi)$, the $W_1$-Wasserstein distance is defined by

$$W_1(\mathbf{m}_1, \mathbf{m}_2) = \sup_{h \in \mathrm{Lip}_1} |\mathrm{E}_{\xi \sim \mathbf{m}_1} h(\xi) - \mathrm{E}_{\zeta \sim \mathbf{m}_2} h(\zeta)|, \quad (16)$$

where $\mathrm{Lip}_1$ is the space of 1-Lipschitz continuous functions $h: \Xi \to \mathbb{R}$.

3. FIRST-ORDER MONOTONE INCLUSIONS

In this section we perform the analysis of the first-order monotone inclusion (9). More precisely, we discuss the existence and uniqueness of solutions as well as the convergence of trajectories. Recall that the dynamic (9) is governed by the operator $F_{\mathbf{m}_x} = A + B_{\mathbf{m}_x}$. We first prove the existence of an equilibrium point $\bar{x}$ which will allow us to reformulate (9) in a suitable form.

The following assumption is essential for the convergence analysis. It describes the sensitivity of the distribution to shifts in the index (here trajectory). It is widely used in the literature (see, e.g., [34, 42, 46, 52]). We will discuss how it closely relates to the so-called coarse Ricci curvature in Section 5.

**Assumption 2 (Lipschitz distributions).** There exists $\tau > 0$ such that

$$W_1(\mathbf{m}_x, \mathbf{m}_y) \le \tau \|x - y\|, \text{ for all } x, y \in \mathcal{H}.$$

The following two assumptions are standard monotonicity assumptions and will be crucial for the well-posedness of the dynamics and to prove existence and uniqueness of equilibria.

**Assumption 3 (Strong monotonicity).** $\exists \mu > 0$ such that for all $x \in \mathcal{H}$, $F_{\mathbf{m}_x}$ is $\mu$-strongly monotone.

The strong monotonicity Assumption 3 can be weakened to uniform monotonicity in the following sense.

**Assumption 4 (Uniform monotonicity).** There exists a function $\phi$ satisfying

$$\phi(t) > \beta\tau t, \forall t > 0,$$ such that for all $x \in \mathcal{H}$, $F_{\mathbf{m}_x}$ is uniformly monotone with a modulus $\phi$.

Finally, define the following parameter $\rho := \frac{\beta\tau}{\mu}$. As we will see, (see also [34, 46]), the parameter regime $\rho < 1$ will play a crucial role in the analysis of the convergence of the trajectories.

3.1. Existence and uniqueness of equilibria.

**Definition 5.** (Equilibrium point) We say that $\bar{x} \in \mathcal{H}$ is at equilibrium with respect to the family of probability measures $(\mathbf{m}_x)_{x \in \mathcal{H}}$ if

$$0 \in F_{\mathbf{m}_x}(\bar{x}). \quad (17)$$

In case $A = \partial g$ and $B = \nabla f$ where $g \in \Gamma_0(\mathcal{H})$ and $f(x, \cdot) \in C_{\beta}^{1,1}(\Xi)$, this definition is to be compared to the one introduced in [46] (see also [34]). Indeed, (17) reduces to:

$$\bar{x} \in \underset{x \in \mathcal{H}}{\operatorname{argmin}} \mathrm{E}_{\xi \sim \mathbf{m}_x} (f(x, \xi)) + g(x). \quad (18)$$

Solutions of (18) are exactly the fixed points of the repeated minimization procedure, that is, starting from some $x_0$, we generate the following sequence for $t \ge 0$

$$x_{t+1} = S(x_t) := \underset{x \in \mathcal{H}}{\operatorname{argmin}} \mathrm{E}_{\xi \sim \mathbf{m}_{x_t}} (f(x, \xi)) + g(x). \quad (19)$$
In [46, Theorem 3.5] it is shown that if $f$ is $C^1$ in both variables, $\xi \mapsto \nabla f(x, \xi)$ is $\beta$-Lipschitz and $E_{\xi \sim m_x}(f(\cdot, \xi))$ is $\mu$-strongly convex for all $x \in \mathcal{H}$ with $\rho < 1$, then, under Assumption 2, the iterates of (19) converge to a unique stable point. Their proof is essentially based on a fixed point argument. In [46, Proposition 4.1], they show the existence of equilibrium points under weaker assumptions on the loss $f$. Specifically, they demonstrate that if $f$ is convex and jointly continuous, then equilibrium points exist provided $\text{dom}(g)$ is compact. However, in this case the equilibrium is not necessarily unique. In the following lemma, we show the existence of equilibrium in the sense of (17).

**Theorem 1 (Existence and uniqueness of equilibrium point).** Under Assumption 4, the map

$$S: x \in \mathcal{H} \mapsto \text{zer}(F_{m_x}) = \{u \in \mathcal{H} : 0 \in F_{m_x}(u)\},$$

is a contraction. In particular, the equilibrium $\bar{x}$ is unique. If moreover, Assumption 3 holds instead, i.e., $F_{m_x}$ is $\mu$-strongly monotone for $\mu > 0$, the mapping $S$ is $\rho$-Lipschitz with $\rho := \frac{\beta\tau}{\mu}$. Thus for $\rho < 1$, there exists a unique equilibrium point $\bar{x}$.

**Proof.** First, we see that $S$ is well defined. Indeed, for any $x \in \mathcal{H}$, $\text{zer}(F_{m_x})$ is nonempty, and is in fact a singleton due to Assumption 4. To see this, we argue as in [23, Proposition 22.11]. Fix $[y_0, u_0] \in \text{gra } F_{m_x}$. We have for any $[y, u] \in \text{gra } F_{m_x}$:

$$ \begin{aligned} \|y - y_0\| \|u\| \ge \langle y - y_0, u \rangle &= \langle y - y_0, u - u_0 \rangle + \langle y - y_0, u_0 \rangle \\ &\ge \|y - y_0\| \phi (\|y - y_0\|) - \|y - y_0\| \|u_0\|. \end{aligned} \quad (20) $$

Since $\lim_{t \to \infty} \phi(t) = \infty$ we infer that $\inf_{u \in \text{gra } F_{m_x}(y)} \|u\| \to \infty$ as $\|y\| \to \infty$ and thus $F_{m_x}$ is surjective (see [23, Corollary 21.25]). Moreover, in view of strict monotonicity of $F_{m_x}$, $\text{zer}(F_{m_x})$ is a singleton (see, e.g., [23, Proposition 23.35]).

Now, pick $x, y \in \mathcal{H}$. We have that $0 \in F_{m_x}(S(x))$ and $0 \in F_{m_y}(S(y))$. In particular, $-B_{m_y}(S(y)) \in A(S(y))$, which gives that $B_{m_x}(S(y)) - B_{m_y}(S(y)) \in F_{m_x}(S(y))$. We get, thanks to Assumption 4

$$ \|S(x) - S(y)\| \phi (\|S(x) - S(y)\|) \le \langle u - v, S(x) - S(y) \rangle, \text{ for any } (u, v) \in F_{m_x}(S(y)) \times F_{m_x}(S(x)). $$

Then, taking $u = B_{m_x}(S(y)) - B_{m_y}(S(y))$ and $v = 0$, we get, using Cauchy-Schwarz inequality

$$ \phi (\|S(x) - S(y)\|) \le \|B_{m_x}(S(y)) - B_{m_y}(S(y))\|. $$

We get, using Corollary 1 below

$$ \phi (\|S(x) - S(y)\|) \le \beta\tau \|x - y\|. \quad (21) $$

Since $\phi$ is strictly increasing, the last inequality gives, thanks to Assumption 4,

$$ \|S(x) - S(y)\| \le \phi^{-1} (\beta\tau \|x - y\|) < \|x - y\|, $$

and by the Banach fixed-point theorem (see Theorem 6), $S$ has a unique fixed point $\bar{x}$.

Now if $F_{m_x}$ is $\mu$-strongly monotone, it is in particular uniformly monotone with modulus $\phi(t) = \mu t$. Equation (21) gives

$$ \|S(x) - S(y)\| \le \rho \|x - y\|, $$

with $\frac{\beta\tau}{\mu} := \rho$. Again, we conclude using Theorem 6. ■

**Remark 3.1.** In the strongly monotone case, Assumption 4 incorporates the parameter regime $\rho < 1$ which appears in particular in [46, Theorem 3.5].
3.2. **Well-posedness.** Before stating the main result of this section, let us fix some extra notation and properties. Assume that Assumption 4 holds, and denote by

2

$$ \varphi(t) = \phi(t) - \beta\tau t, \qquad (22) $$

where $\phi$ is the modulus of uniform monotonicity of $F_{m_x}$.

We prove the following.

**Lemma 1.** Let $a > 0$ and define

5

$$ \theta(z) := \int_z^a \frac{ds}{\varphi(s)}. \qquad (23) $$

Then $\theta$ is nonincreasing and $\lim_{z \to 0^+} \theta(z) = \infty$.

Proof. Indeed, we have $\dot{\theta}(z) = -1/\varphi(z) < 0$ since $\phi$ satisfies Assumption 4. Moreover, since $\phi(t) \ge \varphi(t)$, and $\lim_{z \to 0^+} \int_z^a \frac{ds}{\varphi(s)} = \infty$, the result follows. ■

Remark 3.2. In the literature of ordinary differential equations, the above lemma is related to the fact that $\varphi$ is somehow an Osgood modulus of continuity (see, e.g., [20, Definitions 2.108 and 3.1]). If $\varphi(s) = s$, which corresponds to Lipschitz regularity, and $a = 1$ then $\theta(z) = \log_+(z) = \max\{0, \log(1/z)\}$. If $a = 1/e$ and $\varphi(s) = s \log(1/s)$, which corresponds to log-Lipschitz regularity, then $\theta(z) = \log \log_+(z)$. More generally, $\varphi(s) = s (\log(1/s))^r$ for $r \le 1$ are admissible choices.

10

11

12

13

14

Now let us define the following gap

15

$$ \mathbf{e}_{\bar{x}}(x) = B_{m_x}(x) - B_{m_x}(x). \qquad (24) $$

Using the notation in (10), we may rewrite (9) in the following form

16

$$ \left\{ \begin{array}{l} \dot{x}(t) + F_{m_x}(x(t)) + \mathbf{e}_{\bar{x}}(x(t)) \ni 0, \text{ a.e } t > t_0 \\ x(t_0) = x_0. \end{array} \right. \qquad (\text{p-SMI}) $$

One advantage of this formulation is that the mapping $x \mapsto \mathbf{e}_{\bar{x}}(x)$ exhibits Lipschitz behaviour (see Lemma 4), and now only the operator $F_{m_x}$ appears instead of $F_{m_x(t)}$. This allows us to treat (p-SMI) within the framework of evolution equations governed by Lipschitz perturbations of maximal monotone operators (cf. [26, Chapter III]). This is behind the notation (p-SMI), which stands for perturbed stochastic monotone inclusion. This being said, our aim is to use [26, Proposition 3.13] and show the existence of a unique strong solution (cf. Definition 4) to (p-SMI) (see also [26, Definition 3.1]).

17

18

19

20

21

22

23

We start with the following properties.

24

**Lemma 2.** Under Assumption 1, we have, for any $\mathbf{m}, \nu \in \mathcal{P}(\Xi)$ and $x \in \mathcal{H}$

25

$$ \sup_{x \in \mathcal{H}} \|B_{\mathbf{m}}(x) - B_{\nu}(x)\| \le \beta W_1(\mathbf{m}, \nu). $$

Proof. Let $v \in \mathcal{H}$ be a unit norm vector. By Assumption 1, $\langle v, B(x, \cdot) \rangle$ is $\beta$-Lipschitz continuous for every $x \in \mathcal{H}$. We then have from the definition of the $W_1$-Wasserstein distance (16)

26

27

$$ \begin{aligned} \|B_{\mathbf{m}}(x) - B_{\nu}(x)\| &= \sup_{v \in \mathcal{H}, \|v\|=1} \langle v, B_{\mathbf{m}}(x) - B_{\nu}(x) \rangle \\ &= \sup_{v \in \mathcal{H}, \|v\|=1} \mathbb{E}_{\xi \sim \mathbf{m}} \langle v, B(x, \xi) \rangle - \mathbb{E}_{\zeta \sim \nu} \langle v, B(x, \zeta) \rangle \le \beta W_1(\mathbf{m}, \nu). \end{aligned} $$

Taking the supremum over $x \in \mathcal{H}$ yields the result. ■

29

<figure><img src="image_3.png" /></figure>

Combining Lemma 2 and Assumption 2 we get the following.

30
1 **Corollary 1.** Under Assumption 1 and Assumption 2, for all $y, z \in \mathcal{H}$

$$ \sup_{x \in \mathcal{H}} \|B_{m_y}(x) - B_{m_z}(x)\| \le \beta\tau \|y - z\|. $$

2 **Assumption 5** (Lipschitz continuity of $B_{m_x}$). $\exists L > 0$ such that $B_{m_x}$ is $L$-Lipschitz continuous for every $x \in \mathcal{H}$.

4 This assumption is true if for instance $B(\cdot, \xi)$ is $L$-Lipschitz continuous for $m_x$-almost every $\xi \in \Xi$. Indeed

5 $B_{m_x}(z) - B_{m_x}(y) = E_{\xi \sim m_x} (B(z, \xi) - B(y, \xi)) \le L \|z - y\|$. We are now in position to characterize the properties of $F_{m_x}$.

7 **Lemma 3.** Assume that Assumptions 1 and 5 hold and that $B_{m_x}$ is monotone. Then the operator $F_{m_x}$ is maximally monotone. If, in addition, Assumption 3 (resp. Assumption 4) then $F_{m_x}$ is maximally uniformly (resp. strongly) monotone.

10 **Proof.** Maximal monotonicity of $B_{m_x}$ follows from [23, Corollary 20.28]. Combine this with maximal monotonicity of A stated in Assumption 1, and [26, Lemma 2.4] yields the claim. The proof of the last statement is immediate. ■

13 The key ingredient to use [26, Proposition 3.13] is the Lipschitz continuity of the perturbation $e_{\bar{x}}(\cdot)$. This is the content of the following statement.

15 **Lemma 4.** Suppose that Assumptions 1, 2 and 5 hold, then $e_{\bar{x}}(\cdot)$ is $(2L + \beta\tau)$-Lipschitz continuous.

17 **Proof.** Let $x, z \in \mathcal{H}$. We then have

$$ \begin{aligned} \|e_{\bar{x}}(x) - e_{\bar{x}}(z)\| &= \|(B_{m_x}(x) - B_{m_x}(z)) - (B_{m_z}(z) - B_{m_x}(z))\| \\ &= \|(B_{m_x}(x) - B_{m_x}(z)) + (B_{m_x}(z) - B_{m_x}(x)) + (B_{m_x}(z) - B_{m_z}(z))\| \\ &\le L\|x - z\| + L\|x - z\| + \beta\tau\|x - z\| \\ &= (2L + \beta\tau)\|x - z\|, \end{aligned} \quad (25) $$

18 where we have used Assumption 5 twice and Corollary 1 once in the inequality. ■

19 **Remark 3.3.** Since $\mathrm{zer}(F_{m_x}) = \{\bar{x}\}$, we clearly see that $\bar{x} \in \mathrm{zer}(F_{m_x} + e_{\bar{x}})$. Indeed, taking into account (24), we have $e_{\bar{x}}(\bar{x}) = 0$ so that $0 \in (F_{m_x} + e_{\bar{x}})(\bar{x})$.

21 **Proposition 1.** Suppose that Assumptions 1, 2 and 5 hold, that either Assumption 3 or Assumption 4 hold and let $T > t_0$. Then, given $x_0 \in \mathrm{dom}(F_{m_x})$, the dynamic (p-SMI), and hence (SMI), admits a unique strong solution $x$ in the sense of Definition 4.

24 **Proof.** Thanks to Lemma 3 and Lemma 4 $F_{m_x}$ is maximally monotone and $e_{\bar{x}}$ is Lipschitz continuous. Since $t \mapsto e_{\bar{x}}(x)$ is trivially in $L^\infty(t_0, T; \mathcal{H})$ for any $x \in \mathcal{H}$, we deduce, thanks to [26, Proposition 3.13] the existence of a unique solution $x \in W^{1,1}(t_0, T; \mathcal{H})$ to (SMI). ■

27 **Remark 3.4.** Evolution problems of the form (p-SMI) were also studied in [13], where the operator $e_{\bar{x}}$ is potentially multivalued and depends on both time and space variables. While standard time-variable measurability is assumed, the authors' analysis relied on the upper semicontinuity of $x \mapsto e_{\bar{x}}(t, x)$, rather than Lipschitz continuity, for which weaker results were proved. We also observe that the domain condition in Assumption 1 can be removed in finite dimension (see e.g., [13, Theorem 1.2]) or when A is the subdifferential of a function in $\Gamma_0(\mathcal{H})$.
1. **3.3. Convergence properties.** We now establish the main convergence result of the unique solution trajectory of (SMI).

2. **Theorem 2.** Let $x$ be the unique solution trajectory of (SMI) under Assumptions 1, 2 and 5, and suppose that Assumption 4 also holds. Then, for all $t \ge t_0$, we have

$$ \|x(t) - \bar{x}\| \le \theta^{-1}(2t - \hat{t}), \qquad (26) $$

5. where $\theta$ is defined in Lemma 1 and $\hat{t} = 2t_0 - \theta (\|x_0 - \bar{x}\|)$.

6. Furthermore, if Assumption 3 holds instead and $\rho := \frac{\beta\tau}{\mu} < 1$, we have

$$ \|x(t) - \bar{x}\| \le Ce^{-2\mu(1-\rho)t} \quad \text{for all } t \ge t_0, \qquad (27) $$

7. with $C = \|x_0 - \bar{x}\|e^{2\mu(1-\rho)t_0}$.

8. **Proof.** We observe thanks to Proposition 1 that

$$ \frac{1}{2} \frac{d}{dt} \|x(t) - \bar{x}\|^2 = \langle \dot{x}(t), x(t) - \bar{x} \rangle. $$

9. We then have, for any selection $u(t)$ of $F_{m_x}(x(t))$

$$ \begin{aligned} \frac{1}{2} \frac{d}{dt} \|x(t) - y\|^2 &= \langle \dot{x}(t), x(t) - \bar{x} \rangle = -\langle u(t), x(t) - \bar{x} \rangle - \langle \mathbf{e}_{\bar{x}}(x(t)), x(t) - \bar{x} \rangle \\ &= \langle u(t) - 0, x(t) - \bar{x} \rangle - \langle \mathbf{e}_{\bar{x}}(x(t)), x(t) - \bar{x} \rangle. \end{aligned} \qquad (28) $$

10. We then get by Assumption 4

$$ \langle u(t) - 0, x(t) - \bar{x} \rangle \ge \|x(t) - \bar{x}\| \phi (\|x(t) - \bar{x}\|). $$

11. Then, thanks to Corollary 1, (28) gives

$$ \begin{aligned} \frac{1}{2} \frac{d}{dt} \|x(t) - \bar{x}\|^2 &\le -\|x(t) - \bar{x}\| \phi (\|x(t) - \bar{x}\|) + \|\mathbf{e}_{\bar{x}}(x(t))\| \|x(t) - \bar{x}\| \\ &\le -\|x(t) - \bar{x}\| \phi (\|x(t) - \bar{x}\|) + \beta \tau \|x(t) - \bar{x}\|^2 \\ &= -\|x(t) - \bar{x}\| \varphi (\|x(t) - \bar{x}\|), \end{aligned} \qquad (29) $$

12. where $\varphi(t) = \phi(t) - \beta\tau t$ is as introduced in (22). When $h(t) := \|x(t) - \bar{x}\| = 0$, then there is

13. nothing to prove and we are done. Otherwise, $h(t) \ne 0$. Thus dividing both sides of (29) by

14. $h(t)$, we get

$$ \dot{h}(t) \le -2\varphi(h(t)). \qquad (30) $$

15. From (29) and (30), we infer that

$$ \frac{d}{dt} \theta(h(t)) = \frac{-1}{\varphi(h(t))} \dot{h}(t) \ge 2. \qquad (31) $$

16. Integrating (31) between $t_0$ and $t > 0$, we get

$$ \theta(h(t)) - \theta(h(t_0)) \ge 2(t - t_0). $$

17. Using Lemma 1, we deduce that

$$ h(t) \le \theta^{-1}(2t - \hat{t}), \qquad (32) $$

18. where $\hat{t} = 2t_0 - \theta(h(t_0))$. This proves (26), as desired.

19. To prove the second claim, we observe that if $F_{m_x}$ is $\mu$-strongly monotone, $\phi(t) = \mu t$ so that $\varphi(t) = (\mu - \beta\tau)t$. Since $\mu > \beta\tau$, we get $\theta(s) = \frac{-\log(s)}{\mu - \beta\tau}$. So that $\theta^{-1}(s) = e^{-(\mu - \beta\tau)s}$.

20. Plugging this into (32), we obtain

$$ h(t) \le \|x_0 - \bar{x}\| e^{-2(\mu - \beta\tau)(t - t_0)}, \qquad (33) $$

21. as desired.

<figure><img src="image_4.png" /></figure>
<footnote n="1">1 Remark 3.5. Observe that one can obtain (27) from (29) using Gronwall's Lemma 8. Indeed, taking $\phi(t) = \beta\mu t$, we infer from (29)</footnote>

<footnote n="2">2 $\frac{1}{2}\frac{d}{dt}\|x(t) - \bar{x}\|^2 \le -(\mu - \beta\tau)\|x(t) - \bar{x}\|^2$.</footnote>

<footnote n="3">3 By Assumption 2 and Theorem 2, we obtain a rate of convergence of the measure $m_{x(t)}$ to $m_{\bar{x}}$ in the $W_1$ distance as $t \to \infty$.</footnote>

<footnote n="4">4 Corollary 2. Let $x: [t_0, \infty) \to \mathbb{R}$ be the solution of (SMI) under the assumptions of Proposition 1, and suppose that Assumption 4 also holds. We have for all $t \ge t_0$</footnote>

<footnote n="5">5 $W_1(m_{x(t)}, m_{\bar{x}}) \le \tau\theta^{-1}(2t - \hat{t})$,</footnote>

<footnote n="6">6 where $\theta$ is defined in Lemma 1 and $\hat{t} = 2t_0 - \theta(\|x_0 - \bar{x}\|)$. If Assumption 3 holds instead and $\rho := \frac{\beta\tau}{\mu} < 1$, we have</footnote>

<footnote n="7">7 $W_1(m_{x(t)}, m_{\bar{x}}) \le Ce^{-2\mu(1-\rho)t}$,</footnote>

<footnote n="8">8 with $C = \tau\|x_0 - \bar{x}\|e^{2\mu(1-\rho)t_0}$.</footnote>

<footnote n="9">9 4. INERTIAL SECOND-ORDER SYSTEM WITH VISCOUS AND HESSIAN DAMPING</footnote>

<footnote n="10">10 In this section, to avoid technicalities, we restrict ourselves to the smooth convex optimization case where i.e., $A = \nabla g$ so that (recall (10))</footnote>

<footnote n="11">11 $F_m = \nabla G_m$ with $G_m := g + f_m$ and $f_m := E_{\xi \sim m}(f(\cdot, \xi))$. (34)</footnote>

<footnote n="12">12 In this case, Assumptions 1 and 5 read as follows:</footnote>

<footnote n="13">13 **Assumption 6.**</footnote>

<footnote n="14">14 $g \in C^1(\mathcal{H}) \cap \Gamma_0(\mathcal{H})$;</footnote>

<footnote n="15">15 $f \in C^1(\mathcal{H} \times \Xi)$, $f(x, \cdot)$ is $m_x$-measurable and $C_{\beta}^{1,1}(\Xi)$ for every $x \in \mathcal{H}$, and $f_{m_x} \in C_{L}^{1,1}(\mathcal{H})$ for every $x \in \mathcal{H}$.</footnote>

<footnote n="16">16 On the other hand, Assumption 3 specializes to</footnote>

<footnote n="17">17 **Assumption 7 (Strong convexity).** $\exists \mu > 0$ such that $G_{m_x}$ is $\mu$-strongly convex for every $x \in \mathcal{H}$.</footnote>

<footnote n="18">18 This last assumption is true for instance if $f(\cdot, \xi)$ is $\mu$-strongly convex for $m_x$-almost every $\xi \in \Xi$.</footnote>

<footnote n="19">19 Let us first start with a discussion of the following second-order in time ODE</footnote>

<footnote n="20">20 $\begin{cases} \ddot{x}(t) + \gamma(t)\dot{x}(t) + \nabla G_{m_x(t)}(x(t)) = 0, \text{ a.e } t \in [t_0, T] \\ x(t_0) = x_0 \in \mathcal{H}, \end{cases}$ (35)</footnote>

<footnote n="21">21 where $\gamma: \mathbb{R}^+ \to \mathbb{R}^+$ is a continuous function, usually referred to as the viscous damping coefficient. In a smooth convex optimization deterministic setting, i.e., $F_{m_x} = \nabla f$ where $f$ is a smooth, strongly convex function, systems of the form (35) were first studies by Polyak in [47] with a fixed viscous damping coefficient. Later on, this kind of systems were studied by [15] and then [50] where the authors establish the link between the continuous dynamics with $\gamma(t) = \frac{3}{t}$ and the Nesterov's gradient method [43]. This is a very active research field (see, e.g., [5, 11, 16, 38] and the references therein).</footnote>

<footnote>22</footnote>
Following what we have done Section 3, (35) can be reformulated again as a Lipschitzian perturbation of an inertial gradient system as follows

$$ \ddot{x}(t) + \gamma(t)\dot{x}(t) + \nabla G_{m_x}(x(t)) + \mathbf{e}_{\bar{x}}(x(t)) = 0. \quad (36) $$

Even in absence of perturbations, it is well-known known that the inertial system (36) may suffer from transverse oscillations, and it is desirable to attenuate them. This is precisely the motivation behind the introduction of geometric Hessian-driven damping in [3] (sometimes called Newton-type inertial dynamics). The inertial system we consider then features both viscous and Hessian-driven damping and reads

$$ \ddot{x}(t) + \gamma(t)\dot{x}(t) + \omega(t)\frac{d}{dt}(\nabla G_{m_x}(x(t))) + \nabla G_{m_x(t)}(x(t)) = 0, $$

where $\omega: [0, \infty) \to \mathbb{R}^+$ is a continuous function usually referred to as the Hessian-driven damping coefficient, and which will be taken to be constant, i.e., $\omega(t) \equiv \omega > 0$. The case $\omega = 0$ then recovers (35). Following the reasoning of Section 3, then under Assumptions 2, 6 and 7 this system is also equivalent to

$$ \ddot{x}(t) + \gamma(t)\dot{x}(t) + \nabla G_{m_x}(x(t)) + \omega \frac{d}{dt}(\nabla G_{m_x}(x(t))) + \mathbf{e}_{\bar{x}}(x(t)) + \omega \frac{d}{dt}\mathbf{e}_{\bar{x}}(x(t)) = 0 \quad \text{a.e. } t \in [t_0, T]. \quad (\text{ISEHD}_{\gamma, \omega}) $$

Recall that $\mathbf{e}_{\bar{x}}$ is Lipschitz continuous by Lemma 4, and thus a.e. differentiable so that the last term in system $(\text{ISEHD}_{\gamma, \omega})$ makes sense whenever $x$ is absolutely continuous.

The acronym ISEHD stands for Inertial System with Explicit Hessian Damping. The Hessian damping is said to be explicit since, when $G_{m_x}$ is of class $C^2$ and $x$ is smooth, we have

$$ \frac{d}{dt}(\nabla G_{m_x}(x(t))) = \nabla^2 G_{m_x}(x(t))\dot{x}(t). $$

Variants and generalizations of the inertial systems with Hessian-driven damping were studied by multiple authors (see, e.g., [9, 18, 19]). The study of the effect of perturbations on these systems, i.e., problems of the form $(\text{ISEHD}_{\gamma, \omega})$, was carried out in the recent work [14]. We will take inspiration from this work but our analysis will depart from that of [14] in some important ways. For instance, and most importantly, the perturbations in [14] depend only on time while they depend on the trajectory in this paper. This poses a few challenges that we tackle by exploiting the Lipschitz continuity property enjoyed by our perturbations (see for instance Lemma 4).

## 4.1. Well-posedness.

### 4.1.1. Equivalent first-order formulation.

**Proposition 2.** Suppose that Assumption 6 holds, that $\gamma(t) \ge 0, \omega > 0$ with $\gamma \in C^1([t_0, +\infty])$. For any initial conditions $(x_0, v_0) \in \mathcal{H} \times \mathcal{H}$, the dynamics $(\text{ISEHD}_{\gamma, \omega})$ admits an equivalent formulation of the form

$$ \left\{ \begin{array}{l} \dot{x}(t) + \omega \left( \nabla G_{m_x}(x(t)) + \mathbf{e}_{\bar{x}}(x(t)) \right) - \left( \frac{1}{\omega} - \gamma(t) \right) x(t) + \frac{1}{\omega} y(t) = 0 \\ \\ \dot{y}(t) - \left( \frac{1}{\omega} - \gamma(t) - \dot{\gamma}(t)\omega \right) x(t) + \frac{1}{\omega} y(t) = 0, \end{array} \right. \quad (37) $$

with initial conditions $x(t_0) = x_0, y(t_0) = -\omega(v_0 + \omega\nabla G_{m_x}(x_0)) + (1 - \omega\gamma(t_0))x_0 - \omega^2\mathbf{e}_{\bar{x}}(x_0)$.

*Proof.* Let $(x, y)$ be a solution of (37). By differentiation of the first equation in (37), we get

$$ \ddot{x}(t) + \omega \frac{d}{dt}(\nabla G_{m_x}(x(t))) + \omega \frac{d}{dt}\mathbf{e}_{\bar{x}}(x(t)) + \dot{\gamma}(t)x(t) - \left(\frac{1}{\omega} - \gamma(t)\right)\dot{x}(t) + \frac{1}{\omega}\dot{y}(t) = 0. $$
1. Replacing $\dot{y}$ by its expression from the second equation in (37), we obtain:

$$ 
\begin{aligned}
\ddot{x}(t) + \omega \frac{d}{dt}(\nabla G_{m_{\bar{x}}}(x(t))) + \omega \frac{d}{dt} \mathbf{e}_{\bar{x}}(x(t)) + \dot{\gamma}(t)x(t) \\
- \left(\frac{1}{\omega} - \gamma(t)\right) \dot{x}(t) + \frac{1}{\omega} \left(\left(\frac{1}{\omega} - \gamma(t) - \dot{\gamma}(t)\omega\right)x(t) - \frac{1}{\omega}y(t)\right) = 0,
\end{aligned}
$$

2. and using again the first equation in (37) to eliminate $y(t)$, we get:

$$ 
\begin{aligned}
\ddot{x}(t) + \omega \frac{d}{dt}(\nabla G_{m_{\bar{x}}}(x(t))) + \omega \frac{d}{dt} \mathbf{e}_{\bar{x}}(x(t)) + \dot{\gamma}(t)x(t) - \left(\frac{1}{\omega} - \gamma(t)\right) \dot{x}(t) \\
+ \frac{1}{\omega} \left( \left(\frac{1}{\omega} - \gamma(t) - \dot{\gamma}(t)\omega\right)x(t) + \dot{x}(t) + \omega \left(\nabla G_{m_{\bar{x}}}x(t) + \mathbf{e}_{\bar{x}}(x(t))\right) - \left(\frac{1}{\omega} - \gamma(t)\right)x(t) \right) = 0,
\end{aligned}
$$

3. and after simplifications, we recover $(\text{ISEHD}_{\gamma, \omega})$. Conversely, let $x$ be a trajectory solution

4. to $(\text{ISEHD}_{\gamma, \omega})$ with initial conditions $(x_0, v_0) \in \mathcal{H} \times \mathcal{H}$ and define

$$ y(t) = -\omega \left( \dot{x}(t) + \omega \left( \nabla G_{m_{\bar{x}}} x(t) + \mathbf{e}_{\bar{x}}(x(t)) \right) - \left( \frac{1}{\omega} - \gamma(t) \right) x(t) \right). $$

5. By differentiating the previous formula and using $(\text{ISEHD}_{\gamma, \omega})$, we recover the second equation

6. of (37), as desired.

■

7. In the above, we have taken the derivatives in time as if they exist. We will show shortly

8. that this is actually the case.

9. **4.1.2. Existence and uniqueness of a solution trajectory.** Proposition 2 opens the door to considering even the nonsmooth version of $(\text{ISEHD}_{\gamma, \omega})$ via (37) requiring only that $g \in \Gamma_0(\mathcal{H})$.

10. Thus, writing $\partial G_{m_{\bar{x}}}(x(t)) = \partial g + \mathcal{E}_{\xi \sim m_{\bar{x}}} (\nabla f(\cdot, \xi))$, we can consider the differential inclusion

$$ 
\left\{ 
\begin{array}{l}
\dot{x}(t) + \omega \left( \partial G_{m_{\bar{x}}}(x(t)) + \mathbf{e}_{\bar{x}}(x(t)) \right) - \left( \frac{1}{\omega} - \gamma(t) \right) x(t) + \frac{1}{\omega} y(t) \quad \ni 0 \\
\dot{y}(t) - \left( \frac{1}{\omega} - \gamma(t) - \dot{\gamma}(t)\omega \right) x(t) + \frac{1}{\omega} y(t) \quad = 0,
\end{array}
\right.
\quad (38)
$$

12. with initial conditions $(x(t_0) = x_0, y(t_0) = y_0) \in \text{dom}(G_{m_{\bar{x}}}) \times \mathcal{H}$.

13. One of the main advantages of (38) is that it can be easily recast as a differential inclusion

14. governed by a Lipschitz perturbation of a maximal monotone operator on the product space

15. $\mathcal{H} \times \mathcal{H}$. Indeed, setting $Z(t) = (x(t), y(t))$, $\mathcal{A}(x, y) = (\omega \partial G_{m_{\bar{x}}}(x), 0)$ and

$$ \mathcal{E}(t, x, y) = \left( \omega \mathbf{e}_{\bar{x}}(x(t)) - \left( \frac{1}{\omega} - \gamma(t) \right) x(t) + \frac{1}{\omega} y(t), - \left( \frac{1}{\omega} - \gamma(t) - \dot{\gamma}(t)\omega \right) x(t) + \frac{1}{\omega} y(t) \right), $$

16. we immediately see that (38) can be written as

$$ \dot{Z}(t) + \mathcal{A}(Z(t)) + \mathcal{E}(t, Z(t)) \ni 0_{\mathcal{H} \times \mathcal{H}}, \quad Z(0) = (x_0, y_0). \quad (39) $$

17. Observe that $\mathcal{A}$ is nothing but the subdifferential of the function $H(x, y) = \omega G_{m_{\bar{x}}}(x)$, and

18. thus $\mathcal{A}$ is maximal monotone if the smooth part is monotone (see Lemma 3). Then, (38) can

19. be written as the differential inclusion in the product space $\mathcal{H} \times \mathcal{H}$

$$ 
\left\{ 
\begin{array}{l}
\dot{Z}(t) + \mathcal{A}(Z(t)) + \mathcal{E}(t, Z(t)) \ni 0_{\mathcal{H} \times \mathcal{H}}, \text{ a.e } t \in [t_0, T] \\
Z(0) = (x_0, y_0) \in \text{dom}(g) \times \mathcal{H},
\end{array}
\right.
\quad (\text{MIS})
$$

20. which fits in the framework of Lipschitz perturbations of maximal monotone operators as in

21. Section 3. Notice that this formulation is different from the classical Hamiltonian one.
Before stating the main result, let us recall that we endow the product space with the scalar product $\langle(u, v), (u^*, v^*)\rangle_{H \times H} = \langle u, u^* \rangle + \langle v, v^* \rangle$, and the induced norm $\|(u, v)\|_{H \times H} = \sqrt{\|u\|^2 + \|v\|^2}$. We have the following auxiliary results

**Lemma 5.** Suppose that Assumptions 2 and 6 hold. Consider the operator $\mathcal{E}$ where $\omega > 0$ and $\gamma \in C^1([t_0, +\infty[)$. Then, for any $t \in [t_0, T]$, $\mathcal{E}(t, \cdot, \cdot)$ is Lipschitz continuous on $H \times H$.

The dependence on $t$ of the Lipschitz constant appears only through $|\gamma(t)|$ and $|\dot{\gamma}(t)|$.

With the two most popular choices of $\gamma(t)$, constant (as in heavy ball method) or the asymptotically vanishing viscous damping $\gamma(t) = \alpha/t$, both $|\gamma(t)|$ and $|\dot{\gamma}(t)|$ are uniformly bounded. This makes $\mathcal{E}(t, \cdot, \cdot)$ Lipschitz continuous uniformly in $t > t_0$.

**Proof.** Let $u, v, u^*, v^* \in H$ and set $p = (\frac{1}{\omega} - \gamma(t))$, $q = (\frac{1}{\omega} - \gamma(t) - \dot{\gamma}(t)\omega)$. We have, for $t \in [t_0, T]$

$$ \begin{aligned} & \|\mathcal{E}(t, u, v) - \mathcal{E}(t, u^*, v^*)\|_{H \times H} \\ &= \left\| \left( \omega(\mathbf{e}_{\bar{x}}(u) - \mathbf{e}_{\bar{x}}(u^*)) + p(u^* - u) + \frac{1}{\omega}(v - v^*), \frac{1}{\omega}(v - v^*) + q(u^* - u) \right) \right\|_{H \times H} \\ &= \sqrt{ \left\| \left( \omega(\mathbf{e}_{\bar{x}}(u) - \mathbf{e}_{\bar{x}}(u^*)) + p(u^* - u) + \frac{1}{\omega}(v - v^*) \right)^2 + \left\| \frac{1}{\omega}(v - v^*) + q(u^* - u) \right\|^2 \right\|} \\ &\le \sqrt{ 2\omega^2 \|\mathbf{e}_{\bar{x}}(u) - \mathbf{e}_{\bar{x}}(u^*)\|^2 + (4p^2 + 2q^2) \|u - u^*\|^2 + \frac{6}{\omega^2} \|v - v^*\|^2 } \\ &= \sqrt{ (2\omega^2(2L + \beta\tau)^2 + 4p^2 + 2q^2) \|u - u^*\|^2 + \frac{6}{\omega^2} \|v - v^*\|^2 } \\ &\le \left( \sqrt{2\omega(2L + \beta\tau)} + 2|p| + \sqrt{2}|q| + \frac{\sqrt{6}}{\omega} \right) \|(u, v) - (u^*, v^*)\|_{H \times H} \\ &= K(\omega, \beta, \tau, \gamma, t) \|(u, v) - (u^*, v^*)\|_{H \times H} \end{aligned} $$

where we have used Young's inequality and Lemma 4 for the Lipschitz continuity of the operator $\mathbf{e}_{\bar{x}}$. ■

**Proposition 3.** Assume that Assumptions 2, 6 and 7 hold where we only require $g \in \Gamma_0(H)$ (but possibly nonsmooth). Suppose also that $\omega > 0$ and $\gamma \in C^1([t_0, +\infty[; H)$ such that both $\gamma$ and $\dot{\gamma} \in L^2([t_0, T])$ for all $T > t_0$. Then, for any initial data $x_0 \in \text{dom}(g)$ and $y_0 \in H$, there exists a unique global strong solution $(x, y) : [t_0, +\infty[ \to H \times H$ to (MIS) such that $x(t_0) = x_0$ and $y(t_0) = y_0$. Moreover, the solution $Z = (x, y)$ satisfies the following properties

(i) $y \in C^1([t_0, +\infty[; H)$, and $\dot{y}(t) - (\frac{1}{\omega} - \gamma(t) - \omega\dot{\gamma}(t))x(t) + \frac{1}{\omega}y(t) = 0$, for all $t \ge t_0$;

(ii) $x$ is absolutely continuous on $[t_0, T]$ and $\dot{x} \in L^2(t_0, T; H)$ for all $T > t_0$;

(iii) $x(t) \in \text{dom}(\partial g)$ for all $t > t_0$;

(iv) $x$ is Lipschitz continuous on any compact subinterval of $]t_0, +\infty[$;

(v) the function $t \in [t_0, +\infty[ \mapsto G_{m_{\bar{x}}}(x(t))$ is absolutely continuous on $[t_0, T]$ for all $T > t_0$;

(vi) there exists a function $\xi : [t_0, +\infty[ \to H$ such that

(a) $\xi(t) \in \partial G_{m_{\bar{x}}}(x(t))$ for all $t > t_0$;

(b) $\dot{x}(t) + \omega\xi(t) - (\frac{1}{\omega} - \gamma(t))x(t) + \frac{1}{\omega}y(t) = 0$ for almost every $t > t_0$;

(c) $\xi \in L^2(t_0, T; H)$ for all $T > t_0$;

(d) $\frac{d}{dt}G_{m_{\bar{x}}}(x(t)) = \langle \xi(t), \dot{x}(t) \rangle$ for almost every $t > t_0$.
1
2
3
4

Proof. In view of our assumptions, $\mathcal{A}$ is maximal monotone (see Lemma 3) and $\mathcal{E}$ is Lipschitz continuous thanks to Lemma 5. We deduce the existence of a unique strong global solution $Z = (x, y) : [t_0, T] \to \mathcal{H} \times \mathcal{H}$ of (MIS) via [26, Proposition 3.12]. The verification of items (iii)-(vi) can be done by following the main arguments of [19, Theorem 4.4] ■

Assuming that $g \in C^1(\mathcal{H})$ with a Lipschitz continuous gradient, the conclusions of Proposition 3 can be strengthened as follows.

7 **Proposition 4.** Assume that Assumptions 2, 6 and 7 hold with $\nabla g$ Lipschitz continuous.
8 Suppose also that $\omega > 0$ and $\gamma \in C^1([t_0, +\infty[; \mathcal{H})$ such that both $\gamma$ and $\dot{\gamma} \in L^2([t_0, T])$ for all
9 $T > t_0$. Then, for any $t_0 > 0$, and any Cauchy data $(x_0, \dot{x}_0)$, the system (ISEHD$_{\gamma,\omega}$) admits
10 a unique global solution $x \in C^1([t_0, +\infty[; \mathcal{H})$ satisfying $(x(t_0), \dot{x}(t_0)) = (x_0, \dot{x}_0)$. Moreover,
11 $\dot{x}, \nabla G_{m_x}(x(\cdot))$ and $e_{\bar{x}}(x(\cdot))$ are absolutely continuous.

Proof. Under our regularity assumptions, we get from the conclusion of Proposition 3 and the first equation of (37) that $x$ is a $C^1([t_0, +\infty[; \mathcal{H})$ function. Moreover, $\dot{x}$ is absolutely continuous thanks to Lipschitz continuity of $\nabla G_{m_x}$ and $e_{\bar{x}}$, absolute continuity of $x$ already established in Proposition 3, together with continuity of $\gamma$, $\dot{\gamma}$ and $y$. ■

Remark 4.1. The lack of differentiability of e<sub>x̄</sub> is the main obstacle for showing the existence of a classical solution even if G<sub>m<sub>x̄</sub></sub> is assumed C<sup>2</sup>. Moreover, arguing using the non-autonomous Cauchy-Lipschitz theorem, the convexity assumption on g can be dropped in Proposition 4 at the price of assuming &gamma; to be twice continuously differentiable on [t<sub>0</sub>, +&infin;[. The proof is left to the reader.

4.2. **Convergence properties.** Now let us examine the convergence properties of (ISEHD<sub>γ,ω</sub>). We will work under Assumptions 2, 6 and 7 where $\nabla g$ is Lipschitz continuous. The $\mu$-strong convexity of $G_{m_x}$ allows to tune the viscous damping coefficient to the modulus $\mu$ by taking $\gamma(t) \equiv 2\sqrt{\mu}$ as advocated in [9, 14]. From now on, we focus on the following system

$$ \ddot{x}(t) + 2\sqrt{\mu}\dot{x}(t) + \nabla G_{m_x}(x(t)) + \omega \frac{d}{dt}(\nabla G_{m_x}(x(t))) + e_{\bar{x}}(x(t)) + \omega \frac{d}{dt}e_{\bar{x}}(x(t)) = 0. \quad (\text{ISEHD}_{2\sqrt{\mu},\omega}) $$
Thanks to our assumptions, this system falls within the framework of the previous section so that Proposition 4 applies.

To perform Lyapunov analysis, let us define the following energy function $V : [t_0, \infty[ \to \mathbb{R}^+$ by

$$ V(t) := G_{m_x}(x(t)) - G_{m_x}^* + \frac{1}{2}\|v(t)\|^2, \quad \text{where } v(t) := \sqrt{\mu}(x(t) - \bar{x}) + \dot{x}(t) + \omega \nabla G_{m_x}(x(t)). $$
where $G_{m_x}^* = \min_{x \in \mathcal{H}} G_{m_x}(\mathcal{H})$.

We prove the following result.

**Theorem 3.** Assume that Assumptions 2, 6 and 7 hold and $\nabla g$ is Lipschitz continuous. Let $x : [t_0, +\infty[ \to \mathcal{H}$ be the solution trajectory of (ISEHD<sub>2√μ,ω</sub>). Suppose that $\rho := \frac{\beta\tau}{\mu}$ and the damping coefficient $\omega$ satisfy

$$ 0 \le \omega \le \min \left( \frac{1}{2\sqrt{\mu}}, \frac{\sqrt{\mu}}{2\sqrt{2}(2L + \beta\tau)} \right) \quad \text{and} \quad 16\rho^2 + \omega < 1. \quad (40) $$

Then, we have:
1. (1) for all $t \ge t_0$

$$ \mathcal{V}(t) \le \mathcal{V}(t_0)e^{-\frac{\sqrt{\mu}}{4}(t-t_0)}. $$

In particular

$$ \frac{\mu}{2} \|x(t) - \bar{x}\|^2 \le G_{m_{\bar{x}}}(x(t)) - G_{m_{\bar{x}}}^* \le \mathcal{V}(t_0)e^{-\frac{\sqrt{\mu}}{4}(t-t_0)}. $$

2. (2) There exists $C > 0$ such that,

$$ e^{-\sqrt{\mu}t} \int_{t_0}^{t} e^{\sqrt{\mu}s} \|\nabla G_{m_{\bar{x}}}(x(s))\|^2 ds \le C e^{-\frac{\sqrt{\mu}}{4}t}, \forall t \ge t_0. $$

4. Some remarks are in order before proving this result.

5. Remark 4.2.

6. Compared to the convergence results in Theorem 2 on the first-order system (SMI), Theorem 3 not only provides a fast convergence rate of the trajectory and objective value, but also of the gradient. Observe also the $\sqrt{\mu}$ in the rate in Theorem 3 instead of $\mu$ Theorem 2. This recovers the known result that the system $(\text{ISEHD}_{2\sqrt{\mu},\omega})$ is faster than (SMI) for badly conditioned objectives, approaching the optimal rate of the class of strongly convex functions with Lipschitz continuous gradient.

7. Contrary to [14], no assumption on the integrability of $e_{\bar{x}}(x(\cdot))$ and $\frac{d}{dt}e_{\bar{x}}(x(\cdot))$ is needed in Theorem 3. In fact, thanks to Corollary 1 and Lemma 4, the norms of these error terms will be absorbed in the right hand side of (44).

8. In case $\omega = 0$, i.e., when the inertial dynamic is considered only with the viscous damping coefficient (which can be viewed as a perturbed HBF method), the condition (40) reduces to $\rho < \frac{1}{4}$. This contrasts with the convergence condition for first-order dynamics (cf. Theorem 2), where convergence is ensured in parameter regime $\rho < 1$. One possible explanation for this difference, is the potential occurrence of oscillations typical of inertial system, which may necessitate a stricter compatibility condition between the parameters $\tau, \beta, \mu$.

22. Proof. In view of our assumptions, Proposition 4 applies and we get that $v : [t_0, +\infty[ \to \mathcal{H}$ is absolutely continuous and thus so is $\mathcal{V}$. We then have

23.

$$ \begin{aligned} \dot{\mathcal{V}}(t) &= \langle \nabla G_{m_{\bar{x}}}(x(t)), \dot{x}(t) \rangle + \langle v(t), \sqrt{\mu}\dot{x}(t) + \ddot{x}(t) + \omega \frac{d}{dt}(\nabla G_{m_{\bar{x}}}(x(t))) \rangle \\ &= \langle \nabla G_{m_{\bar{x}}}(x(t)), \dot{x}(t) \rangle + \langle v(t), -\sqrt{\mu}\dot{x}(t) - \nabla G_{m_{\bar{x}}}(x(t)) - e_{\bar{x}}(x(t)) - \omega \frac{d}{dt}e_{\bar{x}}(x(t)) \rangle, \end{aligned} \quad (41) $$

where we used the constitutive equation in $(\text{ISEHD}_{2\sqrt{\mu},\omega})$. Replacing $v$ by its expression and rearranging, we arrive at

$$ \begin{aligned} \dot{\mathcal{V}}(t) + \mu \langle \dot{x}(t), x(t) - \bar{x} \rangle + \sqrt{\mu} \|\dot{x}(t)\|^2 + \sqrt{\mu} \langle \nabla G_{m_{\bar{x}}}(x(t)), x(t) - \bar{x} \rangle \\ + \omega \sqrt{\mu} \langle \nabla G_{m_{\bar{x}}}(x(t)), \dot{x}(t) \rangle + \omega \|\nabla G_{m_{\bar{x}}}(x(t))\|^2 = - \langle v(t), e_{\bar{x}}(x(t)) + \omega \frac{d}{dt} e_{\bar{x}}(x(t)) \rangle. \end{aligned} \quad (42) $$

24. Using $\mu$-strong convexity of $G_{m_{\bar{x}}}$, we have

$$ \langle \nabla G_{m_{\bar{x}}}(x(t)), x(t) - \bar{x} \rangle \ge G_{m_{\bar{x}}}(x(t)) - G_{m_{\bar{x}}}^* + \frac{\mu}{2} \|x(t) - \bar{x}\|^2, \quad (43) $$

25. and plugging this into (42), we get

$$ \dot{\mathcal{V}}(t) + \sqrt{\mu}\Theta(t) \le \|v(t)\| \|e_{\bar{x}}(x(t)) + \omega \frac{d}{dt}e_{\bar{x}}(x(t))\|, \quad (44) $$
1 where

$$
\begin{aligned}
\Theta(t) := G_{m_x}(x(t)) - G_{m_x}^* &+ \frac{\mu}{2} \|x(t) - \bar{x}\|^2 + \sqrt{\mu} \langle \dot{x}(t), x(t) - \bar{x} \rangle + \|\dot{x}(t)\|^2 \\
& + \omega \langle \nabla G_{m_x}(x(t)), \dot{x}(t) \rangle + \frac{\omega}{\sqrt{\mu}} \|\nabla G_{m_x}(x(t))\|^2.
\end{aligned}
$$

2 Using the definition of $\mathcal{V}(t)$, we may rewrite $\Theta(t)$ as

$$
\Theta(t) = \mathcal{V}(t) + \frac{1}{2} \|\dot{x}(t)\|^2 - \omega \sqrt{\mu} \langle \nabla G_{m_x}(x(t)), x(t) - \bar{x} \rangle + \left( \frac{\omega}{\sqrt{\mu}} - \frac{\omega^2}{2} \right) \|\nabla G_{m_x}(x(t))\|^2.
$$

3 Consequently, (44) becomes

$$
\begin{aligned}
& \dot{\mathcal{V}}(t) + \sqrt{\mu}\mathcal{V}(t) + \frac{\sqrt{\mu}}{2}\|\dot{x}(t)\|^2 \\
& \quad + \sqrt{\mu}\left(\left(\frac{\omega}{\sqrt{\mu}} - \frac{\omega^2}{2}\right)\|\nabla G_{m_x}(x(t))\|^2 - \omega\sqrt{\mu}\langle\nabla G_{m_x}(x(t)), x(t) - \bar{x}\rangle\right) \\
& \le \|v(t)\|\|\mathbf{e}_{\bar{x}}(x(t))\| + \omega \frac{d}{dt}\mathbf{e}_{\bar{x}}(x(t))\|. \quad (45)
\end{aligned}
$$

4 Using strong convexity of $G_{m_x}$ again, and discarding the quadratic term in $v(t)$, we have

$$
\mathcal{V}(t) = \frac{1}{2}\mathcal{V}(t) + \frac{1}{2}\mathcal{V}(t) \geq \frac{1}{2}\mathcal{V}(t) + \frac{\mu}{4}\|x(t) - \bar{x}\|^2.
$$

5 Observing that $\frac{\omega}{2\sqrt{\mu}} \le \frac{\omega}{\sqrt{\mu}} - \frac{\omega^2}{2}$ for $0 \le \omega \le \frac{1}{\sqrt{\mu}}$, we end up with

$$
\begin{aligned}
& \dot{\mathcal{V}}(t) + \frac{\sqrt{\mu}}{2}\mathcal{V}(t) + \frac{\sqrt{\mu}}{2}\|\dot{x}(t)\|^2 \\
& \quad + \sqrt{\mu}\left(\frac{\mu}{4}\|x(t) - \bar{x}\|^2 + \frac{\omega}{2\sqrt{\mu}}\|\nabla G_{m_x}(x(t))\|^2 - \omega\sqrt{\mu}\|\nabla G_{m_x}(x(t))\|\|x(t) - \bar{x}\|\right) \\
& \le \|v(t)\|\|\mathbf{e}_{\bar{x}}(x(t))\| + \omega \frac{d}{dt}\mathbf{e}_{\bar{x}}(x(t))\|. \quad (46)
\end{aligned}
$$

6 Now let us treat the right hand side of this inequality. Thanks to Corollary 1, we have

7 $\|\mathbf{e}_{\bar{x}}(x(t))\| = \|\mathbf{e}_{\bar{x}}(x(t)) - \mathbf{e}_{\bar{x}}(\bar{x})\| \le \beta\tau\|x(t) - \bar{x}\|$, and by Lemma 4 $\|\frac{d}{dt}\mathbf{e}_{\bar{x}}(x(t))\| \le (2L +

8 $\beta\tau)\|\dot{x}(t)\|$. Since $\mathcal{V}(t) \ge \frac{1}{2}\|v(t)\|^2$, applying Young's inequality yields

$$
\begin{aligned}
\|v(t)\|\|\mathbf{e}_{\bar{x}}(x(t))\| &+ \omega \frac{d}{dt} \mathbf{e}_{\bar{x}}(x(t)) \| \le \frac{\sqrt{\mu}}{8} \|v(t)\|^2 + \frac{4}{\sqrt{\mu}} \|\mathbf{e}_{\bar{x}}(x(t))\|^2 + \frac{4\omega^2}{\sqrt{\mu}} \|\frac{d}{dt} \mathbf{e}_{\bar{x}}(x(t))\|^2 \\
& \le \frac{\sqrt{\mu}}{4} \mathcal{V}(t) + \frac{4\beta^2\tau^2}{\sqrt{\mu}} \|x(t) - \bar{x}\|^2 + \frac{4\omega^2(2L + \beta\tau)^2}{\sqrt{\mu}} \|\dot{x}(t)\|^2 \quad (47)
\end{aligned}
$$

9 We get after rearranging the terms

$$
\dot{\mathcal{V}}(t) + \frac{\sqrt{\mu}}{4} \mathcal{V}(t) + \left( \frac{\sqrt{\mu}}{2} - \frac{4\omega^2(2L + \beta\tau)^2}{\sqrt{\mu}} \right) \|\dot{x}(t)\|^2 + \sqrt{\mu} \Psi(t) \le 0, \quad (48)
$$

where

$$
\Psi(t) = \left( \frac{\mu}{4} - \frac{4\beta^2\tau^2}{\mu} \right) \|x(t) - \bar{x}\|^2 + \frac{\omega}{2\sqrt{\mu}} \|\nabla G_{m_x}(x(t))\|^2 - \omega\sqrt{\mu} \|\nabla G_{m_x}(x(t))\| \|x(t) - \bar{x}\|.
$$
1 Setting

$$ \mathbf{a} = \frac{\mu}{4} - \frac{4\beta^2\tau^2}{\mu}, \mathbf{b} = \frac{\omega}{2\sqrt{\mu}}, \mathbf{c} = -\frac{\omega\sqrt{\mu}}{2} $$

2 and $X = \|x(t) - \bar{x}\|$ and $Y = \|\nabla G_{m_{\bar{x}}}(x(t))\|$, we see that $\Psi$ can be written as a quadratic

3 form $Q: \mathcal{H} \times \mathcal{H} \to \mathbb{R}$ with $Q(X, Y) = \mathbf{a}\|X\|^2 + 2\mathbf{c}\langle X, Y \rangle + \mathbf{b}\|Y\|^2$. By assumption $\mathbf{a}, \mathbf{b} \ge 0$,

4 and the discriminant $\mathbf{c}^2 - \mathbf{a}\mathbf{b}$ of $Q$ is nonpositive. Indeed, since $0 \le \omega \le \frac{1}{2\sqrt{\mu}}$

$$ \mathbf{c}^2 - \mathbf{a}\mathbf{b} = \frac{\omega^2\mu}{4} - \frac{\omega}{2\sqrt{\mu}}\left(\frac{\mu}{4} - \frac{4\beta^2\tau^2}{\mu}\right) \le \frac{\omega^2\sqrt{\mu}}{8} - \frac{\omega\sqrt{\mu}}{2}\left(\frac{1}{4} - 4\rho^2\right) = \frac{\omega\sqrt{\mu}}{8}(\omega - 1 + 16\rho^2) < 0. $$

5 Hence, $\Psi(t) \ge 0$ and since $\omega \le \frac{\sqrt{\mu}}{2\sqrt{2(2L+\beta\tau)}}$, we get $\left(\frac{\sqrt{\mu}}{2} - \frac{4\omega^2(2L+\beta\tau)^2}{\sqrt{\mu}}\right) \ge 0$, and thus

$$ \dot{V}(t) + \frac{\sqrt{\mu}}{4}V(t) \le 0, $$

6 which gives after integration

$$ V(t) \le V(t_0)e^{-\frac{\sqrt{\mu}}{4}(t-t_0)}. \quad (49) $$

7 Therefore, $\lim_{t \to \infty} V(t) = 0$ and in particular

$$ \lim_{t \to \infty} G_{m_{\bar{x}}}(x(t)) - G_{m_{\bar{x}}}^* = 0 \quad \text{and} \quad \lim_{t \to \infty} \|v(t)\| = 0. \quad (50) $$

8 This implies, using strong convexity of $G_{m_{\bar{x}}}$

$$ \lim_{t \to \infty} \|x(t) - \bar{x}\| = 0, $$

9 which gives that

$$ \lim_{t \to \infty} \|\nabla G_{m_{\bar{x}}}(x(t))\| = 0. $$

10 We deduce from (50) that $\lim_{t \to \infty} \|\dot{x}(t)\| = 0$.

11 Coming back to (49), we have, by definition of $V$ and $\mu$-strong convexity of $G_{m_{\bar{x}}}$, that

$$ \frac{\mu}{2}\|x(t) - \bar{x}\|^2 \le G_{m_{\bar{x}}}(x(t)) - G_{m_{\bar{x}}}^* \le V(t_0)e^{-\frac{\sqrt{\mu}}{4}(t-t_0)} \quad \text{and} \quad \|v(t)\|^2 \le 2V(t_0)e^{-\frac{\sqrt{\mu}}{4}(t-t_0)}. \quad (51) $$

Developing in (51) we have

$$ \mu\|x(t) - \bar{x}\|^2 + \|\dot{x}(t)\|^2 + \omega^2\|\nabla G_{m_{\bar{x}}}(x(t))\|^2 + 2\omega\sqrt{\mu}\langle\nabla G_{m_{\bar{x}}}(x(t)), x(t) - \bar{x}\rangle + 2\omega\langle\nabla G_{m_{\bar{x}}}(x(t)), \dot{x}(t)\rangle + 2\sqrt{\mu}\langle\dot{x}(t), x(t) - \bar{x}\rangle \le Ce^{-\frac{\sqrt{\mu}}{4}t}, \quad (52) $$

12 where $C = 2V(t_0)e^{\frac{\sqrt{\mu}}{4}t_0}$.

13 Since $\langle\nabla G_{m_{\bar{x}}}(x(t)), x(t) - \bar{x}\rangle \ge G_{m_{\bar{x}}}(x(t)) - G_{m_{\bar{x}}}^*$, we deduce from (52) that

$$ \dot{U}(t) + \sqrt{\mu}U(t) + \omega^2 \|\nabla G_{m_{\bar{x}}}(x(t))\|^2 \le Ce^{-\frac{\sqrt{\mu}}{4}t}, \quad (53) $$

14 where $U(t) := \sqrt{\mu}\|x(t) - \bar{x}\|^2 + 2\omega (G_{m_{\bar{x}}}(x(t)) - G_{m_{\bar{x}}}^*)$. Integrating (53), we obtain after

15 elementary computation

$$ e^{-\sqrt{\mu}t} \int_{t_0}^{t} e^{\sqrt{\mu}s} \|\nabla G_{m_{\bar{x}}}(x(s))\|^2 ds \le C_1 e^{-\frac{\sqrt{\mu}}{4}t}, $$

16 as desired.

<figure><img src="image_5.png" /></figure>

17 We end this section with a similar result to Corollary 2, which is a direct consequence of

18 Assumption 2 and Theorem 3(1).

19 Corollary 3. Let $x: [t_0, +\infty[ \to \mathcal{H}$ be the solution of (ISEHD$_{2\sqrt{\mu},\omega}$) where (40) holds. Then

20 $\forall t \ge t_0$

$$ W_1(m_{x(t)}, m_{\bar{x}}) \le \tau \sqrt{\frac{2}{\mu} V(t_0)} e^{-\frac{\sqrt{\mu}}{8}(t-t_0)}. $$
## 5. ON COARSE RICCI CURVATURE

In this section we discuss some dynamical and geometrical properties of the family $(\mathbf{m}_x)_{x \in \mathcal{H}}$, particularly the notion of Ollivier-Ricci curvature and how it tightly related to Assumption 2. In fact, the family of probabilities $\mathbf{m} = (\mathbf{m}_x)_{x \in \mathcal{H}}$ and its Lipschitz behavior with respect to the $W_1$-Wasserstein distance, reveals that a natural setting to address monotone inclusions of the form (9), and thus stochastic optimization problems with decision-dependent distributions is the framework of *metric random walk spaces* (see, e.g., [41, 44]). All definitions of this section can be found in [36, 41].

**5.1. Metric random walk spaces.** Before going further, let us recall the following definitions to introduce a couple of probabilistic notions.

**Definition 6 (Random walks [44]).** Given a Polish space $(X, d)$. A family of probabilities $\mathbf{m} = (\mathbf{m}_x)_{x \in X}$ is a random walk on $X$ if $\mathbf{m}_x \in \mathcal{P}(\Xi)$ for each $x \in X$ and

* $\mathbf{m}_x$ depends measurably on $x \in X$,

* Each $\mathbf{m}_x$ has finite first-order moment, i.e., for some $x^o \in X$, $E_{y \sim \mathbf{m}_x} d(y, x^o) < \infty$.

Then $(X, d)$ equipped with a random walk $\mathbf{m}$ is a metric random walk space (m.r.w.s for short), and we denote it by $[X, d, \mathbf{m}]$.

Let us recall the notion of invariant and ergodic measures.

**Definition 7 (Invariance).** Let $\nu$ be a $\sigma$-finite measure on $X$ and $\mathbf{m}$ a random walk on $(X, \mathcal{B})$. We say that $\nu$ is invariant with respect to $\mathbf{m}$ if $\nu \star \mathbf{m} = \nu$, where $\nu \star \mathbf{m}$ is the convolution of $\nu$ by the random walk $\mathbf{m}$ and is defined by

$$ \nu \star \mathbf{m}(A) = \int_X \mathbf{m}_x(A) d\nu(x) \text{ for all } A \in \mathcal{B}. $$

As pointed out in [44], each measure $\mathbf{m}_x$ can be seen as a replacement of a sphere around $x$. While in a probabilistic framework one think about a Markov chain whose transition kernel from $x$ to $y$ in $n$ steps is defined by

$$ \mathrm{d}\mathbf{m}_x^{*n}(y) = \int_{z \in X} \mathrm{d}\mathbf{m}_x^{*(n-1)}(z) \mathrm{d}\mathbf{m}_z(y), \quad (54) $$

with $\mathbf{m}_x^1 = \mathbf{m}_x$ and $\mathbf{m}_x^0 = \delta_x$.

In the sequel, we assume that $(\mathcal{H}, d)$ is a separable real Hilbert space, and thus a Polish space, where $d(x, y) = (x - y, x - y)^{1/2}$.

**5.2. Feller Property.** Let us recall the following definition.

**Definition 8.** We say that $\mathbf{m} := (\mathbf{m}_x)_x$ has the weak-Feller property if and only if for every sequence $x_n \to x^0 \in \mathcal{H}$ we have $\mathbf{m}_{x_n} \to \mathbf{m}_{x^0}$, i.e., $\int f \mathrm{d}\mathbf{m}_{x_n} \to \int f \mathrm{d}\mathbf{m}_{x^0}$ for any $f \in C_b(\mathcal{H})$.

It turns that Assumption 2 implies directly that the family $\mathbf{m}$ is weak-Feller.

**Proposition 5.** Under Assumption 2, $\mathbf{m}$ has the weak-Feller property. Moreover, for each $x \in \mathcal{H}$, $\mathbf{m}_x$ has finite first-order moment.

*Proof.* Let $x^0 \in \mathcal{H}$ and $(x_n)_n$ a sequence of $\mathcal{H}$ such that $x_n \to x^0$ as $n \to 0$. Then Assumption 2 gives

$$ W_1(\mathbf{m}_{x_n}, \mathbf{m}_{x^0}) \le \tau \|x_n - x^0\|, $$

and thus $\lim_{n \to \infty} W_1(\mathbf{m}_{x_n}, \mathbf{m}_{x^0}) = 0$. Thanks to [4, Proposition 7.1.5], $(\mathbf{m}_{x_n})$ has uniformly integrable $p$-moments with $p \ge 1$ and narrowly converges towards $\mathbf{m}_{x^0}$. In particular $\mathbf{m}$ is weak-Feller and each $\mathbf{m}_x$ has finite first-order moments. ■
<table><caption>**Table 1:** Relation between the values of $\tau$ and $\kappa_m$.</caption><tr><td>Values of τ</td><td>0</td><td>&lt; 1</td><td>≤ 1</td></tr><tr><td>Values of κ<sub>m</sub></td><td>1</td><td>[0, 1]</td><td>[0, 1]</td></tr></table>

<footnote n="1">1 *Remark 5.1.* We already know that $\mathbf{m}_x \in \mathcal{P}(\Xi)$ for each $x \in \mathcal{H}$ and that $x \mapsto \mathbf{m}_x(C)$ is measurable for each $C \in \mathcal{B}$. Moreover, thanks to Proposition 5, we have finiteness of first-order moments of each $\mathbf{m}_x$, so that the family $\mathbf{m}$ satisfies the requirements of Definition 6. This shows that a natural setting to address dynamics of the form (9) is the metric random walk space $[\mathcal{H}, \mathrm{d}, \mathbf{m}]$. Many diffusion and variational problems has been studied within this framework, with allows in particular consider nonlocal continuum problems or problems on weighted graphs (see, e.g., [41] and the references therein).</footnote>

<footnote n="2">2 *Remark 5.2.* Let us point out that if $v$ is an invariant measure with respect to $\mathbf{m}$ then it is also and invariant measure with respect to $\mathbf{m}^{*n}$ for every $n \in \mathbb{N}$, where $\mathbf{m}^{*n}$ is the $n$-step transition probability function given by (54). It turns out that weak-Feller property implies that every weak-$*$ limit $v$ of $(\mathbf{m}^{*n})_n$ is an invariant measure of $\mathbf{m}$ cf. [36, Proposition 7.2.2] (see also [33, Proposition 12.3.4]). However, without assuming at first the existence of an invariant measure with respect to $\mathbf{m}$, the measure $v$ may be trivial. Without further compactness assumptions on the metric space (see, e.g., [36, Theorem 7.2.3]) one needs some Lyapunov like condition to ensure the existence of an invariant measure $v$ of the weak-Feller family $\mathbf{m}$ (see, e.g., [36, Theorem 7.2.4] or [33, Theorem 12.3.3]). As we will see in Corollary 5, another way to obtain the existence of invariant measures is having a positive lower bound on the coarse Ricci curvature of $[\mathcal{H}, \mathrm{d}, \mathbf{m}]$.</footnote>

<footnote n="3">3 *5.3. Ollivier-Ricci curvature.* Let us discuss here the connexion between Assumption 2 and the so-called coarse or Olliver-Ricci curvature (ORC for short). The results can be found in [44] or [45]. A more recent presentation can be found in [41].</footnote>

<footnote n="23">23 **Definition 9** (Ollivier-Ricci curvature [44]). Let $[\mathcal{H}, \mathrm{d}, \mathbf{m}]$ be a m.r.w.s. Then, for any distinct points $x, y \in \mathcal{H}$, the ORC along $(x, y)$ is defined as:</footnote>

<footnote n="24">24 $\kappa_m(x, y) = 1 - \frac{W_1(\mathbf{m}_x, \mathbf{m}_y)}{\mathrm{d}(x, y)}$, (55)</footnote>

<footnote n="25">25 The ORC of $[\mathcal{H}, \mathrm{d}, \mathbf{m}]$ is defined as</footnote>

<footnote n="26">26 $\kappa_m := \inf_{x \neq y} \kappa_m(x, y)$. (56)</footnote>

<footnote n="25">25 We clearly see from (55) that $\kappa_m(x, y) \le 1$. Moreover, rearranging the terms, we have</footnote>

<footnote n="26">26 $W_1(\mathbf{m}_x, \mathbf{m}_y) = (1 - \kappa_m(x, y)) \mathrm{d}(x, y)$. (57)</footnote>

<footnote n="26">26 Consequently, having some lower bound $\kappa_m(x, y) \ge c \in \mathbb{R}$ for any $x, y \in \mathcal{H}$ gives</footnote>

<footnote n="27">27 $W_1(\mathbf{m}_x, \mathbf{m}_y) \le (1-c)\mathrm{d}(x, y)$, (58)</footnote>

<footnote n="27">27 which describes a Lipschitz behavior of the random walk $\mathbf{m}$. This has to be compared to Assumption 2. Indeed, we see from Assumption 2 that, for $x \neq y$</footnote>

<footnote n="28">28 $1 - \tau \le \kappa_m(x, y) \le 1$, (59)</footnote>

<footnote n="29">29 so according to the values of $\tau$ we have different regimes on the ORC $\kappa_m$ (cf. Table 1).</footnote>

<footnote n="29">29 Notice that Assumption 2 excludes both the cases $\kappa_m \equiv 1$ and $\kappa_m < 0$. Typically, $\tau = 0$, would give that $\kappa_m = 1$ in other words $W_1(\mathbf{m}_x, \mathbf{m}_y) = 0$ for any $x, y \in \mathcal{H}$, i.e., the distribution $\mathbf{m}$ is constant.</footnote>
<figure><img src="image_6.png" /><figcaption>Figure 1: Illustration of the ORC.</figcaption></figure>

Moreover, it turns that there is equivalence between the lower bound on $\kappa_m$ in (59) and the Lipschitz behavior (58). This is directly related to a $W_1$-contraction property cf. [44, Proposition 20].

**Proposition 6.** Let $\mathbf{m}$ be a random walk on $(\mathcal{H}, \mathbf{d})$ and assume that $\mathbf{m}_x$ has finite moment for all $x \in \mathcal{H}$. Then

$$ \kappa_m(x, y) \geq c \in \mathbb{R}, \forall x \neq y \iff W_1(\nu_1 \star \mathbf{m}, \nu_2 \star \mathbf{m}) \leq (1-c)W_1(\nu_1, \nu_2) \quad \forall \nu_1, \nu_2 \in \mathcal{P}_1(\Xi). $$

In view of Proposition 6, taking $\nu_1 = \delta_x$ and $\nu_2 = \delta_y$ for $x \neq y$, and $c = 1 - \tau$, we get

$$ W_1(\mathbf{m}_x, \mathbf{m}_y) = W_1(\delta_x \star \mathbf{m}, \delta_y \star \mathbf{m}) \leq \tau W(\delta_x, \delta_y) = \tau d(x, y) $$

which is exactly Assumption 2 since the above inequality is trivial for $x = y$.

In the case of positive curvature, this contraction result implies the existence of a unique invariant measure for the random walk $\mathbf{m}$ when the ORC is positive.

**Corollary 4 ([41]).** Assume that $\kappa_m(x, y) \geq c > 0$ for all $x \neq y$. Then, the random walk $\mathbf{m}$ has a unique invariant measure $v \in \mathcal{P}_1(\Xi)$. Moreover, for any $v \in \mathcal{P}_1(\Xi)$

$$ (1) \quad W_1(\nu \star \mathbf{m}^{*n}, v) \leq (1-c)^n W_1(\nu, v), \forall n \in \mathbb{N}. $$

$$ (2) \quad W_1(\mathbf{m}^{*n}, v) \leq \frac{(1-c)^n}{c} W_1(\delta_x, v), \forall n \in \mathbb{N}, \forall x \in \mathcal{H}. $$

In our setting, the parameter regime $\tau < 1$ would define a positive ORC. Consequently, we have the following

**Corollary 5.** Assume that $\mathbf{m}$ satisfies Assumption 2 with $\tau < 1$. Then, there exists a unique invariant measure $v \in \mathcal{P}_1(\Xi)$ with respect to $\mathbf{m}$.

## 6. APPLICATION: INERTIAL PRIMAL-DUAL ALGORITHM

This section is devoted to the application of the developed results in Section 3 and Section 4 to the following class of saddle-point problems

$$ \inf_{x \in \mathcal{H}} \sup_{y \in \mathcal{K}} \mathbf{f}_{\mathbf{m}_x}^p(x) + \mathbf{g}(x) + \langle y, \mathbf{K}x \rangle - \mathbf{h}(y) - \mathbf{r}_{\mathbf{m}_y^d}(y), \quad (60) $$

where $\mathbf{K}: \mathcal{H} \to \mathcal{K}$ is a bounded linear operator, $\mathcal{K}$ is a real Hilbert space, $(\mathbf{m}_x^p)_{x \in \mathcal{H}}$ and $(\mathbf{m}_y^d)_{y \in \mathcal{K}}$ are two families of probability measures on $\Xi$ and $\mathcal{Z}$ respectively. Recall that $\Xi$ and $\mathcal{Z}$ are two Polish spaces. We equip $\mathcal{H} \times \mathcal{K}$ with the product space inner product and associated norm. In what follows, we make the following assumptions:

**Assumption 8.**

(a) $\mathbf{f}_{\mathbf{m}_x^p}(x) := \mathbb{E}_{\xi \sim \mathbf{m}_x^p}(f(x, \xi)), f \in C^1(\mathcal{H} \times \Xi), f(x, \cdot)$ is $\mathbf{m}_x^p$-measurable and $C_{\beta_f}^{1,1}(\Xi)$ for every $x \in \mathcal{H}$, and $\mathbf{f}_{\mathbf{m}_x^p}(\cdot, \xi) \in C_{L_f}^{1,1}(\mathcal{H})$ for every $\xi \in \Xi$;
(b) $r_{m_y^d}(y) := E_{\zeta \sim m_y^d}(r(y, \zeta))$, with $r \in C^1(\mathcal{K} \times \mathcal{Z})$, $r(\cdot, \cdot)$ is $m_y^d$-measurable and $C_{\beta_r}^{1,1}(\mathcal{Z})$ for every $y \in \mathcal{K}$, and $r_{m_y^d}(\cdot, \zeta) \in C_{L_r}^{1,1}(\mathcal{K})$ for every $\zeta \in \mathcal{Z}$;

(c) $h \in \Gamma_0(\mathcal{K})$ and $g \in \Gamma_0(\mathcal{H})$;

(d) $f_{m_y^p} + g$ is $\mu_p$-strongly convex for all $x \in \mathcal{H}$, and $r_{m_y^d} + h$ is $\mu_d$-strongly convex for all $y \in \mathcal{K}$.

(e) There exists $\tau > 0$ such that
$$W_1(m_x^p \otimes m_y^d, m_{x'}^p \otimes m_{y'}^d) \le \tau \|(x, y) - (x', y')\|_{\mathcal{H} \times \mathcal{K}}$$ for all $x, x' \in \mathcal{H}$, $y, y' \in \mathcal{K}$.

$\beta_f, \beta_r, L_f, L_r, \mu_p, \mu_d$ are all positive constants.

Assumption 8(a)-(c) are to be compared to Assumptions 1 and 5. This will be made clearer when we will cast the saddle point problem (60) and an inclusion problem reminiscent of (9). In the same way, Assumption 8(d) will be sufficient for Assumption 3 to hold. (e) is the version of Assumption 2 on the product space $\mathcal{H} \times \mathcal{K}$.

Problems of the form (60) arise in many fields such as image and signal processing, machine learning and partial differential equations. In the deterministic case where $f$ and/or $r$ do not depend on the distributions $m_x^p$ and $m_y^d$, and where all functions are convex, various forms of such problem were studied in e.g., [29, 30, 48, 51] and many others. In [24], the authors studied a stochastic variant of (60) with non-state dependent measures and convex functions. Specifically, they considered $f(x) = E_{\xi \sim m}(f(x, \xi))$, $g(x) = E_{\xi \sim m}(g(x, \xi))$, $K = E_{\xi \sim m}(K(\xi))$, $r \equiv 0$ and $h(y) = E_{\xi \sim m}(p(x, \xi))$ for some suitable functions $f, g, h, p$ and $K$. Problems of the form (60) in their full generality are difficult to tackle directly because of the presence of the (primal and dual) state-dependent distributions. To the best of our knowledge, general saddle point problems of the form (60) have not been addressed in the literature.

**6.1. Formulation as a monotone inclusion.** As seen in Section 3, the appropriate notion of solutions of (60) is that of equilibria. Thus, our aim is to find an equilibrium point $(\bar{x}, \bar{y})$, i.e., that $(\bar{x}, \bar{y})$ a solution of the static saddle point problem

$$\inf_{x \in \mathcal{H}} \sup_{y \in \mathcal{K}} \mathcal{L}_{m_x^p, m_y^d}(x, y), \qquad (61)$$

where we have defined the Lagrangian for the (primal and dual) probability measures $m^p \in \mathcal{P}(\Xi)$ and $m^d \in \mathcal{P}(\mathcal{Z})$

$$\mathcal{L}_{m^p, m^d}(x, y) := f_{m^p}(x) + g(x) + \langle y, Kx \rangle - h(y) - r_{m^d}(y).$$

If the set of saddle points is nonempty (typically under strong duality, see e.g., [23, Chapter 15]), then $(x^*, y^*)$ is a saddle point for (61) if and only if

$$\mathcal{L}_{m_x^p, m_y^d}(x^*, y) \le \mathcal{L}_{m_x^p, m_y^d}(x^*, y^*) \le \mathcal{L}_{m_x^p, m_y^d}(x, y^*) \quad \forall (x, y) \in \mathcal{H} \times \mathcal{K},$$

or equivalently, if the following optimality condition holds [23, Corollary 19.19]$^1$

$$\left\{ \begin{array}{l} 0 \in \nabla f_{m_x^p}(x^*) + \partial g(x^*) + K^* y^* \\ 0 \in \nabla r_{m_y^d}(y^*) + \partial h(y^*) - K x^*. \end{array} \right. \qquad (62)$$

where $K^*$ is the adjoint operator of $K$. Therefore, $(\bar{x}, \bar{y})$ will be said to be an equilibrium of (60) if it is a solution of (62).

Equivalently, this can be cast as the monotone inclusion problem

$$0_{\mathcal{H} \times \mathcal{K}} \in \mathbf{T}_{m_x^p, m_y^d}(x^*, y^*), \qquad (63)$$

<footnote n="1">¹Here, only convexity is needed rather than Assumption 8(d).</footnote>
1 where

$$
\begin{align}
\mathbf{T}_{m^p, m^d} &= \mathbf{A} + \mathbf{L} + \mathbf{B}_{m^p, m^d} : \mathcal{H} \times \mathcal{K} \Rightarrow \mathcal{H} \times \mathcal{K}, \\
\mathbf{A} &= \begin{pmatrix} \partial g & 0 \\ 0 & \partial h \end{pmatrix}, \quad \mathbf{L} = \begin{pmatrix} 0 & \mathbf{K}^* \\ -\mathbf{K} & 0 \end{pmatrix}, \quad \text{and} \quad \mathbf{B}_{m^p, m^d} = \begin{pmatrix} \nabla f_{m^p} & 0 \\ 0 & \nabla r_{m^d} \end{pmatrix}. \tag{64}
\end{align}
$$

2 Remark 6.1. The authors in [52] studied problems of the form

$$
\min_{x \in X} \max_{y \in Y} \mathbb{E}_{\xi \sim m_{(x, y)}} \phi(x, y, \xi),
$$

3 in finite dimension where X, Y are compact sets and φ is a convex-concave function that

4 plays the role of the Lagrangian in our case. Such problems fall into the scope of (61)

5 where X and Y can be absorbed in the functions g and h through their respective indicator

6 functions l_X and l_Y.

7 The form (64) brings back our problem to the general setting studied in Section 3. The

8 following result shows that the operators **A** and **B**<sub>m<sup>p</sup>,m<sup>d</sup></sub> indeed verify all the appropriate

9 properties required there.

10 **Lemma 6.** Under Assumption 8, the following properties hold for any probability measures

11 m<sup>p</sup> ∈ P(Ξ) and m<sup>d</sup> ∈ P(Z):

12 (i) The operator **A** is maximally monotone.

13 (ii) The operator **B**<sub>m<sup>p</sup>,m<sup>d</sup></sub> is L̃-Lipschitz continuous with L̃ = max(L<sub>f</sub>, L<sub>r</sub>). It is maximally monotone if f(·, ξ) is convex for m<sup>p</sup>-almost every ξ ∈ Ξ and r(·, ζ) is convex

14 for m<sup>d</sup>-almost every ζ ∈ Z. 15

16 (iii) The operator T<sub>m<sup>p</sup>,m<sup>d</sup></sub> is maximally µ-strongly monotone with µ = min(µ<sub>p</sub>, µ<sub>d</sub>).

17 (iv) For any x, x' ∈ H and y, y' ∈ K

$$
\sup_{(x'', y'') \in \mathcal{H} \times \mathcal{K}} \| \mathbf{B}_{m_x^p, m_y^d} (x'', y'') - \mathbf{B}_{m_{x'}^p, m_{y'}^d} (x'', y'') \|_{\mathcal{H} \times \mathcal{K}} \le \tilde{\beta} \tau \|(x, y) - (x', y')\|_{\mathcal{H} \times \mathcal{K}}
$$

18 where $\tilde{\beta} = \max(\beta_f, \beta_r)$.

19 Proof. Thanks to (64) and assumption (c), (i) follows from [23, Theorem 21.2 and Propo-
sition 20.23]. The operator **B**<sub>m<sup>p</sup>,m<sup>d</sup></sub> is the gradient of the separable function **Φ**(x, y) =

20 **f**<sub>m<sup>p</sup></sub>(x) + **r**<sub>m<sup>d</sup></sub>(y). From assumptions (a) and (b), we have **∇**f<sub>m<sup>p</sup></sub> = **E**<sub>ξ~m<sup>p</sup></sub>(**∇**f(·, ξ)) and similarly for **∇**r<sub>m<sup>d</sup></sub>. Thus **Φ** ∈ C<sup>1,1</sup><sub>max(L<sub>f</sub>,L<sub>r</sub>)</sub>(H × K) whence we get claim the first part of

21 (ii). Moreover, under convexity, f<sub>m<sup>p</sup></sub> and r<sub>m<sup>d</sup></sub> are also convex, hence their gradients are maximally monotone, and we conclude using [23, Proposition 20.23]. To show (iii), observe that T<sub>m<sup>p</sup>,m<sup>d</sup></sub>(x, y) = (∂(f<sub>m<sup>p</sup></sub> + g)(x), ∂(r<sub>m<sup>d</sup></sub> + h)(y)) + L(x, y). The first part is maximally strongly monotone by (d) and [23, Example 22.4, Theorem 21.2 and Proposition 20.23].

22 In addition, **L** is a skew-symmetric linear operator, hence maximally monotone (see [23, Example 20.35]). We conclude using [26, Lemma 2.4].
23

24 The proof of (iv) is similar to that of Corollary 1. Indeed, let us fix (x'', y'') ∈ H × K. We
25 then have, for any x, y ∈ H and y, y' ∈ K

$$
\begin{align*}
& \| \mathbf{B}_{m_x^p, m_y^d} (x'', y'') - \mathbf{B}_{m_{x'}^p, m_{y'}^d} (x'', y'') \|_{\mathcal{H} \times \mathcal{K}} \\
&= \sqrt{ \left\| \nabla f_{m_x^p} (x'') - \nabla f_{m_{x'}^p} (x'') \right\|_{\mathcal{H}}^2 + \left\| \nabla r_{m_y^d} (y'') - \nabla r_{m_{y'}^d} (y'') \right\|_{\mathcal{K}}^2 } \\
&\le \sqrt{ \beta_f^2 W_1(m_x^p, m_{x'}^p)^2 + \beta_r^2 W_1(m_y^d, m_{y'}^d)^2 } \\
&\le \max(\beta_f, \beta_r) W_1(m_x^p \otimes m_{x'}^p, m_y^d \otimes m_{y'}^d) \\
&\le \max(\beta_f, \beta_r) \tau \|(x, y) - (x', y')\|_{\mathcal{H} \times \mathcal{K}}.
\end{align*}
$$
1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

**6.2. Existence and uniqueness of equilibrium.** We are now in a position to prove the existence and uniqueness of an equilibrium to the problem (61)

**Theorem 4 (Existence and uniqueness of equilibrium point).** Under Assumption 8, the map

$$S : (x, y) \in \mathcal{H} \times \mathcal{K} \mapsto \operatorname{zer}(\mathbf{T}_{m_x^p, m_y^d}) = \{(u, v) \in \mathcal{H} \times \mathcal{K} : (0,0)_{\mathcal{H} \times \mathcal{K}} \in \mathbf{T}_{m_x^p, m_y^d}(u, v)\}$$

is $\tilde{\rho}$-Lipschitz with

$$\tilde{\rho} := \frac{\tau \tilde{\beta}}{\tilde{\mu}}, \qquad (65)$$

where we recall that $\tilde{\mu} = \min(\mu_p, \mu_d)$ and $\tilde{\beta} = \max(\beta_f, \beta_r)$. In particular, if $\tilde{\rho} < 1$, there is a unique equilibrium point $(\bar{x}, \bar{y})$.

An immediate consequence of this result is that (62) has a unique solution which is $(\bar{x}, \bar{y})$.

*Proof.* Following the same lines as in Theorem 1 and using Lemma 6(iii), we conclude. ■

**Remark 6.2.** Theorem 4 is to be compared to [52, Theorem 2.6]. Notice that the monotone inclusion (63) defining the equilibrium in [52], handled through a variational inequality there to the presence of constraints on $x$ and $y$, can be cast in our setting through normal cones absorbed in $\partial g$ and $\partial h$ under appropriate qualification conditions.

Following the reasoning of Sections 3 and 4, let us define the following gap function

$$\mathbf{E}_{\bar{x}, \bar{y}}(x, y) = \mathbf{B}_{m_x^p, m_y^d}(x, y) - \mathbf{B}_{m_{\bar{x}}, m_{\bar{y}}^d}(x, y). \qquad (66)$$

Mimicking the proof of Lemma 4 using Assumption 8, we prove the following.

**Lemma 7.** Under Assumption 8, $\mathbf{E}_{\bar{x}, \bar{y}}(.)$ is $(2\tilde{L} + \tau\tilde{\beta})$-Lipschitz continuous, where $\tilde{L} = \max(L_f, L_r)$.

**6.3. Related first and second-order dynamics.** To lighten the notation, from now on we set $Z(t) := (x(t), y(t))$, $\bar{Z} := (\bar{x}, \bar{y})$ is the equilibrium of (61) (see Theorem 4), and $m_{\bar{Z}} := m_{\bar{x}}^p \otimes m_{\bar{y}}^d$.

**6.3.1. First order system.** Given an initial data $Z(t_0) = (x(t_0), y(t_0)) \in \operatorname{dom}(g) \times \operatorname{dom}(h)$, we consider the following first-order system associated to the monotone inclusion (63):

$$\dot{Z}(t) + \mathbf{T}_{m_{\bar{Z}}}(Z(t)) + \mathbf{E}_{\bar{Z}}(Z(t)) \ni 0_{\mathcal{H} \times \mathcal{K}}. \qquad (\text{SPDS})$$

Arguing as in Proposition 1² and Theorem 2 we have the following result.

**Proposition 7.** Assume that Assumption 8 holds and $\tilde{\rho} < 1$. Then, for any initial data $Z(t_0) = (x(t_0), y(t_0)) \in \operatorname{dom}(g) \times \operatorname{dom}(h)$, (SPDS) admits a unique strong solution $Z : t \in [t_0, +\infty[ \mapsto (x(t), y(t))$. Moreover,

$$\|Z(t) - \bar{Z}\|_{\mathcal{H} \times \mathcal{K}} \le C e^{-2\tilde{\mu}(1-\tilde{\rho})t}, \quad \forall t \ge t_0,$$

with $C = \|Z_0 - \bar{Z}\|_{\mathcal{H} \times \mathcal{K}} e^{2\tilde{\mu}(1-\tilde{\rho})t_0}$, where $\tilde{\rho}$ and $\tilde{\mu}$ are defined in (65).

<footnote n="2">²Here, we invoke [26, Proposition 3.12] rather than [26, Proposition 3.13] as $\mathbf{T}_{m_{\bar{Z}}}$ is a subdifferential operator perturbed by a bounded linear operator. This explains why we do not need any domain interiority assumption in this case.</footnote>
1. **6.3.2. Second-order system.** As we have done in Section 4, here, we restrict ourselves to the smooth setting where both g and h are smooth with Lipschitz continuous gradient.

2. In this case, the operator A, hence $T_{m_Z}$, is single-valued, Lipschitz continuous and strongly monotone. In explicit forms, for $Z = (x, y)$, we have

$$ T_{m_Z}(Z) = \left( \nabla_x \mathcal{L}_{m_x^p, m_y^d}(x, y), -\nabla_y \mathcal{L}_{m_x^p, m_y^d}(x, y) \right). \quad (67) $$

5. We propose the following inertial system associated to (63)

$$ \ddot{Z}(t) + 2\sqrt{\bar{\mu}}\dot{Z}(t) + \begin{pmatrix} \nabla_x \mathcal{L}_{m_x^p, m_y^d}(x(t), y(t) + \frac{1}{\sqrt{\bar{\mu}}} \dot{y}(t)) \\ -\nabla_y \mathcal{L}_{m_x^p, m_y^d}(x(t) + \frac{1}{\sqrt{\bar{\mu}}} \dot{x}(t), y(t)) \end{pmatrix} + \mathbf{E}_{\bar{Z}}(Z(t)) = 0. \quad (\text{ISPDS}_{2\sqrt{\bar{\mu}}}) $$

6. Unperturbed second order primal-dual dynamical systems for solving saddle point problems have been actively studied in recent years, essentially with vanishing viscous damping; see e.g., [8] and [35] and references therein.

7. Unlike $(ISEHD_{2\sqrt{\bar{\mu}},\omega})$, we did not include geometric damping in $(ISPDS_{2\sqrt{\bar{\mu}}})$ for the sake of simplicity. Moreover, one can see that the system $(ISPDS_{2\sqrt{\bar{\mu}}})$ is driven by the gradient of the Lagrangian but evaluated at extrapolated points, which is standard approach for primal-dual second-order systems. In fact, one has to keep in mind that unlike $(ISEHD_{2\sqrt{\bar{\mu}},\omega})$, the system $(ISPDS_{2\sqrt{\bar{\mu}}})$ is driven by an operator deriving from the Lagrangian which is convex-concave and contains a bilinear form. This is the reason underlying the inclusion of the extrapolation step. This will also necessitate to modify the Lyapunov analysis in the proof Theorem 3 though the main ingredients will remain essentially the same.

10. We first state that $(ISPDS_{2\sqrt{\bar{\mu}}})$ is well posed. This follows from the same arguments as [8, Theorem 5], using the Cauchy-Lipschitz theorem since $\nabla_x \mathcal{L}_{m_x^p, m_y^d}$ and $\nabla_y \mathcal{L}_{m_x^p, m_y^d}$ are globally Lipschitz continuous by assumption, and so is $\mathbf{E}_{\bar{Z}}$ by Lemma 7.

21. **Proposition 8.** Suppose that Assumption 8(a), (c) and (e) hold and, moreover, that g and h are also continuously differentiable with Lipschitz continuous gradient. Then, for any given initial condition $(x(t_0), \dot{x}(t_0)) = (x_0, u_0) \in \mathcal{H} \times \mathcal{H}$ and $(y(t_0), \dot{y}(t_0)) = (y_0, v_0) \in \mathcal{K} \times \mathcal{K}$ the evolution system $(ISPDS_{2\sqrt{\bar{\mu}}})$ has a unique strong global solution $Z(\cdot) = (x(\cdot), y(\cdot))$ with

25. $Z \in C^1([0, +\infty[; \mathcal{H} \times \mathcal{K})$;

26. $Z$ and $\dot{Z}$ are absolutely continuous on every compact subset of the interior of $[t_0, +\infty[$ (hence almost everywhere differentiable);

27. $\dot{Z}(t_0) = (u_0, v_0)$, $\dot{Z}(t_0) = (u_0, v_0)$, $\dot{Z}(t_0) = (u_0, v_0)$.

28. Let us now move to the convergence properties of $(ISPDS_{2\sqrt{\bar{\mu}}})$. We define the following energy function $\mathcal{V}: [t_0, \infty[ \to \mathbb{R}^+$ by

$$ \mathcal{V}(t) := \mathcal{L}_{m_x^p, m_y^d}(x(t), \bar{y}) - \mathcal{L}_{m_x^p, m_y^d}(\bar{x}, y(t)) + \frac{1}{2} \|v(t)\|_{\mathcal{H} \times \mathcal{K}}^2 $$

31. with

$$ v(t) := \sqrt{\bar{\mu}}(Z(t) - \bar{Z}) + \dot{Z}(t), $$

32. Observe that the Lagrangian gap in $\mathcal{V}$ is nonnegative and convex in $Z(t)$.

33. **Theorem 5.** Assume that Assumption 8 holds and, moreover, that g and h are also continuously differentiable with Lipschitz continuous gradient. Let $t \in [t_0, \infty[ \mapsto (x(t), y(t))$ be the solution of $(ISPDS_{2\sqrt{\bar{\mu}}})$. Suppose that $\tilde{\rho} < \frac{\sqrt{2}}{4}$, where $\tilde{\rho}$ is as given in (65). We then have for all $t \ge t_0$:

$$ \frac{\tilde{\mu}}{2} \|Z(t) - \bar{Z}\|_{\mathcal{H} \times \mathcal{K}}^2 \le \mathcal{L}_{m_x^p, m_y^d}(x(t), \bar{y}) - \mathcal{L}_{m_x^p, m_y^d}(\bar{x}, y(t)) \le \mathcal{V}(t_0)e^{-\frac{\sqrt{\bar{\mu}}}{4}(t-t_0)}. \quad (68) $$
1

2 Proof. To lighten notation, we drop dependence of the norm and inner product on the

3 underlying space which is to be understood from the context. At first, we have, using

4 (ISPDS$_{2}$√$_{µ}$$_{)

$$
\begin{align*}
\dot{v}(t) = \sqrt{\tilde{\mu}} \dot{Z}(t) + \ddot{Z}(t) &= -\sqrt{\tilde{\mu}} \dot{Z}(t) - \begin{pmatrix} \nabla_x \mathcal{L}_{m_x^p, m_y^d}(x(t), y(t) + \frac{1}{\sqrt{\tilde{\mu}}} \dot{y}(t)) \\ -\nabla_y \mathcal{L}_{m_x^p, m_y^d}(x(t) + \frac{1}{\sqrt{\tilde{\mu}}} \dot{x}(t), y(t)) \end{pmatrix} - \mathbf{E}_{\bar{Z}}(Z(t)) \\
&= -\sqrt{\tilde{\mu}} \dot{Z}(t) - \begin{pmatrix} \nabla(f_{m_x^p} + g)(x(t)) \\ \nabla(r_{m_y^d}^d + h)(y(t)) \end{pmatrix} - \mathbf{L} \begin{pmatrix} Z(t) + \frac{1}{\sqrt{\tilde{\mu}}} \dot{Z}(t) \\ Z(t) \end{pmatrix} - \mathbf{E}_{\bar{Z}}(Z(t)) \\
&= -\sqrt{\tilde{\mu}} \dot{Z}(t) - \mathbf{T}_{m_z}(Z(t)) - \frac{1}{\sqrt{\tilde{\mu}}} \mathbf{L} \dot{Z}(t) - \mathbf{E}_{\bar{Z}}(Z(t)).
\end{align*}
$$

5 We thus get,

$$
\begin{align*}
\dot{V}(t) &= \langle (\nabla_x \mathcal{L}_{m_x^p, m_y^d}(x(t), \bar{y}), -\nabla_y \mathcal{L}_{m_x^p, m_y^d}(\bar{x}, y(t))) , \dot{Z}(t) \rangle + \langle v(t), \dot{v}(t) \rangle \\
&= \langle (\nabla(f_{m_x^p} + g)(x(t)), \nabla(r_{m_y^d}^d + h)(y(t))) + \mathbf{L}\bar{Z}, \dot{Z}(t) \rangle + \langle v(t), \dot{v}(t) \rangle \\
&= \langle \mathbf{T}_{m_z}(Z(t)) + \mathbf{L}(\bar{Z} - Z(t)), \dot{Z}(t) \rangle \\
&\quad + \langle \sqrt{\tilde{\mu}}(Z(t) - \bar{Z}) + \dot{Z}(t), -\sqrt{\tilde{\mu}}\dot{Z}(t) - \mathbf{T}_{m_z}(Z(t)) - \frac{1}{\sqrt{\tilde{\mu}}}\mathbf{L}\dot{Z}(t) - \mathbf{E}_{\bar{Z}}(Z(t)) \rangle.
\end{align*}
$$

6 After some simplifications and using that L is skew-symmetric, we arrive at

$$
\dot{V}(t) + \sqrt{\tilde{\mu}} \langle \mathbf{T}_{m_z}(Z(t)), Z(t) - \bar{Z} \rangle + \tilde{\mu} \langle \dot{Z}(t), Z(t) - \bar{Z} \rangle + \sqrt{\tilde{\mu}} \|\dot{Z}(t)\|^2 = - \langle v(t), \mathbf{E}_{\bar{Z}}(Z(t)) \rangle. \quad (69)
$$

7 Rather than using $\tilde{\mu}$-strong monotonicity of $\mathbf{T}_{m_z}$ (Lemma 6(iii)), we will now prove an

8 alternative bound on the first inner product in (69) that will be useful for our Lyapunov

9 analysis. We have

$$
\begin{align*}
\langle \mathbf{T}_{m_z}(Z(t)), Z(t) - \bar{Z} \rangle \\
&= \langle \nabla(f_{m_x^p} + g)(x(t)), x(t) - \bar{x} \rangle + \langle \nabla(r_{m_y^d}^d + h)(y(t)), y(t) - \bar{y} \rangle + \langle \mathbf{L}Z(t), Z(t) - \bar{Z} \rangle \\
&= \langle \nabla(f_{m_x^p} + g)(x(t)), x(t) - \bar{x} \rangle + \langle \nabla(r_{m_y^d}^d + h)(y(t)), y(t) - \bar{y} \rangle - \langle K^*y(t), \bar{x} \rangle + \langle Kx(t), \bar{y} \rangle.
\end{align*}
$$

10 Now, by Assumption 8(d), we have

$$
\begin{align*}
\langle \nabla(f_{m_x^p} + g)(x(t)), x(t) - \bar{x} \rangle &\ge (f_{m_x^p} + g)(x(t)) - (f_{m_x^p} + g)(\bar{x}) + \frac{\mu_p}{2} \|x(t) - \bar{x}\|^2 \\
\langle \nabla(r_{m_y^d}^d + h)(y(t)), y(t) - \bar{y} \rangle &\ge (r_{m_y^d}^d + h)(y(t)) - (r_{m_y^d}^d + h)(\bar{y}) + \frac{\mu_d}{2} \|y(t) - \bar{y}\|^2.
\end{align*}
$$

11 Summing we get that

$$
\begin{align*}
& \langle \mathbf{T}_{m_z}(Z(t)), Z(t) - \bar{Z} \rangle \\
&\ge (f_{m_x^p} + g)(x(t)) - (r_{m_y^d}^d + h)(\bar{y}) + (r_{m_y^d}^d + h)(y(t)) - (f_{m_x^p} + g)(\bar{x}) - \langle K^*y(t), \bar{x} \rangle + \langle Kx(t), \bar{y} \rangle + \frac{\tilde{\mu}}{2} \|Z(t) - \bar{Z}\|^2 \\
&= \mathcal{L}_{m_x^p, m_y^d}(x(t), \bar{y}) - \mathcal{L}_{m_x^p, m_y^d}(\bar{x}, y(t)) + \frac{\tilde{\mu}}{2} \|Z(t) - \bar{Z}\|^2.
\end{align*}
$$
1 Plugging this into (69), we get

$$
\begin{align}
& \dot{V}(t) + \sqrt{\tilde{\mu}} \left( \mathcal{L}_{m_x^p, m_y^d}(x(t), \bar{y}) - \mathcal{L}_{m_x^p, m_y^d}(\bar{x}, y(t)) + \frac{\tilde{\mu}}{2} \|Z(t) - \bar{Z}\|^2 + \sqrt{\tilde{\mu}} \langle \dot{Z}(t), Z(t) - \bar{Z} \rangle + \|\dot{Z}(t)\|^2 \right) \nonumber \\
& = \dot{V}(t) + \sqrt{\tilde{\mu}} \left( \mathcal{L}_{m_x^p, m_y^d}(x(t), \bar{y}) - \mathcal{L}_{m_x^p, m_y^d}(\bar{x}, y(t)) + \frac{1}{2} \|\sqrt{\tilde{\mu}}(Z(t) - \bar{Z}) + \dot{Z}(t)\|^2 \right) + \frac{\tilde{\mu}}{2} \|\dot{Z}(t)\|^2 \nonumber \\
& = \dot{V}(t) + \sqrt{\tilde{\mu}} V(t) + \frac{\sqrt{\tilde{\mu}}}{2} \|\dot{Z}(t)\|^2 \le \|v(t)\| \|\mathbf{E}_{\bar{Z}}(Z(t))\|, \tag{70}
\end{align}
$$

2 where we used Cauchy-Schwarz inequality in the last line. Arguing as above, using again

3 Assumption 8(d), we have

$$
\mathcal{V}(t) \geq \mathcal{L}_{m_x^p, m_y^d}(x(t), \bar{y}) - \mathcal{L}_{m_x^p, m_y^d}(\bar{x}, y(t)) \geq \frac{\tilde{\mu}}{2} \|Z(t) - \bar{Z}\|^2. \quad (71)
$$

4 Using this in (70) yields the bound

$$
\dot{V}(t) + \frac{\sqrt{\tilde{\mu}}}{2} \mathcal{V}(t) + \frac{\sqrt{\tilde{\mu}}}{2} \|\dot{Z}(t)\|^2 + \frac{\tilde{\mu}^{3/2}}{4} \|Z(t) - \bar{Z}\|^2 \le \|v(t)\| \|\mathbf{E}_{\bar{Z}}(Z(t))\|. \quad (72)
$$

5 We have, using Young's inequality, Lemma 6(iv) and the fact that $\mathcal{V}(t) \ge \frac{1}{2} \|v(t)\|^2$

$$
\|v(t)\| \|\mathbf{E}_{\bar{Z}}(Z(t))\| \le \frac{\sqrt{\tilde{\mu}}}{8} \|v(t)\|^2 + \frac{2}{\sqrt{\tilde{\mu}}} \|\mathbf{E}_{\bar{Z}}(Z(t))\|^2 \le \frac{\sqrt{\tilde{\mu}}}{4} \mathcal{V}(t) + \frac{2\tilde{\beta}^2\tau^2}{\sqrt{\tilde{\mu}}} \|Z(t) - \bar{Z}\|^2.
$$

6 Inserting into (72) and rearranging gives

$$
\dot{V}(t) + \frac{\sqrt{\tilde{\mu}}}{4} \mathcal{V}(t) + \frac{\sqrt{\tilde{\mu}}}{2} \|\dot{Z}(t)\|^2 + \sqrt{\tilde{\mu}} \left( \frac{\tilde{\mu}}{4} - \frac{2\tilde{\beta}^2\tau^2}{\tilde{\mu}} \right) \|Z(t) - \bar{Z}\|^2 \le 0.
$$

7 Since $\frac{\tilde{\beta}\tau}{\tilde{\mu}} < \frac{\sqrt{2}}{4}$ by assumption, we have

$$
\dot{V}(t) + \frac{\sqrt{\tilde{\mu}}}{4} \mathcal{V}(t) \le 0,
$$

8 which gives after integration yields

$$
\mathcal{V}(t) \le \mathcal{V}(t_0)e^{-\frac{\sqrt{\tilde{\mu}}}{4}(t-t_0)}, \quad \forall t \ge t_0. \qquad (73)
$$

9 Plugging this into (71), we prove the claim.

10

■

11

7. COMMENTS, EXTENSIONS AND FUTURE WORK

12 In this paper, we adopted a dynamical system approach to study some stochastic opti-
mization problems with state-dependent distributions. We investigated the existence and
13 uniqueness of equilibrium points, well-posedness as well as convergence properties of the
14 trajectories, for both first and second-order dynamics. We highlighted some dynamical
15 and geometrical properties of the state-dependent distributions suggesting that the natural
16 framework to study problems of the form (5) is the one of *metric random walk spaces*. More
17 particularly, the notion of coarse Ricci curvature gives a new insight on the geometrical
18 hidden structure of this kind of problems. Finally, we discussed as an application the inertial
19 primal-dual algorithm. We present here some ongoing works, possible extensions as well as
20 some open problems.
21

20

21
**Inertial algorithms.** Relying on the discretization of the dynamics studied in Section 4 and Section 6, more specifically, $(ISEHD_{2\sqrt{\mu},\omega})$ and $(ISPDS_{2\sqrt{\mu}})$, we obtain new inertial algorithms with Hessian-driven damping for stochastic optimization problems with decision-dependent distributions. These algorithms exhibit rapid convergence properties and can also be adapted to the nonsmooth case. This is being addressed in ongoing work.

**Implicit Hessian damping.** We focused in Section 4 on the explicit Hessian damping in the smooth case. Yet, it is possible to consider implicit damping as in [9, 14]

$$ \ddot{x}(t) + \gamma(t)\dot{x}(t) + \nabla G_{m_{\bar{x}}} (x(t) + \omega\dot{x}(t)) + e_{\bar{x}}(x(t)) = 0, \quad (\text{ISIHD}_{m,\gamma}) $$

The dynamics $(\text{ISIHD}_{m,\gamma})$ is referred to as an Inertial System with Implicit Hessian damping, since one can observe, using Taylor expansion

$$ \nabla G_{m_{\bar{x}}}(x(t) + \omega \dot{x}(t)) \approx \nabla G_{m_{\bar{x}}} x(t) + \nabla^2 G_{m_{\bar{x}}}(x(t)) \dot{x}(t). $$

As it was observed in [14], higher-order moments of the perturbation $e_{\bar{x}}$ are required to get fast convergence guarantees in the implicit case compared to the explicit one. Since in our analysis (see Theorem 3) no integrability assumption on $e_{\bar{x}}$ is needed, it is interesting to investigate the effect of implicit Hessian damping, both in the smooth and nonsmooth cases.

**Tikhonov-regularization.** In Section 4 we restricted ourselves, for sake of simplicity, the analysis to the case where the operator $F_{m_{\bar{x}}}$ is smooth. However, it is possible to consider second-order dynamics for general (and possibly nonpotential) operators, by considering, for $\lambda > 0$ the following dynamic

$$ \ddot{x}(t) + \gamma(t)\dot{x}(t) + F_{m_{\bar{x}}(t)}^{\lambda(t)}(x(t)) + \omega \frac{d}{dt} \left( F_{m_{\bar{x}}(t)}^{\lambda(t)}(x(t)) \right) + e_{\bar{x}}(x(t)) = 0. \quad (\text{ISEHD}_{\lambda,\gamma}) $$

where $F_{m_{\bar{x}}}^{\lambda}$ is the so-called Yosida approximation of $F_{m_{\bar{x}}}$ defined by $F_{m_{\bar{x}}}^{\lambda} = \frac{1}{\lambda}(\text{id} - J_{\lambda F_{m_{\bar{x}}}})$ and $J_{\lambda F_{m_{\bar{x}}}} = (\text{id} + JF_{m_{\bar{x}}})^{-1}$ is the resolvent of $F_{m_{\bar{x}}}$. This approach comes with several advantages. First, the Yosida approximation is single-valued so that there is no nonsmoothness to take care of. In addition one can exploit the $\lambda$-cocoercivity of $F_{m_{\bar{x}}}^{\lambda}$ and the fact that zer $F_{m_{\bar{x}}} = \text{zer } F_{m_{\bar{x}}}^{\lambda}$. The approach was used in [17] for $\omega = 0$ and in the recent work [12] for Newton-like dynamics. We are exploring the adaptation of this techniques to stochastic monotone inclusions with state-dependent distribution in an ongoing work.

**Weaker Assumptions.** We have seen that one of the crucial assumptions in the analysis is Assumption 2, which concerns the Lipschitz behavior of the distribution $(m_x)_x$. A natural question that arises is what happens under a weaker assumption. For example, when $x \mapsto m_x$ is Hölder continuous. We are not aware of any existing results in this direction. We plan to investigate this in future work.

### ACKNOWLEDGEMENT

The work of H.E was supported by the ANR grant, reference: ANR-20-CE92-0037.

### A. GRONWALL INEQUALITIES

In this section we list several auxiliary results that we make use of in the paper.

**Lemma 8 (Gronwall's lemma: differential form).** Let $u, v$ be two $C^0$ (resp. $C^1$) nonnegative function on $[0, T]$ and let $w$ be a continuous function on $[0, T]$. We assume that

$$ \frac{1}{2} \frac{d}{dt} u^2(t) \le w(t) u^2(t) + u(t) v(t) \quad \text{on } (0, T), \quad (74) $$
1 then, for any $t \in [0, T]$

$$u(t) \le u(0)e^{K(t)} + \int_{0}^{t} v(s)e^{K(t)-K(s)}ds, \qquad (75)$$

2 where $K(t) = \int_{0}^{t} w(s)ds$.

3 **Lemma 9.** [26, Lemma A.5] Let $v \in L^1(t_0, T; \mathbb{R}^+)$ and $u \in C^0(t_0, T)$ such that

$$\frac{1}{2}u^2(t) \le \frac{1}{2}c^2 + \int_{t_0}^{t} u(s)v(s)ds,$$

4 for some $c \ge 0$ for all $t \in [t_0, T]$. Then

$$|u(t)| \le c + \int_{t_0}^{t} v(s)ds.$$

## B. BANACH FIXED POINT THEOREM & PICARD ITERATIVE METHOD

6 **Theorem 6** (see, e.g., [6, 27]). Let $(\mathcal{X}, d)$ be a complete metric space and $S: \mathcal{X} \to \mathcal{X}$ be a strict contraction, i.e., there exists a constant $\rho < 1$ such that

$$d(S(x), S(y)) \le \rho d(x, y), \forall x, y \in \mathcal{X}.$$

8 Then, there exists a unique $\bar{x} \in \mathcal{X}$ such that $S(\bar{x}) = \bar{x}$. Moreover, for any $x_0 \in \mathcal{X}$, the sequence starting from $x_0$ with $x_{n+1} = S(x_n)$ for all $n \in \mathbb{N}$ converges to $\bar{x}$ as $n$ goes to $\infty$.

## REFERENCES

11 [1] F. Alvarez. On the minimizing property of a second order dissipative system in Hilbert spaces. *SIAM J. Control Optim.*, 38(4):1102–1119, 2000.

12 [2] F. Alvarez and H. Attouch. An inertial proximal method for maximal monotone operators via discretization of a nonlinear oscillator with damping. *Set-Valued Anal.*, 9(1-2):3–11, 2001.

13 [3] F. Alvarez, H. Attouch, J. Bolte, and P. Redont. A second-order gradient-like dissipative dynamical system with hessian-driven damping: Application to optimization and mechanics. *Journal de mathématiques pures et appliquées*, 81(8):747–779, 2002.

14 [4] L. Ambrosio, N. Gigli, and G. Savaré. Gradient flows in metric spaces and in the space of probability measures. Basel: Birkhäuser, 2nd ed. edition, 2008.

15 [5] V. Apidopoulos, J.-F. Aujol, and C. Dossal. Convergence rate of inertial forward-backward algorithm beyond Nesterov’s rule. *Mathematical Programming*, 180(1):137–156, 2020.

16 [6] H. Attouch, G. Buttazzo, and G. Michaille. Variational analysis in Sobolev and BV spaces, volume 17 of MOS-SIAM Series on Optimization. Society for Industrial and Applied Mathematics (SIAM), Philadelphia, PA; Mathematical Optimization Society, Philadelphia, PA, second edition, 2014. Applications to PDEs and optimization.

17 [7] H. Attouch, A. Cabot, and P. Redont. The dynamics of elastic shocks via epigraphical regularization of a differential inclusion. *Barrier and penalty approximations. Adv. Math. Sci. Appl.*, 12(1):273–306, 2002.

18 [8] H. Attouch, Z. Chbani, J. Fadili, and H. Riahi. Fast convergence of dynamical ADMM via time scaling of damped inertial dynamics. *J. Optim. Theory Appl.*, 193(1-3):704–736, 2022.

19 [9] H. Attouch, Z. Chbani, J. Fadili, and H. Riahi. First-order optimization algorithms via inertial systems with Hessian driven damping. *Math. Program.*, 193(1, Ser. A):113–155, 2022.

20 [10] H. Attouch, Z. Chbani, J. Peypouquet, and P. Redont. Fast convergence of inertial dynamics and algorithms with asymptotic vanishing viscosity. *Mathematical Programming*, 168:123–175, 2018.

21 [11] H. Attouch, Z. Chbani, and H. Riahi. Rate of convergence of the nesterov accelerated gradient method in the subcritical case $\alpha \le 3$. *ESAIM: Control, Optimisation and Calculus of Variations*, 25:2, 2019.

22 [12] H. Attouch and S. Csaba László. Continuous Newton-like inertial dynamics for monotone inclusions. *Set-Valued Var. Anal.*, 29(3):555–581, 2021.

23 [13] H. Attouch and A. Damlamian. On multivalued evolution equations in hilbert spaces. *Israel Journal of Mathematics*, 12:373–390, 1972.
1

[14] H. Attouch, J. Fadili, and V. Kungurtsev. On the effect of perturbations in first-order optimization methods with inertia and Hessian driven damping. *Evol. Equ. Control Theory*, 12(1):71–117, 2023.

2

3

[15] H. Attouch, X. Goudou, and P. Redont. The heavy ball with friction method, i. the continuous dynaamical system: global exploration of the local minima of a real-valued function by asymptotic analysis of a dissipative dynamical system. *Communications in Contemporary Mathematics*, 2(01):1–34, 2000.

4

5

6

7

[16] H. Attouch and J. Peypouquet. The rate of convergence of nesterov's accelerated forward-backward method is actually faster than 1/k². *SIAM Journal on Optimization*, 26(3):1824–1834, 2016.

8

10

[17] H. Attouch and J. Peypouquet. Convergence of inertial dynamics and proximal algorithms governed by maximally monotone operators. *Math. Program.*, 174(1-2 (B)):391–432, 2019.

11

[18] H. Attouch, J. Peypouquet, and P. Redont. A dynamical approach to an inertial forward-backward algorithm for convex minimization. *SIAM Journal on Optimization*, 24(1):232–256, 2014.

12

14

[19] H. Attouch, J. Peypouquet, and P. Redont. Fast convex optimization via inertial dynamics with hessian driven damping. *Journal of Differential Equations*, 261(10):5734–5783, 2016.

15

[20] H. Bahouri, J.-Y. Chemin, and R. Danchin. *Fourier analysis and nonlinear partial differential equations*, volume 343 of *Grundelehren der mathematischen Wissenschaften* [Fundamental Principles of Mathematical Sciences]. Springer, Heidelberg, 2011.

16

19

[21] J. Baillon and H. Brézis. Une remarque sur le comportement asymptotique des semigroupes non linéaires. *Houston J. Math.*, 2:5–7, 1976.

20

[22] V. Barbu. *Nonlinear differential equations of monotone types in Banach spaces*. Springer Monographs in Mathematics. Springer, New York, 2010.

21

22

[23] H. H. Bauschke and P. L. Combettes. *Convex analysis and monotone operator theory in Hilbert spaces*. CMS Books in Mathematics/Ouvrages de Mathématiques de la SMC. Springer, New York, 2011. With a foreword by Hédy Attouch.

23

24

25

[24] P. Bianchi, W. Hachem, and A. Salim. A fully stochastic primal-dual algorithm. *Optimization Letters*, 15(2):701–710, 2021.

26

27

[25] R. I. Bołt, E. R. Csetnek, and S. C. László. Tikhonov regularization of a second order dynamical system with Hessian driven damping. *Math. Program.*, 189(1-2 (B)):151–186, 2021.

28

29

[26] H. Brézis. *Opérateurs maximaux monotones et semi-groupes de contractions dans les espaces de Hilbert*. North-Holland Mathematics Studies, No. 5. North-Holland Publishing Co., Amsterdam-London; American Elsevier Publishing Co., Inc., New York, 1973.

30

31

32

[27] H. Brézis. *Functional analysis, Sobolev spaces and partial differential equations*. Universitext. Springer, New York, 2011.

33

34

[28] R. E. Bruck Jr. Asymptotic convergence of nonlinear contraction semigroups in hilbert space. *Journal of Functional Analysis*, 18(1):15–26, 1975.

35

36

[29] A. Chambolle and T. Pock. A first-order primal-dual algorithm for convex problems with applications to imaging. *Journal of mathematical imaging and vision*, 40:120–145, 2011.

37

38

[30] L. Condat. A primal-dual splitting method for convex optimization involving lipschitzian, proximable and linear composite terms. *Journal of optimization theory and applications*, 158(2):460–479, 2013.

39

40

[31] J. Cutler, M. Díaz, and D. Drusvyatskiy. Stochastic approximation with decision-dependent distributions: asymptotic normality and optimality. *arXiv preprint arXiv:2207.04173*, 2022.

41

42

[32] A. Derrow-Pinion, J. She, D. Wong, O. Lange, T. Hester, L. Perez, M. Nunkesser, S. Lee, X. Guo, B. Wiltshire, P. W. Battaglia, V. Gupta, A. Li, Z. Xu, A. Sanchez-Gonzalez, Y. Li, and P. Velickovic. Eta prediction with graph neural networks in google maps. In *Proceedings of the 30th ACM International Conference on Information & Knowledge Management*, CIKM '21, pages 3767–3776, New York, NY, USA, 2021. Association for Computing Machinery.

43

44

45

46

47

[33] R. Douc, E. Moulines, P. Priouret, and P. Soulier. *Markov chains*. Springer Series in Operations Research and Financial Engineering. Springer, Cham, 2018.

48

49

[34] D. Drusvyatskiy and L. Xiao. Stochastic optimization with decision-dependent distributions. *Mathematics of Operations Research*, 2022.

50

51

[35] X. He, R. Hu, and Y. Fang. A second order primal-dual dynamical system for a convex-concave bilinear saddle point problem. *Applied Mathematics & Optimization*, 89(2):30, 2024.

52

53

[36] O. Hernández-Lerma and J. B. Lasserre. *Markov chains and invariant probabilities*, volume 211 of *Prog. Math*. Basel: Birkhäuser, 2003.

54
1

[37] D. Kim. Accelerated proximal point method for maximally monotone operators. *Mathematical Programming*, 190(1):57-87, 2021.

2

3

[38] S. C. László. Convergence rates for an inertial algorithm of gradient type associated to a smooth non-convex minimization. *Mathematical Programming*, 190(1):285-329, 2021.

4

5

[39] T. Lin and M. I. Jordan. A control-theoretic perspective on optimal high-order optimization. *Mathematical Programming*, pages 1-47, 2022.

6

7

[40] J. Macfarlane. Your navigation app is making traffic unmanageable. *IEEE Spectrum*, 19, 2019.

8

9

[41] J. Mazón, M. Solera-Diana, and J. Toledo-Melero. *Variational and Diffusion Problems in Random Walk Spaces*. Progress in Nonlinear Differential Equations and Their Applications. Springer Nature Switzerland, 2023.

10

11

12

13

[42] C. Mendler-Dünner, J. C. Perdomo, T. Zrnic, and M. Hardt. Stochastic optimization for performative prediction. In *Proceedings of the 34th International Conference on Neural Information Processing Systems*, NIPS’20, Red Hook, NY, USA, 2020. Curran Associates Inc.

14

15

[43] Y. Nesterov. A method for solving the convex programming problem with convergence rate $o(1/k^2)$. *Proceedings of the USSR Academy of Sciences*, 269:543-547, 1983.

16

[44] Y. Ollivier. Ricci curvature of Markov chains on metric spaces. *J. Funct. Anal.*, 256(3):810-864, 2009.

17

18

[45] Y. Ollivier. A survey of Ricci curvature for metric spaces and markov chains. In *Probabilistic approach to geometry*, volume 57, pages 343-382. Mathematical Society of Japan, 2010.

19

20

21

[46] J. Perdomo, T. Zrnic, C. Mendler-Dünner, and M. Hardt. Performative prediction. In H. D. III and A. Singh, editors, *Proceedings of the 37th International Conference on Machine Learning*, volume 119 of *Proceedings of Machine Learning Research*, pages 7599-7609. PMLR, 13-18 Jul 2020.

22

23

[47] B. T. Polyak. Some methods of speeding up the convergence of iteration methods. *Ussr computational mathematics and mathematical physics*, 4(5):1-17, 1964.

24

25

[48] H. Raguet, J. Fadili, and G. Peyré. A generalized forward-backward splitting. *SIAM Journal on Imaging Sciences*, 6(3):1199-1226, 2013.

26

27

[49] B. Shi, S. S. Du, M. I. Jordan, and W. J. Su. Understanding the acceleration phenomenon via high-resolution differential equations. *Mathematical Programming*, pages 1-70, 2022.

28

29

[50] W. Su, S. Boyd, and E. J. Candès. A differential equation for modeling Nesterov’s accelerated gradient method: theory and insights. *J. Mach. Learn. Res.*, 17:43, 2016. Id/No 153.

30

31

[51] B. C. Vū. A splitting algorithm for dual monotone inclusions involving cocoercive operators. *Advances in Computational Mathematics*, 38(3):667-681, 2013.

32

33

[52] K. Wood and E. Dall’Anese. Stochastic saddle point problems with decision-dependent distributions. *SIAM Journal on Optimization*, 33(3):1943-1967, 2023.