# Strong existence and uniqueness of a calibrated local stochastic volatility
model

Scander Mustapha  Program in Applied & Computational Mathematics, Princeton
University, Princeton, NJ 08544<EMAIL_ADDRESS>
(Date: May 2024)

###### Abstract.

We study a two-dimensional McKean-Vlasov stochastic differential equation,
whose volatility coefficient depends on the conditional distribution of the
second component with respect to the first component. We prove the strong
existence and uniqueness of the solution, establishing the well-posedness of a
two-factor local stochastic volatility (LSV) model calibrated to the market
prices of European call options. In the spirit of [10], we assume that the
factor driving the volatility of the log-price takes finitely many values.
Additionally, the propagation of chaos of the particle system is established,
giving theoretical justification for the algorithm [5].

soft open fences

##  1\. Introduction

For the family of local stochastic volatility (LSV) models

(1.1) $\begin{cases}dS_{t}=\sigma_{\textrm{lsv}}(t,S_{t})a_{t}S_{t}dB_{t},\
t\geq 0,\\\ S_{t=0}=S_{0},\\\ \end{cases}$

to be calibrated to the market prices of European call options
$(C(t,K))_{t>0,K>0}$ it is enough to take the leverage function
$\sigma_{\textrm{lsv}}$ as

(1.2)
$\sigma_{\textrm{lsv}}(t,K)=\frac{\sigma_{\textrm{loc}}(t,K)}{\sqrt{\mathbb{E}\left[a_{t}^{2}|S_{t}=K\right]}},\
t>0,K>0,$

where $(S_{t})_{t\geq 0}$ is the price process, $(B_{t})_{t\geq 0}$ is a real
standard Brownian motion, $S_{0}$ is a positive random variable independent of
$(B_{t})_{t\geq 0}$,
$\sigma_{\textrm{loc}}(t,K):=\sqrt{\frac{2\partial_{t}C(t,K)}{\partial_{K}^{2}C(t,K)}},t>0,K>0$
is the so-called Dupire volatility [4] and $(a_{t})_{t\geq 0}$ is any
stochastic process. This is justified by Gyöngy’s theorem [7] that asserts
that a stochastic process $(S_{t})_{t\geq 0}$ solving (1.1) and (1.2) has
fixed-time marginal distributions given by the marginal distributions of the
local volatility model
$dS^{\textrm{loc}}_{t}=\sigma_{\textrm{loc}}(t,S^{\textrm{loc}}_{t})S^{\textrm{loc}}_{t}dB_{t}$
and starting at $S_{0}$. In order to make the model tractable, one commonly
assumes that $(a_{t})_{t\geq 0}$ is given by $(\sqrt{f(Y_{t})})_{t\geq 0}$
where $f:\mathbb{R}\to\mathbb{R}_{+}^{*}$ is a bounded smooth function and
$(Y_{t})_{t\geq 0}$ is another real Itô-Lévy process eventually correlated to
$(B_{t})_{t\geq 0}$.

Formally, the joint density $p(t,s,y)$ of
$(S_{t},Y_{t}),(t,s,y)\in\mathbb{R}_{+}\times\mathbb{R}_{+}^{*}\times\mathbb{R}$,
solves a quasi-linear and non-local Fokker-Planck partial differential
equation (PDE), whose coefficients depend upon the non-linear term

(1.3) $\frac{\int p(\cdot,\cdot,z)dz}{\int f(z)p(\cdot,\cdot,z)dz}.$

While the applications are important in calibration to market implied
volatility surfaces ([8], [6, Chapter 11], [1], [15]) the existence and
uniqueness of solutions to the stochastic differential equation (SDE) and the
partial differential equation (PDE) problems is still an open problem. Only
partial results have been obtained in particular cases.

Abergel and Tachet showed in [2, Theorem 3.1], the existence of a classical
solution to the PDE problem in a bounded domain of $\mathbb{R}^{2}$ and with
additional Dirichlet boundary condition - for short time, and for
$\sup|f^{\prime\prime}|$ small enough. In [12, Theorem 1.4], the authors
established the existence and uniqueness of a stationary solution to a similar
SDE (a drift needs to be added to the dynamic of $(S_{t})_{t\geq 0}$ to allow
the possibility of a stationary measure). Jourdain and Zhou proved in [10,
Theorems 2.4, 2.5], the existence of a weak solution to the SDE. Under the
assumption that $(Y_{t})_{t\geq 0}$ is a jump process taking finitely many
values and after writing the PDE problem as a system of parabolic equations in
non-divergence form, the authors are able to write a variational formulation
of the problem, for which existence of a solution can be proved - by the
Galerkin’s method - provided that the range of $f$ is small enough. They are
then able to prove the existence of a weak solution to the SDE. Uniqueness and
propagation of chaos are out of reach in the approach of [10] because higher
regularity of the weak solution is needed. Global regularity of weak solutions
to parabolic systems is still an open problem.

The difficulty of the analysis of the McKean-Vlasov SDE and the related non-
linear PDE stems from the singularity of the denominator $\int
f(z)p(\cdot,\cdot,z)dz$. However, the well-posedness of McKean-Vlasov problems
with coefficients depending smoothly on the density $p$ of the unknown process
has been studied in various papers [9], [14], [11], [3].

The main contribution of our paper is the proof of well-posedness of a
regularized version of SDE (1.1) - (1.2) when $(a_{t})_{t\geq 0}$ is time
independent and given by $f(Y)$, where $Y$ is a fixed random variable
independent of $(B_{t})_{t\geq 0}$ and takes finitely many values in
$\\{1,\dots,N\\}$, $N\geq 2$ and $f:\\{1,\dots,N\\}\to\mathbb{R}_{+}^{*}$. The
regularization is chosen so that the calibration property is conserved. In
other words, the fixed-time marginals of a solution are still given by the
marginals of
$S^{\textrm{loc}}_{t}=S_{0}+\int_{0}^{t}\sigma_{\textrm{loc}}(s,S_{s}^{\textrm{loc}})S_{s}^{\textrm{loc}}dB_{s}$,
$t\geq 0$.

Let $(X_{t})_{t\geq 0}=(\log{S_{t}})_{t\geq 0}$ be the log price and set
$\sigma(t,x):=\sigma_{\textrm{loc}}(t,e^{x}),(t,x)\in\mathbb{R}_{+}\times\mathbb{R}$.
Herein, we assume that $Y$ is a random variable taking finitely many values
$\\{1,\cdots,N\\}$, $N\geq 2$. Let $\varepsilon>0$ and consider the McKean-
Vlasov SDE

(1.4)
$\begin{cases}dX_{t}=-\frac{1}{2}\frac{\varepsilon+f(Y)p(t,X_{t})}{\varepsilon+\mathbb{E}\left[f(Y)|X_{t}\right]p(t,X_{t})}\sigma(t,X_{t})^{2}dt+\sqrt{\frac{\varepsilon+f(Y)p(t,X_{t})}{\varepsilon+\mathbb{E}\left[f(Y)|X_{t}\right]p(t,X_{t})}}\sigma(t,X_{t})dB_{t},\\\
\mathbb{P}\left(X_{t}\in dx\right)=p(t,x)dx,\\\ X_{t=0}=X_{0},\end{cases}$

where $X_{0}$ is a real random variable independent of $(B_{t})_{t\geq 0}$.

We notice immediately

(1.5)
$\mathbb{E}\left[\left(\sqrt{\frac{\varepsilon+f(Y)p(t,X_{t})}{\varepsilon+\mathbb{E}\left[f(Y)|X_{t}\right]p(t,X_{t})}}\sigma(t,X_{t})\right)^{2}\Big{|}X_{t}\right]=\sigma(t,X_{t})^{2},$

and consequently the model is calibrated, meaning $\exp(X_{t})$ and
$S_{t}^{\textrm{loc}}$ have the same law for all $t\geq 0$.

In the general case where $(a_{t})_{t\geq 0}=(f(Y_{t}))_{t\geq 0}$ and $Y$ is
an Itô process, the uniform ellipticity of the non-linear PDE - written in
divergence form - solved by the joint density $p(t,x,y)$ of $(X_{t},Y_{t})$,
$t\geq 0$, does not hold a priori. This property is a key element in
establishing uniqueness in [9, Equation 1.8]. Therefore, we made the
restrictive assumption that $Y$ is a time-independent discrete random variable
and takes finitely many values. Assuming that the range of $f$ is small
enough, we are able to derive a uniform ellipticity property for the PDE
problem, in the spirit of what is done in [10].

We assume in the whole paper that $\sigma$ is smooth, has bounded derivatives
and there exist $0<\sigma_{0}<\sigma_{1}$ such that
$\sigma_{0}\leq\sigma(t,x)\leq\sigma_{1}$ for all
$(t,x)\in\mathbb{R}_{+}\times\mathbb{R}$. Moreover, we assume that the measure
$\mathbb{P}\left(X_{0}\in dx\cap Y=n\right)$ admits a density
$P_{n}:\mathbb{R}\to\mathbb{R}_{+}$ \- of total mass
$\mathbb{P}\left(Y=n\right)$ \- with bounded first and second derivatives and
such that for some $\alpha\in(0,1)$, the second derivative $P^{(2)}_{n}$ is
$\alpha$-Hölder and

(1.6) $\|P_{n}\|_{C^{2+\alpha}}:=\sum_{0\leq k\leq
2}\|P^{(k)}_{n}\|_{\infty}+\sup_{x,y\in\mathbb{R},x\neq
y}\frac{|P_{n}^{(k)}(x)-P_{n}^{(k)}(y)|}{|x-y|^{\alpha}}<\infty.$

Denote $f_{\max}=\max_{1\leq n\leq N}f(n)$, $f_{\min}=\min_{1\leq n\leq
N}f(n)$ and

(1.7) $\bar{f}=\frac{1}{N}\sum_{n=1}^{N}f(n).$

We introduce the following small range condition on $f$.

######  Condition 1.1.

(1.8) $\begin{split}\frac{1}{2}\left[N+1-\max_{1\leq k\leq
N}\sqrt{\sum_{n=1,n\neq k}^{N}f(n)\sum_{n=1,n\neq
k}^{N}\frac{1}{f(n)}}\right]\wedge
1>&\frac{1}{f_{\min}}\sqrt{\sum_{n=1}^{N}(f(n)-\bar{f})^{2}}\\\
&+\frac{f_{\max}-f_{\min}}{f_{\min}}.\\\ \end{split}$

We are now ready to state our main result.

######  Theorem 1.2.

Let $T>0$. If Condition 1.1 is satisfied, $P$ belongs to
$(C^{2+\alpha})^{N},\alpha\in(0,1)$ and the norm
$\sum_{n=1}^{N}\|P_{n}\|_{C^{2+\alpha}}$ is small enough, then there exists a
unique strong solution to SDE (1.4) on $[0,T]$. Moreover, $X$ admits a smooth
density $p\in\left(C^{1+\alpha/4,2+\alpha/2}\right)^{N}$.

We study as well the question of propagation of chaos for SDE (1.4). Let
$M\geq 1$, $(\delta_{M})_{M\geq 1}\in(0,1)^{\mathbb{N}^{*}}$ be a sequence
converging to $0$ and $W_{1}:\mathbb{R}\to\mathbb{R}_{+}^{*}$ be a bounded
kernel function with a bounded derivative and satisfying $\int W_{1}(x)dx=1$
and $\int xW_{1}(x)dx=0$. Denote

(1.9)
$W_{\delta_{M}}:=\frac{1}{\delta_{M}}W_{1}\left(\frac{\cdot}{\delta_{M}}\right).$

Let $(B_{t}^{i})_{t\geq 0,i\geq 1}$ be a collection of i.i.d standard Brownian
motions and $(X_{0}^{i},Y^{i})_{i\geq 1}$ be a collection of i.i.d random
variables of law $\mathbb{P}\left(Y^{i}=n\cap X_{0}^{i}\in
dx\right)=P_{n}(x)dx,\ 1\leq n\leq N$.

For $M\geq 1$, we introduce the approximating system of $M$ particles, already
considered in [5] (for $\varepsilon=0$)

(1.10)
$\begin{split}dX^{i,M}_{t}=&-\frac{1}{2}\frac{\varepsilon+f(Y^{i})\frac{1}{M}\sum_{j=1}^{M}W_{\delta_{M}}(X^{i,M}_{t}-X_{t}^{j,M})}{\varepsilon+\frac{1}{M}\sum_{j=1}^{M}f(Y^{j})W_{\delta_{M}}(X^{i,M}_{t}-X^{j,M})}\sigma(t,X^{i,M}_{t})^{2}dt\\\
&+\sqrt{\frac{\varepsilon+f(Y^{i})\frac{1}{M}\sum_{j=1}^{M}W_{\delta_{M}}(X^{i,M}_{t}-X_{t}^{j,M})}{\varepsilon+\frac{1}{M}\sum_{j=1}^{M}f(Y^{j})W_{\delta_{M}}(X^{i,M}_{t}-X^{j,M})}}\sigma(t,X^{i,M}_{t})dB^{i}_{t},\\\
\end{split}$

initialized at $X^{i,M}_{t=0}=X_{0}^{i}$, $1\leq i\leq N$.

Assume that the conditions of Theorem 1.3 are satisfied. For each $1\leq i\leq
M$, let $(\hat{X}^{i})_{t\geq 0}$ be the particle starting at $X^{i,M}_{0}$
solving (1.4) and driven by the Brownian motion $(B^{i}_{t})_{t\geq 0}$ and
$Y^{i}$

(1.11) $\begin{cases}\begin{aligned}
d\hat{X}^{i}_{t}=&-\frac{1}{2}\frac{\varepsilon+f(Y^{i})\sum_{n=1}^{N}\hat{p}_{n}(t,\hat{X}^{i}_{t})}{\varepsilon+\sum_{n=1}^{N}f(n)\hat{p}_{n}(t,\hat{X}^{i}_{t})}\sigma(t,\hat{X}_{t}^{i})^{2}dt\\\
&+\sqrt{\frac{\varepsilon+f(Y^{i})\sum_{n=1}^{N}\hat{p}_{n}(t,\hat{X}^{i}_{t})}{\varepsilon+\sum_{n=1}^{N}f(n)\hat{p}_{n}(t,\hat{X}^{i}_{t})}}\sigma(t,\hat{X}_{t}^{i})dB^{i}_{t},\\\
\end{aligned}\\\ \mathbb{P}\left(\hat{X}^{i}_{t}\in dx\cap
Y^{i}=n\right)=\hat{p}_{n}(t,x)dx,\\\ 0\leq t\leq T.\end{cases}.$

The existence and uniqueness of the process $(X^{i,M}_{t})_{0\leq t\leq
T,1\leq i\leq N}$ is ensured by the regularity of the drift and diffusion
coefficients, given by Theorem 1.2.

Under Condition 1.1, propagation of chaos holds.

######  Theorem 1.3.

