# Swarm-based gradient descent
meets simulated annealing

Zhiyan Ding  Department of Mathematics, University of California, Berkeley
<EMAIL_ADDRESS>,  Martin Guerra  Department of Mathematics, University
of Wisconsin-Madison<EMAIL_ADDRESS>,  Qin Li  Department of Mathematics,
University of Wisconsin-Madison<EMAIL_ADDRESS>and  Eitan Tadmor
Department of Mathematics and Institute for Physical Science & Technology
University of Maryland, College Park<EMAIL_ADDRESS>
###### Abstract.

We introduce a novel method for non-convex optimization which is at the
interface between the swarm-based gradient-descent (SBGD) [21, 28], and
Simulated Annealing (SA) [5, 18, 11]. We follow the methodology of SBGD in
which a swarm of agents, each identified with a position, $\mathbf{x}$ _and_
mass $m$, explores the ambient space. The agents proceed in gradient descent
direction, and are subject to Brownian motion with annealing-rate dictated by
a decreasing function of their mass. Thus, instead of the SA protocol for
time-decreasing temperature, we let the swarm decide how to ‘cool down’
agents, depending on their accumulated mass over time. The dynamics of masses
is coupled with the dynamics of positions: agents at higher ground transfer
(part of) their mass to those at lower ground. Consequently, the swarm is
dynamically divided between heavier, cooler agents viewed as ‘leaders’ and
lighter, warmer agents viewed as ‘explorers’. Mean-field convergence analysis
and benchmark optimizations demonstrate the effectiveness of the swarm-based
method as a multi-dimensional global optimizer.

###### Key words and phrases:

Optimization, stochastic gradient descent, swarm-based gradient descent,
provisional minimum, simulated annealing.

###### 1991 Mathematics Subject Classification:

90C26,65K10,92D25

Acknowledgments. Research was supported by NSF grant DMS-2012292 (MG), by NSF
grant DMS-2308440 and ONR grant N00014-21-1-2140 (QL) and by ONR grant
N00014-2112773 (ET). ET thanks the Fondations Sciences Mathématiques des Paris
(FSMP) and LJLL at the Sorbonne University for the support and hospitality.

###### Contents

  1. 1 Introduction
  2. 2 Statement of main results
  3. 3 From mean-field to macroscopic description
  4. 4 Pitfalls
  5. 5 Proofs of the main results
  6. 6 Numerical experiments
  7. A Interacting vs. non-interacting agents
  8. B Technical lemmas

##  1\. Introduction

We introduce a new swarm-based optimization method to compute the global
minimum of non-convex objective functions,
$\displaystyle\mathbf{x}_{*}:=\operatorname*{\textnormal{argmin}}_{\Omega\subset\mathbb{R}^{d}}F(\mathbf{x})$.
The swarm consists of $N$ agents — enumerated $j=1,2,\ldots,N$, each of which
is identified with a time-dependent position,
$\mathbf{x}^{j}_{t}\in{\mathbb{R}}^{d}$, and time-dependent mass,
$m^{j}_{t}\in{\mathbb{R}}_{+}$. The positions are governed by overdamped
Langevin process111 A word about notation. Throughout the paper the time-
dependent quantities in deterministic models are denoted $\square(t,\cdots)$
while time-dependent quantities which arise from stochastic descriptions are
denoted $\square_{t}(\cdots)$.,

$\textnormal{d}\mathbf{x}^{j}_{t}=-\nabla
F\left(\mathbf{x}^{j}_{t}\right)\textnormal{d}t+\sqrt{2\sigma^{j}_{t}}\textnormal{d}W^{j}_{t},\qquad
j=1,2,\ldots,N,$

where $\left\\{W^{j}_{t}\right\\}^{N}_{j=1}$ are independent Brownian motions
with amplitude $\sqrt{2\sigma^{j}_{t}}$, where $\sigma^{j}_{t}$ is the scaled
“temperature” (dating back to Langevin [20]), or is referred to as the
_annealing-rate_ , [18, 5, 12, 11, 7]. The intricate aspect of the dynamics is
a proper tuning of these annealing-rates $\sigma^{j}_{t}$. This is where the
masses, $\\{m^{j}_{t}\\}$, come into play: they are driven by

$\textnormal{d}m^{j}_{t}=-m^{j}_{t}\left(F\left(\mathbf{x}^{j}_{t}\right)-\overline{F}^{N}_{t}\right)\textnormal{d}t,\qquad
j=1,2,\ldots,N,$

where $\overline{F}^{N}_{t}$ is the mass-weighted average of the $N$ agents,
$\displaystyle\overline{F}^{N}_{t}:=\frac{\sum_{j=1}^{N}m^{j}_{t}F(\mathbf{x}^{j}_{t})}{\sum_{j=1}^{N}m^{j}_{t}}$.
We tune the annealing rates as a decreasing function of the mass,
$\sigma^{j}_{t}=\sigma(m^{j}_{t})$. Thus, our protocol for annealing-rate
prevents lighter agents from being trapped in basins of _local_ minima while
‘cools down’ heavier agents at the basin of attraction of _global_ minimum of
$F$. This is where the swarm-based gradient descent [21] meets Simulated
Annealing (SA), [18, 10, 5], except that the protocol for “cooling process” is
different: instead of an explicit recipe for decreasing the temperature as
done in SA, our swarm-based approach let different agents at a lower ground
cool down, depending on their increasing mass
$\sigma^{j}_{t}=\sigma(m^{j}_{t})$.

Why provisional minimum? We refer to the mass-weighted average
$\overline{F}^{N}_{t}:=\frac{\sum_{j=1}^{N}m^{j}_{t}F(\mathbf{x}^{j}_{t})}{\sum_{j=1}^{N}m^{j}_{t}}$
as the _provisional minimum_. To clarify why provisional minimum is needed, we
note that according to the mass equation above, agents that are above the
provisional minimum will shed a fraction of their mass — masses that are
transferred to agents that are below the provisional minimum. Since the total
mass remains constant, $\sum_{j}\textnormal{d}m^{j}_{t}/\textnormal{d}t=0$,
masses tend to concentrate with agents near or below the provisional minimum,
and it therefore makes sense to adjust the annealing rate, $\sigma^{j}_{t}$,
as a _decreasing_ function of the mass,

$\sigma^{j}_{t}=\sigma(m^{j}_{t}):(0,\infty)\mapsto{\mathbb{R}}_{+},\qquad
j=1,2,\ldots,N.$

Given this adjustment of mass-dependent annealing, _it is expected_ that for
large crowd of agents, $N\gg 1$, the provisional minimum will approach the
global minimum,

$\lim_{N\rightarrow\infty}\overline{F}^{N}_{t}\stackrel{{\scriptstyle
t\rightarrow\infty}}{{\longrightarrow}}F_{*}:=F(\mathbf{x}_{*}).$

Indeed, our main result states that, under the appropriate assumptions
outlined in Theorem 2.3 below, this is true.

The method described above lies at the interface of the deterministic swarm-
based gradient descent (SBGD) approach and the stochastic-based SA. We recall
that the SBGD is governed by a deterministic system which governs the swarm of
$N$-agents