Let $T>0$, $(X^{i,M}_{t})_{t\geq 0}$ and $(\hat{X}^{i}_{t})_{t\geq 0}$ be
given by (1.10) and (1.11) respectively. Assume that $(\delta_{M})_{M\geq 1}$
converges slowly enough towards zero such that

(1.12)
$\lim_{M\to\infty}\frac{1}{M}\frac{\exp{(CT\delta_{M}^{-4})}}{\delta_{M}^{2}}=0.$

for all constants $C>0$.

Then for all $1\leq i\leq N$

(1.13) $\lim_{M\to\infty}\mathbb{E}\left[\sup_{0\leq t\leq
T}|X^{i}_{t}-\hat{X}^{i}_{t}|^{2}\right]=0.$

Consequently, for all $k\geq 1$ the law of the process
$(X^{1,M},X^{2,M},\dots,X^{k,M})$ converges weakly towards
$\mu\otimes\mu\otimes\dots\otimes\mu$, where $\mu$ is the law of
$\hat{X}^{1}$.

The main application of this result is the justification of the particular
method for calibrating the LSV model described in [5], [6, Section 11.6.1].

The rest of the article is structured as follows. In Section 2 we introduce
notation and prove Theorem 1.2. Existence is proved in Subsection 2.1 and
uniqueness in Subsection 2.2. Section 3 is devoted to establishing propagation
of chaos.

##  2\. Existence and uniqueness

Throughout the rest of the paper $T>0$ is fixed, and

  * •

For any $d\geq 1$ and $X\in\mathbb{R}^{d}$, we denote
$\|X\|_{2}=\sqrt{\sum_{i=1}^{d}X_{i}^{2}}$ and $\|X\|_{\infty}=\max_{1\leq
i\leq d}|X_{i}|$.

  * •

The scalar product is denoted $\langle\cdot,\cdot\rangle$

(2.1) $\langle X,Y\rangle=\sum_{i=1}^{d}X_{i}Y_{i},\
X\in\mathbb{R}^{d},Y\in\mathbb{R}^{d}.$

  * •

We denote by $L^{p},p\geq 1$ the space of measurable real-valued functions
$\phi$ for which $\|\phi\|_{L^{p}}=\left(\int|\phi(x)|^{p}dx\right)^{1/p}$ is
finite.

  * •

We denote by $C(0,T,L^{2})$ the space of measurable functions
$\phi:(t,x)\in[0,T]\times\mathbb{R}\to\mathbb{R}$ such that the mapping
$t\in[0,T]\to\phi(t,\cdot)$ takes values in $(L^{2},\|\|_{L^{2}})$ and is
continuous.

  * •

For $\alpha\in(0,1)$ we denote by $C^{2+\alpha}$ the space of real functions
$\phi$ on $\mathbb{R}$ for which $\|\phi\|_{C^{2+\alpha}}<\infty$. We denote
by $C^{(k+\alpha)/2,k+\alpha}$ the space of real functions $\phi$ on
$[0,T]\times\mathbb{R}$ which are continuous together with their derivatives
$\partial_{t}^{r}\partial_{x}^{l}\phi$, $2r+l\leq k$ and admit a finite norm

(2.2) $\begin{split}\|\phi\|_{C^{(k+\alpha)/2,k+\alpha}}&=\sum_{2r+l\leq
k}\sup_{(t,x)\in[0,T]\times\mathbb{R}}|\partial_{t}^{r}\partial_{x}^{l}\phi(t,x)|\\\
&+\sum_{k-1\leq 2r+l\leq
k}\sup_{x\in\mathbb{R},t,s\in[0,T]}\frac{|\partial_{t}^{r}\partial_{x}^{l}\phi(t,x)-\partial^{r}_{t}\partial^{l}_{x}\phi(s,x)|}{|t-s|^{(k-2r-l+\alpha)/2}}\\\
&+\sum_{2r+l=k}\sup_{t\in[0,T],x,y\in\mathbb{R}}\frac{|\partial_{t}^{r}\partial_{x}^{l}\phi(t,x)-\partial^{r}_{t}\partial^{l}_{x}\phi(t,y)|}{|x-y|^{\alpha}},\end{split}$

where $\partial_{t}$ and $\partial_{x}$ are the partial derivatives with
respect to $t$ and $x$.

  * •

$S^{++}$ is the set of symmetric and positive $N\times N$ real matrices.

For each $1\leq n\leq N$ and $t\geq 0$, denote by $p_{n}(t,x)$ the conditional
density of $X_{t}$ given $Y=n$, multiplied by $\mathbb{P}\left(Y=n\right)$. In
other words, $p_{n}(t,x)dx=\mathbb{P}\left(X_{t}\in dx\cap
Y=n\right),x\in\mathbb{R},1\leq n\leq N$. Set
$P_{n}(x)dx=\mathbb{P}\left(X_{0}\in dx\cap Y=n\right)$ for $1\leq n\leq N$.

The PDE problem solved by $(p_{n})_{1\leq n\leq N}$ can be formulated as a
system of $N$ parabolic PDEs

(2.3)
$\begin{cases}\partial_{t}p=\frac{1}{2}\partial_{xx}[\sigma^{2}B^{\varepsilon}(p)p]+\frac{1}{2}\partial_{x}\left[\sigma^{2}B^{\varepsilon}(p)p\right],\
(t,x)\in[0,T]\times\mathbb{R},\\\ p_{n}/\mathbb{P}\left(Y=n\right)\textrm{is a
probability density,}\\\ p_{n}(0,x)=P_{n}(x),\ x\in\mathbb{R},1\leq n\leq
N,\\\ \end{cases}$

where $B^{\varepsilon}(u),u\in(\mathbb{R}_{+})^{N}$ is a diagonal matrix whose
diagonal elements are

(2.4)
$B^{\varepsilon}_{nn}(u)=B^{\varepsilon}_{n}(u):=\frac{\varepsilon+f(n)\sum_{k=1}^{N}u_{k}}{\varepsilon+\sum_{k=1}^{N}f(k)u_{k}},\quad
u\in\mathbb{R}_{+}^{N},1\leq n\leq N.$

The proof of Theorem 1.2 is a direct application of Propositions 2.1 and 2.4
below.

###  2.1. Existence

Under the condition that the norm $\sum_{n=1}^{N}\|P_{n}\|_{C^{2+\alpha}}$ is
small enough, we prove the existence, in the classical sense, of a solution to
Problem (2.3) and the existence of a stochastic process $(X_{t})_{t\geq 0}$
solving SDE (1.4). The proof follows the approach of [9, Proof of Proposition
2.2] adapted to our specific McKean-Vlasov equation, whose coefficients depend
on the marginal $\sum_{n=1}^{N}p_{n}$ and the quantity
$\sum_{n=1}^{N}f(n)p_{n}$.

######  Proposition 2.1.

There exists a constant $C$, depending on $T$, $f_{\max}$, $f_{\min}$ and
$\varepsilon$, such that if $\sum_{n=1}^{N}\|P_{n}\|_{C^{2+\alpha}}\leq C$,
then there exist a solution $(p_{n})_{1\leq n\leq
N}\in(C^{1+\alpha/2,2+\alpha})^{N}$ to Problem (2.3), and a solution
$(X_{t})_{t\in[0,T]}$ to SDE (1.4).

###### Proof.

For $1\leq n\leq N$ and non-negative $u\in(C^{1+\alpha/2,2+\alpha})^{N}$,
define the operator

(2.5)
$\partial_{t}v-L_{n,u}v:=\partial_{t}v-\frac{1}{2}\partial_{xx}[\sigma^{2}B^{\varepsilon}_{n}(u)v]-\frac{1}{2}\partial_{x}[\sigma^{2}B^{\varepsilon}_{n}(u)v].$

$L_{n,u}$ is a uniformly parabolic operator of second order with coefficients
in $C^{\alpha/2,\alpha}$. According to [9, Proposition 1.1](or [13, Chapter
IV, Theorem 5.1]), there exists a unique $v_{n}\in C^{1+\alpha/2,2+\alpha}$
such that $v_{n}/\mathbb{P}\left(Y_{t}=n\right)$ is a probability density and
solves in $[0,T]\times\mathbb{R}$

(2.6) $\begin{cases}\partial_{t}v_{n}=L_{n,u}v_{n},\\\ v_{n}(0,x)=P_{n}(x),\
x\in\mathbb{R}.\end{cases}$

Moreover, there exists a constant $C^{\prime}$ depending only on the
regularity of the coefficients of $L_{n,u}$ (namely
$f_{\max},f_{\min},\varepsilon$ and
$\sum_{l=1}^{N}\|u_{l}\|_{1+\alpha/2,2+\alpha}$) such that

(2.7) $\|v_{n}\|_{C^{1+\alpha/2,2+\alpha}}\leq
C^{\prime}\left(f_{\max},f_{\min},\varepsilon,\sum_{l=1}^{N}\|u_{l}\|_{C^{1+\alpha/2,2+\alpha}}\right)\|P_{n}\|_{C^{2+\alpha}}.$

A useful fact about the constant $C^{\prime}$ is that it is non-decreasing in
its last argument.

We construct by induction a sequence of solutions as follows. Start by
$p^{0}_{n}(t,x)=P_{n}(x)$, $t\in[0,T],1\leq n\leq N$ and $x\in\mathbb{R}$. For
$m\geq 1$, let $p^{m}_{n}\in C^{1+\alpha/2,2+\alpha}$, for each $1\leq n\leq
N$, be the solution of

(2.8) $\begin{cases}\partial_{t}p^{m}_{n}=L_{n,p^{m-1}}p^{m}_{n},\\\
p^{m}_{n}/\mathbb{P}\left(Y=n\right)\textrm{is a probability density,}\\\
p^{m}_{n}(0,x)=P_{n}(x),\ x\in\mathbb{R}.\end{cases}$

Using estimate (2.7), we have for each $m\geq 1$

(2.9) $\sum_{n=1}^{N}\|p^{m}_{n}\|_{C^{1+\alpha/2,2+\alpha}}\leq
C^{\prime}\left(f_{\max},f_{\min},\varepsilon,\sum_{l=1}^{N}\|p^{m-1}_{l}\|_{C^{1+\alpha/2,2+\alpha}}\right)\sum_{n=1}^{N}\|P_{n}\|_{C^{2+\alpha}}.$

Under the assumption $\sum_{n=1}^{N}\|P_{n}\|_{2+\alpha}\leq
C:=1\wedge\frac{1}{C^{\prime}(f_{\max},f_{\min},\varepsilon,1)}$, we have by
induction that $\sum_{n=1}^{N}\|p^{m}_{n}\|_{C^{1+\alpha/2,2+\alpha}}\leq 1$
for each $m\geq 1$. Indeed, the result holds for $m=0$. If
$\sum_{n=1}^{N}\|p^{m-1}_{n}\|_{C^{1+\alpha/2,2+\alpha}}\leq 1$ then

(2.10) $\sum_{n=1}^{N}\|p^{m}_{n}\|_{C^{1+\alpha/2,2+\alpha}}\leq
C^{\prime}\left(f_{\max},f_{\min},\varepsilon,1\right)\sum_{n=1}^{N}\|P\|_{C^{2+\alpha}}\leq
1.$

We proved that $(p^{m})_{m\geq 0}$ is a bounded sequence of elements of
$\left(C^{1+\alpha/2,2+\alpha}\right)^{N}$, and similarly to what is done in
[9, Proposition 2.2], we can extract a limit point $(p_{n})_{1\leq n\leq
N}\in(C^{1+\alpha/4,2+\alpha/2})^{N}$, such that
$p_{n}/\mathbb{P}\left(Y=n\right)$ is a density of probability and
$(p_{n})_{1\leq n\leq N}$ solves Problem (2.3).

Now, for each $1\leq n\leq N$, $(p_{m})_{1\leq m\leq N}$ being regular, there
exists a strong solution $(X_{t},Y_{t})_{t\in[0,T]}$ to the SDE

(2.11)
$\begin{split}X_{t}=X_{0}&-\frac{1}{2}\int_{0}^{t}\frac{\varepsilon+f(Y)\sum_{m=1}^{N}p_{m}(s,X^{n}_{s})}{\varepsilon+\sum_{m=1}^{N}f(m)p_{m}(s,X^{n}_{s})}\sigma(t,X_{s})^{2}ds\\\
&+\int_{0}^{t}\sqrt{\frac{\varepsilon+f(Y)\sum_{m=1}^{N}p_{m}(s,X^{n}_{s})}{\varepsilon+\sum_{m=1}^{N}f(m)p_{m}(s,X^{n}_{s})}}\sigma(t,X_{s})dB_{s}.\end{split}$

Moreover, for each $1\leq n\leq N$ and $t\in[0,T]$, the law of $X_{t}1_{Y=n}$
admits a density $q_{n}(t,\cdot)$ solving $\partial_{t}q_{n}=L_{n,p}q_{n}$,
with initial data $P_{n}$. The solution of this problem being unique and given
by $p_{n}$, we conclude that $q_{n}=p_{n}$ for all $1\leq n\leq N$. ∎

###  2.2. Uniqueness

Uniqueness is proved by first writing system (2.3) in divergence form and then
proving the uniform ellipticity of the differential operator. For all $u\in
C^{1}(\mathbb{R}\times\mathbb{R}_{+})^{N}$, we can rewrite

(2.12) $\partial_{x}[B^{\varepsilon}(u)u]=A^{\varepsilon}(u)\partial_{x}u,$

where $(A^{\varepsilon}_{nk})_{1\leq n,k\leq
N}:\mathbb{R}_{+}^{N}\to\mathbb{R}^{N\times N}$ is defined by

(2.13)
$A^{\varepsilon}_{nn}(u):=B^{\varepsilon}_{n}(u)+u_{n}\frac{f(n)\sum_{l=1}^{N}(f(l)-f(n))u_{l}}{(\varepsilon+\sum_{l=1}^{N}f(l)u_{l})^{2}},\
1\leq n\leq N,$ (2.14)
$A^{\varepsilon}_{nk}(u):=u_{n}\frac{\varepsilon(f(n)-f(k))+f(n)\sum_{l=1}^{N}(f(l)-f(k))u_{l}}{(\varepsilon+\sum_{l=1}^{N}f(l)u_{l})^{2}},\
1\leq n\neq k\leq N.$

In these terms, PDEs (2.3) can be rewritten

(2.15)
$\begin{cases}\partial_{t}p=\frac{1}{2}\partial_{x}[\sigma^{2}A^{\varepsilon}(p)\partial_{x}p]+\partial_{x}[\left(\sigma\partial_{x}\sigma+\frac{1}{2}\sigma^{2}\right)B^{\varepsilon}(p)p],\
(t,x)\in[0,T]\times\mathbb{R},\\\ p_{n}/\mathbb{P}\left(Y=n\right)\textrm{is a
probability density,}\\\ p_{n}(0,x)=P_{n}(x),\ x\in\mathbb{R},1\leq n\leq
N.\end{cases}$

We recall the following result of [10, Proof of Proposition 2.3, Corollary
B.3].

######  Proposition 2.2.

Let

(2.16) $\kappa_{0}:=\frac{1}{2}\left[N+1-\max_{1\leq k\leq
N}\sqrt{\sum_{n=1,n\neq k}^{N}f(i)\sum_{n=1,n\neq
k}^{N}\frac{1}{f(n)}}\right].$

Then for all $\delta>0$, $\eta>0$, $u\in(\mathbb{R}_{+}^{*})^{N}$,
$U\in\textbf{1}\mathbb{R}$ and $V\in\textbf{1}^{\perp}$,
$\textbf{1}=(1)_{1\leq n\leq N}$,

(2.17) $\begin{split}\langle U+V,(J+\delta
I)A^{0}(u)(U+V)\rangle&\geq\frac{N}{2}\left[N-\delta\left(1+\frac{f_{\max}}{f_{\min}}\right)\left(1+\frac{1}{2\eta}\right)\right]\|U\|_{2}^{2}\\\
&+\delta\left(\kappa_{0}-N^{2}\left(1+\frac{f_{\max}}{f_{\min}}\right)\eta\right)\|V\|_{2}^{2},\end{split}$

where $I=(1_{i=j})_{1\leq i,j\leq N}$ and $J=(1)_{1\leq i,j\leq N}$.

We prove the following result.

######  Proposition 2.3.

Under Condition 1.1, there exist a matrix $\Gamma\in S^{++}$ and $\kappa>0$
such that for all $u\in\mathbb{R}_{+}^{N}$ and $X\in\mathbb{R}^{N}$

(2.18) $\langle X,\Gamma A^{\varepsilon}(u)X\rangle\geq\kappa\|X\|_{2}^{2}.$

###### Proof.

If $\sum_{n=1}^{N}u_{n}=0$ then $A^{\varepsilon}(u)=I$ and the result holds.
Assume that $\sum_{n=1}^{N}u_{n}>0$ and set

(2.19)
$\rho:=\frac{\sum_{l=1}^{N}f(l)u_{l}}{\varepsilon+\sum_{l=1}^{N}f(l)u_{l}}\in[0,1].$

Rewrite

(2.20)
$\begin{split}A^{\varepsilon}_{nn}(u)&=B^{\varepsilon}_{n}(u)+u_{n}\frac{f(n)\sum_{l=1}^{N}(f(l)-f(n))u_{l}}{(\varepsilon+\sum_{l=1}^{N}f(l)u_{l})^{2}}\\\
&=1-\rho+\rho
B_{n}^{0}(u)+\rho^{2}u_{n}\frac{f(n)\sum_{l=1}^{N}(f(l)-f(n))u_{l}}{(\sum_{l=1}^{N}f(l)u_{l})^{2}}\\\
&=(1-\rho)(1+\rho B^{0}_{n}(u))+\rho^{2}A^{0}_{nn}(u),\end{split}$

and

(2.21)
$\begin{split}A^{\varepsilon}_{nk}(u)&=u_{n}\frac{\varepsilon(f(n)-f(k))+f(n)\sum_{l=1}^{N}(f(l)-f(k))u_{l}}{(\varepsilon+\sum_{l=1}^{N}f(l)u_{l})^{2}}\\\
&=\rho(1-\rho)\frac{(f(n)-f(k))u_{n}}{\sum_{l=1}^{N}f(l)u_{l}}+\rho^{2}A^{0}_{nk}(u).\end{split}$

Define the matrix $D\in\mathbb{R}^{N\times N}$ by

(2.22) $D_{nk}:=\frac{(f(n)-f(k))u_{n}}{\sum_{l=1}^{N}f(l)u_{l}},\ 1\leq n\neq
k\leq N,$

and

(2.23) $D_{nn}:=-\sum_{m\neq
n}D_{mn}=\frac{\sum_{m=1}^{N}(f(n)-f(m))u_{m}}{\sum_{l=1}^{N}f(l)u_{l}}=B^{0}_{n}(u)-1,\
1\leq n\leq N.$

We see immediately that $\sum_{n=1}^{N}D_{nk}=0$, $JD=0$, and that

(2.24) $\max_{1\leq n,k\leq
N}|D_{nk}(u)|\leq\frac{f_{\max}-f_{\min}}{f_{\min}}.$

In view of

(2.25) $(1-\rho)(1+\rho B^{0}_{n})-D_{nn}=(1-\rho)(1+\rho
B^{0}_{n}(u))-\rho(1-\rho)(B^{0}_{m}(u)-1)=1-\rho^{2},$

we find that

(2.26) $A^{\varepsilon}(u)=\rho^{2}A^{0}(u)+(1-\rho^{2})I+\rho(1-\rho)D.$

Set $M:=A^{0}(u)-I$. For all $1\leq k\leq N$, $\sum_{n=1}^{N}M_{nk}=0$
according to [10, Proof of Lemma 3.13] and therefore $JM=0$. The previous
equation can be written

(2.27) $A^{\varepsilon}(u)=I+\rho^{2}M+\rho(1-\rho)D.$

Exploiting the idea of [10, Proof of Proposition 2.3], we are going to show
that for $\delta\in(0,1)$ small enough, there exists $\kappa>0$ such that for
all $X\in\mathbb{R}^{N}$ and $u\in(\mathbb{R}^{*}_{+})^{N}$

(2.28) $\langle X,(J+\delta
I)A^{\varepsilon}(u)X\rangle\geq\kappa\|X\|^{2}_{2}.$

Let $\delta\in(0,1)$ and $X\in\mathbb{R}^{N}$. We decompose

(2.29) $(J+\delta I)A^{\varepsilon}(u)=J+\delta I+\rho^{2}\delta
M+\rho(1-\rho)\delta D.$

If $\langle X,MX\rangle\leq 0$, then using $\rho\in[0,1]$,

(2.30) $\begin{split}\langle X,(J+\delta I)A^{\varepsilon}(u)X\rangle&=\langle
X,(J+\delta I)X\rangle+\rho^{2}\delta\langle
X,MX\rangle+\rho(1-\rho)\delta\langle X,DX\rangle\\\ &\geq\langle X,(J+\delta
I)X\rangle+\delta\langle X,MX\rangle+\rho(1-\rho)\delta\langle X,DX\rangle\\\
&=\langle X,(J+\delta I)A^{0}X\rangle+\rho(1-\rho)\delta\langle
X,DX\rangle.\end{split}$

If $\langle X,MX\rangle>0$ then

(2.31) $\langle X,(J+\delta I)A^{\varepsilon}(u)X\rangle\geq\langle
X,(J+\delta I)X\rangle+\rho(1-\rho)\delta\langle X,DX\rangle.$

Write $X=U+V$ where $U\in\textbf{1}\mathbb{R}$ and $V\in\textbf{1}^{\perp}$.
Then $DU=0$ and

(2.32) $\langle X,DX\rangle=\langle V,DV\rangle+\langle U,DV\rangle.$

On the one hand,

(2.33) $\begin{split}\langle
V,DV\rangle&=\sum_{n=1}^{N}D_{nn}V^{2}_{n}+\sum^{N}_{n,k=1}\frac{u_{n}(f(n)-f(k))}{\sum_{l=1}^{N}f(l)u_{l}}V_{n}V_{k}\\\
&=\sum_{n=1}^{N}D_{nn}V^{2}_{n}+\sum^{N}_{n,k=1}\frac{u_{n}((f(n)-\bar{f})-(f(k)-\bar{f}))}{\sum_{l=1}^{N}f(l)u_{l}}V_{n}V_{k}\\\
&=\sum_{n=1}^{N}D_{nn}V^{2}_{n}-\frac{\sum_{n,k=1}^{N}u_{n}(f(k)-\bar{f})V_{n}V_{k}}{\sum_{l=1}^{N}f(l)u_{l}},\end{split}$

where we used that $\sum_{k=1}^{N}V_{k}=0$ between the second and third line.

We estimate

(2.34)
$\begin{split}\left|\sum_{n,k=1}^{N}u_{n}(f(k)-\bar{f})V_{n}V_{k}\right|&=\left|\sum_{n=1}^{N}u_{n}V_{n}\sum_{k=1}^{N}(f(k)-\bar{f})V_{k}\right|\\\
&\leq\sqrt{\sum_{n=1}^{N}u_{n}^{2}}\sqrt{\sum_{n=1}^{N}(f(n)-\bar{f})^{2}}\|V\|_{2}^{2}.\end{split}$

Using that $u_{n}>0$ for $1\leq n\leq N$

(2.35)
$\frac{\sqrt{\sum_{n=1}^{N}u_{n}^{2}}}{{\sum_{l=1}^{N}f(l)u_{l}}}\leq\frac{1}{f_{\min}}\sqrt{\frac{\sum_{n=1}^{N}u_{n}^{2}}{\left(\sum_{l=1}^{N}u_{l}\right)^{2}}}\leq\frac{1}{f_{\min}},$

we readily conclude that

(2.36) $\rho(1-\rho)\delta\langle
V,DV\rangle\geq-\delta\beta(f)\|V\|_{2}^{2}.$

where
$\beta(f):=\left(\frac{f_{\max}-f_{\min}}{f_{\min}}+\frac{1}{f_{\min}}\sqrt{\sum_{n=1}^{N}(f(n)-\bar{f})^{2}}\right)$,

On the second hand, for all $\eta>0$

(2.37) $\begin{split}\rho(1-\rho)\delta\langle
U,DV\rangle&\geq-\delta\max_{1\leq n,k\leq N}|D_{nk}|N\|U\|_{2}\|V\|_{2}\\\
&\geq-\delta\beta(f)\sqrt{N}|\|U\|_{2}\|V\|_{2}\\\
&\geq-\delta\frac{\beta(f)N}{\eta}\|U\|^{2}_{2}-\delta\beta(f)N^{2}\eta\|V\|_{2}^{2}.\end{split}$

Thus combining (2.36), (2.37) and Proposition 2.2, we infer in the case
$\langle X,MX\rangle\leq 0$

(2.38) $\begin{split}\langle X,(J+\delta
I)A^{\varepsilon}(u)X\rangle&\geq\frac{N}{2}\left[N-\delta\left[\left(1+\frac{f_{\max}}{f_{\min}}\right)\left(1+\frac{1}{2\eta}\right)+\frac{2\beta(f)}{\eta}\right]\right]\|U\|_{2}^{2}\\\
&+\delta\left(\kappa_{0}-\beta(f)-N^{2}\left(1+\frac{f_{\max}}{f_{\min}}+\beta(f)\right)\eta\right)\|V\|_{2}^{2}.\end{split}$

Condition 1.1 yields $\kappa_{0}>\beta(f)$ and we can choose

(2.39)
$\eta<\eta_{-}:=\frac{1}{2}\frac{(\kappa_{0}-\beta(f))}{N^{2}}\left(1+\frac{f_{\max}}{f_{\min}}+\beta(f)\right)^{-1},$

and

(2.40)
$\delta<\delta_{-}(\delta):=\frac{N}{\left[\left(1+\frac{f_{\max}}{f_{\min}}\right)\left(1+\frac{1}{2\eta}\right)+\frac{2\beta(f)}{\eta}\right]}.$

We check that with such choice, if $\langle X,MX\rangle\leq 0$

(2.41) $\langle X,(J+\delta
I)A^{\varepsilon}(u)X\rangle\geq\kappa_{-}\|X\|_{2}^{2},$

where

(2.42)
$\kappa_{-}:=\min\left(\frac{N}{2}\left[N-\delta\left[\left(1+\frac{f_{\max}}{f_{\min}}\right)\left(1+\frac{1}{2\eta}\right)+\frac{2\beta(f)}{\eta}\right]\right],\frac{\delta}{2}(\kappa_{0}-\beta(f))\right)>0.$

Now in the case where $\langle X,MX\rangle>0$

(2.43) $\begin{split}\langle X,(J+\delta I)A^{\varepsilon}(u)X\rangle&\geq
N\left[1-\beta(f)\frac{\delta}{\eta}\right]\|U\|_{2}^{2}\\\
&+\delta\left[1-\beta(f)-\beta(f)N^{2}\eta\right]\|V\|_{2}^{2}.\end{split}$

Therefore, under Condition 1.1, $1>\beta(f)$, and we can choose

(2.44) $\eta<\eta_{+}:=\frac{1}{2}\frac{1-\beta(f)}{\beta(f)N^{2}},$

and

(2.45) $\delta<\delta_{+}(\eta):=\eta\beta(f).$

Thus

(2.46) $\langle U+V,(J+\delta
I)A^{\varepsilon}(u)(U+V)\rangle\geq\kappa_{+}\|X\|_{2}^{2},$

where

(2.47)
$\kappa_{+}:=\min\left(N\left[1-\beta(f)\frac{\delta}{\eta}\right],\frac{\delta}{2}\left[1-\beta(f)\right]\right)>0.$

To conclude, take $\eta<\min(\eta_{-},\eta_{+})$ and
$\delta<\min(\delta_{-}(\eta),\delta_{+}(\eta))$, to obtain

(2.48) $\langle X,(J+\delta
I)A^{\varepsilon}(u)X\rangle\geq\kappa\|X\|_{2}^{2},$

where $\kappa:=\min(\kappa_{-},\kappa_{+})>0$. ∎

We are now ready to prove uniqueness.

######  Proposition 2.4.

Under Condition 1.1, there exists at most one solution in
$C^{1+\alpha/4,2+\alpha/2}$ to Problem (2.3).

###### Proof.

Let $p$ and $q$ be two solutions in $(C^{1+\alpha/4,2+\alpha/2})^{N}$. By
integrating by parts, we have the estimate

(2.49)
$\|\partial_{x}p_{n}\|_{L^{2}}^{2}\leq\|\partial_{xx}p_{n}\|_{L^{\infty}}\|p_{n}\|_{L^{1}},\
1\leq n\leq N$

we see that $\partial_{x}p,\partial_{x}q\in C(0,T,L^{2})^{N}$.

Let $\Gamma\in S^{++}$ and $\kappa>0$ be given by Proposition 2.3, such that
for all $u\in(\mathbb{R}_{+})^{N}$

(2.50) $\langle X,\Gamma A^{\varepsilon}(u)X\rangle\geq\kappa\|X\|^{2}_{2}.$

Set $p^{\prime}=\sqrt{\Gamma}p$ and $q^{\prime}=\sqrt{\Gamma}q$. $p^{\prime}$
and $q^{\prime}$ are in $C^{1+\alpha/4,2+\alpha/2}$, and their spatial
gradient in $C(0,T,L^{2})$, and they solve respectively the system of PDEs