(1.1)
$\left\\{\begin{split}\frac{\textnormal{d}}{\textnormal{d}t}\mathbf{x}^{j}(t)&=-\alpha(m^{j}(t))\nabla
F(\mathbf{x}^{j}(t)),\\\
\frac{\textnormal{d}}{\textnormal{d}t}m^{j}(t)&=-\big{(}F(\mathbf{x}^{j}(t))-\overline{F}^{N}(t)\big{)}m^{j}(t),\end{split}\right.\qquad
j=1,2,\ldots,N.$

The key feature advocated in (1.1), is swarm dynamics of both — positions and
masses embedded in ${\mathbb{R}}^{d}\times{\mathbb{R}}_{+}$. The additional
dimension of mass serves as a platform for _communication among agents_ ,
encoded in the provisional minimum, which in [21] was taken as the minimum of
the crowd at the given time,
$\overline{F}^{N}(t)=\min_{j}F(\mathbf{x}^{j}(t))$. The same swarm-based
methodology is used in this work, except two distinct features. First, the
provisional minimum is given here by the weighted average
$\overline{F}^{N}_{t}:=\frac{\sum_{j}m^{j}_{t}F(\mathbf{x}^{j}_{t})}{\sum_{j}m^{j}_{t}}$
which eventually is expected to approach the global minimum,
$\min_{j}F(\mathbf{x}^{j}(t))$, used in [21]. The second and more essential
distinction, is how masses are being used to adjust the dynamics of positions
along the gradient descent: in SBGD, it is the _time step_ which is adjusted
as a decreasing function of the mass, $\alpha=\alpha(m^{j}(t))$. In the swarm-
based optimization proposed in this work, masses are used to adjust the
annealing rate, $\sigma=\sigma(m^{j}_{t})$, as decreasing function of mass,

(1.2) $\left\\{\begin{split}\textnormal{d}\mathbf{x}^{j}_{t}&=-\nabla
F(\mathbf{x}^{j}_{t})\textnormal{d}t+\sqrt{2\sigma(m^{j}_{t})}\textnormal{d}W^{j}_{t},\\\
\textnormal{d}m^{j}_{t}&=-\big{(}F(\mathbf{x}^{j}_{t})-\overline{F}^{N}_{t}\big{)}m^{j}_{t}\textnormal{d}t,\quad\overline{F}^{N}_{t}:=\frac{\sum_{j=1}^{N}m^{j}_{t}F(\mathbf{x}^{j}_{t})}{\sum_{j=1}^{N}m^{j}_{t}},\end{split}\right.\qquad
j=1,2,\ldots,N.$

######  Remark 1.1 (On the choice of provisional minimum).

The SBGD method (1.1) employs the actual minimum,
$\overline{F}(t)=\min_{j}F(\mathbf{x}^{j}(t))$ as its provisional minimum [21,
28]. This, however, is not amenable to mean-field limit that we will pursue in
Section 2.1. The mean-field limit requires us to tour the whole landscape,
ending with

$\lim_{N\rightarrow\infty}\overline{F}^{N}(t)=\lim_{N\rightarrow\infty}\min_{j=1,\ldots,N}F(\mathbf{x}^{j}(t)),$

which is already the global minimum we are looking for. Instead, the use of
mass-weighted average $\overline{F}^{N}_{t}$ in (1.2) as a provisional
minimum, admits a mean-field interpretation which is shown to be driven
towards the same desired global minimum.

###  1.1. In-swarm communication combined with stochastic search

Application of deterministic swarm-based methods reveals that such methods
_succeed in the accuracy sense but fail in the probability sense_ , namely —
when these methods succeed, then they find the global optimizer, yet there is
non-zero probability for the methods to fail to do so. In contrast, stochastic
methods _succeed in the probability sense but fail in the accuracy sense_ ,
namely — such methods can always find the global basin, but they do not
necessarily return a true global optimizer. In this context, our swarm-based
optimization with adjusted annealing rate takes advantage of randomness in
(1.2)1 combined with communication-based swarming in (1.2)2, with the aim of
succeeding in both the probability and the accuracy sense. This is achieved by
a dynamic process that combines swarm communication and randomness, as
described below.

The essential role of communication in SBGD was already emphasized in [21,
28]. This is amplified in the present context, upon setting
$\sigma(\cdot)\equiv 0$ in (1.2)1. One ends up with a crowd of $N$ independent
agents driven by gradient descent,

$\dot{\mathbf{x}}^{j}(t)=-\nabla F(\mathbf{x}^{j}(t)),\qquad j=1,2,\ldots,N,.$

The success of such a non-communicating crowd in exploring the ambient
landscape is _significantly worse_ than that of a communication-based crowd.
In Section 4, we show that the method fails in the probability sense.

Communication is encoded in the provisional minimum, where agents
‘communicate’ their height relative to $\overline{F}^{N}_{t}$. In our swarm-
based optimization with adjusted annealing rate (1.2), $\overline{F}^{N}_{t}$
is taken as the weighted average height, weighted by the different masses.
Consequently, there is a dynamic distinction between agents above the
provisional minimum and agents below or at the level of the provisional
minimum. Agents above the provisional minimum carry lighter mass and therefore
explore the ambient landscape with a relatively large Brownian motion; these
are viewed as the _explorers_ of the crowd. The heavier agents — those that
are below or at the level of the provisional minimum, evolve in the gradient
direction with a smaller Brownian motion and are viewed as _leaders_ of the
crowd. Of course, once an exploring light agent ‘hits’ a new lower ground at,
say, $\mathbf{x}^{j}_{t}$ where
$F(\mathbf{x}^{j}_{t})\ll-\overline{F}^{N}_{t}$, it is expected to accumulate
more mass from the higher agents, to cool down, and eventually to become a
leader, driving the crowd to a new lower point.

We note that by fixing a constant annealing rate, (1.2)1 with
$\sigma(m^{j}_{t})\equiv\sigma>0$, we also encounter a problem. In this case,
we recover Langevin Monte Carlo

$\textnormal{d}\mathbf{x}^{j}_{t}=-\nabla
F(\mathbf{x}^{j}_{t})\textnormal{d}t+\sqrt{2\sigma}\textnormal{d}W^{j}_{t},\qquad
j=1,2,\ldots,N,$

whose long time invariant measure suggests $\mathbf{x}\propto
e^{-F(\mathbf{x})/\sigma}$. It admits the proper convergence in the zero-
variance limit, $\mathbf{x}\stackrel{{\scriptstyle\sigma\rightarrow
0}}{{\longrightarrow}}\operatorname*{\textnormal{argmin}}_{\mathbf{x}}F(\mathbf{x})$,
but, as we show in Section 4 below, the resulting method fails in the accuracy
sense for any finite $\sigma$. Thus, the dynamic process of relabeling
explorers and leaders, coupled through a mass-dependent annealing rate,
$\sigma^{j}_{t}=\sigma(m^{j}_{t})$, is at the heart of matter.

In Section 5 we prove global convergence of (1.2) in both — sense of
probability and sense of accuracy. This is achieved in a two-step argument:
first, large-crowd, mean-field analysis summarized in Theorem 2.1 below, which
allows us to pass from finitely many Stochastic-ODEs governing the particle
system (1.2), to the limiting dynamics for the probability measure,
$\mu_{t}(\mathbf{x},m)$, governed by the limiting PDE (2.4); and second, the
large-time behavior of the latter, summarized in Theorem 2.2, implies that in
this mean-field limit, the PDE in long time returns the the provisional
minimum as the global minimum. Combining these results we conclude in Theorem
2.3 the claimed convergence
$\displaystyle\mathop{\liminf}_{t\rightarrow\infty}\lim_{N\rightarrow\infty}\overline{F}^{N}_{t}=F_{*}$
with a precise convergence rate of polynomial type. Further, in Section 3 we
derive the corresponding macroscopic description of the system.

###  1.2. The swarm-based optimization with adjusted annealing rate

The exchange of mass encoded in (1.2) assigns smaller masses to agents with
higher values, while the decreasing property of $\sigma(m)$ ensures that these
agents experience larger random perturbations. As a result, lighter agents
serve as explorers, while heavier agents serve as leaders. As noted in [21],
upon normalization of the total mass,
$\sum_{j}m^{j}_{t}\equiv\sum_{j}m^{j}_{0}=1$, masses can be interpreted as
‘probabilities of finding a global minimum’; heavier agents have larger
probability in doing so. The computation proceeds by discretization of (1.2)
using Euler-Maruyama’s formula with a proper time step $h>0$,

(1.3) $\mathbf{x}_{n+1}^{j}=\mathbf{x}_{n}^{j}-h\nabla
F(\mathbf{x}_{n}^{j})+\sqrt{2h\sigma\left(m_{n}^{j}\right)}\boldsymbol{\xi}^{j}$

for $j=1,2,...,N$ and $\boldsymbol{\xi}^{j}\sim\mathcal{N}(0,\mathbf{I}_{d})$
i.i.d for all $j$. The detailed algorithm is summarized in Algorithm 1.

The performance of the algorithm is presented in Section 6. This follows the
analytical study in Sections 2 and 5 in which we state, and respectively,
prove the main analytical results of the paper, regarding the large-crowd,
large-time dynamics associated with (1.2).

Algorithm 1 Swarm-based stochastic gradient descent

1:Input: Initial distribution $\mu_{0}(x,m)$, $h$ as step size, and $n_{T}$ as
$\\#$ of iterations

2:$\\{(\mathbf{x}_{0}^{j},m_{0}^{j})\\}_{j=1}^{N}$ are i.i.d. drawn from
$\mu_{0}$ for $j=1,2,...,N$

3:$\displaystyle\overline{F}_{0}\leftarrow\frac{\sum_{j=1}^{N}m_{0}^{j}F(\mathbf{x}_{0}^{j})}{\sum^{N}_{j=1}m^{j}_{0}}$
%
$\displaystyle\overline{F}_{n}\equiv\overline{F}_{n}^{N}=\frac{\sum^{N}_{j=1}m^{j}_{n}F(\mathbf{x}^{j}_{n})}{\sum^{N}_{j=1}m^{j}_{n}}$
the provisional minimum at $t^{n}$

4:for $n=0,1,2,...,n_{T}-1$ do

5: for $j=1,2,...,N$ do

6: $m_{n+1}^{j}\leftarrow
m_{n}^{j}-m_{n}^{j}h\big{(}F(\mathbf{x}_{n}^{j})-\overline{F}_{n}\big{)}$ for
$j=1,2,...,N$

7: end for

8: Generate $N$ samples i.i.d such that
$\boldsymbol{\xi}^{j}\sim\mathcal{N}(0,\mathbf{I}_{d})$

9: $\mathbf{x}_{n+1}^{j}\leftarrow\mathbf{x}_{n}^{j}-h\nabla
F(\mathbf{x}_{n}^{j})+\sqrt{2h\sigma(m_{n}^{j})}\boldsymbol{\xi}^{j}$ for
$j=1,2,...,N$

10:
$\displaystyle\overline{F}_{n+1}\leftarrow\sum_{j=1}^{N}m_{n+1}^{j}F(\mathbf{x}_{n+1}^{j})/\sum_{j=1}^{N}m_{n+1}^{j}$

11:end for

12:$j_{\text{opt}}\leftarrow\operatorname*{\textnormal{argmin}}_{j}F(\mathbf{x}^{j}_{n+1})$

13:Output: $\mathbf{x}_{n_{T}}^{j_{\text{opt}}}$ and $\overline{F}_{n_{T}}$

###  1.3. Related work

There are several well-known non-convex optimization algorithms which combine
swarm-based and/or stochastic effects. We mention biologically inspired
methods of ant colony optimization [23], artificial bee colony optimization
[16], and firefly optimization [30], and physically inspired methods of wind-
driven optimization, [1]. We also mention particle swarm optimization methods
which explore the state space with randomized drifts toward the best global
position, [17, 26] and its stochastic version in [14].

The stochastic part of our method is motivated by simulated annealing driven
by stochastic noise that is ‘cooled down’ as time evolves, [15, 24]. The key
feature in one-particle SA dynamics is the protocol for ‘cooling down’, i.e.,
setting a properly tuned decreasing-in-time annealing rate $\sigma_{t}\propto
c/\sqrt{\log{t}}$, [13, 12, 11, 19, 7]. This should be compared with our
swarm-based approach, in which the annealing rate of each agent is dictated as
a decreasing function of its mass222in fact – we allow a larger class of mass
scaling, covered in Assumptions 2.1 and 2.2.,
$\sigma_{t}\mapsto\sigma(m^{j}_{t})$; essentially, we let the swarm decide how
to ‘cool down’ the agents as they accumulate more mass.

The idea of gaining from the stochastic interaction of more than just one
agent can be found in the more recent work [6]. Here Langevin diffusion swaps
between two agents — a “global explorer” and a “local explorer”, with two
annealing rates which correspond to high and low temperatures. The dynamics
takes place in ${\mathbb{R}}^{d}\times{\mathbb{R}}^{d}$. In our method, there
is a crowd of $N$ agents which evolve their mass-dependent annealing rates in
$\big{(}\mathbb{R}^{d}\times\mathbb{R}_{+}\big{)}^{\otimes N}$.

In the present context, it is instructive to calibrate our method with the
Consensus Based Optimization (CBO), [25, 2, 29, 3, 4, 8]. The CBO is governed
by a swarm-based dynamics of SDEs

(1.4)
$\textnormal{d}\mathbf{x}^{j}_{t}=-\lambda\left(\mathbf{x}^{j}_{t}-\overline{F}^{\lambda}_{t}\right)\textnormal{d}t+\sigma\left|\mathbf{x}^{j}_{t}-\overline{F}^{\lambda}_{t}\right|\textnormal{d}W^{j}_{t},\qquad\overline{F}^{\lambda}_{t}:=\sum_{j=1}^{N}\mathbf{x}^{j}_{t}\left(\frac{\exp(-\lambda
F(\mathbf{x}^{j}_{t}))}{\sum_{i=1}^{N}\exp(-\lambda
F(\mathbf{x}^{i}_{t}))}\right),$

which are steering towards the weighted minimum $\overline{F}^{\lambda}_{t}$.
Indeed, the weighted min $\overline{F}^{\lambda}_{t}$ can be viewed as the
provisional minimum driving the CBO (1.4), in view of the Laplace principle,

$\lim_{\lambda\rightarrow\infty}\left(-\frac{1}{\lambda}\log\left(\int\omega^{\lambda}_{F}(\mathbf{x})\textnormal{d}\rho\right)\right)=\min_{\textnormal{supp}(\rho)}F,\qquad\omega^{\lambda}_{F}(\mathbf{x}):=\exp(-\lambda
F(\mathbf{x})).$

It is here that we observe the main difference between our swarm-based
optimization with adjusted annealing rate (1.2) and the CBO (1.4): while the
latter requires a user tuning of parameter $\lambda$ in
$\overline{F}^{\lambda}$, the swarm-based framework of [21] provides us with
an _adaptive_ mechanism which dynamically adjusts the provisional minimum
based on the mass distribution in different parts of the landscape,
$\overline{F}^{N}_{t}$, and the corresponding mean-field
$\overline{\mathscr{F}}^{\mu}_{t}$, outlined in (2.4) below.

##  2\. Statement of main results

###  2.1. From empirical distribution to mean-field

We define the _empirical distribution_ $\mu^{N}_{t}$, which records the
_ensemble_ distribution of $\\{\mathbf{x}^{j}_{t},m^{j}_{t}\\}_{j=1}^{N}$ that
obey (1.2), as a measure on the feature space $(\mathbf{x},m)$,

(2.1)
$\mu^{N}_{t}(\mathbf{x},m)=\frac{1}{N}\sum_{j=1}^{N}\delta_{\mathbf{x}^{j}_{t}}(\mathbf{x})\otimes\delta_{m^{j}_{t}}(m)\,.$

Recall that the provisional minimum is given by the weighted average in (1.2),
which is expressed as

(2.2)
$\overline{F}^{N}_{t}=\frac{\displaystyle\sum_{j=1}^{N}m^{j}_{t}F(\mathbf{x}^{j}_{t})}{\displaystyle\sum_{j=1}^{N}m^{j}_{t}}=\frac{\displaystyle\iint
mF(\mathbf{x})\hskip
0.56917pt\textnormal{d}\mu^{N}_{t}(\mathbf{x},m)}{\displaystyle\iint m\hskip
0.56917pt\textnormal{d}\mu^{N}_{t}(\mathbf{x},m)}.$

To formally derive the mean-field limit, we use a smooth test function
$\phi(x,m)$ that is compactly supported on $\mathbb{R}^{d}\times\mathbb{R}$.
By definition,

$\iint\phi(\mathbf{x},m)\mu^{N}_{t}(\mathbf{x},m)\textnormal{d}x\textnormal{d}m=\frac{1}{N}\sum^{N}_{j=1}\phi(\mathbf{x}^{j}_{t},m^{j}_{t})\,.$

Differentiating, calling Itô’s formula, and using (1.2) yield

$\displaystyle\frac{1}{N}\sum_{j=1}^{N}$
$\displaystyle\mathbb{E}\left[\textnormal{d}\phi(\mathbf{x}^{j}_{t},m^{j}_{t})\right]$
$\displaystyle=\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\left[\nabla_{\mathbf{x}}\phi(\mathbf{x}^{j}_{t},m^{j}_{t})\cdot\textnormal{d}\mathbf{x}^{j}_{t}+\frac{1}{2}\langle
D^{2}_{\mathbf{x}}\phi,2\sigma(m^{j}_{t})\rangle_{\text{Tr}}\textnormal{d}t+\partial_{m}\phi(\mathbf{x}^{j}_{t},m^{j}_{t})\textnormal{d}m^{j}_{t}\right]$
$\displaystyle=\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\underbrace{\left[\nabla_{\mathbf{x}}\phi\cdot\big{(}-\nabla
F(\mathbf{x}^{j}_{t})\textnormal{d}t+\sqrt{2\sigma(m^{j}_{t})}\textnormal{d}W^{j}_{t}\big{)}\right]}_{\text{Term
I}}+\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\underbrace{\left[\langle
D^{2}_{\mathbf{x}}\phi,\sigma(m^{j}_{t})\rangle_{\text{Tr}}\textnormal{d}t\right]}_{\text{Term
II}}$
$\displaystyle\quad+\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\underbrace{\left[-\partial_{m}\phi(\mathbf{x}^{j}_{t},m^{j}_{t})(F(\mathbf{x}^{j}_{t})-\overline{F}^{N}_{t})m^{j}_{t}\textnormal{d}t\right]}_{\text{Term
III}},$

where $D^{2}_{\mathbf{x}}\phi$ is the Hessian. The first term,
$\mathrm{I}=\mathrm{I}_{1}+\mathrm{I}_{2}$, where

$\begin{split}\mathrm{I}_{1}&:=\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\left[\nabla_{\mathbf{x}}\phi(\mathbf{x}^{j}_{t},m^{j}_{t})\cdot-\nabla
F(\mathbf{x}^{j}_{t})\right]=-\iint\nabla_{\mathbf{x}}\phi(\mathbf{x},m)\cdot\nabla
F(\mathbf{x})\mu^{N}_{t}(\mathbf{x},m)\textnormal{d}\mathbf{x}\textnormal{d}m\\\
&\ =\iint\phi(\mathbf{x},m)\nabla_{\mathbf{x}}\cdot(\mu^{N}_{t}\nabla
F(\mathbf{x}))\textnormal{d}\mathbf{x}\textnormal{d}m,\end{split}$

and $\mathrm{I}_{2}:=\mathbb{E}[\textnormal{d}W^{j}_{t}]=0$. For the second
term, we have

$\begin{split}\mathrm{II}&:=\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}[\langle
D^{2}_{\mathbf{x}}\phi(\mathbf{x}^{j}_{t},m^{j}_{t}),\sigma(m^{j}_{t})\rangle_{\text{Tr}}\textnormal{d}t]=\iint\Delta_{\mathbf{x}}\phi(\mathbf{x},m)\sigma(m)\mu^{N}_{t}(\mathbf{x},m)\textnormal{d}\mathbf{x}\textnormal{d}m\\\
&\
=\iint\phi(\mathbf{x},m)\sigma(m)\Delta_{\mathbf{x}}\mu^{N}_{t}(\mathbf{x},m)\textnormal{d}\mathbf{x}\textnormal{d}m,\end{split}$

and finally, the third term

$\displaystyle\mathrm{III}$
$\displaystyle:=\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\left[-\partial_{m}\phi(\mathbf{x}^{j}_{t}m^{j}_{t})m^{j}_{t}(F(\mathbf{x}^{j}_{t})-\overline{F}^{N}_{t})\right]$
$\displaystyle=-\iint\partial_{m}\phi(\mathbf{x},m)(m(F(\mathbf{x})-\overline{F}^{N}_{t}))\mu^{N}_{t}(\mathbf{x},m)\textnormal{d}\mathbf{x}\textnormal{d}m$
$\displaystyle=\iint\phi(\mathbf{x},m)\partial_{m}(m\mu^{N}_{t}(\mathbf{x},m))(F(\mathbf{x})-\overline{F}^{N}_{t})\textnormal{d}\mathbf{x}\textnormal{d}m\,.$

Adding all these terms together, we conclude that the dynamics of the
empirical distribution, $\mu^{N}_{t}$, is governed by the Vlasov equation

(2.3) $\partial_{t}\mu^{N}_{t}=\nabla_{\mathbf{x}}\cdot(\mu^{N}_{t}\nabla
F)+\big{(}F(\mathbf{x})-\overline{F}^{N}_{t}\big{)}\partial_{m}(m\mu^{N}_{t})+\sigma(m)\Delta_{\mathbf{x}}\mu^{N}_{t}\,.$

This formally provides the mean-field equation. In the mean-field limit,
$\mu=\mu_{t}(\mathbf{x},m)$ is governed by:

(2.4a) $\partial_{t}\mu=\nabla_{\mathbf{x}}\cdot(\mu\nabla
F)+\big{(}F(\mathbf{x})-\overline{\mathscr{F}}^{\mu}_{t})\big{)}\partial_{m}(m\mu)+\sigma(m)\Delta_{\mathbf{x}}\mu\,,$
driven by the corresponding provisional minimum (2.4b)
$\overline{\mathscr{F}}^{\mu}_{t}:=\frac{\mathbb{E}_{\mu}[mF(\mathbf{x})]}{\mathbb{E}_{\mu}[m]}=\frac{\displaystyle\iint
mF(\mathbf{x})\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)}{\displaystyle\iint m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)}\,.$

Observe that $\overline{\mathscr{F}}^{\mu}_{t}$ denotes the provisional
minimum associated with mean-field $\mu$, which is to be distinguished from
the discrete provisional minimum in (2.2), denoted $\overline{F}^{N}_{t}$. We
note that, similar to the discrete mass conservation in (1.2),
$\frac{1}{N}\sum_{j}m^{j}_{t}\equiv 1$, the total mass in (2.4a) is conserved
and equals one

$\frac{\textnormal{d}}{\textnormal{d}t}\int m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)=0,\quad\int m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)=1\,.$

Our first result quantifies the convergence of the empirical distribution,
$\mu^{N}_{t}$, to its limiting mean-field $\mu_{t}$. To this end, we make the
following assumption.

######  Assumption 2.1.

  * •

_Lipschitz bound_ — $F$, $\nabla F$ and $\sqrt{\sigma(\cdot)}$ are all
Lipschitz continuous functions with the Lipschitz constant denoted by $L$;

  * •

_Boundedness_ — There exists $D>0$ such that $|F|+\sqrt{\sigma}\leqslant D$.
Without loss of generality, we set $D\geqslant 1$ and can further assume
$F\geqslant 0$;

  * •

_Bounded support uniformly in $m$_ — there exists $M_{\max}<\infty$ such that

$\text{Supp}_{\mathbf{x}}(\mu_{0}(\mathbf{x},m))\subset[0,M_{\max}],\qquad\forall
m\in(0,\infty).$

The next theorem shows that the empirical distribution $\mu^{N}_{t}$ and the
mean-field $\mu$ stay close. The proof is postponed to Section 5.1.

######  Theorem 2.1 (Mean-field limit).

Assume that Assumption 2.1 holds. Let $\mu_{t}=\mu_{t}(\mathbf{x},m)$ be the
mean-field solution of (2.4a) and let
$\displaystyle\mu^{N}_{t}=\frac{1}{N}\sum_{j=1}^{N}\delta_{\mathbf{x}^{j}_{t}}(\mathbf{x})\otimes\delta_{m^{j}_{t}}(m)$
be the empirical distribution associated with the ensemble of swarm-based
solutions (1.2) subject to compatible initial data so that
$\\{\mathbf{x}^{j}_{0}\,,m^{j}_{0}\\}$ are i.i.d. samples of
$\mu_{t=0}(\mathbf{x},m)$. Then $\mu_{t}$ and $\mu^{N}_{t}$ are close in the
Wasserstein sense:

(2.5)
$\lim_{N\rightarrow\infty}\mathbb{E}\left[W_{2}\left(\mu_{t},\mu^{N}_{t}\right)\right]=0\,,$

and the corresponding provisional minimum (2.2), $\overline{F}^{N}_{t}$,
converges to $\overline{\mathscr{F}}^{\mu}_{t}$ in (2.4b), with the law of
large numbers rate: there exists a constant $C_{t}=C_{t}(\mu_{0},F,d)>0$
independent of $N$ such that

(2.6)
$\mathbb{E}\left[\left|\overline{F}^{N}_{t}-\overline{\mathscr{F}}^{\mu}_{t}\right|\right]<\frac{C_{t}}{\sqrt{N}}\,.$

The bound (2.6) provides the crucial convergence bound, which translates the
provisional minimum of the many-agent system to that of the limiting mean-
field equation.

###  2.2. Large-time convergence— from mean-field to global minimum

Upon translating the convergence of $\overline{F}^{N}_{t}$ to that of
$\overline{\mathscr{F}}^{\mu}_{t}$, we now switch gear to study the
convergence of the mean-field limiting PDE. To this end, we need another set
of assumptions for $F,\sigma$, and $\mu_{0}$:

######  Assumption 2.2.

We assume that the following conditions hold true:

  * •

$\sigma(m)$ has a bounded support: there exists an $m_{c}>0$ such that
$\sigma(m)=0$ when $m>m_{c}$.

  * •

$(m_{c},\infty)\cap\mathrm{supp}_{m}(\mu_{0})\neq\emptyset$ and $\displaystyle
M(0)=\iint m\hskip 0.56917pt\textnormal{d}\mu_{0}(\mathbf{x},m)=1$.

  * •

$F$ admits a unique global minimum
$\mathbf{x}_{*}:=\operatorname*{\textnormal{argmin}}_{\mathbf{x}\in\Omega}F(\mathbf{x})$,
and thus $\nabla F(\mathbf{x}_{*})=0$.

The next theorem shows that the provisional minimum of the mean-field,
$\overline{\mathscr{F}}^{\mu}_{t}$, is converging towards the global minimum
of $F$. The proof is postponed to Section 5.2 below.

######  Theorem 2.2 (Large time behavior).

Assume that Assumptions 2.1, 2.2 hold. Let $\mu$ be the mean-field with
provisional minimum $\overline{\mathscr{F}}^{\mu}_{t}$. We have

(2.7)
$\liminf_{t\rightarrow\infty}\overline{\mathscr{F}}^{\mu}_{t}=\min_{\mathbf{x}}F(\mathbf{x})=F_{*}\,.$

Specifically, — for any $\epsilon>0$, we let ${\mathcal{N}}_{\epsilon}$ denote
the following $\epsilon$ neighborhood:

(2.8)
${\mathcal{N}}_{\epsilon}:=\left\\{\mathbf{x}\,|\,F(\mathbf{x})<F_{*}+\epsilon\right\\}\,.$

Then, there exists a time $t_{\epsilon}(m_{c})<\infty$ defined by333Note that
$\int_{{\mathcal{N}}_{\epsilon}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)\leqslant\iint m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)=1$ for mass conservation. The
term in the logarithm is always bigger than $1$, so that the expression on the
right hand side of (2.9) quantifies positive time.:

(2.9)
$t_{\epsilon}(m_{c})=\frac{2}{\epsilon}\log\left(\left(\sup_{t\in[0,1]}\int_{{\mathcal{N}}_{\epsilon}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)\right)^{-1}\right)+1,$

and $t\leqslant t_{\epsilon}(m_{c})$ such that

(2.10) $\overline{\mathscr{F}}^{\mu}_{t}<F_{*}+\epsilon.$

Observe that even if
${\mathcal{N}}_{\epsilon}\cap\mathrm{supp}_{\mathbf{x}}\mu_{0}(\cdot,m)=\emptyset$,
then the presence of the Brownian motion implies that there exists a time,
$t\in(0,1]$, such that
${\mathcal{N}}_{\epsilon}\cap\mathrm{supp}_{\mathbf{x}}\mu_{t}(\cdot,m)\neq\emptyset$,
and consequently, (2.9) always yields a finite number threshold,
$t_{\epsilon}$. However, we do not claim that the ‘convergence rate’ worked
out below, $t_{\epsilon}\lesssim\epsilon^{-1}$, is optimal: our numerical
findings reported in Section 6 indicate the possibility of a much better
convergence rate in the sense that $t_{\epsilon}=\log(1/\epsilon)$ is already
sufficient to ensure (2.10).

Combing Theorems 2.1 and 2.2 we conclude the following main result of this
paper, namely — the convergence of the empirical distribution to the mean-
field.

######  Theorem 2.3.

Assume that Assumptions 2.1, 2.2 hold. Let $\mu^{N}_{t}(\mathbf{x},m)$ be the
empirical distribution associated with an ensemble of solutions of the swarm-
based optimization with adjusted annealing rate (1.2),
$(\mathbf{x}^{j}_{t},m^{j}_{t})$, subject to the initial data drawn i.i.d.
from $\mu_{t=0}$. Fix an arbitrarily $\epsilon>0$. Then, for a large enough
crowd spelled out in (2.6), depending on time, $t_{\epsilon/2}$ defined in
(2.9),

$N>\frac{4C_{t_{\epsilon/2}}^{2}}{\epsilon^{2}}\ \ \textnormal{where}\ \
t_{\epsilon/2}<\frac{4}{\epsilon}\log\left(\left(\sup_{t\in[0,1]}\int_{{\mathcal{N}}_{\epsilon/2}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)\right)^{-1}\right)+1,$

so that there is a $t<t_{\epsilon/2}$ and the provisional minimum
$\overline{F}^{N}_{t}$ is within $\epsilon$ of
$\min_{\mathbf{x}}F(\mathbf{x})$,

$\mathbb{E}\left[\left|\overline{F}^{N}_{t}-\min_{\mathbf{x}}F(\mathbf{x})\right|\right]<\epsilon.$

_Proof of Theorem 2.3._ By Theorem 2.2, there is a $t<t_{\epsilon/2}$ so that

$\left|\overline{\mathscr{F}}^{\mu}_{t}-\min_{\mathbf{x}}F(\mathbf{x})\right|<\frac{\epsilon}{2},$

and according to Theorem 2.1, there exists a constant, $C_{t}$, such that

$\mathbb{E}\left[\left|\overline{F}^{N}_{t}-\overline{\mathscr{F}}^{\mu}_{t}\right|\right]<\frac{C_{t}}{\sqrt{N}}.$

Hence, for large enough $N$ such that
$\displaystyle\frac{C_{t}}{\sqrt{N}}<\frac{\epsilon}{2}$, or equivalently
$N>\frac{4C^{2}_{t_{\epsilon/2}}}{\epsilon^{2}}$ (note
$C_{t_{\epsilon/2}}>C_{t}$), we conclude

$\begin{split}\mathbb{E}\left[\left|\overline{F}^{N}_{t}-\min_{\mathbf{x}}F(\mathbf{x})\right|\right]\leqslant\mathbb{E}\left[\left|\overline{F}^{N}_{t}-\overline{\mathscr{F}}^{\mu}_{t}\right|\right]+\left|\overline{\mathscr{F}}^{\mu}_{t}-\min_{\mathbf{x}}F(\mathbf{x})\right|\leqslant\epsilon\,.\quad\square\end{split}$

##  3\. From mean-field to macroscopic description

We shall briefly comment on the macroscopic description of the swarm-based
optimization with adjusted annealing rate, furnished in terms of the firsts
two moments — the density,
$\rho:{\mathbb{R}}_{+}\times{\mathbb{R}}^{d}\mapsto{\mathbb{R}}_{+}$, and
momentum, $\rho{\sf
v}:{\mathbb{R}}_{+}\times{\mathbb{R}}^{d}\mapsto{\mathbb{R}}$, (observe that
time-dependence in macroscopic quantities is denoted as usual
$\rho=\rho(t,\mathbf{x}),{\sf v}={\sf v}(t,\mathbf{x})$ etc.)

$\rho(t,\mathbf{x}):=\int\mu_{t}(\mathbf{x},m)\textnormal{d}m,\qquad\rho{\sf
v}(t,\mathbf{x}):=\int m\mu_{t}(\mathbf{x},m)\textnormal{d}m.$

Taking the first two moments of the mean-field (2.4), we find