(2.51)
$\begin{cases}\partial_{t}p^{\prime}=\frac{1}{2}\partial_{x}[\sigma^{2}A^{\prime}(p)\partial_{x}p^{\prime}]++\partial_{x}[(\sigma\partial_{x}\sigma+\frac{1}{2}\sigma^{2})B^{\prime}(p)p^{\prime}],\
(t,x)\in[0,T]\times\mathbb{R},\\\ p^{\prime}_{n}(0,x)=\sqrt{\Gamma}P_{n}(x),\
x\in\mathbb{R},1\leq n\leq N,\end{cases}$

and

(2.52)
$\begin{cases}\partial_{t}q^{\prime}=\frac{1}{2}\partial_{x}[\sigma^{2}A^{\prime}(q)\partial_{x}q^{\prime}]+\partial_{x}[(\sigma\partial_{x}\sigma+\frac{1}{2}\sigma^{2})B^{\prime}(q)q^{\prime}],\
(t,x)\in[0,T]\times\mathbb{R},\\\ q^{\prime}_{n}(0,x)=\sqrt{\Gamma}P_{n}(x),\
x\in\mathbb{R},1\leq n\leq N,\end{cases}$

where the operators $A^{\prime}$ and $B^{\prime}$ are respectively defined by
$A^{\prime}(u):=\sqrt{\Gamma}A^{\varepsilon}(u)\sqrt{\Gamma}^{-1}$ and
$B^{\prime}(u):=\sqrt{\Gamma}B^{\varepsilon}(u)\sqrt{\Gamma}^{-1}$,
$u\in(\mathbb{R}_{+})^{N}$. $A^{\prime}$ satisfies the same coercivity
property as $A$. Indeed, for all $X\in\mathbb{R}$ and
$u\in(\mathbb{R}_{+})^{N}$

(2.53)
$\langle\sqrt{\Gamma}X,\sqrt{\Gamma}A^{\varepsilon}(u)\sqrt{\Gamma}^{-1}\sqrt{\Gamma}X\rangle=\langle
X,\Gamma
A^{\varepsilon}(u)X\rangle\geq\kappa\|X\|_{2}^{2}\geq\tilde{\kappa}\|\sqrt{\Gamma}X\|_{2}^{2},$

where $\tilde{\kappa}:=\frac{\kappa}{N+\delta}$, $(N+\delta)^{-1}$ being the
smallest eigenvalue of $\Gamma^{-1}$.

Multiplying (2.51) by $p^{\prime}-q^{\prime}$ and integrating gives for all
$t\in[0,T]$

(2.54) $\begin{split}\int\langle
p^{\prime}-q^{\prime},\partial_{t}p^{\prime}\rangle&=-\frac{1}{2}\int\langle\partial_{x}(p^{\prime}-q^{\prime}),A^{\prime}(p)\partial_{x}p^{\prime}\rangle\\\
&-\int\left(\sigma\partial_{x}\sigma+\frac{1}{2}\sigma^{2}\right)\langle\partial_{x}(p^{\prime}-q^{\prime}),B^{\prime}(p)p^{\prime}\rangle.\end{split}$

A similar equation holds for $q^{\prime}$. Taking the difference with the
previous equation and integrating over $[0,t]$ gives

(2.55)
$\begin{split}\frac{1}{2}\int\|p^{\prime}-q^{\prime}\|_{2}^{2}&=-\frac{1}{2}\int_{0}^{t}\int\langle\partial_{x}(p^{\prime}-q^{\prime}),A^{\prime}(p)\partial_{x}p^{\prime}-A^{\prime}(q)\partial_{x}q^{\prime}\rangle
dt\\\
&-\int_{0}^{t}\int\left(\sigma\partial_{x}\sigma+\frac{1}{2}\sigma^{2}\right)\langle\partial_{x}(p^{\prime}-q^{\prime}),B^{\prime}(p)p^{\prime}-B^{\prime}(q)q^{\prime}\rangle
dt.\end{split}$

Rewrite the integrand of the first integral of the right-hand side as

(2.56)
$\begin{split}\langle\partial_{x}(p^{\prime}-q^{\prime}),A^{\prime}(p)\partial_{x}p^{\prime}-A^{\prime}(q)\partial_{x}q^{\prime}\rangle&=\langle\partial_{x}(p^{\prime}-q^{\prime}),A^{\prime}(p)\partial_{x}(p^{\prime}-q^{\prime})\rangle\\\
&+\langle\partial_{x}(p^{\prime}-q^{\prime}),(A^{\prime}(p)-A^{\prime}(q))\partial_{x}q^{\prime}\rangle\\\
&\geq-\tilde{\kappa}\|\partial_{x}(p^{\prime}-q^{\prime})\|_{2}^{2}\\\
&-C\|q\|_{(C^{1+\alpha/2,2+\alpha})^{N}}\|\partial_{x}(p^{\prime}-q^{\prime})\|_{2}\|A^{\prime}(p)-A^{\prime}(q)\|_{2}.\end{split}$

Likewise for the second integral

(2.57)
$\begin{split}\sigma\partial_{x}\sigma\langle\partial_{x}(p^{\prime}-q^{\prime}),B^{\prime}(p)p^{\prime}-B^{\prime}(q)q^{\prime}\rangle&=\sigma\partial_{x}\sigma\langle\partial_{x}(p^{\prime}-q^{\prime}),B^{\prime}(p)(p^{\prime}-q^{\prime})\rangle\\\
&+\sigma\partial_{x}\sigma\langle\partial_{x}(p^{\prime}-q^{\prime}),(B^{\prime}(p)-B^{\prime}(q))q^{\prime}\rangle\\\
&\geq-C\|\partial_{x}(p^{\prime}-q^{\prime})\|_{2}\|p^{\prime}-q^{\prime}\|_{2}\\\
&-C\|q\|_{(C^{1+\alpha/2,2+\alpha})^{N}}\|\partial_{x}(p^{\prime}-q^{\prime})\|_{2}\|B^{\prime}(p)-B^{\prime}(q)\|_{2},\end{split}$

where we used the boundedness of $B^{\prime}(p)$, $\sigma$ and
$\partial_{x}\sigma$.

We check easily the existence of a constant $C$ depending on $\Gamma$,
$f_{\max}$, $f_{\min}$ and $\varepsilon$ such that

(2.58)
$\|A^{\prime}(p)-A^{\prime}(q)\|_{2}^{2}+\|B^{\prime}(p)-B^{\prime}(q)\|_{2}^{2}\leq
C\|p^{\prime}-q^{\prime}\|_{2}^{2}.$

In the view of the upper-bounds

(2.59)
$\|\partial_{x}(p^{\prime}-q^{\prime})\|_{2}\|A^{\prime}(p)-A^{\prime}(q)\|_{2}\leq\frac{\kappa}{4}\|\partial_{x}(p^{\prime}-q^{\prime})\|_{2}^{2}+\frac{4}{\kappa}\|A^{\prime}(p)-A^{\prime}(q)\|_{2}^{2},$
(2.60)
$\|\partial_{x}(p^{\prime}-q^{\prime})\|_{2}\|p^{\prime}-q^{\prime}\|_{2}\leq\frac{\kappa}{4}\|\partial_{x}(p^{\prime}-q^{\prime})\|_{2}^{2}+\frac{4}{\kappa}\|p^{\prime}-q^{\prime}\|_{2}^{2},$

and

(2.61)
$\|\partial_{x}(p^{\prime}-q^{\prime})\|_{2}\|B^{\prime}(p)-B^{\prime}(q)\|_{2}\leq\frac{\kappa}{4}\|\partial_{x}(p^{\prime}-q^{\prime})\|_{2}^{2}+\frac{4}{\kappa}\|B^{\prime}(p)-B^{\prime}(q)\|^{2}_{2},$

we conclude that for all $t\in[0,T]$ and some constant $C>0$

(2.62) $\int\|p^{\prime}-q^{\prime}\|^{2}_{2}\leq
C\|q\|_{(C^{1+\alpha/4,2+\alpha/2})^{N}}^{2}\int_{0}^{t}\int\|p^{\prime}-q^{\prime}\|^{2}_{2}dt.$

Gronwall’s inequality yields $p^{\prime}=q^{\prime}$ a.e. on
$[0,T]\times\mathbb{R}$ and by continuity of $p$ and $q$, we have immediately
that $p=q$.

∎

##  3\. Propagation of chaos

In this section, we prove the propagation of chaos and introduce an
intermediate mollified version of SDE (1.4). For each $M\geq 1$, let
$(\widetilde{X}^{i,M}_{t})_{i\geq 1,t\geq 0}$ be a particle system such that
the i-th particle starts at $X^{i}_{0}$ and evolves according to the dynamics

(3.1) $\begin{cases}\begin{aligned}
d\widetilde{X}^{i,M}_{t}=&-\frac{1}{2}\frac{\varepsilon+f(Y^{i})\sum_{n=1}^{N}W_{\delta_{M}}*\widetilde{p}_{n}(t,\widetilde{X}^{i,M}_{t})}{\varepsilon+\sum_{n=1}^{N}f(n)W_{\delta_{M}}*\widetilde{p}_{n}(t,\widetilde{X}^{i,M}_{t})}\sigma(t,\widetilde{X}^{i,M}_{t})^{2}dt\\\
&+\sqrt{\frac{\varepsilon+f(Y^{i})\sum_{n=1}^{N}W_{\delta_{M}}*\widetilde{p}_{n}(t,\widetilde{X}^{i,M}_{t})}{\varepsilon+\sum_{n=1}^{N}f(n)W_{\delta_{M}}*\widetilde{p}_{n}(t,\widetilde{X}^{i,M}_{t})}}\sigma(t,\widetilde{X}^{i,M}_{t})dB^{i}_{t},\end{aligned}\\\
\mathbb{P}\left(\widetilde{X}^{i,M}_{t}\in dx\cap
Y^{i}=n\right)=\widetilde{p}_{n}(t,x)dx,\end{cases},$

where we denote

(3.2) $W_{\delta_{M}}*\phi(x)=\int W_{\delta_{M}}(x-y)\phi(y)dy,\quad\phi\in
L^{2}(\mathbb{R}),\ x\in\mathbb{R}.$

We introduce the PDE system associated to SDE (3.1)

(3.3)
$\begin{cases}\partial_{t}\widetilde{p}=\frac{1}{2}\partial_{xx}\left[\sigma^{2}B^{\varepsilon}(W*\widetilde{p})\widetilde{p}\right]+\frac{1}{2}\partial_{x}\left[\sigma^{2}B^{\varepsilon}(W*\widetilde{p})\widetilde{p}\right],\
(t,x)\in[0,T]\times\mathbb{R},\\\ \tilde{p}_{n}(0,x)=P_{n}(x),\
x\in\mathbb{R},\\\ 1\leq n\leq N.\end{cases}.$

The existence of $(\widetilde{X}^{i,M}_{t})_{1\leq i\leq M,t\geq 0}$ is
ensured by the following proposition.

######  Proposition 3.1.

If $(P_{n})_{1\leq n\leq N}$ satisfies the assumption of Theorem 1.2, then
there exists a strong solution
$(\widetilde{X}^{i,M}_{t},\widetilde{p}(t,\cdot))_{t\geq 0}$ to SDE (1.11).
Moreover, $(\widetilde{p}_{n})_{1\leq n\leq N}$ belongs to the space
$(C^{1+\alpha/4,2+\alpha/2})^{N}\cap C(0,T,L^{2})^{N}$ and there exists a
constant $C>0$ independent of $M$ and $\delta_{M}$ such that

(3.4) $\sum_{n=1}^{N}\|\tilde{p}_{n}\|_{C^{1+\alpha/4,2+\alpha/2}}\leq C,$

and

(3.5) $\sup_{0\leq t\leq T}\int\|\widetilde{p}-\hat{p}\|_{2}^{2}\leq
C\delta_{M},\ \forall t\geq 0.$

###### Proof.

The existence of $(\tilde{p}_{n})_{1\leq n\leq N}$ and
$(\tilde{X}^{i,M}_{t})_{t\geq 0},1\leq i\leq M$, can easily be established by
noticing that the mapping $\phi\in C^{1+\alpha/4,2+\alpha/2}\mapsto
W_{\delta_{M}}*u\in C^{1+\alpha/4,2+\alpha/2}$ is non-expansive for the norm
$\|\cdot\|_{C^{1+\alpha/4,2+\alpha/2}}$ and by using the techniques of the
proof of Proposition 2.1 (see as well [9, Proposition 2.2]). The non-
expansivity of the mentioned mapping makes the estimate (2.7) still valid and
independent of $\delta_{M}$, and therefore (3.4) is justified.

In order to establish (3.5), we shall show that the proof of [9, Lemma 2.6]
carries to parabolic systems. Mimicking the computation of [9, Equation
(2.8)], we can write that the difference $\hat{p}-\widetilde{p}$ solves the
system

(3.6)
$\partial_{t}(\hat{p}-\widetilde{p})=\frac{1}{2}\partial_{x}[A^{\varepsilon}(\hat{p})\partial_{x}(\hat{p}-\widetilde{p})]+a_{1}\partial_{x}(\hat{p}-\widetilde{p})+a_{2}(\hat{p}-\widetilde{p})+\phi.$

where the matrix valued functions $a_{1}(t,x)$ and $a_{2}(t,x)$ are bounded,
depend on $\hat{p}$ and $\tilde{p}$, and the vector valued function
$\phi(t,x)$ depends on $\widetilde{p}-W*\widetilde{p}$ and its first and
second derivatives. Moreover, $\phi$ satisfies

(3.7) $\sum_{n=1}^{N}\sup_{x\in\mathbb{R}}|\phi_{n}(x)|\leq C\delta_{M},$

for some constant $C>0$ independent of $\delta_{M}$.

Next, we observe that under Condition 1.1, there exist $\Gamma\in S^{++}$ and
$\kappa>0$ such that for all $X\in\mathbb{R}^{N}$

(3.8) $\langle X,\Gamma
A^{\varepsilon}(\hat{p})X\rangle\geq\kappa\|X\|_{2}^{2}.$

Like in the proof of Proposition 2.4, set
$\hat{p}^{\prime}=\sqrt{\Gamma}\hat{p}$,
$\widetilde{p}^{\prime}=\sqrt{\Gamma}\widetilde{p}$,
$A^{\prime}(\hat{p})=\sqrt{\Gamma}A^{\varepsilon}(\hat{p})\sqrt{\Gamma}^{-1}$,
$a_{1}^{\prime}=\sqrt{\Gamma}a_{1}\sqrt{\Gamma}^{-1}$,
$a_{2}^{\prime}=\sqrt{\Gamma}a_{2}\sqrt{\Gamma}^{-1}$ and
$\phi^{\prime}=\sqrt{\Gamma}\phi$. With (2.53), there exists
$\kappa^{\prime}>0$ such that for all $X\in\mathbb{R}^{N}$

(3.9) $\langle X,A^{\prime}X\rangle\geq\kappa^{\prime}\|X\|^{2}_{2}.$