(3.1) $\left\\{\begin{split}\rho_{t}+\nabla_{\mathbf{x}}\cdot(\rho\nabla
F)&=\Delta_{\mathbf{x}}\int\sigma(m)\mu_{t}(\mathbf{x},m)\textnormal{d}m,\\\
(\rho{\sf v})_{t}+\nabla_{\mathbf{x}}\cdot(\rho{\sf v}\nabla
F)&=\big{(}\overline{F}(t)-F(\mathbf{x})\big{)}\rho{\sf
v}+\Delta_{\mathbf{x}}\int
m\sigma(m)\mu_{t}(\mathbf{x},m)\textnormal{d}m.\end{split}\right.$

We normalize the total mass and momentum,
$\displaystyle\int\rho(t,\mathbf{x})\textnormal{d}\mathbf{x}\equiv 1$ and
$\displaystyle\int\rho{\sf v}(t,\mathbf{x})\textnormal{d}\mathbf{x}\equiv 1$,
which leaves us with a provisional minimum,
$\overline{F}(t)=\overline{\mathscr{F}}^{\mu}_{t}$, given by

$\overline{F}(t)=\int F(\mathbf{x})\rho{\sf
v}(t,\mathbf{x})\textnormal{d}\mathbf{x}.$

To unlock a closed form of the diffusion on the right, we assume the simplest
closure based on a mono-kinetic pseudo-Maxwellian ,
$\mu_{t}(\mathbf{x},m)=\rho(t,\mathbf{x})\delta(m-{\sf v}(t,\mathbf{x}))$.
This yields

(3.2a) $\displaystyle\rho_{t}+\nabla_{\mathbf{x}}\cdot(\rho\nabla F)$
$\displaystyle=\Delta_{\mathbf{x}}(\sigma(t,\mathbf{x})\rho(t,\mathbf{x})),\qquad\sigma(t,\mathbf{x}):=\sigma({\sf
v}(t,\mathbf{x}))$ (3.2b) $\displaystyle(\rho{\sf
v})_{t}+\nabla_{\mathbf{x}}\cdot(\rho{\sf v}\nabla F)$
$\displaystyle=\big{(}\overline{F}(t)-F(\mathbf{x})\big{)}\rho{\sf
v}+\Delta_{\mathbf{x}}\big{(}\sigma(t,\mathbf{x})\rho{\sf
v}(t,\mathbf{x})\big{)}.$

The momentum (3.2b) can be converted into a drift-diffusion equation for the
velocity ${\sf v}$,

(3.3a) $\displaystyle\rho_{t}+\nabla_{\mathbf{x}}\cdot(\rho\nabla F)$
$\displaystyle=\Delta_{\mathbf{x}}(\sigma(t,\mathbf{x})\rho(t,\mathbf{x})),$
(3.3b) $\displaystyle{\sf v}_{t}+\nabla F\cdot\nabla_{\mathbf{x}}{\sf v}$
$\displaystyle=\big{(}\overline{F}(t)-F(\mathbf{x})\big{)}{\sf
v}+\frac{1}{\rho}\Delta_{\mathbf{x}}\big{(}\sigma(t,\mathbf{x})\rho{\sf
v}(t,\mathbf{x})\big{)},\quad\mathbf{x}\in\textnormal{supp}\\{\rho(t,\cdot)\\}.$

We shall not dwell on a detailed study of (3.3), but note the decent estimate

$\frac{\textnormal{d}}{\textnormal{d}t}\int(F(\mathbf{x})-F_{*})\rho(t,\mathbf{x})\textnormal{d}\mathbf{x}=-\int|\nabla
F(\mathbf{x})|^{2}\rho(t,\mathbf{x})\textnormal{d}\mathbf{x}+\int\Delta_{\mathbf{x}}F(\mathbf{x})\sigma(t,\mathbf{x})\rho(t,\mathbf{x})\textnormal{d}\mathbf{x},$

with _expected_ long-time behavior $\rho(t,\mathbf{x})\stackrel{{\scriptstyle
t\rightarrow\infty}}{{\longrightarrow}}\delta(\mathbf{x}-\mathbf{x}_{*})$ and
${\sf v}(t,\mathbf{x})\stackrel{{\scriptstyle
t\rightarrow\infty}}{{\longrightarrow}}{\mathds{1}}(\mathbf{x}_{*})$. Observe
that $\sigma(t,\mathbf{x})$ should vanish as
$\mathbf{x}\rightarrow\mathbf{x}_{*}$.

##  4\. Pitfalls

We examine the swarm-based optimization with adjusted annealing rate in two
extreme cases — when $\sigma(m)\equiv 0$ and $\sigma(m)\equiv 1$. The former
case turns off stochasticity and amounts to a deterministic system of $N$
independent gradient descending agents; we show that their convergence fails
in probability sense. The latter turns on stochasticity uniformly for all
agents, which amounts to a stochastic system of $N$ independent gradient
descending agents; we show that their convergence fails in accuracy sense.

We note that in both cases, the resulting systems are passive systems in the
sense that the weights of agents do _not_ affect their trajectories. As a
consequence, agents roam the landscape with equal randomness with no
distinction between leaders and explorers. The convergence of the provisional
minimum $\overline{\mathscr{F}}^{\mu}_{t}$ is naturally worse than in the
communicating swarm-based dynamics.

###  4.1. Swarming with no communication

What is the effect of Brownian motion in our swarm-based optimization with
adjusted annealing rate? If we turn off the amplitude, $\sigma(m)\equiv 0$,
then in (1.2) is reduced to a deterministic system of _non-interacting_ agents
(and again we switch notations of time-dependence in deterministic quantities,
$\mathbf{x}^{j}(t),m^{j}(t),\overline{F}^{N}(t)$ etc.)

(4.1)
$\left\\{\begin{split}\frac{\textnormal{d}}{\textnormal{d}t}\mathbf{x}^{j}(t)&=-\nabla
F(\mathbf{x}^{j}(t))\\\
\frac{\textnormal{d}}{\textnormal{d}t}m^{j}(t)&=-m^{j}(t)\big{(}F(\mathbf{x}^{j}(t))-\overline{F}^{N}(t)\big{)},\quad\end{split}\right.j=1,...,N.$

We still assume random initial configurations: the initial drawing is
independent, with $\\{\mathbf{x}^{j}(0)\\}^{N}_{j=1}$ drawn from an initial
distribution $\rho_{0}$, and $m^{j}(0)=\frac{1}{N}$. Since the total mass is
conserved, we can be $M(t)=1$, so that provisional minimum is reduced to the
usual average
$\overline{F}^{N}(t)=\sum_{j=1}^{N}m^{j}(t)F(\mathbf{x}^{j}(t))$. The main
feature _missing_ in (4.1) is communication. Different agents at different
positions proceed along gradient decent, independently of each other: since
communication through provisional minimum is missing, the dynamics lacks the
collective engagement of the agents as a self-organized swarm.

The essential role of communication was already highlighted in the context of
SBGD method, [21]: once we decouple the step-size and mass, setting a
_uniform_ step-size in (1.1), that is, $\alpha(m)\equiv 1$, then the success
rate of the resulting non-communicating swarm-based method decreases
dramatically. The system (4.1) provides yet another realization of the
essential role of communication. We note that since all particles are
“descending” it is also straightforward to show that the collective mean
objective $\overline{F}^{N}(t)$ is always decreasing:
$\frac{\textnormal{d}}{\textnormal{d}t}\overline{F}^{N}(t)\leqslant 0$. See
Section A.1 for the proof. However, even though $\overline{F}^{N}(t)$
continues to decay in time, it does not guarantee the desired convergence to
the global minimum,
$\displaystyle\liminf_{t\rightarrow\infty}\overline{F}^{N}(t)=\min_{\mathbf{x}}F(\mathbf{x})$.
This is due to the passive nature of the gradient descent of the agents. In
particular, it is likely that all agents will become trapped at local minima
or saddle points, especially when the number of particles is not sufficiently
large. Specifically, we have the following.

######  Theorem 4.1 (Lack of communication and failure in probability).

Assume $F$ has a unique global minimum $\mathbf{x}_{*}$ inside an open set
$\Omega$ so that $F_{*}=F(\mathbf{x}_{*})$. Moreover, $\nabla F(\mathbf{x})=0$
on $\partial\Omega$. If we let $\\{\mathbf{x}^{j}(0)\\}$ drawn from
$\rho_{0}$, $m^{j}(0)=\frac{1}{N}$ and evolve according to (4.1), then the
success rate is bounded:

(4.2)
$\mathbb{P}\left(\lim_{t\rightarrow\infty}\overline{F}^{N}(t)=F_{*}\right)\leqslant
1-\left(1-\int_{\Omega}\textnormal{d}\rho_{0}\right)^{N}\,.$

In particular, if $\mathrm{supp}_{\mathbf{x}}(\rho_{0})\cap\Omega=\emptyset$,
the success rate is $0$, meaning
$\mathbb{P}\left(\lim_{t\rightarrow\infty}\overline{F}^{N}(t)=F_{*}\right)=0$
for all $N$.

_Proof_. To prove (4.2), we first notice that without the Brownian motion, the
deterministic system drives the samples to the local minima of the local basin
in which the samples are located at $t=0$, so the event of
$\left\\{\lim_{t\rightarrow\infty}\overline{F}^{N}(t)=\min
F\right\\}\subset\left\\{\exists\mathbf{x}^{i}(0)\in\Omega\right\\}$. This
implies

$\mathbb{P}\left(\lim_{t\rightarrow\infty}\overline{F}^{N}(t)=\min
F\right)\leqslant\mathbb{P}\left(\left\\{\exists\mathbf{x}^{j}(0)\in\Omega\right\\}\right)=1-\left(1-\int_{\Omega}\textnormal{d}\rho_{0}\right)^{N}\,.$

When $\mathrm{supp}_{\mathbf{x}}(\rho_{0})\cap\Omega=\emptyset$, we have
$\int_{\Omega}\textnormal{d}\rho_{0}=0$, so using this formula, we have for
all $N$,

$\mathbb{P}\left(\lim_{t\rightarrow\infty}\overline{F}^{N}(t)=\min
F\right)=0\,,$

Since $\overline{F}^{N}(t)$ is a decreasing in $t$, the above equality leads
to $\mathbb{P}\left(\lim_{t\rightarrow\infty}\overline{F}^{N}(t)>\min
F\right)=1.\hfill\square$

This theorem explicitly spells out the “probability failure”, so there is a
nonzero probability that the mean objectives obtained by the particles cannot
ever achieve the global optimizer.

###  4.2. Stochastic system

The limitation of (4.1) clearly originates from the fact that particle
evolution closely resembles classical gradient descent, lacking a mechanism
for domain exploration and avoiding entrapment in local minima. The other end
of the spectrum is to introduce “too much” exploration. Indeed, setting
$\sigma=1$, we have:

(4.3) $\left\\{\begin{split}\textnormal{d}\mathbf{x}^{j}_{t}&=-\nabla
F(\mathbf{x}^{j}_{t})\textnormal{d}t+\sqrt{2}\textnormal{d}W^{j}_{t}\\\
\textnormal{d}m^{j}_{t}&=-m^{j}_{t}\big{(}F(\mathbf{x}^{j}_{t})-\overline{F}^{N}_{t}\big{)}\textnormal{d}t\end{split}\right.,\quad
j=1,...,N$

where $W^{j}_{t}$ is the Brownian motion at time $t$ associated with particle
$j$. The only difference compared to (4.1) is the addition of a constant
Brownian motion term. This Brownian motion term allows more effective
exploration of the domain by preventing agents from becoming trapped at any
local minima.

The trajectory of $\\{\mathbf{x}^{j}_{t}\\}$ resembles that of Langevin Monte
Carlo, and under appropriate assumptions, it is well known that the
distribution of $\mathbf{x}^{j}_{t}$ converges to
$\rho_{\infty}\propto\exp(-F(\mathbf{x}))$ as $t\rightarrow\infty$.
Incorporating mass exchange in (4.3) does not eliminate this issue. Carrying
out the same mean-field analysis as was done in Section 2.1, we find that the
limiting PDE is:

(4.4) $\begin{split}\partial_{t}\mu=\nabla_{\mathbf{x}}\cdot(\mu\nabla
F)+\big{(}F(\mathbf{x})-\overline{\mathscr{F}}^{\mu}_{t}\big{)}\partial_{m}(m\mu)&+\Delta_{\mathbf{x}}\mu,\\\
\displaystyle\overline{\mathscr{F}}^{\mu}_{t}=\frac{\displaystyle\iint
mF(\mathbf{x})\mu_{t}(\mathbf{x},m)\textnormal{d}\mathbf{x}\textnormal{d}m}{\displaystyle\iint
m\mu_{t}(\mathbf{x},m)\textnormal{d}\mathbf{x}\textnormal{d}m}.\end{split}$

Analyzing this mean-field PDE, we could show that there is always a positive
gap between $\overline{\mathscr{F}}^{\mu}_{t}$ and $F_{*}$.

######  Theorem 4.2.

Consider a function $F$ satisfying the conditions outlined in Assumptions 2.1
and 2.2, and further assume it is uniformly strongly convex in a ball around
the global minimum $\mathbf{x}_{*}$, meaning $\Delta_{\mathbf{x}}F\geqslant\xi
d$ for $\xi>0$ and $\mathbf{x}\in B_{R}(\mathbf{x}_{*})$. Then there exists
$\epsilon>0$ (depending on $L,D,\xi,d$ and $R$) such that the provisional
minimum associated with the mean field $\mu$ in (4.4) satisfies

(4.5)
$\overline{\mathscr{F}}^{\mu}_{t}\geqslant\min\\{F_{*}+\epsilon,\overline{\mathscr{F}}^{\mu}_{0}\\},\qquad\forall
t>0.$

For brevity, the proof of this theorem is in Appendix A.2. This theorem
explicitly spells out the “accuracy failure.” The global optimizer can not be
achieved by the mean objective.

##  5\. Proofs of the main results

In this section we prove Theorem 2.1 on the convergence to the mean-field
limit, and Theorem 2.2 on the large-time behavior of the limiting equation.

###  5.1. Convergence to the mean-field limit

Theorem 2.1 states that the empirical distribution $\mu^{N}_{t}$ that
assembles all samples that follow the coupled SDE system (1.2), is close to
the limiting distribution $\mu$ that solves the PDE (2.4a).

There are a few different but related approaches for justifying mean-field
limits, including the coupling method, the energy/entropy estimates, the
utilization of BBGKY hierarchies and running the tightness argument to confine
a sequence of measures for the convergence; see references [22, 27]. Here we
employ the coupling method, introduced in [27]. The coupling method places a
stronger requirement on the field, but the machinery is rather simple to use,
and it has the advantage of providing precise convergence rate in terms of the
number of particles, $N$.

To apply the coupling method, one usually designs an auxiliary system that is
pushed forward by the underlying field. This auxiliary system serves as a
bridge linking the discrete SDE system and the limiting mean-field PDE. On the
one hand, it is driven by the underlying field and thus its ensemble closely
resembles $\rho$, and on the other hand, this system adopts the SDE structure
and can be easily compared with the original self-consistent SDE system. In
our context, we define such auxiliary system as
$\left\\{\widetilde{\mathbf{x}}^{j}_{t}\,,\widetilde{m}^{j}_{t}\right\\}_{j=1}^{N}$,
that are governed by:

(5.1)
$\left\\{\begin{array}[]{rcl}\textnormal{d}\widetilde{\mathbf{x}}^{j}_{t}&=&-\nabla
F\left(\widetilde{\mathbf{x}}^{j}_{t}\right)\textnormal{d}t+\sqrt{2\sigma\left(\widetilde{m}^{j}_{t}\right)}\textnormal{d}W^{j}_{t}\\\
\textnormal{d}\widetilde{m}^{j}_{t}&=&-\widetilde{m}^{j}_{t}(F\left(\widetilde{\mathbf{x}}^{j}_{t}\right)-\overline{\mathscr{F}}^{\mu}_{t})\textnormal{d}t,\end{array}\right.$

with the initial data set to be
$\\{\widetilde{\mathbf{x}}^{j}_{0}=\mathbf{x}^{j}_{0}\,,\widetilde{m}^{j}_{0}=m^{j}_{0}\\}$,
and
$\overline{\mathscr{F}}^{\mu}=\mathbb{E}_{\mu}[mF(\mathbf{x})]/\mathbb{E}_{\mu}[m]$
was defined in (2.4b). For this system, we define the corresponding empirical
distribution and the provisional minimum:

(5.2)
$\widetilde{\mu}^{N}_{t}=\frac{1}{N}\sum_{i}\delta_{\widetilde{\mathbf{x}}^{j}_{t}}(\mathbf{x})\otimes\delta_{\widetilde{m}^{j}_{t}}(m)\,,\quad\text{and}\quad\widetilde{F}^{N}_{t}=\frac{\sum_{j}\widetilde{m}^{j}_{t}F\left(\widetilde{\mathbf{x}}^{j}_{t}\right)}{\sum_{j}\widetilde{m}^{j}_{t}}.$

We note that this system has a clear similarity with the original SDE (1.2)
but is passive, in the sense that the dynamics of $\widetilde{m}^{j}_{t}$ is
driven by $\overline{\mathscr{F}}^{\mu}$ instead of $\widetilde{F}^{N}_{t}$.
While the dynamics of $\widetilde{\mathbf{x}}_{t}$ agrees with that of
$\mathbf{x}_{t}$, the mass changes call for different rates. In particular, in
(1.2), we compare $F(\mathbf{x}_{t})$ with $\overline{F}^{N}_{t}$ that is
self-generated in (2.2), and here, $F(\widetilde{\mathbf{x}}_{t})$ is compared
with $\overline{\mathscr{F}}^{\mu}$ defined from the underlying PDE system
(2.4a). As a consequence, this new
$\\{\widetilde{\mathbf{x}}^{j}_{0}=\mathbf{x}^{j}_{0}\,,\widetilde{m}^{j}_{0}=m^{j}_{0}\\}$
system gets passively pushed forward by the underlying $\mu$.

This design of the new auxiliary system prompts two lemmas.

######  Lemma 5.1.

Suppose Assumption 2.1 holds. Let
$\\{\widetilde{\mathbf{x}}^{j}_{t}\,,\widetilde{m}^{j}_{t}\\}_{j=1}^{N}$ be
the system defined in (5.1) and let $\mu_{t}$ solve (2.4a). If the initial
condition
$\\{\widetilde{\mathbf{x}}^{j}_{0}=\mathbf{x}^{j}_{0}\,,\widetilde{m}^{j}_{0}=m^{j}_{0}\\}$
is drawn from $\mu_{t=0}$, then

  * •

$\mu_{t}$ and $\widetilde{\mu}^{N}_{t}$ are close in Wasserstein sense:

$\lim_{N\to\infty}\mathbb{E}(W_{2}(\mu_{t},\widetilde{\mu}^{N}_{t}))=0.$

  * •

There exists a constant $C=C(t,d,\mu_{0},F)$ independent of $N$ such that
$\widetilde{F}^{N}_{t}$ and $\overline{\mathscr{F}}^{\mu}_{t}$ are close

$\mathbb{E}(|\widetilde{F}^{N}_{t}-\overline{\mathscr{F}}^{\mu}_{t}|)<\frac{C}{\sqrt{N}}\,,$

where $\widetilde{F}^{N}_{t}$ and $\overline{\mathscr{F}}^{\mu}_{t}$ are
defined in (5.2) and (2.4b) respectively.

and

######  Lemma 5.2.

Suppose Assumption 2.1 holds. Let
$\\{\widetilde{\mathbf{x}}^{j}_{t}\,,\widetilde{m}^{j}_{t}\\}_{j=1}^{N}$ be
the system defined in (5.1), and let
$\\{{\mathbf{x}}^{j}_{t}\,,{m}^{j}_{t}\\}_{j=1}^{N}$ solve (1.2). Suppose the
two systems have the same initial conditions, then

  * •

$\mu^{N}_{t}$ and $\widetilde{\mu}^{N}_{t}$ are close in Wasserstein sense:

$\lim_{N\to\infty}\mathbb{E}(W_{2}(\mu^{N}_{t},\widetilde{\mu}^{N}_{t}))=0$

  * •

There exists a constant $C=C(t,d,\mu_{0},F)$ independent of $N$ such that
$\widetilde{F}^{N}_{t}$ converges to $\overline{F}^{N}_{t}$ with the rate of

$\mathbb{E}(|\widetilde{F}^{N}_{t}-\overline{F}^{N}_{t}|)<\frac{C}{\sqrt{N}}\,.$

Here $\widetilde{F}^{N}_{t}$ and $\overline{F}^{N}_{t}$ are defined in (5.2)
and (2.2).

With these two lemmas, Theorem 2.1 is a straightforward corollary, calling
triangle inequality:

$\mathbb{E}\left(W_{2}(\mu_{t},\mu^{N}_{t})\right)\leqslant\mathbb{E}\left(W_{2}(\mu_{t},\widetilde{\mu}^{N}_{t})\right)+\mathbb{E}\left(W_{2}(\widetilde{\mu}^{N}_{t},\mu^{N}_{t})\right)\,.$

We omit the details here and move on to provide the proof of the two lemmas.

###### Proof for Lemma 5.1.

Calling [9],

$\mathbb{E}\left(W_{2}(\mu_{t},\widetilde{\mu}^{N}_{t})\right)=\mathcal{O}(N^{-\theta})\to
0\,,\quad\text{as}\quad N\to\infty\,,$

with $\theta=\max\\{1/2,1/d\\}$.

For weak convergence, we note that $(\widetilde{\mathbf{x}},\widetilde{m})$
are i.i.d. samples from $\mu$ and thus:

$\mathbb{E}\left(\left|\frac{1}{N}\sum_{j=1}^{N}\widetilde{m}^{j}_{t}F\left(\widetilde{\mathbf{x}}^{j}_{t}\right)-\int
mF(\mathbf{x})\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)\right|\right)<\frac{\mathrm{Var}^{1/2}_{\mu}(mF)}{\sqrt{N}}\leqslant\frac{C(t)M}{\sqrt{N}}\,.$

∎

The proof for Lemma 5.2 is the standard coupling method. This amounts to
compare the original SDE (1.2) with the auxiliary SDE (5.1). In particular, we
note that the $\mathbf{x}$ component of the two systems are almost identical:

$\mathbf{x}^{j}_{t}=-\nabla
F(\mathbf{x}^{j}_{t})\textnormal{d}t+\sqrt{2\sigma(m^{j}_{t})}\textnormal{d}W^{j}_{t}\,,\quad\text{vs}\quad\widetilde{\mathbf{x}}^{j}_{t}=-\nabla
F(\widetilde{\mathbf{x}}^{j}_{t})\textnormal{d}t+\sqrt{2\sigma(\widetilde{m}^{j}_{t})}\textnormal{d}W^{j}_{t}\,.$

Then we can subtract one from the other and trace the growth of the error. The
same strategy is applied to tracking $m^{j}_{t}-\widetilde{m}^{j}_{t}$. The
combination of the two pieces of error gives the control of the Wasserstein
distance of $W_{2}(\mu^{N}_{t},\widetilde{\mu}^{N}_{t})$.

###### Proof for Lemma 5.2.

Noting that
$\mu^{N}_{t}=\frac{1}{N}\sum_{i}\delta_{\mathbf{x}^{j}_{t}}(\mathbf{x})\otimes\delta_{m^{j}_{t}}(m)$
and
$\widetilde{\mu}^{N}_{t}=\frac{1}{N}\sum_{i}\delta_{\widetilde{\mathbf{x}}^{j}_{t}}(\mathbf{x})\otimes\delta_{\widetilde{m}^{j}_{t}}(m)$,
then recall the definition of the Wasserstein metric, we have:

(5.3)
$W_{2}(\widetilde{\mu}^{N}_{t},\mu^{N}_{t})^{2}\leqslant\left(\frac{1}{N}\sum^{N}_{j=1}\left[|\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}|^{2}+\left|m^{j}_{t}-\widetilde{m}^{j}_{t}\right|^{2}\right]\right)\,.$

Similarly, for the mean objective, we have, according to the definition:

(5.4) $\displaystyle\mathbb{E}(|\widetilde{F}^{N}_{t}-\overline{F}^{N}_{t}|)$
$\displaystyle=$
$\displaystyle\mathbb{E}\left(\left|\frac{1}{N}\sum^{N}_{j=1}\left(m^{j}_{t}F(\mathbf{x}^{j}_{t})-\widetilde{m}^{j}_{t}F\left(\widetilde{\mathbf{x}}^{j}_{t}\right)\right)\right|\right)$
$\displaystyle=$
$\displaystyle\mathbb{E}\left(\left|\frac{1}{N}\sum^{N}_{j=1}\left(m^{j}_{t}-\widetilde{m}^{j}_{t}\right)F(\mathbf{x}^{j}_{t})-\frac{1}{N}\sum^{N}_{j=1}\widetilde{m}^{j}_{t}\left(F(\mathbf{x}^{j}_{t})-F\left(\widetilde{\mathbf{x}}^{j}_{t}\right)\right)\right|\right)$
$\displaystyle\leqslant$
$\displaystyle\mathbb{E}\sqrt{\frac{C}{N}\sum^{N}_{j=1}\left[|\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}|^{2}+\left|m^{j}_{t}-\widetilde{m}^{j}_{t}\right|^{2}\right]}\,,$

where in the last inequality, we used the triangle inequality, Cauchy-Schwarz
inequality, Lipschitz continuity of $F$ and upper boundedness of both
$m$444Firstly we notice that, according to the Boundedness in $m$ condition in
Assumption 2.1, the initial data for $m^{j}_{t}$ is bounded by $M_{\max}$,
meaning $|m^{j}_{0}|<M_{\max}$ for all $j$ at the initial time. Given that
$|F|\leqslant D$, then according to the $m^{j}_{t}$ equation in (1.2), we call
Grönwall inequality to obtain that  (5.5) $|m^{j}_{t}|\leqslant
M_{\max}e^{2Dt}$ for all time. The same argument applies to bound
$\widetilde{m}^{j}_{t}$. We denote the upper bound $c$ that depends on
$M_{\max}$, $D$ and $t$.  and $F$. The constant $C$ depends on $D$, $L$ and
$M_{\max}$ (from Assumption 2.1).

Then controlling (5.3) and (5.4) gets translated to controlling the following
term:

(5.6)
$G(t)=\mathbb{E}\left(\frac{1}{N}\sum^{N}_{j=1}\left[|\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}|^{2}+\left|m^{j}_{t}-\widetilde{m}^{j}_{t}\right|^{2}\right]\right)\,.$