Moreover

(3.10)
$\begin{split}\frac{1}{2}\int\|\hat{p}^{\prime}-\widetilde{p}^{\prime}\|_{2}^{2}=&-\frac{1}{2}\int_{0}^{t}\int\langle\partial_{x}(\hat{p}^{\prime}-\widetilde{p}^{\prime}),A^{\prime}(\hat{p})\partial_{x}(\hat{p}^{\prime}-\widetilde{p}^{\prime})\rangle
dt\\\
&+\int_{0}^{t}\int\langle\hat{p}^{\prime}-\widetilde{p}^{\prime},a_{1}^{\prime}\partial_{x}(\hat{p}^{\prime}-\widetilde{p}^{\prime})\rangle
dt+\int_{0}^{t}\int\langle\hat{p}^{\prime}-\widetilde{p}^{\prime},a_{2}(\hat{p}^{\prime}-\widetilde{p}^{\prime})\rangle
dt\\\
&+\int_{0}^{t}\int\langle\hat{p}^{\prime}-\widetilde{p}^{\prime},\phi^{\prime}\rangle
dt.\end{split}$

Using the boundedness of $a_{1}^{\prime}$, $a_{2}^{\prime}$ and standard
techniques, it is easy to see that for some appropriate constant $C>0$ and for
all $t\in[0,T]$

(3.11)
$\begin{split}&\int\langle\hat{p}^{\prime}-\widetilde{p}^{\prime},a_{1}^{\prime}\partial_{x}(\hat{p}^{\prime}-\widetilde{p}^{\prime})\rangle+\int\langle\hat{p}^{\prime}-\widetilde{p}^{\prime},a_{2}(\hat{p}^{\prime}-\widetilde{p}^{\prime})\rangle\\\
&\leq\frac{\kappa^{\prime}}{2}\int\|\partial_{x}(\hat{p}^{\prime}-\tilde{p}^{\prime})\|^{2}_{2}+C\int\|\hat{p}^{\prime}-\tilde{p}^{\prime}\|^{2}_{2}.\end{split}$

In the view of the uniform coercivity of $A^{\prime}$, and $\phi^{\prime}$,
and the estimate

(3.12)
$\int\langle\hat{p}^{\prime}-\widetilde{p}^{\prime},\phi^{\prime}\rangle\leq
C\delta_{M}\sum_{n=1}^{N}\int(\hat{p}_{n}+\widetilde{p}_{n})=C\delta_{M},$

we conclude that there exists some constant $C>0$ such that for all
$t\in[0,T]$

(3.13) $\frac{1}{2}\int\|\hat{p}^{\prime}-\widetilde{p}^{\prime}\|_{2}^{2}\leq
C\delta_{M}+C\int_{0}^{t}\frac{1}{2}\int\|\hat{p}^{\prime}-\widetilde{p}^{\prime}\|_{2}^{2}dt.$

The result is proved by applying Gronwall’s inequality and using the obvious
inequality
$\|\hat{p}-\widetilde{p}\|_{2}\leq\|\hat{p}^{\prime}-\widetilde{p}^{\prime}\|_{2}$,
for some constant $C>0$ depending only on $\Gamma$ and independent of
$\delta_{M}$. ∎

Let $1\leq i\leq M$. To simplify the notation, we will drop the subscripts and
superscripts in $i$ and $M$, and we denote $\delta=\delta_{M}$,
$W=W_{\delta_{M}}$, $X=X^{i}$, $\widetilde{X}=\widetilde{X}^{i}$,
$\hat{X}=\hat{X}^{i}$, and for $j\neq i$, $X^{j}=X^{j,M}$ and
$\widetilde{X}^{j}=\widetilde{X}^{j,M}$. Denote
$W_{\delta_{M}}*X_{t}(x)=\frac{1}{M}\sum_{j=1}^{M}W_{\delta_{M}}(x-X^{j}_{t})$
and
$W_{\delta_{M}}^{f}*X_{t}(x)=\frac{1}{M}\sum_{j=1}^{M}f(Y^{j})W_{\delta_{M}}(x-X^{j}_{t})$.
Under this notation, the volatility coefficients of $X$, $\widetilde{X}$ and
$\hat{X}$ are respectively given by

(3.14)
$\Sigma_{t}:=f(Y)\sqrt{\frac{\varepsilon+W*X_{t}(X_{t})}{\varepsilon+W^{f}*X_{t}(X_{t})}}\sigma(t,X_{t}),$
(3.15)
$\widetilde{\Sigma}_{t}:=\sqrt{\frac{\varepsilon+f(Y)\sum_{n=1}^{N}W*\widetilde{p}_{n}(t,\widetilde{X}_{t})}{\varepsilon+\sum_{n=1}^{N}f(n)W*\widetilde{p}_{n}(t,\widetilde{X}_{t})}}\sigma(t,\widetilde{X}_{t}),$

and

(3.16)
$\hat{\Sigma}_{t}:=\sqrt{\frac{\varepsilon+f(Y)\sum_{n=1}^{N}\hat{p}_{n}(t,\hat{X}_{t})}{\varepsilon+\sum_{n=1}^{N}f(n)\hat{p}_{n}(t,\widetilde{X}_{t})}}\sigma(t,\hat{X}_{t}).$

Notice immediately that

(3.17) $\sigma_{0}^{2}\frac{\inf f}{\sup
f}\leq\Sigma^{2},\widetilde{\Sigma}^{2},\hat{\Sigma}^{2}\leq\frac{\sup f}{\inf
f}\sigma_{1}^{2}.$

In the rest of the section, we denote, for notational simplicity, by $C$ any
constant depending on $N$, $\varepsilon$, $f$, $\sigma_{0}$, $\sigma_{1}$, and
$|W_{1}|_{\infty}$.

The proof of Theorem 1.3 requires the Propositions 3.2 and 3.3 below.

######  Proposition 3.2.

Let $1\leq i\leq M$. Then for all $T>0$.

(3.18) $\mathbb{E}\left[\sup_{0\leq t\leq
T}(X^{i,M}_{t}-\widetilde{X}^{i,M}_{t})^{2}\right]\leq\frac{C}{M}\frac{\exp{(CT\delta_{M}^{-4})}}{\delta_{M}^{2}}.$

###### Proof.

According to the Burkholder–Davis–Gundy inequality, for any $s\in[0,T]$

(3.19) $\begin{split}\mathbb{E}\left[\sup_{0\leq t\leq
s}(X_{t}-\widetilde{X}_{t})^{2}\right]&\leq C\mathbb{E}\left[\langle
X-\widetilde{X}\rangle_{s}\right]\\\
&=C\int_{0}^{s}\mathbb{E}\left[(\Sigma_{t}-\widetilde{\Sigma}_{t})^{2}\right]dt.\end{split}$

For $x\in\mathbb{R}$ and $(U^{j})_{1\leq j\leq M}\in\mathbb{R}^{M}$, we
abbreviate $\frac{1}{M}\sum_{j=1}^{M}W(x-U^{j})$ by $W*U(x)$ and
$\frac{1}{M}\sum_{j=1}^{M}f(Y^{j})W(x-U^{j})$ by $W^{f}*U(x)$.

The boundedness $\Sigma_{t}$, $\widetilde{\Sigma}_{t}$, $\sigma$ and
$\partial_{x}\sigma$ gives for all $t\in[0,T]$

(3.20) $\begin{split}|\Sigma_{t}-\widetilde{\Sigma}_{t}|&\leq
C|\Sigma_{t}^{2}-\widetilde{\Sigma}_{t}^{2}|\\\ &\leq
C\left|\frac{\varepsilon+W^{f}*X_{t}(X_{t})}{\varepsilon+W*X_{t}(X_{t})}-\frac{\varepsilon+\sum_{n=1}^{N}f(n)W*\widetilde{p}_{n}(t,\widetilde{X}_{t})}{\varepsilon+\sum_{n=1}^{N}W*\widetilde{p}_{n}(t,\widetilde{X}_{t})}\right|\\\
&+C|\sigma(t,X_{t})-\sigma(t,\widetilde{X}_{t})|\\\ &\leq
C\left|\frac{\varepsilon+W^{f}*X(X_{t})}{\varepsilon+W*X_{t}(X_{t})}-\frac{\varepsilon+W^{f}*X_{t}(X_{t})}{\varepsilon+\sum_{n=1}^{N}W*\widetilde{p}_{n}(t,\widetilde{X}_{t})}\right|\\\
&+C\left|\frac{\varepsilon+W^{f}*X(X_{t})}{\varepsilon+\sum_{n=1}^{N}W*\widetilde{p}_{n}(t,\widetilde{X}_{t})}-\frac{\varepsilon+\sum_{n=1}^{N}f(n)W*\widetilde{p}_{n}(t,\widetilde{X}_{t})}{\varepsilon+\sum_{n=1}^{N}W*\widetilde{p}_{n}(t,\widetilde{X}_{t})}\right|\\\
&+C|X_{t}-\widetilde{X}_{t}|\\\ &\leq
C\left|W*X_{t}(X_{t})-\sum_{n=1}^{N}W*\widetilde{p}_{n}(t,\widetilde{X}_{t})\right|\\\
&+C\left|W^{f}*X_{t}(X_{t})-\sum_{n=1}^{N}f(n)W*\widetilde{p}_{n}(t,\widetilde{X}_{t})\right|+C|X_{t}-\widetilde{X}_{t}|.\end{split}$

Let $t\in[0,T]$. After estimating the term
$\mathbb{E}\left[\left|W^{f}*X_{t}(X_{t})-\sum_{n=1}^{N}f(n)W*\widetilde{p}_{n}(t,\widetilde{X}_{t})\right|^{2}\right]$,
we will apply the estimate with $f=1$ to treat the term
$\mathbb{E}\left[\left|W*X_{t}(X_{t})-\sum_{n=1}^{N}W*\widetilde{p}_{n}(t,\widetilde{X}_{t})\right|^{2}\right]$.

Upper-bound

(3.21)
$\begin{split}\left|W^{f}*X_{t}(X_{t})-\sum_{n=1}^{N}f(n)W*p_{n}(t,\widetilde{X}_{t})\right|&\leq(\textrm{I})+(\textrm{II})+(\textrm{III}),\end{split}$

where each term $(\textrm{I}),(\textrm{II})$ and $(\textrm{III})$ is
respectively defined by

(3.22)
$(\textrm{I}):=\left|W^{f}*X_{t}(X_{t})-W^{f}*X_{t}(\widetilde{X}_{t})\right|,$
(3.23)
$(\textrm{II}):=\left|W^{f}*X_{t}(\widetilde{X}_{t})-W^{f}*\widetilde{X}_{t}(\widetilde{X}_{t})\right|,$

and

(3.24)
$(\textrm{III}):=\left|W^{f}*\widetilde{X}_{t}(\widetilde{X}_{t})-\sum_{n=1}^{N}f(n)W*\widetilde{p}_{n}(t,\widetilde{X}_{t})\right|.$

By using $\sup_{x\in\mathbb{R}}|W^{\prime}(x)|\leq\frac{C}{\delta^{2}}$, we
find

(3.25)
$\begin{split}\mathbb{E}\left[(\textrm{I})^{2}\right]&=\mathbb{E}\left[\left(\frac{1}{M}\sum_{j=1}^{M}\mathbb{E}\left[f(Y^{j})(W(X_{t}-X^{j}_{t})-W(\widetilde{X}_{t}-X^{j}_{t}))\right]\right)^{2}\right]\\\
&\leq\frac{C}{\delta^{4}}\mathbb{E}\left[|X_{t}-\widetilde{X}_{t}|^{2}\right].\end{split}$

Likewise

(3.26)
$\begin{split}\mathbb{E}\left[(\textrm{II})^{2}\right]&=\mathbb{E}\left[\left(\frac{1}{M}\sum_{j=1}^{M}f(Y^{j})(W(\widetilde{X}_{t}-X^{j}_{t})-W(\widetilde{X}_{t}-\widetilde{X}^{j}_{t}))\right)^{2}\right]\\\
&\leq\frac{C}{\delta^{4}M}\mathbb{E}\left[\left(\sum_{j=1}^{M}|X^{j}_{t}-\widetilde{X}^{j}_{t}|^{2}\right)\right]\\\
&\leq\frac{C}{\delta^{4}}\mathbb{E}\left[|X_{t}-\widetilde{X}_{t}|^{2}\right],\end{split}$

where we used that by symmetry
$\mathbb{E}\left[|X^{j}_{t}-\widetilde{X}^{j}_{t}|^{2}\right]=\mathbb{E}\left[|X_{t}-\widetilde{X}_{t}|^{2}\right]$
for all $j=1,\dots,M$.

By independence of the $(\widetilde{X}^{j}_{t})_{1\leq j\leq M}$ and using
that

(3.27)
$\sum_{n=1}^{N}f(n)W*\widetilde{p}_{n}(t,x)=\mathbb{E}\left[f(Y^{j})W(x-\widetilde{X}^{j}_{t})\right],\
1\leq j\neq i\leq N,$

we get

(3.28)
$\begin{array}[]{c}\mathbb{E}\left[(\textrm{III})^{2}\right]=\displaystyle\int\mathbb{E}\left[\left(\frac{1}{M}\sum_{j\neq
i}f(Y^{j})W(x-\widetilde{X}^{j}_{t})-\mathbb{E}\left[f(Y^{j})W(x-\widetilde{X}^{j}_{t})\right]\right)^{2}\right]\\\
\times\mathbb{P}\left(\widetilde{X}_{t}\in dx\right).\end{array}$

Therefore

(3.29) $\begin{split}\mathbb{E}\left[(\textrm{III})^{2}\right]&\leq
C\int\mathbb{E}\left[\left(\frac{1}{M}\sum_{j\neq
i}\left(f(Y^{j})W(x-\widetilde{X}^{j}_{t})-\mathbb{E}\left[f(Y^{j})W(x-\widetilde{X}^{j}_{t})\right]\right)\right)^{2}\right]\mathbb{P}\left(\widetilde{X}_{t}\in
dx\right)\\\
&+\frac{C}{M^{2}}\int\mathbb{E}\left[f(Y)W(x-\widetilde{X}_{t})\right]^{2}\mathbb{P}\left(\widetilde{X}_{t}\in
dx\right).\end{split}$

Making use of the bound $\sup_{x\in\mathbb{R}}|W(x)|\leq\frac{C}{\delta}$, we
see that the second integral in the right-hand side of (3.29) is dominated by

(3.30)
$\begin{split}\int\mathbb{E}\left[f(Y)W(x-\widetilde{X}_{t})\right]^{2}\mathbb{P}\left(\widetilde{X}_{t}\in
dx\right)&\leq C\int W(x-y)^{2}\mathbb{P}\left(\widetilde{X}_{t}\in
dx\right)\mathbb{P}\left(\widetilde{X}_{t}\in dy\right)\\\
&\leq\frac{C}{\delta^{2}}.\end{split}$

In estimating the first integral of (3.29), we utilize the independence of the
$(\tilde{X}^{j}_{t})_{1\leq j\leq N}$, for all $x\in\mathbb{R}$, which gives

(3.31) $\begin{split}&\mathbb{E}\left[\left(\frac{1}{M}\sum_{j\neq
i}\left(f(Y^{j})W(x-\widetilde{X}^{j}_{t})-\mathbb{E}\left[f(Y^{j})W(x-\widetilde{X}^{j}_{t})\right]\right)\right)^{2}\right]\\\
&=\textrm{Var}\left(\frac{1}{M}\sum_{j\neq
i}f(Y^{j})W(x-\widetilde{X}^{j}_{t})\right)\\\
&=\frac{M-1}{M^{2}}\textrm{Var}\left(f(Y)W(x-\widetilde{X}_{t})\right)\\\
&\leq\frac{C}{M\delta^{2}}.\end{split}$

Putting together (3.25), (3.26) and (3.29)-(3.31), we infer

(3.32)
$\mathbb{E}\left[\left|W^{f}*X_{t}(X_{t})-\sum_{n=1}^{N}f(n)W*\widetilde{p}_{n}(t,\widetilde{X}_{t})\right|^{2}\right]\leq\frac{C}{M\delta^{2}}+\frac{C}{\delta^{4}}\mathbb{E}\left[|X_{t}-\widetilde{X}_{t}|^{2}\right].$

Taking $f=1$ leads to the same estimate for
$\mathbb{E}\left[\left|W*X_{t}(X_{t})-\sum_{n=1}^{N}W*\widetilde{p}_{n}(t,\widetilde{X}_{t})\right|^{2}\right]$.

If we insert this in (3.20) and use (3.19), we obtain

(3.33) $\mathbb{E}\left[\sup_{0\leq s\leq
t}(X_{s}-\widetilde{X}_{s})^{2}\right]\leq\frac{1}{M}\frac{C}{\delta^{2}}+\frac{C}{\delta^{4}}\int_{0}^{t}\mathbb{E}\left[\sup_{0\leq
s\leq r}(X_{s}-\widetilde{X}_{s})^{2}\right]dr.$

Applying Gronwall’s inequality proves (3.18) and concludes the proof of the
proposition. ∎

######  Proposition 3.3.

Let $1\leq i\leq M$. Then for some constant $C>0$ independent of $\delta$

(3.34) $\mathbb{E}\left[\sup_{0\leq t\leq
T}|\widetilde{X}^{i,M}_{t}-\hat{X}^{i}_{t}|^{2}\right]\leq C\delta_{M}.$

###### Proof.

In the same manner of (3.20), for some constant $C>0$

(3.35)
$|\widetilde{\Sigma}_{t}-\hat{\Sigma}_{t}|\leq\frac{C}{\delta}\sum_{n=1}^{N}\left|W*\widetilde{p}_{n}(t,\widetilde{X}_{t})-\hat{p}_{n}(t,\hat{X}_{t})\right|+C|\widetilde{X}_{t}-\hat{X}_{t}|.$

For all $t\in[0,T]$

(3.36)
$\begin{split}\left|W*\widetilde{p}_{n}(t,\widetilde{X}_{t})-\hat{p}_{n}(t,\hat{X}_{t})\right|&\leq\left|W*\hat{p}_{n}(t,\hat{X}_{t})-\hat{p}_{n}(t,\hat{X}_{t})\right|\\\
&+\left|W*\hat{p}_{n}(t,\widetilde{X}_{t})-W*\hat{p}_{n}(t,\hat{X}_{t})\right|\\\
&+\left|W*\widetilde{p}_{n}(t,\widetilde{X}_{t})-W*\hat{p}_{n}(t,\widetilde{X}_{t})\right|.\end{split}$

Next, using the boundedness of the gradient of $\hat{p}_{n}$, we estimate the
first two terms of the left-hand side of (3.36) as follows

(3.37)
$\begin{split}\left|W*\hat{p}_{n}(t,\hat{X}_{t})-\hat{p}_{n}(t,\hat{X}_{t})\right|&\leq\sup_{x\in\mathbb{R}}\int
W_{1}(y)|\hat{p}_{n}(t,x-\delta y)-\hat{p}_{n}(t,x)|dy\\\
&\leq\delta\|\hat{p}_{n}\|_{C^{1+\alpha/2,2+\alpha}}\int|y|W_{1}(y)dy,\end{split}$

and

(3.38)
$\left|W*\hat{p}_{n}(t,\widetilde{X}_{t})-W*\hat{p}_{n}(t,\hat{X}_{t})\right|\leq\|\hat{p}_{n}\|_{C^{1+\alpha/2,2+\alpha}}|\widetilde{X}_{t}-\hat{X}_{t}|.$

Therefore, there exists a constant $C>0$ independent of $\delta$ such that

(3.39)
$\left|W*\hat{p}_{n}(t,\hat{X}_{t})-\hat{p}_{n}(t,\hat{X}_{t})\right|+\left|W*\hat{p}_{n}(t,\widetilde{X}_{t})-W*\hat{p}_{n}(t,\hat{X}_{t})\right|\leq
C(\delta+|\widetilde{X}_{t}-\hat{X}_{t}|).$

For the last term of (3.36), write

(3.40)
$\begin{split}\mathbb{E}\left[\left|W*\widetilde{p}_{n}(t,\widetilde{X}_{t})-W*\hat{p}_{n}(t,\widetilde{X}_{t})\right|^{2}\right]&=\sum_{m=1}^{N}\int|W*(\widetilde{p}_{n}-\hat{p}_{n})|^{2}p_{m}\\\
&\leq\|p\|_{(C^{1+\alpha/4,2+\alpha/2})^{N}}\int|W*(\widetilde{p}_{n}-\hat{p}_{n})|^{2}\\\
&\leq\|p\|_{(C^{1+\alpha/4,2+\alpha/2})^{N}}\int|\widetilde{p}_{n}-\hat{p}_{n}|^{2}\\\
&\leq\|p\|_{(C^{1+\alpha/4,2+\alpha/2})^{N}}\delta,\end{split}$

where we used Proposition 3.1.

Combining the previous estimate and (3.36) we deduce the existence of $C>0$
independent of $\delta$ such that

(3.41)
$\mathbb{E}\left[|\widetilde{\sigma}_{t}-\hat{\sigma}_{t}|^{2}\right]\leq
C\left(\delta^{2}+\mathbb{E}\left[|\widetilde{X}_{t}-\hat{X}_{t}|^{2}\right]\right).$

Burkholder–Davis–Gundy inequality and Gronwall’s inequality yield

(3.42) $\mathbb{E}\left[\sup_{0\leq t\leq
T}|\widetilde{X}_{t}-\hat{X}_{t}|^{2}\right]\leq C\delta.$

∎

We are now ready for the proof of Theorem 1.3.

###### Proof of Theorem 1.3.

According to Propositions 3.2 and 3.3 above, there exists a constant $C>0$
such that for all $1\leq i\leq N$

(3.43) $\mathbb{E}\left[\sup_{0\leq t\leq
T}\left(X^{i,M}_{t}-\hat{X}^{i,M}_{t}\right)^{2}\right]\leq\frac{C}{M}\frac{\exp{(CT\delta_{M}^{-4})}}{\delta_{M}^{2}}+C\delta_{M}.$

It is enough to choose $(\delta_{M})_{M\geq 1}$ such that
$\lim_{M\to\infty}\delta_{M}=0$ and

(3.44)
$\lim_{M\to\infty}\frac{1}{M}\frac{\exp{(CT\delta_{M}^{-4})}}{\delta_{M}^{2}}=0.$

∎

## References

  * [1] Frédéric Abergel, Rémy Tachet Des Combes, Riadh Zaatour, et al.  Nonparametric model calibration for derivatives.  Journal of Mathematical Finance, 7(03):571, 2017.
  * [2] Frédéric Abergel and Rémi Tachet.  A nonlinear partial integro-differential equation from mathematical finance.  Discrete and Continuous Dynamical Systems-Series A, 27(3):907–917, 2010.
  * [3] Mireille Bossy and Jean-François Jabir.  On the wellposedness of some mckean models with moderated or singular diffusion coefficient.  In Frontiers in Stochastic Analysis–BSDEs, SPDEs and their Applications: Edinburgh, July 2017 Selected, Revised and Extended Contributions 8, pages 43–87. Springer, 2019.
  * [4] Bruno Dupire et al.  Pricing with a smile.  Risk, 7(1):18–20, 1994.
  * [5] Julien Guyon and Pierre Henry-Labordère.  Being particular about calibration.  Risk, 25(1):88, 2012.
  * [6] Julien Guyon and Pierre Henry-Labordere.  Nonlinear option pricing.  CRC Press, 2013.
  * [7] István Gyöngy.  Mimicking the one-dimensional marginal distributions of processes having an itô differential.  Probability theory and related fields, 71(4):501–516, 1986.
  * [8] Pierre Henry-Labordère.  Calibration of local stochastic volatility models to market smiles: A monte-carlo approach.  Risk Magazine, September, 2009.
  * [9] Benjamin Jourdain and Sylvie Méléard.  Propagation of chaos and fluctuations for a moderate model with smooth initial data.  In Annales de l’Institut Henri Poincare (B) Probability and Statistics, volume 34, pages 727–766. Elsevier, 1998.
  * [10] Benjamin Jourdain and Alexandre Zhou.  Existence of a calibrated regime switching local volatility model.  Mathematical Finance, 30(2):501–546, 2020.
  * [11] Arturo Kohatsu-Higa and Shigeyoshi Ogawa.  Weak rate of convergence for an euler scheme of nonlinear sde’s.  1997\.
  * [12] Daniel Lacker, Mykhaylo Shkolnikov, and Jiacheng Zhang.  Inverting the markovian projection, with an application to local stochastic volatility models.  2020\.
  * [13] Olga Aleksandrovna Ladyzhenskaya.  Linear and quasilinear elliptic equations.  translated by Scripta Technica, 1968.
  * [14] Karl Oelschläger.  A law of large numbers for moderately interacting diffusion processes.  Zeitschrift für Wahrscheinlichkeitstheorie und verwandte Gebiete, 69(2):279–322, 1985.
  * [15] Yuri F Saporito, Xu Yang, and Jorge P Zubelli.  The calibration of stochastic local-volatility models: An inverse problem perspective.  Computers & Mathematics with Applications, 77(12):3054–3067, 2019.


# Optical conductivity of a topological system driven using a realistic pulse

Ranjani Seshadri1,2<EMAIL_ADDRESS>T. Pereg-Barnea2
<EMAIL_ADDRESS>1Department of Physics, Ben-Gurion University of the
Negev, Beer-Sheva 84105, Israel  2Department of Physics, McGill University,
Montréal,Québec H3A 2T8, Canada

###### Abstract

The effect of a time-periodic perturbation, such as radiation, on a system
otherwise at equilibrium has been studied in the context of Floquet theory
with stationary states replaced by Floquet states and the energy replaced by
quasienergy. These quasienergy bands in general differ from the energy bands
in their dispersion and, especially in the presence of spin-orbit coupling, in
their states. This may, in some cases, alter the topology when the quasienergy
bands exhibit different topological invariants than their stationary
counterparts. In this work, motivated by advances in pump-probe techniques, we
consider the optical response of driven topological systems when the drive is
not purely periodic but is instead multiplied by a pulse shape/envelope
function. We use real time-evolved states to calculate the optical
conductivity and compare it to the response calculated using Floquet theory.
We find that the conductivity bears a memory of the initial equilibrium state
even when the pump is turned on slowly and the measurement is taken well after
the ramp. The response of the time-evolved system is interpreted as coming
from Floquet bands whose population has been determined by their overlap with
the initial equilibrium state. In particular, at band inversion points in the
Brillouin zone the population of the Floquet bands is inverted as well.

Introduction- The theoretical prediction and experimental realization of
topological insulators (TIs) [1, 2, 3, 4, 5, 6, 7, 8, 9] has been one of the
greatest developments in condensed matter physics in the last decade. Not only
do topological insulators represent a paradigm shift in condensed matter
physics, they are also predicted to have a variety of applications [10, 11,
12, 13, 14].

While spin-orbit coupling is a key ingredient, it need not always lead to non-
trivial topology as band inversion may not always occur or the presence of a
Fermi surface may not be avoided. It has therefore been proposed to use a
time-periodic perturbation in order to control the topology [15, 16, 17, 18,
19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37,
38]. When a time-periodic perturbation is added to a Hamiltonian the system is
no longer invariant under an arbitrary translation in time. However a reduced
discrete time-translation symmetry still exists. This allows finding solutions
to the time-dependent Schrödinger equation using Floquet theory. These
solutions or Floquet states are eigenstates of the time evolution operator
over a single drive cycle. In other words, Floquet states are periodic up to a
phase which is interpreted as $-\epsilon T$ where $\epsilon$ is the
quasienergy and $T$ is the drive period. The quasienergies and Floquet states
(say, spinors) in general differ from the equilibrium energies and eigenstates
respectively. Therefore, along with other properties, topological invariants
can change as a result of irradiation leading to driven topological phase
transitions. For example, graphene may be driven into a topological phase [31]
where gaps appear in the spectrum. In case of Weyl semimetals, Weyl nodes may
split into Dirac points or gap out and give rise to Chern bands [39] while
spin-orbit coupled insulators are predicted to become topological upon diving
[15].

While there are several theoretical predictions, the experimental realization
of such Floquet-driven topological phase transitions seems to be challenging.
Notably, analogue photonic systems were the first to realize some of these
predictions [40] and recently driven graphene has shown signs of topology [41]
while the general idea of Floquet bands has been demonstrated by time-resolved
ARPES [42, 43].

Several obstacles occur while trying to realize such Floquet-driven
topological transitions. These include sample heating, damping and disorder.
But perhaps the most elementary deviation from the pure Floquet drive is the
unavoidable pulse shape. The drive can not be turned on at time $t=-\infty$
and therefore the state of the system is always connected to that of the
equilibrium state. One might expect that at long times after the turning on of
the drive, the state will resemble a Floquet state. However, as will be shown
below, the notion of adiabaticity does not hold at relevant drive frequencies.
In particular, as will be discussed here, the system does not forget its
initial conditions. While Floquet states may be a good approximation for the
single particle states at long time after the perturbation has been turned on,
their population is highly dependent on the initial state. Therefore, one
should not expect completely filled or empty Floquet bands at low temperature
meaning that the full potential of topological invariance may not be realized.
Moreover, it seems that relaxation effects do not necessarily lead to the
desired population as the quasienergy is periodic and energy may not
necessarily relax to one band [44]. Similarly, when connecting a Floquet spin-
Hall insulator to leads, one can not measure quantized conductivity due to
mismatch between the equilibrium states of the leads and the driven system
[45, 46].