If we can show $G(t)\sim 1/N$, the lemma is complete. The rest of the proof is
dedicated to bounding $G(t)$. To do so, we will trace the evolution of
$|\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}|^{2}$ and
$|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}$ respectively.

  * •

To trace $|\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}|^{2}$, we compare
the $\mathbf{x}^{j}_{t}$-equation in (1.2) and the
$\widetilde{\mathbf{x}}^{j}_{t}$-equation in (5.1) to have:

$g(\mathbf{x})=|\mathbf{x}|^{2},\quad\nabla_{\mathbf{x}}g(\mathbf{x})=2\mathbf{x},\quad
D^{2}_{\mathbf{x}}g=2{\mathbb{I}}\,.$

Then
$g(\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t})=|\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}|^{2}=\langle\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t},\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}\rangle$,
and with Itô’s formula we get

$\displaystyle\textnormal{d}g(\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t})$
$\displaystyle=\left(\langle\nabla g,-(\nabla F(\mathbf{x}^{j}_{t})-\nabla
F\left(\widetilde{\mathbf{x}}^{j}_{t}\right))\rangle+\text{Tr}\left(\left(\sqrt{\sigma(m^{j}_{t})}-\sqrt{\sigma\left(\widetilde{m}^{j}_{t}\right)}\right)^{2}D^{2}_{\mathbf{x}}g\right)\right)\textnormal{d}t$
$\displaystyle\quad+\left\langle\nabla
g,\sqrt{2\sigma(m^{j}_{t})}-\sqrt{2\sigma\left(\widetilde{m}^{j}_{t}\right)}\right\rangle\textnormal{d}W^{j}_{t}$
$\displaystyle=\left(-2\langle\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t},\nabla
F(\mathbf{x}^{j}_{t})-\nabla
F\left(\widetilde{\mathbf{x}}^{j}_{t}\right)\rangle+2d\left|\sqrt{\sigma(m^{j}_{t})}-\sqrt{\sigma\left(\widetilde{m}^{j}_{t}\right)}\right|^{2}\right)\textnormal{d}t$
$\displaystyle\quad+2\left\langle\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t},\sqrt{2\sigma(m^{j}_{t})}-\sqrt{2\sigma\left(\widetilde{m}^{j}_{t}\right)}\right\rangle\textnormal{d}W^{j}_{t}\,.$

Taking the expectation of the whole formula, the mean of the Brownian motion
gets dropped out and we have:

(5.7)
$\displaystyle\frac{\textnormal{d}}{\textnormal{d}t}\mathbb{E}[g(\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t})]$
$\displaystyle=\frac{\textnormal{d}}{\textnormal{d}t}\mathbb{E}[|\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}|^{2}]\leqslant
2(L+L^{2})\left(\mathbb{E}[|\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}|^{2}]+d\mathbb{E}[|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}]\right)\,,$

where we use Assumption 2.1 that both $\nabla F$ and $\sqrt{\sigma}$ are
$L$-Lipschitz.

  * •

To trace $|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}$, we compare the $m$-equation
in (1.2) and the $m$-equation in (5.1) to have:

(5.8)
$\displaystyle\frac{\textnormal{d}}{\textnormal{d}t}|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}$
$\displaystyle=-2(m^{j}_{t}-\widetilde{m}^{j}_{t})(m^{j}_{t}F(\mathbf{x}^{j}_{t})-\widetilde{m}^{j}_{t}F\left(\widetilde{\mathbf{x}}^{j}_{t}\right))$
$\displaystyle\quad+2(m^{j}_{t}-\widetilde{m}^{j}_{t})\left(\frac{1}{N}\sum_{i=1}^{N}m^{i}_{t}F(\mathbf{x}^{i}_{t})-\mathbb{E}_{\mu}[mF(\mathbf{x})]\right)$
$\displaystyle=-2(m^{j}_{t}-\widetilde{m}^{j}_{t})^{2}F(\mathbf{x}^{j}_{t})-2(m^{j}_{t}-\widetilde{m}^{j}_{t})\widetilde{m}^{j}_{t}(F(\mathbf{x}^{j}_{t})-F\left(\widetilde{\mathbf{x}}^{j}_{t}\right))$
$\displaystyle\quad+2(m^{j}_{t}-\widetilde{m}^{j}_{t})\left(\frac{1}{N}\sum_{i=1}^{N}m^{i}_{t}F(\mathbf{x}^{i}_{t})-\frac{1}{N}\sum_{i=1}^{N}\widetilde{m}^{i}_{t}F(\widetilde{\mathbf{x}}^{i}_{t})\right)$
$\displaystyle\quad+2(m^{j}_{t}-\widetilde{m}^{j}_{t})\left(\frac{1}{N}\sum_{i=1}^{N}\widetilde{m}^{i}_{t}F(\widetilde{\mathbf{x}}^{i}_{t})-\mathbb{E}_{\mu}[mF(\mathbf{x})]\right)$

Re-order all terms, we have

(5.9)
$\displaystyle\frac{\textnormal{d}}{\textnormal{d}t}|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}$
$\displaystyle=\underbrace{-2(m^{j}_{t}-\widetilde{m}^{j}_{t})^{2}F(\mathbf{x}^{j}_{t})}_{\text{(i)}}\underbrace{-2(m^{j}_{t}-\widetilde{m}^{j}_{t})\widetilde{m}^{j}_{t}(F(\mathbf{x}^{j}_{t})-F\left(\widetilde{\mathbf{x}}^{j}_{t}\right))}_{\text{(ii)}}$
$\displaystyle\quad\underbrace{+2(m^{j}_{t}-\widetilde{m}^{j}_{t})\left(\frac{1}{N}\sum_{i=1}^{N}(m^{i}_{t}-\widetilde{m}^{i}_{t})F(\mathbf{x}^{i}_{t})\right)}_{\text{(iii)}}$
$\displaystyle\quad+\underbrace{2(m^{j}_{t}-\widetilde{m}^{j}_{t})\left(\frac{1}{N}\sum_{i=1}^{N}\widetilde{m}^{i}_{t}(F(\mathbf{x}^{i}_{t})-F(\widetilde{\mathbf{x}}^{i}_{t}))\right)}_{\text{(iv)}}$
$\displaystyle\quad\underbrace{+2(m^{j}_{t}-\widetilde{m}^{j}_{t})\left(\frac{1}{N}\sum_{i=1}^{N}\widetilde{m}^{i}_{t}F(\widetilde{\mathbf{x}}^{i}_{t})-\mathbb{E}_{\mu}[mF(\mathbf{x})]\right)}_{\text{(v)}}\,.$

We analyze each term separately,

    1. (i)

According to Assumption 2.1, $F(\mathbf{x}^{j}_{t})\geqslant 0$ so,
$\displaystyle-2(m^{j}_{t}-\widetilde{m}^{j}_{t})^{2}F(\mathbf{x}^{j}_{t})\leqslant
0$.

    2. (ii)

Using Young’s inequality,

$\begin{split}-2(m^{j}_{t}-\widetilde{m}^{j}_{t})\widetilde{m}^{j}_{t}(F(\mathbf{x}^{j}_{t})-F\left(\widetilde{\mathbf{x}}^{j}_{t}\right))&\leqslant\widetilde{m}^{j}_{t}(|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}+|F(\mathbf{x}^{j}_{t})-F\left(\widetilde{\mathbf{x}}^{j}_{t}\right)|^{2})\\\
&\leqslant
c[|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}+L^{2}|\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}|^{2}]\,,\end{split}$

where $c$ comes from the boundedness of $m$, seen in Footnote 4 and $L$ is the
Lipschitz constant from Assumption 2.1.

    3. (iii)

Similarly:

$\begin{split}2(m^{j}_{t}-\widetilde{m}^{j}_{t})\left(\frac{1}{N}\sum_{i=1}^{N}(m^{i}_{t}-\widetilde{m}^{i}_{t})F(\mathbf{x}^{i}_{t})\right)&\leqslant|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}+\left(\frac{1}{N}\sum_{i=1}^{N}(m^{i}_{t}-\widetilde{m}^{i}_{t})F(\mathbf{x}^{i}_{t})\right)^{2}\\\
&\leqslant
D^{2}\left[|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}+\left(\frac{1}{N}\sum_{i=1}^{N}(m^{i}_{t}-\widetilde{m}^{i}_{t})\right)^{2}\right]\\\
&\leqslant
D^{2}\left[|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}+\frac{1}{N}\sum_{i=1}^{N}|m^{i}_{t}-\widetilde{m}^{i}_{t}|^{2}\right]\,,\end{split}$

where we used Young’s in the first inequality, boundedness of $F$ in the
second (following Assumption 2.1), and the Hölder’s inequality in the last.

    4. (iv)

Using the same argument as in (iii) we get

$\begin{split}2(m^{j}_{t}-\widetilde{m}^{j}_{t})\Big{(}\frac{1}{N}\sum_{i=1}^{N}&\widetilde{m}^{i}_{t}(F(\mathbf{x}^{i}_{t})-F(\widetilde{\mathbf{x}}^{i}_{t}))\Big{)}\\\
\leqslant&|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}+\Big{(}\frac{1}{N}\sum_{i=1}^{N}\widetilde{m}^{i}_{t}(F(\mathbf{x}^{i}_{t})-F(\widetilde{\mathbf{x}}^{i}_{t}))\Big{)}^{2}\\\
\leqslant&c\left[|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}+\Big{(}\frac{1}{N}\sum_{i=1}^{N}(F(\mathbf{x}^{i}_{t})-F(\widetilde{\mathbf{x}}^{i}_{t}))\Big{)}^{2}\right]\\\
\leqslant&c\left[|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}+\frac{1}{N}\sum_{i=1}^{N}|F(\mathbf{x}^{i}_{t})-F(\widetilde{\mathbf{x}}^{i}_{t})|^{2}\right]\\\
\leqslant&c\left[|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}+\frac{L^{2}}{N}\sum_{i=1}^{N}|\mathbf{x}^{i}_{t}-\widetilde{\mathbf{x}}^{i}_{t}|^{2}\right]\end{split}$

where we used Young’s in the first inequality, boundedness of
$\widetilde{m}^{j}_{t}$ in the second (see Footnote 4), Hölder’s in the third,
and finally the Lipschitz continuity of $F$ in the last.

    5. (v)

Using the same argument,

$\displaystyle 2(m^{j}_{t}-\widetilde{m}^{j}_{t})$
$\displaystyle\left(\frac{1}{N}\sum_{i=1}^{N}\widetilde{m}^{i}_{t}F(\widetilde{\mathbf{x}}^{i}_{t})-\mathbb{E}_{\mu}[mF(\mathbf{x})]\right)$
$\displaystyle\leqslant$
$\displaystyle|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}+\left|\frac{1}{N}\sum_{i=1}^{N}\widetilde{m}^{i}_{t}F(\widetilde{\mathbf{x}}^{i}_{t})-\mathbb{E}_{\mu}[mF(\mathbf{x})]\right|^{2}\,.$

Collecting these terms together to feed in (5.9), and sum up against $j$, we
have:

(5.10)
$\displaystyle\frac{\textnormal{d}}{\textnormal{d}t}\left(\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\left[|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}\right]\right)\leqslant\widetilde{C}\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\left[|\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}|^{2}+|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}\right]$
$\displaystyle\hskip
144.54pt+\mathbb{E}\left[\left|\frac{1}{N}\sum_{j=1}^{N}\widetilde{m}^{j}_{t}F\left(\widetilde{\mathbf{x}}^{j}_{t}\right)-\mathbb{E}_{\mu}[mF(\mathbf{x})]\right|^{2}\right]$

where $\widetilde{C}$ depends on $D$, $M_{\max}$, $t$ and $L$.

Recall the definition of $G$ in (5.6), we collect (5.7) and (5.10) to have:

(5.11) $\displaystyle\frac{\textnormal{d}}{\textnormal{d}t}G(t)$
$\displaystyle\leqslant{C}(t)\mathbb{E}\left[\frac{1}{N}\sum_{j=1}^{N}|\mathbf{x}^{j}_{t}-\widetilde{\mathbf{x}}^{j}_{t}|^{2}+\frac{1}{N}\sum_{j=1}^{N}|m^{j}_{t}-\widetilde{m}^{j}_{t}|^{2}\right]$
$\displaystyle\qquad+\underbrace{\mathbb{E}\left[\left|\frac{1}{N}\sum_{j=1}^{N}\widetilde{m}^{j}_{t}F\left(\widetilde{\mathbf{x}}^{j}_{t}\right)-\mathbb{E}_{\mu}[mF(\mathbf{x})]\right|^{2}\right]}_{\text{(vi)}}$
$\displaystyle={C}(t)G(t)+\text{\text{(vi)}},$

where ${C}(t)$ depends on $D$, $M_{\max}$, $t$, $L$ and $d$. Noting that term
(vi) can be further controlled using the law of large numbers: let
$X=mF(\mathbf{x})$ be a random variable, then
$\widetilde{m}^{j}_{t}F\left(\widetilde{\mathbf{x}}^{j}_{t}\right)$ is a
realization of $X$ and we have

$\overline{X}_{N}=\frac{1}{N}\sum_{j=1}^{N}\widetilde{m}^{j}_{t}F(\widetilde{\mathbf{x}}^{j}_{t})$

then, since
$\\{(\widetilde{\mathbf{x}}^{j}_{t},\widetilde{m}^{j}_{t})\\}_{j=1}^{N}$ are
i.i.d. samples of $\mu$ we have

$\text{\text{(vi)}}=\mathbb{E}[|\overline{X}_{N}-\mathbb{E}[X]|^{2}]=\text{Var}[\overline{X}_{N}]=\frac{1}{N}\text{Var}[X]=\mathcal{O}\left(\frac{1}{N}\right)\,.$

Combining with (5.11), we have:

$\frac{\textnormal{d}}{\textnormal{d}t}G(t)\leqslant
C(t)G(t)+\mathcal{O}\left(\frac{1}{N}\right)\quad\text{with}\quad G(0)=0\,,$

then, calling Grönwall’s inequality, we have
$G(t)\leqslant\mathcal{O}\left(\frac{1}{N}\right),$ which completes the proof.
∎

###  5.2. Large time behavior

This section is dedicated to the proof of Theorem 2.2. This theorem quantifies
the large time behavior of the PDE (2.4a). We note that (2.7) and (2.10) are
asymptotic and rate-specific respectively. Both proofs are by contradiction.
To start off, we first present two lemmas.