The task at hand is therefore to accurately describe a system driven by a
pulse of light whose width is in the range of a few to many time periods, as
appropriate for pump-probe techniques and understand the relation between the
optical conductivity and Floquet band population. In this work we are
interested in the physics of such pump-probe measurements and how topological
systems respond to a perturbation that breaks time periodicity. We look at the
behavior of the Bernevig-Hughes-Zhang (BHZ) model of a two-dimensional TI [6]
in presence of a perturbation in the form of a short pulse and compare it with
the response of an exactly periodic (Floquet) drive.

Driven BHZ Model-  The equilibrium Hamiltonian in momentum-space is written as

$H_{0}({\bf k})={\bf d({\bf k})}\cdot{\boldsymbol{\sigma}}$ (1)

with ${\bf d}({\bf k})=(A\sin k_{x},A\sin k_{y},M-2B(2-\cos k_{x}-\cos
k_{y}))$ and ${\boldsymbol{\sigma}}=(\sigma^{x},\sigma^{y},\sigma^{z})$ are
the $2\times 2$ Pauli matrices. The mass $M$ and hopping amplitude $B$ are
expressed in units of the spin-orbit coupling strength $A$. The spectrum, in
general, is insulating in the bulk with a finite band gap. We work in a
parameter regime ($M=0.2A$ and $B=0.2A$) where the system is topological. The
equilibrium Chern numbers are calculated numerically as explained in the
Supplemental material [47] following the method prescribed in Ref [48] and are
found to be $C^{Eq}_{\pm}=\pm 1$ for the top and bottom band respectively.

A time-dependence $H(t)=H_{0}+V(t)\sigma^{z}$ effectively makes the mass time-
dependent. In the ideal case, $V(t)$ is perfectly periodic with a frequency
$\Omega$. However, in reality, this perfect periodicity cannot exist forever
and is instead approximated by realistic cases of a slow quench or a Gaussian
pulse,

$\displaystyle V(t)=\begin{cases}V_{{Floq}}(t)=V_{0}\sin(\Omega t)\\\
V_{{Pump}}(t)=V_{0}\sin(\Omega t)e^{-\frac{t^{2}}{2\Delta^{2}}}\\\
V_{{Quench}}(t)=V_{0}\sin(\Omega t)\frac{1+\tanh(\beta t)}{2}\end{cases}$ (2)

where $V_{{Floq}},V_{{Pump}}$ and $V_{{Quench}}$ have different envelope
functions - constant, Gaussian and a smooth ramp respectively. Here $\beta$
and $\Delta$ are the rate of the quench and the width of the pulse,
respectively. $V_{0}$ ($=0.35A$) is the peak amplitude of the perturbation in
all cases. The $V_{Floq}(t)$ with frequency $\Omega=2A$ drives the system to a
topological phase with Chern numbers $C^{Fl}_{\pm}=\pm 2$.

Figure 1: (a) The gapped equilibrium spectrum with Chern number $\pm 1$. (b)
With a periodic drive of frequency $\Omega=2A$ the ideal Floquet spectrum is
also gapped with Chern numbers $\pm 2$. The transitions marked in (b)
correspond to features marked as ${*}$s in Fig. 2. (c) The drive protocols of
in Eq. (2) - the ideal Floquet case (red dotted), the slow quench (black dot
dash) and the Gaussian pump (green solid) .

In the perfectly periodic case we employ Floquet theory to find the wave
functions which are then used to calculate the response functions, assuming
that one of the Floquet bands is completely filled while the other is
completely empty. In the cases of quench and pump, the analysis requires an
actual time evolution over several drive cycles since the perfect periodicity
is lost due to the envelope. The response to a non-periodic drive (Gaussian or
quench) is calculated using the time-evolved states starting from the
equilibrium states of the lower band of the undriven system. We use these
states in the Kubo formula as described below.

Linear Response Theory -  According to Kubo formula the susceptibility, which
is in our case is the response of the driven system to a small probe field is
given by,

$\displaystyle\chi_{AB}(t,t^{\prime})=\lim_{\eta\rightarrow 0^{+}}e^{\eta
t^{\prime}}\Big{(}$ $\displaystyle i\Theta(t-t^{\prime})$ $\displaystyle
Tr\Big{\\{}g_{0}[A^{I}(t^{\prime}),B^{I}(t)]\Big{\\}}$ (3)
$\displaystyle+\delta(t-t^{\prime})$ $\displaystyle
Tr\Big{\\{}g_{0}M^{I}(t)\Big{\\}}\Big{)}$

where $A^{I}$, $B^{I}$ and $M^{I}$ are operators in the interaction
representation. The density matrix $g_{0}$ determines the initial state of the
system and $\eta>0$ is a small parameter used to smoothen the response
function. The Heaviside step-function $\Theta(t-t^{\prime})$ ensures that
causality is not violated. The diamagnetic term $M^{I}(t)$ contributes only to
the DC conductivity in the limit $\omega\rightarrow 0$.

For computing electrical conductivity, both $A$ and $B$ are current operators.
As explained in the Supplemental material [47], following reference [49], Eq.
(3) becomes,

$\displaystyle\chi_{uv}(t,t^{\prime})$ $\displaystyle=\lim_{\eta\rightarrow
0^{+}}e^{\eta
t^{\prime}}\sum_{\alpha\gamma}g_{0\alpha}\Big{[}2i\Theta(t-t^{\prime})$ (4)
$\displaystyle\times$
$\displaystyle\Big{(}\langle\psi_{\alpha}(t^{\prime})|J_{u}|\psi_{\gamma}(t^{\prime})\rangle\langle\psi_{\gamma}(t)|J_{v}|\psi_{\alpha}(t)\rangle-{u\leftrightarrow
v}\Big{)}$ $\displaystyle+$
$\displaystyle\delta(t-t^{\prime})\langle\psi_{\alpha}(t)|M_{uv}|\psi_{\alpha}(t)\rangle\Big{]}$

Here $|\psi_{\alpha}(t))\rangle$ is the state corresponding to band $\alpha$
at time $t$, $g_{0\alpha}$ gives the occupation of states at the initial time;
the current operator $J_{u}=\partial_{k_{u}}H$ and the inverse mass
$M_{uv}=\partial_{k_{u}}\partial_{k_{v}}H$. The subscripts $u$ and $v$ are the
in-plane directions with $u=(\neq)v$ being longitudinal (transverse)
conductivity. and the ${\bf k}$-dependence has been skipped for brevity. In
the specific case of the model we have considered, the diamagnetic term
contributes only to the longitudinal conductivity, since $M_{uv}=0$
identically when $u\neq v$. To obtain the frequency response, we Fourier
transform Eq. (4) with respect to the time difference $\tau=t-t^{\prime}$,

$\chi_{uv}(\omega,t)=\int_{\tau=-\infty}^{\tau=0}d\tau~{}~{}\chi_{uv}(t,t+\tau)e^{-i\omega\tau}.$
(5)

In general, this depends on the probe time $t$. This is especially important
when the perturbation breaks time-periodicity as in the case of a pump-probe
experiment and the results are sensitive to the time of measurement.
Additionally, we average this over one cycle around $t$ to take into account
the small but finite width of the probe,

$\bar{\chi}_{uv}(\omega,t)=\frac{1}{T}\int_{t}^{t+T}dt^{\prime}\chi_{uv}(\omega,t^{\prime}).$
(6)

The electrical conductivity is then expressed as

$\sigma_{uv}(\omega)={\bar{\chi}_{uv}(\omega)}/{\omega}.$ (7)

Response of a system at equilibrium -  For the special case of an unperturbed
system, by noting that stationary states evolve as
$|\psi_{\alpha}(t)\rangle=e^{-iE_{\alpha}t}|\psi_{\alpha}(0)\rangle$ with
$E_{\alpha}$ being the energy of the $\alpha$th band, Eq. (5) becomes,

$\displaystyle\chi^{Eq}_{uv}(\omega)=i\sum_{\alpha\gamma{\bf k}}$
$\displaystyle g_{0\alpha}$
$\displaystyle\Big{[}\frac{\langle\psi_{\alpha}|H_{u}|\psi_{\gamma}\rangle\langle\psi_{\gamma}|H_{v}|\psi_{\alpha}\rangle}{\omega+(E_{\alpha}-E_{\gamma})+i\eta}$
$\displaystyle-$
$\displaystyle\frac{\langle\psi_{\gamma}|H_{u}|\psi_{\alpha}\rangle\langle\psi_{\alpha}|H_{v}|\psi_{\gamma}\rangle}{\omega-(E_{\alpha}-E_{\gamma})+i\eta}\Big{]}$
$\displaystyle+\sum_{\alpha{\bf k}}$ $\displaystyle g_{0\alpha}$
$\displaystyle\langle\psi_{\alpha}|M_{uv}|\psi_{\alpha}\rangle.$ (8)

Note that in the absence of a drive, there is no dependence on the final time
$t$ as the system is actually time independent and the averaging in Eq. (6) is
skipped.

Response of a perfect Floquet drive - Similarly, we derive a simpler
expression for a Floquet system by noting that the Floquet states can be
written in terms of the Fourier components $|\phi_{\alpha}^{(n)}\rangle$ i.e.

$\displaystyle|\Psi_{\alpha}(t)\rangle=e^{-i\epsilon_{\alpha}t}|\phi_{\alpha}(t)\rangle=\sum_{n}e^{-i(\epsilon_{\alpha}-\Omega
n)t}|\phi_{\alpha}^{(n)}\rangle$
$\displaystyle|\phi_{\alpha}^{(n)}\rangle=\frac{1}{T}\int_{0}^{T}dt~{}e^{-in\Omega
t}|\phi_{\alpha}(t)\rangle.$ (9)

The details of calculating these quasi-mode wavefunctions are given in the
Supplemental Material [47]. The expression for the homodyne [49]
susceptibility is then modified to,

$\displaystyle\chi^{Fl}_{uv}(\omega)$ $\displaystyle=$ $\displaystyle
i\sum_{\alpha\gamma m{\bf
k}}g_{0\alpha}\Bigg{[}\frac{\sum_{l}\langle\phi_{\alpha}^{(l)}|H_{u}|\phi_{\gamma}^{(l+m)}\rangle{\sum_{l^{\prime}}\langle\phi_{\gamma}^{(l^{\prime}+m)}|H_{v}|\phi_{\alpha}^{(l^{\prime})}\rangle}}{\omega+(\epsilon_{\alpha}-\epsilon_{\gamma}-m\Omega)+i\eta}$
(10) $\displaystyle-$
$\displaystyle\frac{\sum_{l}\langle\phi_{\alpha}^{(l)}|H_{v}|\phi_{\gamma}^{(l+m)}\rangle{\sum_{l^{\prime}}\langle\phi_{\gamma}^{(l^{\prime}+m)}|H_{u}|\phi_{\alpha}^{(l^{\prime})}\rangle}}{\omega-(\epsilon_{\alpha}-\epsilon_{\gamma}+m\Omega)+i\eta}\Bigg{]}$
$\displaystyle+$ $\displaystyle\sum_{\alpha{\bf
k}l}g_{0\alpha}\langle\phi_{\alpha}^{(l)}|M_{uv}|\phi_{\alpha}^{(l)}\rangle.$

The terms in the above expression correspond to optical transitions from the
$\alpha$ band of the $l$th Floquet zone to the $\gamma$ band of the $(l+m)$th
Floquet zone. While in principle the Fourier indices being summed over should
span all integers from $-\infty$ to $+\infty$, in practice only a few Fourier
components of each state are significant. This can be seen by solving a
simpler case of a driven single band system where the weight of the $n$th
Fourier mode is proportional to the Bessel function
$\mathcal{J}_{n}\left(\frac{V_{0}}{\Omega}\right)$ [45]. For a small ratio
$V_{0}/\Omega$ this drops rapidly with $|n|$. In our case we find that
$|\phi_{\alpha}^{(n)}\rangle$ is negligible beyond $n=\pm 3$ for the drive
parameters that we are working with. Therefore, in order to numerically
evaluate the conductivity, we terminate the sums in Eq. (10) at $l=\pm 3$ and
$m=\pm 3$. For lower drive frequencies or higher amplitudes a higher cut off
may be required.

Importantly, since we assume a perfectly periodic drive we take the population
of the levels to have the simple form,

$\displaystyle g_{0\alpha}=\begin{cases}1&\mbox{for the lower
band~{}~{}}(\alpha=-)\\\ 0&\mbox{for the upper
band~{}~{}}(\alpha=+).\end{cases}$ (11)

While this is never the case for a driven system, since all drives are turned
on at some finite time, many authors resort to this population as it is the
simplest.

Figure 2: Conductivities for the undriven case (grey solid), pure Floquet
drive (red dotted), a Gaussian pump of width $\Delta=20$ cycles (green solid)
and a slow quench (black dash-dot) with $\beta=0.05$. As the response to a
Gaussian and quench is almost identical, it is safe to infer that the response
is almost entirely dependent on the instantaneous drive amplitude.

Gaussian and quench pumps -  We now turn our attention to a realistic scenario
where the drive is a Gaussian pulse. We first compute the conductivity from
Eq. (4) using the real time evolution for a Gaussian pump as well as a quench.
We compare that to the response of a perfectly periodic drive as well as the
unperturbed (equilibrium) case which are obtained from Eq. (10) and Eq. (8)
respectively. For the Floquet response we have used the form of $g_{0\alpha}$
given in Eq. (11).

This comparison is shown in Fig. 2 for $\Omega=2A$ for the real and imaginary
parts of longitudinal and transverse conductivity, i.e. (a)
$\sigma^{Re}_{xx}$, (b) $\sigma^{Im}_{xx}$, (c) $\sigma^{Re}_{xy}$ and (d)
$\sigma^{Im}_{xy}$. The width of the Gaussian is $\Delta=20$ cycles and the
quench ramp time is $1/\beta=20$ cycles. The conductivities are shown at the
peak of the Gaussian (green solid line) and after the quench has reached
saturation (black dot-dash line). Although some features seem to agree, there
is a significant difference between the ideal Floquet response and the actual
response with a Gaussian drive or a quench. This difference is more pronounced
when the probe frequency is higher than the drive frequency, i.e.,
$\omega>\Omega$ where the sign of certain features is inverted. Moreover, we
see that some features which are very strong in the Floquet response are
suppressed in the Gaussian/quench response.

Figure 3: Comparison of the response of a Gaussian evaluated using a real time
evolution (green) and the modulated Floquet response (blue) using the
$g_{0\alpha}$ in Eq. (12) at the peak of the Gaussian pictured in Fig. 1, when
the response is expected to have the most similarity to the ideal Floquet
case.