######  Lemma 5.3.

Under the condition of Theorem 2.2, assume
$\liminf_{t\rightarrow\infty}\overline{\mathscr{F}}^{\mu}_{t}>F_{*}$, then:
there exists $\epsilon_{0}>0$ such that for any $t\geqslant 0$, we have
$\overline{\mathscr{F}}^{\mu}_{t}>F_{*}+\epsilon_{0}$.

At the same time,

######  Lemma 5.4.

Under the condition of Theorem 2.2, if there is an $\epsilon>0$ so that
$\overline{\mathscr{F}}^{\mu}_{t}>F_{*}+\epsilon$ for all $t\geqslant 0$, then
the total mass $\int\int m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)\to\infty$ as $t\to\infty$ and is
not conserved.

Putting these two lemmas together, we can prove the theorem.

###### Proof for (2.7) in Theorem 2.2.

We deploy the proof by contradiction. This is to assume (2.7) is not true,
meaning $\liminf_{t\rightarrow\infty}\overline{\mathscr{F}}^{\mu}_{t}>F_{*}$.
Then according to Lemma 5.3, there exists $\epsilon_{0}>0$ such that
$\overline{\mathscr{F}}^{\mu}_{t}>F_{*}+\epsilon_{0}$ for all $t\geqslant 0$.
Setting $\epsilon=\epsilon_{0}$ in Lemma 5.4, mass is not conserved,
contradicting the property of the equation. Hence (2.7) has to be true,
namely:

$\liminf_{t\to\infty}\overline{\mathscr{F}}^{\mu}_{t}=F_{*}\,.$

∎

We now spell out details for the two lemmas.

###### Proof for Lemma 5.3.

This is also completed by contradiction. Assuming the statement is not true,
meaning there is no such $\epsilon_{0}$ to make
$\overline{\mathscr{F}}^{\mu}_{t}>F_{*}+\epsilon_{0}$, then $\inf_{t\geqslant
0}\overline{\mathscr{F}}^{\mu}_{t}=F_{*}$. Considering
$\liminf_{t\rightarrow\infty}\overline{\mathscr{F}}^{\mu}_{t}>F_{*}$, there
must exist a finite time, denoted as $t_{*}\geqslant 0$ such that
$\overline{\mathscr{F}}^{\mu}_{t_{*}}=F_{*}$. According to Assumption 2.2, $F$
has a unique global minimum at $\mathbf{x}_{*}$, so at $t_{*}$,
$\mu_{t_{*}}(\mathbf{x},m)$ must have the form of

(5.12)
$\mu_{t_{*}}(\mathbf{x},m)=\delta_{\mathbf{x}_{*}}(\mathbf{x})\otimes\eta(m)\,.$

where $\eta$ is a probability measure over $m$-dimension supported on
$\mathbb{R}_{+}$.

We will show below that $\mu_{t_{*}}(\mathbf{x},m)$ in the form of (5.12) is
an equilibrium. If so, $\mu_{t}(\mathbf{x},m)=\mu_{t_{*}}(\mathbf{x},m)$ for
all $t>t_{*}$, and thus $\overline{\mathscr{F}}^{\mu}_{t}=F_{*}$ for all time
after $t_{*}$. This contradicts the assumption that
$\liminf_{t\rightarrow\infty}\overline{\mathscr{F}}^{\mu}_{t}>F_{*}$,
completing the proof.

To show $\mu_{t_{*}}(\mathbf{x},m)$ in the form of (5.12) is an equilibrium is
amount to plugging the form into the equation (2.4a) and show the right hand
side vanishes. Since this is a probability-measured solution with a Dirac
delta on $\mathbf{x}$-domain, the proof has to be conducted in the weak sense.
Specifically, given any smooth test function $\phi(\mathbf{x},m)$, we show the
right hand side of (2.4a), when tested by $\phi(\mathbf{x},m)$, vanishes:

(5.13) $\displaystyle\iint$
$\displaystyle\mathrm{R.H.S}~{}\eqref{eqn:mean_field_pde}\left.\right|_{\mu=\delta_{\mathbf{x}_{*}}(\mathbf{x})\otimes\eta(m)}\times\phi(\mathbf{x},m)\textnormal{d}{\mathbf{x}}\textnormal{d}{m}$
$\displaystyle=$
$\displaystyle\underbrace{-\iint(F(\mathbf{x})-\overline{\mathscr{F}}^{\mu}_{t})m\partial_{m}\phi(\mathbf{x},m)\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)}_{\text{Term I}}$
$\displaystyle-\iint\left(\underbrace{\nabla_{\mathbf{x}}\phi(\mathbf{x},m)\cdot\nabla
F(\mathbf{x})}_{\text{Term
II}}+\underbrace{\sigma(m)\Delta_{\mathbf{x}}\phi(\mathbf{x},m)}_{\text{Term
III}}\right)\hskip 0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)\,.$

We will show all three terms are zero. Indeed,

  * •

Since $\mu_{t}$ takes on the form of (5.12),

(5.14) $\displaystyle\mathrm{Term\ I}=$
$\displaystyle\iint(F(\mathbf{x})-\overline{\mathscr{F}}^{\mu}_{t}m\partial_{m}\phi(\mathbf{x},m)\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)$ $\displaystyle=$
$\displaystyle\iint(F_{*}-F_{*})m\partial_{m}\phi(\mathbf{x},m)\textnormal{d}\eta(m)=0\,.$

  * •

Similarly,

(5.15) $\displaystyle\mathrm{Term\ II}=$
$\displaystyle\iint\nabla_{\mathbf{x}}\phi(\mathbf{x},m)\cdot\nabla
F(\mathbf{x})\delta_{x_{*}}(\mathbf{x})\textnormal{d}{\mathbf{x}}\textnormal{d}\eta(m)$
$\displaystyle=$ $\displaystyle\nabla
F(\mathbf{x}_{*})\cdot\int\nabla_{\mathbf{x}}\phi(\mathbf{x}_{*},m)\textnormal{d}\eta(m)=0\,.$

  * •

To show $\mathrm{Term\ III}=0$, we note that

(5.16) $\displaystyle\mathrm{Term\ III}=$
$\displaystyle\iint\sigma(m)\Delta_{\mathbf{x}}\phi(\mathbf{x},m)\delta_{x_{*}}(\mathbf{x})\textnormal{d}{\mathbf{x}}\textnormal{d}\eta(m)$
$\displaystyle=$
$\displaystyle\int\sigma(m)\Delta_{x}\phi(\mathbf{x}_{*},m)\textnormal{d}\eta(m)=0\,,$

while in the last equation, we used Lemma B.1 and that
$\mathrm{Supp}(\eta(m))\subset[m_{c},\infty)$.

Plugging (5.14)-(5.16) into (5.13), it is clear that
$\mu_{t_{*}}(\mathbf{x},m)=\delta_{\mathbf{x}_{*}}(\mathbf{x})\otimes\eta(m)$
is an equilibrium, and thus $\mu_{t}(\mathbf{x},m)=\mu_{t_{*}}(\mathbf{x},m)$
for all $t>t_{*}$, making $\overline{\mathscr{F}}^{\mu}_{t}=F_{*}$ and
contradicting the assumption that
$\liminf_{t\rightarrow\infty}\overline{\mathscr{F}}^{\mu}_{t}>F_{*}$. ∎

###### Proof for Lemma 5.4.

According to Assumption 2.2, the initial distribution satisfies
$(m_{c},\infty)\cap\mathrm{supp}_{m}(\mu_{0}(\mathbf{x},m))\neq\emptyset$, so
using the diffusivity property in $\mathbf{x}$-direction of the Fokker-Planck
equation, for ${\mathcal{N}}_{\epsilon/2}$ as defined in (2.8), there must
exist $T_{c}>0$ such that

$\int_{{\mathcal{N}}_{\epsilon/2}}\int^{\infty}_{m_{c}}\hskip
0.56917pt\textnormal{d}\mu_{T_{c}}(\mathbf{x},m)>0\,,$

which means there is no trivial probability at $T_{c}$ to find a sample in the
set of ${\mathcal{N}}_{\epsilon/2}\times(m_{c},\infty)$. Denote one such
sample $(\mathbf{x}_{s}(T_{c}),m_{s}(T_{c}))$, we are to show that

  1. (a)

$(\mathbf{x}_{s},m_{s})(t)\in{\mathcal{N}}_{\epsilon/2}\times(m_{c},\infty)$
for all $t>T_{c}$;

  2. (b)

$m_{s}(t)>\exp(\frac{\epsilon}{2}(t-T_{c}))m_{s}(T_{c})>\exp(\frac{\epsilon}{2}(t-T_{c}))m_{c}$
for any $t>T_{c}$.

Indeed, we recall its trajectory governed by the SDE:

$\left\\{\begin{array}[]{rcl}\textnormal{d}\mathbf{x}_{s}(t)&=&-\nabla
F(\mathbf{x}_{s}(t))\textnormal{d}t+\sqrt{2\sigma(m_{s}(t))}\textnormal{d}W_{t}\\\
\frac{\textnormal{d}}{\textnormal{d}t}m_{s}(t)&=&-m_{s}(t)\left(F(\mathbf{x}_{s}(t))-\overline{\mathscr{F}}^{\mu}_{t}\right)\,.\end{array}\right.$

Since $\mathbf{x}_{s}(T_{c})\in{\mathcal{N}}_{\epsilon/2}$ and
$m_{s}(T_{c})>m_{c}$, then noting $\sigma(m_{s}(T_{c}))=0$, we obtain

$\left.\frac{\textnormal{d}\mathbf{x}_{s}(t)}{\textnormal{d}t}\right|_{t=T_{c}}=-\nabla
F(\mathbf{x}_{s}(t)),\quad\left.\frac{\textnormal{d}m_{s}(t)}{\textnormal{d}t}\right|_{t=T_{c}}>\frac{\epsilon}{2}m_{s}(t)$

where we used the fact that $\overline{\mathscr{F}}^{\mu}_{t}>F_{*}+\epsilon$
while $F(\mathbf{x}_{s}(T_{c}))<F_{*}+\epsilon/2$ to obtain the equation for
$m_{s}$. This dynamics suggests that $F(\mathbf{x}_{s}(t))$ keeps decreasing
so $\mathbf{x}_{s}\in{\mathcal{N}}_{\epsilon/2}$, and $m_{s}(t)$ keeps
increasing with an exponential rate:

$m_{s}(t)>\exp\left(\frac{\epsilon}{2}(t-T_{c})\right)m_{s}(T_{c})>m_{c}\,,$

confirming (a)-(b). Accordingly,

(5.17) $\displaystyle\int\limits_{{\mathcal{N}}_{\epsilon/2}}$
$\displaystyle\int\limits^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)=\mathbb{E}(m(t)|\mathbf{x}(t)\in{\mathcal{N}}_{\epsilon/2},m(t)>m_{c})\mathbb{P}\left(\left\\{\mathbf{x}(t)\in{\mathcal{N}}_{\epsilon/2},m(t)>m_{c}\right\\}\right)$
$\displaystyle\ \geqslant$
$\displaystyle\mathbb{E}(m(t)|\mathbf{x}(T)\in{\mathcal{N}}_{\epsilon/2},m(T)>m_{c})\mathbb{P}\left(\left\\{\mathbf{x}(T)\in{\mathcal{N}}_{\epsilon/2},m(T)>m_{c}\right\\}\right)$
$\displaystyle\ \geqslant$
$\displaystyle\mathbb{E}(\exp(\frac{\epsilon}{2}(t-T))m(T)|\mathbf{x}(T)\in{\mathcal{N}}_{\epsilon/2},m(T)>m_{c})\mathbb{P}\left(\left\\{\mathbf{x}(T)\in{\mathcal{N}}_{\epsilon/2},m(T)>m_{c}\right\\}\right)$
$\displaystyle\ =$
$\displaystyle\exp(\frac{\epsilon}{2}(t-T))\int_{{\mathcal{N}}_{\epsilon/2}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{T}(\mathbf{x},m)\,,$

where we use (a) in the first inequality and (b) in the second inequality.
This leads to the fact that

$\limsup_{t\rightarrow\infty}\int_{\mathbb{R}^{d}}\int^{\infty}_{0}m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)\geqslant\limsup_{t\rightarrow\infty}\int_{{\mathcal{N}}_{\epsilon}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)=\infty\,,$

completing the proof. ∎

###### Proof for (2.10).

We again deploy a proof by contradiction and the argument is rather similar to
that of (2.7). Assuming that it is not true, then for all $t\leqslant
t_{\epsilon}(m_{c})$,

$\overline{\mathscr{F}}^{\mu}_{t}\geqslant F_{*}+\epsilon\,.$

We will show this leads to a contradiction. Define $\eta_{\epsilon}$ to be the
time so that