Memory of initial state \- The comparison between the response of the
periodically driven systems and the ones with a pulse shape leads us to
speculate that the initial state is not forgotten even after several cycles of
the drive. To illustrate this we devise an approximate expression for the
time-evolved conductivity as follows. For a measurement of the Gaussian-driven
system at time $\tilde{t}$ we calculate the Floquet states of a system driven
by an ideal sinusoidal drive whose amplitude is $V_{{Pump}}(\tilde{t})$. We
then use these states to calculate the response using Eq. (10), albeit with
one important difference. We replace the simple population $g_{0\alpha}$ of
Eq. (11) by the overlap of the Floquet state with the equilibrium state,

$g_{0\alpha}(\tilde{t})=|\langle\phi_{\alpha}^{\tilde{t}}|\Psi_{0}\rangle|^{2}$
(12)

where $|\Psi_{0}\rangle$ is the initial state, which in our case we have taken
to be lying in the lower band of the equilibrium spectrum. The Floquet state
$|\phi_{\alpha}^{\tilde{t}}\rangle$ is the eigenstate of the Floquet operator
corresponding to the drive frequency and the instantaneous amplitude, i.e. the
magnitude of the envelope function at the probe time $\tilde{t}$. We call the
response thus obtained as the “Modulated Floquet response”. The inset in Fig.
3 shows $g_{0\alpha}$ at the peak, with the colors showing the population of
the lower band of a given Floquet zone. Here red represents the regions where
the lower band is populated and the response in these regions is very close to
the Floquet-like response. On the other hand blue is where the upper band is
populated and the response is exactly inverted from the Floquet response.
Moreover, at the momenta where the original bands have folded to create the
Floquet bands, both bands are partially populated. This intermediate regime is
responsible for the behaviour around the drive frequency and optical
transitions around there are suppressed. This is marked with an asterisk (*)
in Figs. 2 and 3 and corresponds to the transitions shown in the Floquet
spectrum of Fig. 1(b).

We compare the above “modulated Floquet response” to that of the time evolved
system for the case of a Gaussian drive at various probe times (See
Supplemental material [47]). The response well before the peak/ramp is found
to be close to the equilibrium response, as naturally expected. However, even
at the peak of the Gaussian pump $(t=0)$, i.e. where the drive amplitude is
the highest, where there is significant mismatch between the time-evolved
response and either the Floquet or equilibrium response, the modulated Floquet
response gives a good fit. The same holds for the quench scenario, where even
though the drive amplitude is kept on for a significantly long time, the
response saturates and does not replicate the case of a pure Floquet drive.
This means that even when an external driving is switched on very slowly, the
system never forgets its initial state and never goes into a pure Floquet
regime where only one of the Floquet bands is fully populated.

Figure 4: Comparison of the response of a Gaussian evaluated using a real time
evolution (green) and the modulated Floquet response (blue) using the
$g_{0\alpha}$ in Eq. (12) at $t=-2\Delta$, i.e. well before the peak of the
gaussian.

We also plot the response obtained from the real time evolution and the
modulated Floquet response for probe time $t=-2\Delta$ in Fig. 4, where, again
the two behaviours are in agreement. Similar plots for intermediate probe
times and for a higher drive frequency are shown in the Supplemental material
[47].

Conclusions - The agreement between the real time evolved conductivity and the
modulated Floquet response is a clear sign of the importance of the initial
state at any time of probe, even at the center of a wide Gaussian shaped pulse
or late after the ramp time of a quench. While the quasi-modes which
contribute to the relevant optical transitions can be approximated by an
instantaneous Floquet theory, one must keep in mind that band inversion may
invert energies but does not invert the population of the bands. This
unfortunately means that, unlike at equilibrium, a situation in which
quasienergy bands with topological character are not, in general completely
full or empty and quantized DC conductivity is not likely.

Acknowledgements -  The authors thank Babak Seradjeh and Martin Rodriguez-Vega
for useful discussions. We acknowledge financial support from the Natural
Sciences and Engineering Research Council of Canada (NSERC) and Fonds de
recherche du Québec – Nature et technologies (FRQNT). RS acknowledges Dganit
Meidan for financial support from Israel Science Foundation (ISF) grant. The
computations presented here were conducted in the computing resources provided
by the Digital Research Alliance of Canada and Calcul Quebec.

## References

  * Moore [2010] J. E. Moore, The birth of topological insulators, Nature 464, 194 (2010).
  * Moore and Balents [2007] J. E. Moore and L. Balents, Topological invariants of time-reversal-invariant band structures, Phys. Rev. B 75, 121306 (2007).
  * Fu  _et al._ [2007] L. Fu, C. L. Kane, and E. J. Mele, Topological insulators in three dimensions, Phys. Rev. Lett. 98, 106803 (2007).
  * Hasan and Kane [2010] M. Z. Hasan and C. L. Kane, Colloquium: Topological insulators, Rev. Mod. Phys. 82, 3045 (2010).
  * Qi and Zhang [2011] X.-L. Qi and S.-C. Zhang, Topological insulators and superconductors, Rev. Mod. Phys. 83, 1057 (2011).
  * Bernevig  _et al._ [2006] B. A. Bernevig, T. L. Hughes, and S.-C. Zhang, Quantum spin Hall effect and topological phase transition in HgTe quantum wells, Science 314, 1757 (2006).
  * Hsieh  _et al._ [2008] D. Hsieh, D. Qian, L. Wray, Y. Xia, Y. S. Hor, R. J. Cava, and M. Z. Hasan, A topological Dirac insulator in a quantum spin hall phase, Nature 452, 970 (2008).
  * König  _et al._ [2007] M. König, S. Wiedmann, C. Brüne, A. Roth, H. Buhmann, L. W. Molenkamp, X.-L. Qi, and S.-C. Zhang, Quantum spin Hall insulator state in hgte quantum wells, Science 318, 766 (2007).
  * Roth  _et al._ [2009] A. Roth, C. Brüne, H. Buhmann, L. W. Molenkamp, J. Maciejko, X.-L. Qi, and S.-C. Zhang, Nonlocal transport in the quantum spin Hall state, Science 325, 294 (2009).
  * Wang  _et al._ [2019] M. Wang, Q. Fu, L. Yan, W. Pi, G. Wang, Z. Zheng, and W. Luo, A $Bi_{2}Te_{3}$ topological insulator as a new and outstanding counter electrode material for high-efficiency and endurable flexible perovskite solar cells, ACS Applied Materials & Interfaces 11, 47868 (2019), pMID: 31799822.
  * Legg  _et al._ [2022] H. F. Legg, M. Rößler, F. Münning, D. Fan, O. Breunig, A. Bliesener, G. Lippertz, A. Uday, A. A. Taskin, D. Loss, J. Klinovaja, and Y. Ando, Giant magnetochiral anisotropy from quantum-confined surface states of topological insulator nanowires, Nature Nanotechnology 17, 696 (2022).
  * Ali Shameli and Yousefi [2022] M. Ali Shameli and L. Yousefi, Light trapping in thin film crystalline silicon solar cells using multi-scale photonic topological insulators, Optics and Laser Technology 145, 107457 (2022).
  * Yao  _et al._ [2012] X.-C. Yao, T.-X. Wang, H.-Z. Chen, W.-B. Gao, A. G. Fowler, R. Raussendorf, Z.-B. Chen, N.-L. Liu, C.-Y. Lu, Y.-J. Deng, Y.-A. Chen, and J.-W. Pan, Experimental demonstration of topological error correction, Nature 482, 489 (2012).
  * Yue  _et al._ [2016] Z. Yue, B. Cai, L. Wang, X. Wang, and M. Gu, Intrinsically core-shell plasmonic dielectric nanostructures with ultrahigh refractive index, Science Advances 2, e1501536 (2016).
  * Lindner  _et al._ [2011] N. H. Lindner, G. Refael, and V. Galitski, Floquet topological insulator in semiconductor quantum wells, Nature Physics 7, 490 (2011).
  * Holthaus [2015] M. Holthaus, Floquet engineering with quasienergy bands of periodically driven optical lattices, Journal of Physics B: Atomic, Molecular and Optical Physics 49, 013001 (2015).
  * Kundu  _et al._ [2014] A. Kundu, H. A. Fertig, and B. Seradjeh, Effective theory of Floquet topological transitions, Phys. Rev. Lett. 113, 236803 (2014).
  * Calvo  _et al._ [2011] H. L. Calvo, H. M. Pastawski, S. Roche, and L. E. F. F. Torres, Tuning laser-induced band gaps in graphene, Applied Physics Letters 98, 232103 (2011).
  * Oka and Aoki [2009] T. Oka and H. Aoki, Photovoltaic Hall effect in graphene, Phys. Rev. B 79, 081406 (2009).
  * Kitagawa  _et al._ [2010] T. Kitagawa, E. Berg, M. Rudner, and E. Demler, Topological characterization of periodically driven quantum systems, Phys. Rev. B 82, 235114 (2010).
  * Rudner  _et al._ [2013] M. S. Rudner, N. H. Lindner, E. Berg, and M. Levin, Anomalous edge states and the bulk-edge correspondence for periodically driven two-dimensional systems, Phys. Rev. X 3, 031005 (2013).
  * Nathan and Rudner [2015] F. Nathan and M. S. Rudner, Topological singularities and the general classification of Floquet–Bloch systems, New Journal of Physics 17, 125014 (2015).
  * Kitagawa  _et al._ [2011] T. Kitagawa, T. Oka, A. Brataas, L. Fu, and E. Demler, Transport properties of nonequilibrium systems under the application of light: Photoinduced quantum Hall insulators without landau levels, Phys. Rev. B 84, 235108 (2011).
  * Tenenbaum Katan and Podolsky [2013] Y. Tenenbaum Katan and D. Podolsky, Generation and manipulation of localized modes in floquet topological insulators, Phys. Rev. B 88, 224106 (2013).
  * Titum  _et al._ [2015] P. Titum, N. H. Lindner, M. C. Rechtsman, and G. Refael, Disorder-induced floquet topological insulators, Phys. Rev. Lett. 114, 056801 (2015).
  * Fregoso  _et al._ [2013] B. M. Fregoso, Y. H. Wang, N. Gedik, and V. Galitski, Driven electronic states at the surface of a topological insulator, Phys. Rev. B 88, 155129 (2013).
  * Foa Torres  _et al._ [2014] L. E. F. Foa Torres, P. M. Perez-Piskunow, C. A. Balseiro, and G. Usaj, Multiterminal conductance of a floquet topological insulator, Phys. Rev. Lett. 113, 266801 (2014).
  * Saha [2016] K. Saha, Photoinduced chern insulating states in semi-dirac materials, Phys. Rev. B 94, 081103 (2016).
  * Seshadri and Sen [2022] R. Seshadri and D. Sen, Engineering Floquet topological phases using elliptically polarized light, Phys. Rev. B 106, 245401 (2022).
  * Seshadri [2023] R. Seshadri, Floquet topological phases on a honeycomb lattice using elliptically polarized light, Materials Research Express 10, 024002 (2023).
  * Gu  _et al._ [2011] Z. Gu, H. A. Fertig, D. P. Arovas, and A. Auerbach, Floquet spectrum and transport through an irradiated graphene ribbon, Phys. Rev. Lett. 107, 216601 (2011).
  * Harper  _et al._ [2020] F. Harper, R. Roy, M. S. Rudner, and S. Sondhi, Topology and broken symmetry in Floquet systems, Annual Review of Condensed Matter Physics 11, 345 (2020).
  * Rudner and Lindner [2020] M. S. Rudner and N. H. Lindner, Band structure engineering and non-equilibrium dynamics in floquet topological insulators, Nature Reviews Physics 2, 229 (2020).
  * Lindner  _et al._ [2013] N. H. Lindner, D. L. Bergman, G. Refael, and V. Galitski, Topological Floquet spectrum in three dimensions via a two-photon resonance, Phys. Rev. B 87, 235131 (2013).
  * Dóra  _et al._ [2012] B. Dóra, J. Cayssol, F. Simon, and R. Moessner, Optically engineering the topological properties of a spin hall insulator, Phys. Rev. Lett. 108, 056602 (2012).
  * Thakurathi  _et al._ [2013] M. Thakurathi, A. A. Patel, D. Sen, and A. Dutta, Floquet generation of majorana end modes and topological invariants, Phys. Rev. B 88, 155133 (2013).
  * Thakurathi  _et al._ [2014] M. Thakurathi, K. Sengupta, and D. Sen, Majorana edge modes in the kitaev model, Phys. Rev. B 89, 235434 (2014).
  * Katan and Podolsky [2013] Y. T. Katan and D. Podolsky, Modulated Floquet topological insulators, Phys. Rev. Lett. 110, 016802 (2013).
  * Hübener  _et al._ [2017] H. Hübener, M. A. Sentef, U. De Giovannini, A. F. Kemper, and A. Rubio, Creating stable floquet–weyl semimetals by laser-driving of 3D Dirac materials, Nature Communications 8, 13940 (2017).
  * Rechtsman  _et al._ [2013] M. C. Rechtsman, J. M. Zeuner, Y. Plotnik, Y. Lumer, D. Podolsky, F. Dreisow, S. Nolte, M. Segev, and A. Szameit, Photonic floquet topological insulators, Nature 496, 196 (2013).
  * McIver  _et al._ [2020] J. W. McIver, B. Schulte, F.-U. Stein, T. Matsuyama, G. Jotzu, G. Meier, and A. Cavalleri, Light-induced anomalous hall effect in graphene, Nature Physics 16, 38 (2020).
  * Farrell  _et al._ [2016] A. Farrell, A. Arsenault, and T. Pereg-Barnea, Dirac cones, floquet side bands, and theory of time-resolved angle-resolved photoemission, Phys. Rev. B 94, 155304 (2016).
  * Wang  _et al._ [2013] Y. H. Wang, H. Steinberg, P. Jarillo-Herrero, and N. Gedik, Observation of floquet-bloch states on the surface of a topological insulator, Science 342, 453 (2013).
  * Dehghani  _et al._ [2014] H. Dehghani, T. Oka, and A. Mitra, Dissipative Floquet topological systems, Phys. Rev. B 90, 195429 (2014).
  * Farrell and Pereg-Barnea [2015] A. Farrell and T. Pereg-Barnea, Photon-inhibited topological transport in quantum well heterostructures, Phys. Rev. Lett. 115, 106403 (2015).
  * Farrell and Pereg-Barnea [2016] A. Farrell and T. Pereg-Barnea, Edge-state transport in Floquet topological insulators, Phys. Rev. B 93, 045121 (2016).
  * [47] Supplemental material, .
  * Fukui  _et al._ [2005] T. Fukui, Y. Hatsugai, and H. Suzuki, Chern numbers in discretized Brillouin zone: Efficient method of computing (spin) hall conductances, Journal of the Physical Society of Japan 74, 1674 (2005).
  * Kumar  _et al._ [2020] A. Kumar, M. Rodriguez-Vega, T. Pereg-Barnea, and B. Seradjeh, Linear response theory and optical conductivity of Floquet topological insulators, Phys. Rev. B 101, 174314 (2020).