$\eta_{\epsilon}=\arg\sup_{t\in[0,1]}\int_{{\mathcal{N}}_{\epsilon}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)\,\quad\Rightarrow\quad\int_{{\mathcal{N}}_{\epsilon}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{\eta_{\epsilon}}(\mathbf{x},m)>0\,.$

Since the probability is non-trivial, one can sample a point in the set of
${\mathcal{N}}_{\epsilon}\times(m_{c},\infty)$. Starting with
$\mathbf{x}(\eta_{\epsilon})\in{\mathcal{N}}_{\epsilon}$ and
$m(\eta_{\epsilon})>m_{c}$, with the same calculation as in (5.17), we get:

$\mathbf{x}(t)\in{\mathcal{N}}_{\epsilon},\quad
m(t)>e^{\epsilon(t-\eta_{\epsilon})/2}m(\eta_{\epsilon})\geqslant
e^{\epsilon(t-1)/2}m(\eta_{\epsilon})$

for any $t>1\geqslant\eta_{\epsilon}$, which makes:

$\int_{{\mathcal{N}}_{\epsilon}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)>e^{\epsilon(t-1)/2}\int_{{\mathcal{N}}_{\epsilon}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{\eta_{\epsilon}}(\mathbf{x},m)\,.$

In particular, we choose $t=t_{\epsilon}(m_{c})$, then according to the
definition (2.9):

$\int_{{\mathcal{N}}_{\epsilon}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{t_{\epsilon}(m_{c})}(\mathbf{x},m)>e^{\epsilon(t_{\epsilon}(m_{c})-1)/2}\int_{{\mathcal{N}}_{\epsilon}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{\eta_{\epsilon}}(\mathbf{x},m)\geqslant 1\,.$

This contradicts the mass conservation:

$\int_{{\mathcal{N}}_{\epsilon}}\int^{\infty}_{m_{c}}m\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)\leqslant 1\,,\quad\forall
t\geqslant 1\,.$

Hence there is a $t\leqslant t_{\epsilon}(m_{c})$ so that
$\overline{\mathscr{F}}^{\mu}_{t}<F_{*}+\epsilon$. ∎

##  6\. Numerical experiments

This section is dedicated to numerical results. We will showcase numerical
evidence that validate theoretical statements on deriving the mean-field limit
and the long time behavior of $\overline{F}^{N}$. Furthermore, we will also
present some numerical findings that are yet to be understood theoretically.

All numerical experiments are ran on two classical non-convex functions: the
Ackley type functions (in 1D, 2D and 10D) and the Rastrigin type functions in
1D, defined as follows

$\displaystyle F_{A}(\mathbf{x})$
$\displaystyle=-A\exp\left(-C\sqrt{\frac{1}{d}\sum_{i=1}^{d}x_{i}^{2}}\right)-\exp\left(\frac{1}{d}\sum_{i=1}^{d}\cos(2\pi
x_{i})\right)+A+e+D,$ $\displaystyle F_{R}(\mathbf{x})$
$\displaystyle=10d+\sum_{i=1}^{d}[x_{i}^{2}-10\cos(2\pi x_{i})]+D\,.$

Throughout the numerical section, we chose the following parameters $A=20$,
$C=0.2$ and $D=0$. The profiles of the two functions in 1D are presented in
Figure 6.1. They are both highly oscillatory, presenting abundant local
minima. Moreover, the local minima for the Rastrigin function give very close
function values to the global minimum, and can easily confuse samples on
convergence.

(a) Ackley function.

(b) Rastrigin function.

Figure 6.1. Non-convex functions that will be used for the numerical examples
when $d=1$.

For the numerical evidence we will use the following $\sigma$ function in the
application of Algorithm 1:

(6.1)
$\sigma_{\lambda,\beta}(m)=\begin{cases}\lambda\exp\left(\frac{m}{m-\beta}\right)&m<\beta\\\
0&m\geqslant\beta\end{cases}\,.$

The smooth decay of the function (see Figure 2(a)) gives different magnitudes
of noise to each particle depending on their mass, but eventually for
sufficiently high mass the noise will disappear.

For the numerical study in Section 6.2 we will use the following $\sigma$
function so we have a more strict change in the application of the noise after
certain amount of mass:

(6.2)
$\sigma_{\lambda,\beta}(m)=\lambda\left(-\frac{1}{2}\tanh(1000(m-\beta))+\frac{1}{2}\right)\,.$

The hyperbolic tangent form of the function is chosen to closely resemble a
discontinuous cut-off function for samples whose weight is below $\beta$ (see
Figure 2(b)): this means particles with smaller weights are encouraged to
explore the landscape of the objective with the same rate $\lambda$ while all
other heavier particles are given zero stochasticity. Experimentally we
observe a very neat relation between $\beta$ and $N$ for a consistent
convergence rate, and we leave that to Section 6.2.

(a) Smooth $\sigma$.

(b) Step-like $\sigma$.

Figure 6.2. $\sigma$ function for different values of the parameter $\beta$
and $\lambda=1$.

To demonstrate the performance of the algorithm, we choose the experiments
where the initial swarm does not cover the global minimum basin, to rule out
the easy-to-converge situations. The convergence of the algorithm purely
depends on the balancing between Brownian motion and local gradient descent
through mass-transfer.

###  6.1. Validation of theoretical results

The first part of our simulation is to showcase the theoretical results we had
for Algorithm 1. For this, we will apply the algorithm on the Ackley function
for $d=1,\,2$, then Rastrigin for $d=1$ and, finally to Ackley function with
$d=10$. In all cases we will show the performance of the algorithm and check
its dependence on $T$ and $N$.

####  6.1.1. Ackley function in lower dimensions

We experiment with the Ackley function for $d=1$ and $d=2$. First, for $d=1$,
we pick $N=8$ particles, $\beta=1/8$, $\lambda=1$ and run the algorithm for
step size $h=10^{-4}$ and $20000$ iterations for a total time of $T=2$. Figure
6.3 demonstrates a few time instances of the particles: It is clear that they
gradually move towards the global minimum. We then tested the convergence with
$N=8,\,16,\,32,\,64$. Each run of the algorithm provides one trajectory of
$\overline{F}^{N}_{t}$, and in simulation we ran the algorithm $20$ times and
compute its expectation, first and third quartiles. The solid lines in Figure
6.4 demonstrate the convergence of $\mathbb{E}[\overline{F}^{N}_{t}]$ in time,
and the shaded areas shows the difference between the first and third
quartiles.

Figure 6.3. Swarm movement (in red) for Ackley function in $d=1$. Figure 6.4.
Expectation of $\overline{F}^{N}_{t}$ for different values of $N$ for Ackley
function in $d=1$ in log-scale. The shaded area shows the difference between
the lower quartile and the higher quartile for $\overline{F}^{N}_{t}$.

For $d=2$, we pick $N=16$ particles, $\beta=1/16$, $\lambda=1$ and we run the
algorithm for step size $h=10^{-4}$ and $20000$ iterations for a total time of
$T=2$. In Figure 6.5 we visualize the convergence of particles to the global
minimum in a few time instances. We repeat the experiment for $d=1$ case and
compute the mean, first and third quartiles of $\overline{F}^{N}$ for
different choices of $N$, and they are plotted in Figure 6.6.

Figure 6.5. Swarm movement (in black) for Ackley function in $d=2$. Figure
6.6. Expectation of $\overline{F}^{N}_{t}$ for different values of $N$ for
Ackley function in $d=2$ in log-scale. The shaded area shows the difference
between the lower quartile and the higher quartile for $\overline{F}^{N}_{t}$.

####  6.1.2. Rastrigin function in low dimension

The third test is conducted for the Rastrigin function for $d=1$. The choice
of parameters are $N=8$, $\beta=1/8$ and $\lambda=2$. The step-size is
$h=10^{-4}$ and we run it up to $T=500$.

In Figure 6.7 we can see that the particles eventually find the global minimum
after a long time. The convergence of $\overline{F}^{N}_{t}$, however, is
significantly faster than finding the global minimum. We also evaluate the
mean-field convergence. Running the experiment for $20$ runs, we obtain 20
trajectories of $\overline{F}^{N}_{t}$. We compute the mean, the first and
third quartiles, and plot them in Figure 6.8. It is clear that bigger $N$
gives a higher chance of achieving a lower value for $F$ in time.

Figure 6.7. Swarm movement (in red) for Rastrigin function in $d=1$. Figure
6.8. Expectation of $\overline{F}^{N}_{t}$ for different values of $N$ for
Rastrigin function in $d=1$ in log-scale. The shaded area shows the difference
between the lower quartile and the higher quartile for $\overline{F}^{N}_{t}$.

####  6.1.3. Ackley function in higher dimensions

The fourth experiment is to evaluate the performance of the algorithm on a
high dimensional function (Ackley function for $d=10$). The hope is the
Brownian motion in the dynamics can help overcoming the curse of
dimensionality. We set the parameters to be $\lambda=1$, $\beta=1/N$. In
Figure 6.9, we show the convergence in time of $20$ runs of experiments using
Algorithm 1, and compare it to the solution to the deterministic system shown
in (4.1). The observation is that the deterministic system converges very
quickly at the initial steps, but it saturates at local minima, and the global
optimization is never achieved. On the contrary, Algorithm 1, using either
$2000$ or $200$ samples both achieve the global minimum, albeit with slower
convergence rate. We attribute this to the randomness in the method and this
observation confirms Theorem 4.1. We also emphasize that the performance of
the algorithm using smaller number of samples ($200$) is compatible to that
using $2000$, indicating the overcome of dimensionality.

Figure 6.9. Long time behavior of $\overline{F}^{N}_{t}$ for different values
of $N$ for Ackley function in $d=10$, in log-scale, in the absence and
presence of noise.

###  6.2. Numerical study

Though not the focus of the paper, we observe, through running numerical
simulations, some interesting algorithmic performance. We collect these
findings in this subsection. We first present studies on parameter-tuning. As
shown earlier, the performance of the algorithm depends on the choice of the
$\sigma$ function that relies on the configuration of the following three
parameters: $\lambda$, $\beta$ and $N$. While $\lambda$ determines the
strength of the Brownian motion, the cut-off point $\beta$ and its relation to
$N$ fully characterize the percentage of particles that are allowed to “roam.”
To fully understand the relation, we run Algorithm 1 on the Ackley functions
for $d=1$ up to $T=1$ with $\sigma$ defined as in (6.2) for $\lambda=1$ using
many choices of $(N,\beta)$. To be more specific, we choose
$N=8,\,16,\,32,\,64$ and $\beta=1/8,\,1/16,\,1/32,\,1/64$ and for each
$(N,\beta)$ configuration, we run the experiment $40$ times for a collection
of $40$ trajectories of $\overline{F}^{N}_{t}$. The mean and variance of these
$\overline{F}^{N}_{t}$ for different $(N,\beta)$ choices are presented, on
log-scale, in Figure 6.10. We note that experimental-wise, it seems that
$N=2/\beta$ seems to give the most ideal convergence. Another observation is
that the performance of the algorithm is consistent across the cases when
$\beta N$ are set to be constants. Namely, the convergence rate seem to agree
along diagonals in Figure 6.10.

Figure 6.10. Convergence rates for the Ackley functions for different values
of $N$ and $\beta$.

We then present the long time behavior of the algorithm performed on a
variation of the Rastrigin function from Section 6.1.2:

(6.3) $F_{R}(\mathbf{x})=10d+\sum_{i=1}^{d}[2x_{i}^{2}-10\cos(2\pi
x_{i})]+D\,,$

with parameters chosen to be $d=1$ and $D=0$, and the same choices of
$\\{x_{i}\\}$. In this case, we run the simulation for a very long time till
$t\sim 700$. As seen in Figure 6.11, the error seems to saturate at around
$O(1)$ after $t=10$, and the second dip occurs at $t=690$, dragging the error
down to $0$. The convergence takes place in a staggered pattern. The
conjecture is that the value of the local minimum of the Rastrigin function is
very close to that of the global one, and the samples stuck at the local
minimum for a very long stretch of time, accumulating all the mass. The
lighter “explorer" uses this time to explore and find the global minimum.
Before the very small basin gets found, local minimum takes the lead and the
error stays at $O(1)$, and the second dip occurs when the global basin finally
was found by one explorer, and the mass gets transferred to it. We find this
staggered pattern enlightening, and future research is needed for spelling out
the explicit dependence of this behavior on different parameters.

(a) $t\in[0,700]$.

(b) Zoom-in at $t\in[0,20]$.

(c) Zoom-in at $t\in[680,720]$.

Figure 6.11. Long time behavior of $\overline{F}^{N}_{t}$ for (6.3).

## References

  * [1] Zikri Bayraktar, Muge Komurcu, Jeremy A Bossard, and Douglas H Werner.  The wind driven optimization technique and its application in electromagnetics.  IEEE transactions on antennas and propagation, 61(5):2745–2757, 2013.
  * [2] José A Carrillo, Young-Pil Choi, Claudia Totzeck, and Oliver Tse.  An analytical framework for consensus-based global optimization method.  Mathematical Models and Methods in Applied Sciences, 28(06):1037–1066, 2018.
  * [3] José A. Carrillo, Shi Jin, Lei Li, and Yuhua Zhu.  A consensus-based global optimization method for high dimensional machine learning problems.  ESAIM: Control, Optimisation and Calculus of Variations, 27, 2021\.
  * [4] José A. Carrillo, Claudia Totzeck, and Urbain Vaes.  Consensus-based optimization and ensemble Kalman inversion for global optimization problems with constraints.  Modeling and Simulation for Collective Dynamics, pages 195–230, 2023.
  * [5] Vladimír Černỳ.  Thermodynamical approach to the traveling salesman problem: An efficient simulation algorithm.  Journal of optimization theory and applications, 45:41–51, 1985\.
  * [6] Yi Chen, Jinglin Chen, Jing Dong, Jian Peng, and Zhaoran Wang.  Accelerating nonconvex learning via replica exchange Langevin diffusion.  ICLR 2019, 2019.
  * [7] Tzuu-Shuh Chiang, Chii-Ruey Hwang, and Shuenn Jyi Sheu.  Diffusion for global optimization in $\mathbb{R}^{n}$.  SIAM Journal on Control and Optimization, 25(3):737–753, 1987.
  * [8] Massimo Fornasier, Timo Klock, and Konstantin Riedl.  Consensus-based optimization methods converge globally.  ArXiv:2103.15130v4, 2021.
  * [9] Nicolas Fournier and Arnaud Guillin.  On the rate of convergence in Wasserstein distance of the empirical measure.  Probability Theory and Related Fields, 162(3):707–738, 2015.
  * [10] Stuart Geman and Donald Geman.  Stochastic relaxation, Gibbs distributions, and the Bayesian restoration of images.  IEEE Transactions on pattern analysis and machine intelligence, (6):721–741, 1984.
  * [11] Stuart Geman and Chii-Ruey Hwang.  Diffusions for global optimization.  SIAM Journal on Control and Optimization, 24(5):1031–1043, 1986\.
  * [12] Basilis Gidas.  Global optimization via the Langevin equation.  In 1985 24th IEEE Conference on Decision and Control, pages 774–778. IEEE, 1985.
  * [13] Basilis Gidas.  Nonstationary Markov chains and convergence of the annealing algorithm.  Journal of Statistical Physics, 39:73–131, 1985.
  * [14] Sara Grassi, Hui Huang, Lorenzo Pareschi, and Jinnao Qiu.  Mean-field particle swarm optimization.  Modeling and Simulation for Collective Dynamics (W. Bao, P. Markovich, B. Perthame and E. Tadmor, eds.), 27:127–193, 2023.
  * [15] Johnson AW Henderson D, Jacobson SH.  The theory and practice of simulated annealing, volume 57.  Springer Science & Business Media, 2003.
  * [16] Dervis Karaboga, Beyza Gorkemli, Celal Ozturk, and Nurhan Karaboga.  A comprehensive survey: artificial bee colony (abc) algorithm and applications.  Artificial intelligence review, 42:21–57, 2014.
  * [17] James Kennedy and Russell Eberhart.  Particle swarm optimization.  In Proceedings of ICNN’95-international conference on neural networks, volume 4, pages 1942–1948. IEEE, 1995.
  * [18] S. Kirkpatrick, C. D. Gelatt, and M. P. Vecchi.  Optimization by simulated annealing.  Science, 220(4598):671–680, 1983.
  * [19] Peter J. M. Laarhoven and Emile H. L. Aarts.  Simulated Annealing: Theory and Applications, chapter Simulated Annealing, pages 7–15.  Springer Dordrecht, Oxford, 1987.
  * [20] Don S Lemons and Anthony Gythiel.  Paul Langevin’s 1908 paper" on the theory of brownian motion"[" sur la thiorie du mouvement brownien," cr acad. sci.(paris) 146, 530-533 (1908)].  American Journal of Physics, 65:1079–1081, 1997.
  * [21] Jingcheng Lu, Eitan Tadmor, and Anil Zenginoglu.  Swarm-based gradient descent method for non-convex optimization.  ArXiv:2211.17157, 2022.
  * [22] Sylvie Méléard.  Asymptotic behaviour of some interacting particle systems; McKean-Vlasov and Boltzmann models, pages 42–95.  Springer Berlin Heidelberg, Berlin, Heidelberg, 1996.
  * [23] B Chandra Mohan and R Baskaran.  A survey: Ant colony optimization based recent research and implementation on several engineering domain.  Expert Systems with Applications, 39(4):4618–4627, 2012.
  * [24] Pierre Monmarché.  Hypocoercivity in metastable settings and kinetic simulated annealing.  Probability Theory and Related Fields, 172(3-4):1215–1248, 2018\.
  * [25] René Pinnau, Claudia Totzeck, Oliver Tse, and Stephan Martin.  A consensus-based model for global optimization and its mean-field limit.  Mathematical Models and Methods in Applied Sciences, 27(01):183–204, 2017.
  * [26] Riccardo Poli, James Kennedy, and Tim Blackwell.  Particle swarm optimization: An overview.  Swarm intelligence, 1:33–57, 2007.
  * [27] Alain-Sol Sznitman.  Topics in propagation of chaos.  In Paul-Louis Hennequin, editor, Ecole d’Eté de Probabilités de Saint-Flour XIX — 1989, pages 165–251, Berlin, Heidelberg, 1991. Springer Berlin Heidelberg.
  * [28] Eitan Tadmor and Anil Zenginoglu.  Swarm-based optimization with random descent.  Acta Applicandae Mathematicae, 190, 2024.
  * [29] Claudia Totzeck.  Trends in consensus-based optimization.  In Active Particles, Volume 3: Advances in Theory, Models, and Applications (Bellomo, N., Carrillo, JA and Tadmor, E., eds), pages 201–226. Springer, 2021.
  * [30] Xin-She Yang.  Firefly algorithms for multimodal optimization.  In International symposium on stochastic algorithms, pages 169–178. Springer, 2009.

##  Appendix A Interacting vs. non-interacting agents

In this appendix, we investigate the effect of communication, and study the
behavior of two systems that produce pitfalls. In subsection A.1 we study
basic properties of the deterministic system (4.1) which governs non-
interacting agents. In subsection A.2, we switch gears to interacting agents,
governed by the stochastic system (4.3) with mass communication passively
adjusted according to $\mathbf{x}$-values, and prove Theorem 4.2.

###  A.1. Properties of the deterministic system (non interacting agents).

######  Proposition A.1.

Given the $N-$particle system defined in (4.1), we have the provisional
minimum $\overline{F}^{N}(t)$ always decay in time555Note that the time
dependent quantities in the deterministic case are denotes
$\square(t,\cdots)$, e.g.,
$\mathbf{x}^{j}(t),\overline{F}^{N}(t),\nu^{N}(t,\mathbf{x})$ etc.:

(A.1) $\frac{\textnormal{d}}{\textnormal{d}t}\overline{F}^{N}(t)\leqslant
0\,,\quad\forall t\geqslant
0\,\qquad\overline{F}^{N}(t)=\frac{1}{M(0)}\sum_{j=1}^{N}m^{j}(t)F(\mathbf{x}^{j}(t)).$

###### Proof.

This comes from direct computation:

$\displaystyle\frac{\textnormal{d}}{\textnormal{d}t}\overline{F}^{N}(t)$
$\displaystyle=\frac{1}{M(0)}\left(\sum^{N}_{j=1}\left(\frac{\textnormal{d}}{\textnormal{d}t}m^{j}(t)\right)F(\mathbf{x}^{j}(t))+m^{j}(t)\frac{\textnormal{d}}{\textnormal{d}t}F(\mathbf{x}^{j}(t))\right)$
$\displaystyle=\frac{1}{M(0)}\left(-\sum^{N}_{j=1}m^{j}_{t}\left(F(\mathbf{x}^{j}(t))-\overline{F}^{N}(t)\right)F(\mathbf{x}^{j}(t))-m^{j}(t)\left|\nabla
F(\mathbf{x}^{j}(t))\right|^{2}\right)$
$\displaystyle=-\frac{1}{M(0)}\sum^{N}_{j=1}m^{j}_{t}\left(F(\mathbf{x}^{j}(t))-\overline{F}^{N}(t)\right)^{2}-\frac{1}{M(0)}\sum^{N}_{j=1}m^{j}_{t}\left|\nabla
F(\mathbf{x}^{j}(t))\right|^{2}\leqslant 0\,.$

This concludes the proof. ∎

Similar to its stochastic counterpart, we can also derive the mean-field
limiting equation in the deterministic case. This is expressed in terms of the
limiting distribution $\nu(t,\mathbf{x})$. To be specific, we have the
following proposition:

######  Proposition A.2.

Let $\\{\mathbf{x}^{j}(t),m^{j}(t)\\}_{j=1}^{N}$ denote the crowd of agents
satisfying the deterministic system (4.1) subject to initial data
$\\{\mathbf{x}^{j}(0),m^{j}(0)=\frac{1}{N}\\}^{N}_{j=1}$, independently drawn
from random distribution $\rho_{0}$, and let $\nu^{N}(t,\mathbf{x})$ denote
its empirical distribution

(A.2)
$\nu^{N}(t,\mathbf{x})=\sum_{j=1}^{N}m^{j}(t)\delta_{\mathbf{x}^{j}(t)}(\mathbf{x}).$

Then its mean-field limit of the ensemble distribution,
$\nu^{N}(t,\mathbf{x})\rightarrow\nu(t,\mathbf{x})$, satisfies

(A.3)
$\partial_{t}\nu(t,\mathbf{x})=\nabla_{\mathbf{x}}\cdot(\nu(t,\mathbf{x})\nabla
F(\mathbf{x}))-(F(\mathbf{x})-\overline{\mathscr{F}}^{\nu}(t))\nu(t,\mathbf{x}),$

subject to $\nu(0,\cdot)=\rho_{0}(\cdot)$. Here
$\overline{\mathscr{F}}^{\nu}(t)$ is the weighted average is:

(A.4) $\overline{\mathscr{F}}^{\nu}(t)=\frac{\displaystyle\int
F(\mathbf{x})\textnormal{d}\nu(t,\mathbf{x})}{\displaystyle\int\textnormal{d}\nu(t,\mathbf{x})}\,.$

###### Proof.

To derive the mean-field limit, we test the system with a test smooth function
$\phi$:

$\mathbb{E}_{\nu(t)}[\phi]=\int\phi(\mathbf{x})\textnormal{d}\nu(t,\mathbf{x})=\sum_{j=1}^{N}m^{j}(t)\phi(\mathbf{x}^{j}(t))\,.$

We take the time derivative on the two sides of the equation. The left hand
side provides:

(A.5)
$\frac{\textnormal{d}}{\textnormal{d}t}\mathbb{E}_{\nu(t)}[\phi]=\int\partial_{t}\nu(t,\mathbf{x})\phi(\mathbf{x})\mathrm{d}\mathbf{x}.$

Here expectation is with respect to the randomly drawn initial data. The time
derivative of the right hand side yields

(A.6)
$\begin{split}\frac{\textnormal{d}}{\textnormal{d}t}\left(\sum_{j=1}^{N}\right.&\left.m^{j}(t)\phi(\mathbf{x}^{j}(t))\right)=\\\
&-\sum_{j=1}^{N}m^{j}_{t}\underbrace{\nabla_{\mathbf{x}}\phi(\mathbf{x}^{j}(t))\cdot\nabla
F(\mathbf{x}^{j}(t))}_{\dot{\mathbf{x}}=-\nabla
F(\mathbf{x}(t))}-\sum_{j=1}^{N}\phi(\mathbf{x}^{j}(t))\underbrace{m^{j}(t)(F(\mathbf{x}^{j}(t))-\overline{F}(t))}_{\dot{m}=-m(t)(F(\mathbf{x}(t)-\overline{\mathscr{F}}^{\nu}(t))},\end{split}$

where we integrated the ODE (4.1). Note that the first term in (A.6) can be
written as

$-\sum_{j=1}^{N}m^{j}(t)\nabla_{\mathbf{x}}\phi(\mathbf{x}^{j}(t))\cdot\nabla
F(\mathbf{x}^{j}(t))=-\mathbb{E}_{\nu(t)}[\nabla_{\mathbf{x}}\phi\cdot\nabla
F]=\int\phi(\mathbf{x})\nabla_{\mathbf{x}}\cdot(\nu(t,\mathbf{x})\nabla
F(\mathbf{x}))\textnormal{d}\mathbf{x}$

using integration-by-parts, and the second term of (A.6) can be written as

$\displaystyle\sum_{j=1}^{N}\phi(\mathbf{x}^{j}(t))m^{j}(t)(F(\mathbf{x}^{j}(t))-\overline{\mathscr{F}}^{\nu}(t))$
$\displaystyle=\mathbb{E}_{\nu(t)}[\phi
F]-\mathbb{E}_{\nu(t)}[\phi]\mathbb{E}_{\nu(t)}[F]$
$\displaystyle=\int\nu(t,\mathbf{x})\phi(\mathbf{x})F(\mathbf{x})\textnormal{d}\mathbf{x}-\overline{\mathscr{F}}^{\nu}(t)\int\phi(\mathbf{x})\nu(\mathbf{x})\textnormal{d}\mathbf{x}.$

Plugging the latter into (A.6) and combining it with (A.5) we recover (A.3). ∎

Some properties of this continuous limit are straightforward. In particular,
obvious observations made earlier in our discussion yield the following.

######  Proposition A.3.

Let $\nu$ be the mean-field limit satisfying (A.3). There holds:
$\bullet$ Mass conservation: the total mass $\displaystyle
M(t):=\int\nu(t,\mathbf{x})\textnormal{d}\mathbf{x}$, is conserved in time,
$\frac{\textnormal{d}}{\textnormal{d}t}M(t)=0$;
$\bullet$ Decrease of provisional mean: the weighted mean
$\overline{\mathscr{F}}^{\nu}(t)$ decreases in time,
$\frac{\textnormal{d}}{\textnormal{d}t}\overline{\mathscr{F}}^{\nu}(t)\leqslant
0$.

###  A.2. Proprieties of the stochastic system (interacting agents)

We turn to the proof of Theorem 4.2.

###### Proof.

We may assume without loss of generality that $F_{*}=0$, and that the total
mass is normalized

$\iint m\hskip 0.56917pt\textnormal{d}\mu_{0}(\mathbf{x},m)=1,$

so that $\displaystyle\rho_{t}(\cdot)=\int m\hskip
0.56917pt\textnormal{d}\mu_{t}(\cdot,m)$ can be viewed as a probability
density function in $\mathbf{x}$. Multiplying (4.4) by $mF(\mathbf{x})$ and
integrate, we obtain that, with a few integration by parts:

$\frac{\textnormal{d}}{\textnormal{d}t}\overline{\mathscr{F}}^{\mu}_{t}=-\iint
m|\nabla F|^{2}\hskip 0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)-\iint
m(F-\overline{\mathscr{F}}^{\mu}_{t})^{2}\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)+\iint m\Delta F\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)\,,$

and rewriting it using $\rho$,

$\frac{\textnormal{d}}{\textnormal{d}t}\overline{\mathscr{F}}^{\mu}_{t}=-\mathbb{E}_{\rho}(|\nabla
F|^{2})-\mathrm{Var}_{\rho}(F)+\mathbb{E}_{\rho}(\Delta F)\,.$

Because $\nabla_{\mathbf{x}}F$ is $L$-Lipschitz according to the first
condition in Assumption 2.1 and $F_{*}=0$, we obtain

$|\nabla F|^{2}\leqslant 2LF(\mathbf{x})\,.$

This implies the first term can be lower bounded by the expectation of
$F(\mathbf{x})$:

(A.7) $-\mathbb{E}_{\rho}(|\nabla
F|^{2})\geqslant-2L\overline{\mathscr{F}}^{\mu}_{t}\,.$

Since $\mathbf{x}_{\ast}$ is the unique global minimal, there exists
$\epsilon^{\prime}>0$ such that

$F^{-1}([0,\epsilon^{\prime}])\subset B_{R}(\mathbf{x}_{\ast})\,.$

Now, we lower bound the last two terms:

  * •

Because $F\leqslant D$ according to the second condition in Assumption 2.1 and
$F_{*}=0$, we obtain

$\mathrm{Var}_{\rho}(F)=\mathbb{E}_{\rho}(F^{2})-(\overline{\mathscr{F}}^{\mu}_{t})^{2}\leqslant
D\overline{\mathscr{F}}^{\mu}_{t}\left(1-\frac{\overline{\mathscr{F}}^{\mu}_{t}}{D}\right)\leqslant
D\overline{\mathscr{F}}^{\mu}_{t}\,,$

where we use Assumption 2.1, $0\leqslant F\leqslant D$ in the inequalities.
This implies

(A.8) $-\mathrm{Var}_{\rho}(F)\geqslant-D\overline{\mathscr{F}}^{\mu}_{t}\,.$

  * •

Because $F(x)$ is $\xi$-strongly convex in the ball of
$B_{R}(\mathbf{x}_{*})$, we obtain that

$\Delta F\geqslant\xi d,\quad\forall\mathbf{x}\in B_{R}(\mathbf{x}_{*})$

and because $\nabla F$ is $L$-Lipschitz, we have:

$\left|\Delta F\right|\leqslant dL\,,\forall\mathbf{x}\in\mathbb{R}^{d}\,.$

To lower bound the last term, we first deploy Markov inequality:

$\mathbb{P}(\mathbf{x}\in
B_{R}(\mathbf{x}_{\ast}))\geqslant\mathbb{P}(\mathbf{x}\in
F^{-1}([0,\epsilon^{\prime}]))=1-\mathbb{P}(F(\mathbf{x})\geqslant\epsilon^{\prime})\geqslant
1-\frac{\overline{\mathscr{F}}^{\mu}_{t}}{\epsilon^{\prime}}\,.$

Then, we have

(A.9) $\mathbb{E}_{\rho}(\Delta F)\geqslant
d\left(\xi\underbrace{\left(1-\frac{\overline{\mathscr{F}}^{\mu}_{t}}{\epsilon^{\prime}}\right)}_{\text{contribution
from
}B_{R}(\mathbf{x}_{\ast})}-\underbrace{L\left(\frac{\overline{\mathscr{F}}^{\mu}_{t}}{\epsilon^{\prime}}\right)}_{\text{contribution
outside }B_{R}(\mathbf{x}_{\ast})}\right)\,.$

Combining (A.7), (A.8), and (A.9), we obtain

$\frac{\textnormal{d}}{\textnormal{d}t}\overline{\mathscr{F}}^{\mu}_{t}\geqslant
d\xi-\overline{\mathscr{F}}^{\mu}_{t}\left(2L+D+\frac{d\xi}{\epsilon^{\prime}}+\frac{dL}{\epsilon^{\prime}}\right)\,.$

Set
$\epsilon=d\xi/\left(2L+D+\frac{d\xi}{\epsilon^{\prime}}+\frac{dL}{\epsilon^{\prime}}\right)$,
we obtain that
$\frac{\textnormal{d}}{\textnormal{d}t}\overline{\mathscr{F}}^{\mu}_{t}\geqslant
0$ if $\overline{\mathscr{F}}^{\mu}_{t}\leqslant\epsilon$. This proves (4.5).
∎

##  Appendix B Technical lemmas

The following lemma is deployed in the proof for the long time behavior Lemma
5.4.

######  Lemma B.1.

Under the condition of Theorem 2.2. Let $\mu_{t}(\mathbf{x},m)$ solve (2.4a).
If at a time $t_{*}$ so that $\mu_{t_{*}}(\mathbf{x},m)$ takes on the form of
$\delta_{\mathbf{x}_{*}}\otimes\eta(m)$ for some probability measure $\eta$,
then $\mathrm{Supp}_{m}(\cdot)\subset[m_{c},\infty)$.

###### Proof.

This can be seen by deploying a special test function. Set the test function
$\phi_{*}(\mathbf{x},m)=|\mathbf{x}-\mathbf{x}_{*}|^{2}\xi(\mathbf{x})$, where
$\xi(\mathbf{x})$ is a compactly supported, non-negative, smooth function such
that $\xi(\mathbf{x})=1$ in a neighborhood of $\mathbf{x}_{*}$. Naturally by
definition, $\partial_{m}\phi_{*}(\mathbf{x},m)=0$. Testing it on the mean-
field PDE (2.4a), we obtain that

(B.1)
$\displaystyle\underbrace{\frac{\textnormal{d}}{\textnormal{d}t}\iint\phi_{*}(\mathbf{x},m)\hskip
0.56917pt\textnormal{d}\mu_{t_{*}}(\mathbf{x},m)}_{\mathrm{Term\ I}}$
$\displaystyle\qquad=-\iint\left(\underbrace{\nabla\phi_{*}(\mathbf{x},m)\cdot\nabla
F(\mathbf{x})}_{\mathrm{Term\
II}}+\sigma(m)\Delta\phi_{*}(\mathbf{x},m)\right)\hskip
0.56917pt\textnormal{d}\mu_{t_{*}}(\mathbf{x},m)\,,$

where the
$\big{(}F(\mathbf{x})-\overline{\mathscr{F}}^{\mu}_{t}\big{)}\partial_{m}(m\mu)\left.\right|_{\mu_{t_{*}}(\mathbf{x},m)}=0$
and is dropped out upon testing, using the same argument that showed (5.14).
We will show that both $\mathrm{Term\ I}$ and $\mathrm{Term\ II}$ are zero.

Indeed, denote $A(t):=\iint\phi_{*}(\mathbf{x},m)\hskip
0.56917pt\textnormal{d}\mu_{t}(\mathbf{x},m)$. According to the definition of
$\phi_{*}$, $A(t)\geqslant 0$ for all $t$. Meanwhile,
$A(t_{*})=\int\phi_{*}(\mathbf{x}_{*},m)\textnormal{d}\eta(m)=0$. As a
consequence,

$\mathrm{Term\
I}=\frac{\textnormal{d}}{\textnormal{d}t}\iint\phi_{*}(\mathbf{x},m)\hskip
0.56917pt\textnormal{d}\mu_{t_{*}}(\mathbf{x},m)=\frac{\textnormal{d}}{\textnormal{d}t}A(t_{*})=0\,.$

Similarly $\mathrm{Term\ II}=0$ is due to the fact that $\nabla
F(\mathbf{x}_{*})=0$.

These together suggests $\iint\sigma(m)\Delta\phi_{*}(\mathbf{x},m)\hskip
0.56917pt\textnormal{d}\mu_{t_{*}}(\mathbf{x},m)=0$. Considering
$\mu_{t_{*}}(\mathbf{x},m)=\delta_{\mathbf{x}_{*}}\otimes\eta(m)$ and that
$\Delta\phi_{*}(\mathbf{x}_{*})=2d$, we have
$\int\sigma(m)\textnormal{d}\eta(m)=0$, and thus
$\mathrm{Supp}(\eta(\cdot))\subset[m_{c},\infty)$. ∎