$\left\\{\begin{aligned}
dQ_{t}^{i}&=b_{t}\,\Lambda\big{(}\delta^{*}(-Y_{t}^{i}-\beta_{t}^{\tilde{1},b})\big{)}\,dt-
a_{t}\,\Lambda\big{(}\delta^{*}(Y_{t}^{i}-\beta_{t}^{\tilde{4},a})\big{)}\,dt,\\\
\,dY_{t}^{i}&=2\phi_{t}\,Q_{t}^{i}\,dt+dM_{t}^{i},\\\
Q_{0}^{i}&=q_{0}^{i},\quad Y_{T}^{i}=-2A\,Q_{T}^{i}.\end{aligned}\right.$ (51)

Similar to the techniques in [26], to solve the non-Lipschitz FBSDE (51) we
first impose the regularizer $\xi>0$ to the forward equation so that it
becomes

$dQ_{t}^{i}=b_{t}\,\Lambda\big{(}\delta^{*}(-Y_{t}^{i}-\beta_{t}^{\tilde{1},b})\vee(-\xi)\wedge\xi\big{)}\,dt-
a_{t}\,\Lambda\big{(}\delta^{*}(Y_{t}^{i}-\beta_{t}^{\tilde{4},a})\vee(-\xi)\wedge\xi\big{)}\,dt.$

Since the equation is then Lipschitz, we prove the well-posedness the
regularized FBSDE. Consequently, the solution $Y^{i}$ turns out to be bounded
by some constant independent of $\xi$, which helps us remove the regularizer
$\xi$ and obtain a solution solving the original equation. Finally, one can
see the solution
$(Q^{i},Y^{i},M^{i})\in\mathbb{S}^{2}\times\mathbb{S}^{2}\times\mathbb{M}$
obtained for FBSDE (51) is also unique by a continuation argument as in [35].
Note that the control problem (50) and (51) are also solved by artificial
agent $\tilde{2}$ and $\tilde{3}$; for illustration, we have

$\beta_{t}^{\tilde{2},a}=\beta_{t}^{\tilde{4},a}+\delta^{*}(Y_{t}^{\tilde{2}}-\beta_{t}^{\tilde{4},a})\text{
\; and \;
}\beta_{t}^{\tilde{2},b}=\beta_{t}^{\tilde{1},b}+\delta^{*}(-Y_{t}^{\tilde{2}}-\beta_{t}^{\tilde{1},b}),$

where $Y^{\tilde{2}}$ solves

$\left\\{\begin{aligned}
dQ_{t}^{\tilde{2}}&=b_{t}\,\Lambda\big{(}\delta^{*}(-Y_{t}^{\tilde{2}}-\beta_{t}^{\tilde{1},b})\big{)}\,dt-
a_{t}\,\Lambda\big{(}\delta^{*}(Y_{t}^{\tilde{2}}-\beta_{t}^{\tilde{4},a})\big{)}\,dt,\\\
\,dY_{t}^{\tilde{2}}&=2\phi_{t}\,Q_{t}^{\tilde{2}}\,dt+dM_{t}^{\tilde{2}},\\\
Q_{0}^{\tilde{2}}&=q_{0}^{\tilde{2}},\quad
Y_{T}^{\tilde{2}}=-2A\,Q_{T}^{\tilde{2}}.\end{aligned}\right.$

For any $i\in\mathcal{I}$, referring to [26], the monotonicity of the control
with respect to the initial inventory yields

  * (1)

$\hat{\delta}_{t}^{i,a}\geq\beta_{t}^{\tilde{3},a}=\beta_{t}^{\tilde{4},a}$,

  * (2)

$\hat{\delta}_{t}^{i,b}\geq\beta_{t}^{\tilde{2},b}=\beta_{t}^{\tilde{1},a}$,

  * (3)

$\hat{\delta}_{t}^{i,a}-\beta_{t}^{\tilde{3},a}\leq
C\,(q_{0}^{\tilde{3}}-q_{0}^{i})$ and
$\hat{\delta}_{t}^{i,b}-\beta_{t}^{\tilde{2},b}\leq
C\,(q_{0}^{i}-q_{0}^{\tilde{2}})$

for all $t$, where $C>0$ is some constant. Considering the assumption on
sequences $(k_{n})_{n\in\mathbb{N}}$ and $(l_{n})_{n\in\mathbb{N}}$, we
conclude that

$\inf_{i\in\mathcal{I}}\hat{\delta}_{t}^{i,b}=\beta_{t}^{\tilde{2},b}=\lim_{n\to\infty}\hat{\delta}_{t}^{k_{n},b}\text{
\; and \;
}\inf_{i\in\mathcal{I}}\hat{\delta}_{t}^{i,a}=\beta_{t}^{\tilde{3},a}=\lim_{n\to\infty}\hat{\delta}_{t}^{l_{n},a}.$

Finally, recalling that both $(k_{n})_{n\in\mathbb{N}}$ and
$(l_{n})_{n\in\mathbb{N}}$ are non-duplicate, the strategy profile
$(Q^{i},Y^{i},M^{i})_{i\in\mathcal{I}}$ is a Nash equilibrium since

$\displaystyle dQ_{t}^{i}$
$\displaystyle=b_{t}\,\Lambda\big{(}\delta^{*}(-Y_{t}^{i}-\beta_{t}^{\tilde{1},b})\big{)}\,dt-
a_{t}\,\Lambda\big{(}\delta^{*}(Y_{t}^{i}-\beta_{t}^{\tilde{4},a})\big{)}\,dt$
$\displaystyle=b_{t}\,\Lambda\big{(}\delta^{*}(-Y_{t}^{i}-\inf_{\mathcal{I}\ni
j\neq i}\hat{\delta}_{t}^{j,b})\big{)}\,dt-
a_{t}\,\Lambda\big{(}\delta^{*}(Y_{t}^{i}-\inf_{\mathcal{I}\ni j\neq
i}\hat{\delta}_{t}^{j,a})\big{)}\,dt.$

In other words, processes $\beta^{\tilde{4},a}$ and $\beta^{\tilde{1},b}$ are
the genuine best ask and bid strategies. ∎

The remainder of this section is dedicated to the well-posedness of BSRE (43).
When the coefficients $a,b,\phi$, and $A$ are all deterministic, then
$\mathcal{B}$ is also deterministic and thus the term $\mathcal{Z}$ becomes
zero, simplifying BSRE (43) into a matrix Riccati equation. Due to its
extensive applications, the well-posedness of the matrix Riccati equation has
attracted significant attention, with an early exploration credited to [38].
However, existing literature has predominantly focused on equations where the
coefficients are either positive definite or symmetric. When $a,b,\phi$ and
$A$ are random, the term $\mathcal{Z}$ must be non-zero, rendering (43)
genuinely stochastic. While [33] introduced a stochastic adaptation of
Bellman’s quasi-linearization method in [38], similar to the deterministic
scenario, cases when coefficients are neither symmetric nor positive definite
are rarely studied. To circumvent such difficulties, we first study the
deterministic (43) using the well-known Radon’s lemma, followed by an analysis
of the stochastic version, building upon insights gained in the previous step.

Assume $a,b,\phi$, and $A$ are all deterministic and consider the following
linear matrix differential equations:

$\begin{pmatrix}V^{\prime}(t)\\\
U^{\prime}(t)\end{pmatrix}=\begin{pmatrix}0&\mathcal{B}_{t}\\\
2\phi_{t}\,I&0\end{pmatrix}\begin{pmatrix}V(t)\\\
U(t)\end{pmatrix};\quad\begin{pmatrix}V(T)\\\
U(T)\end{pmatrix}=\begin{pmatrix}I\\\ -2A\,I\end{pmatrix}.$ (52)

The purpose of Radon’s lemma is to establish a connection between the matrix
Riccati equation (43) and the linear equation (52):

######  Lemma 6.2 ([18]).

(1) Let $\mathcal{X}$ be a solution of the deterministic Riccati equation (43)
on some interval $\mathcal{J}\subseteq[0,T]$ such that $T\in\mathcal{J}$.
Denote by $V$ the unique solution of the linear equation

$V^{\prime}(t)=\mathcal{B}_{t}\,\mathcal{X}_{t}\,V(t);\quad V(T)=I,$

for $t\in\mathcal{J}$ and set $U(t)=\mathcal{X}_{t}\,V(t)$. Then, the matrix
$\begin{pmatrix}V(t)\\\ U(t)\end{pmatrix}$ defines for $t\in\mathcal{J}$ the
solution of the linear differential equation (52);

(2) Suppose $\begin{pmatrix}V(t)\\\ U(t)\end{pmatrix}$ is on some interval
$\mathcal{J}\subseteq[0,T]$ a solution of the linear differential equation
(52) such that $\det V(t)\neq 0$ for all $t\in\mathcal{J}$, then

$\mathcal{X}:\mathcal{J}\to\mathbb{R}^{N\times N},\quad t\mapsto
U(t)\,V(t)^{-1}=:\mathcal{X}_{t}$ (53)

is a solution to the deterministic Riccati equation (43) on $\mathcal{J}$.

Radon’s lemma reveals that the matrix Riccati equation (43) is locally
equivalent to the linear differential equation (52). This equivalence persists
until a potentially finite blow-up time, as indicated by (53). Furthermore, it
is evident from (53) that the solution $\mathcal{X}$ of the matrix Riccati
equation experiences blow-ups at moments when $\det V(t)$ vanishes. Hence, we
shift our focus to the linear differential equation (52) and verify the non-
singularity of $V$.

######  Theorem 6.3.

Assume the coefficients $a,b,\phi$, and $A$ are all deterministic. Then the
matrix Riccati equation (43) accepts a unique bounded solution $\mathcal{X}$
for cases:

  * (i)

when $\phi_{t}=0$ for all $t$;

  * (ii)

when $N=2$.

Moreover, in both cases, the solution $\mathcal{X}_{t}$ has row sum
$-2A-2\int_{t}^{T}\phi_{s}\,ds$ for all $t$.

###### PROOF.

If $\phi\equiv 0$ as in case (i), then $U$ and $V$ can be solved one-by-one in
the linear system (52). This allows us to write the solution explicitly as

$U(t)=-2A\,I\text{ \; and \; }V(t)=I+2A\,\int_{t}^{T}\mathcal{B}_{s}\,ds.$

Recalling that $\mathcal{B}_{s}$ is an $M_{0}$-matrix for all $s$, then
$\int_{t}^{T}\mathcal{B}_{s}\,ds$ is also an $M_{0}$-matrix and $V(t)$ is then
a $Z_{+}$-matrix with row sum $1$. By Lemma 4.4, we know $V(t)$ is an
$M$-matrix, the non-singularity of which can be deduced from its strict
diagonal dominance. According to Lemma 6.2, the unique solution $\mathcal{X}$
of matrix Riccati equation (43) can be represented by
$\mathcal{X}_{t}=U(t)\,V(t)^{-1}$. While $V(t)$ is a strictly diagonally
dominant $M$-matrix, it further implies $V(t)^{-1}$ is a non-negative matrix
and is strictly diagonally dominant of column entries. Defining
$\vec{1}=(1,\dots,1)$, the fact $V(t)\cdot\vec{1}=\vec{1}$ yields

$\vec{1}=I\cdot\vec{1}=V(t)^{-1}\cdot V(t)\cdot\vec{1}=V(t)^{-1}\cdot\vec{1}.$

This tells us that $V(t)^{-1}$ has row sum $1$. Denote by
$C([0,T];\mathbb{R}^{2\times 2})$ the set of matrix-valued functions on
$[0,T]$ that is element-wise continuous. Let $C_{1}([0,T];\mathbb{R}^{2\times
2})\subset C([0,T];\mathbb{R}^{2\times 2})$ be the class such that, for any
fixed $t$, the associated matrix is a $Z_{+}$-matrix with row sum $1$. In the
set $C([0,T];\mathbb{R}^{2\times 2})$, consider the norm

$\|P\|_{\ell}:=\sup_{t\in[0,T]}e^{-\ell\,(T-t)}\cdot\|P(t)\|_{F},$

where $\ell$ is a positive constant to be specified, and $\|\cdot\|_{F}$
represents the Frobenius norm. The supremum in the $\|\cdot\|_{\ell}$ ensures
that both $C([0,T];\mathbb{R}^{2\times 2})$ and
$C_{1}([0,T];\mathbb{R}^{2\times 2})$ are complete with respect to
$\|\cdot\|_{\ell}$. Picking any $P\in C_{1}([0,T];\mathbb{R}^{2\times 2})$,
let us solve the linear equation

$U^{\prime}(t)=2\phi_{t}\,P(t);\quad U(T)=-2A\,I,$

and denote by $U$ its solution. This yields

$-U(t)=2A\,I+2\int_{t}^{T}\phi_{s}\,P(s)\,ds$

and consequently $-U(t)$ is a $Z_{+}$-matrix with row sum
$2A+2\int_{t}^{T}\phi_{s}\,ds$ for any $t$. Given such $U$, we move on to the
other linear equation

$V^{\prime}(t)=\mathcal{B}_{t}\,U(t);\quad V(T)=I$

and its solution $V$ reads

$V(t)=I+\int_{t}^{T}\mathcal{B}_{s}\big{[}-U(s)\big{]}\,ds.$

Considering the following facts: (1) $\mathcal{B}_{s}$ is an $M_{0}$-matrix;
(2) $-U(s)$ is $Z_{+}$-matrix with a uniform row sum
$2A+2\int_{t}^{T}\phi_{s}\,ds$; (3) both of them are $2$-by-$2$ matrices, some
direct calculations infer $\mathcal{B}_{s}[-U(s)]$ is an $M_{0}$-matrix.
Indeed, we apply the $2$-dimensional condition only for the property that
$M$-matrices are closed under multiplication in $2$-dimension. The matrix
$V(t)$ is thus a $Z_{+}$-matrix with row sum $1$. In summary, the above
procedure defines a map

$C_{1}([0,T];\mathbb{R}^{2\times 2})\ni P\hookrightarrow\mathfrak{T}(P)=V\in
C_{1}([0,T];\mathbb{R}^{2\times 2}).$

We proceed to show $\mathfrak{T}$ is a contraction with respect to
$\|\cdot\|_{\ell}$ for some large $\ell$.

For any $P,\tilde{P}\in C_{1}([0,T];\mathbb{R}^{2\times 2})$, write
$V=\mathfrak{T}(P)$ and $\tilde{V}=\mathfrak{T}(\tilde{P})$. Utilizing the
Fubini’s theorem, the map $\mathfrak{T}$ can be briefly written as

$\displaystyle V(t)$
$\displaystyle=I+2\int_{t}^{T}\mathcal{B}_{s}\Big{(}A\,I+\int_{s}^{T}\phi_{u}\,P(u)\,du\Big{)}\,ds$
$\displaystyle=I+2A\int_{t}^{T}\mathcal{B}_{s}\,ds+2\int_{t}^{T}\Big{(}\int_{t}^{u}\mathcal{B}_{s}\,ds\Big{)}\,\phi_{u}\,P(u)\,du.$

We can then see it holds for all $t$ that

$\displaystyle\|\tilde{V}(t)-V(t)\|_{F}$ $\displaystyle\leq
2\int_{t}^{T}\big{\|}\phi_{u}\,\Big{(}\int_{t}^{u}\mathcal{B}_{s}\,ds\Big{)}\big{\|}_{F}\cdot\|\tilde{P}_{u}-P_{u}\|_{F}\,du$
$\displaystyle\leq C\int_{t}^{T}e^{\ell\,(T-u)}\cdot
e^{-\ell\,(T-u)}\,\|\tilde{P}_{u}-P_{u}\|_{F}\,du$ $\displaystyle\leq
C\,\frac{e^{\ell\,(T-t)}-1}{\ell}\,\|\tilde{P}-P\|_{\ell},$

which immediately yields

$\displaystyle e^{-\ell\,(T-t)}\,\|\tilde{V}(t)-V(t)\|_{F}$ $\displaystyle\leq
C\,\frac{e^{\ell\,(T-t)}-1}{\ell\,e^{\ell\,(T-t)}}\,\|\tilde{P}-P\|_{\ell},$
$\displaystyle\|\tilde{V}-V\|_{\ell}$
$\displaystyle\leq\frac{C}{\ell}\,\|\tilde{P}-P\|_{\ell}.$

It suffices to pick $\ell>C$ for $\mathfrak{T}$ to be a contraction. By the
Banach fixed-point theorem, the map $\mathfrak{T}$ accepts a unique fixed
point, corresponding to the unique solution of the linear differential
equation (52). Therefore, we can conclude for all $t$ that:

  * (1)

$V(t)$ as the solution of (52) is a $Z_{+}$-matrix with row sum $1$ and thus
is non-singular;

  * (2)

$-U(t)$ is a $Z_{+}$-matrix with row sum $2A+2\int_{t}^{T}\phi_{s}\,ds$.

Discussed above, the matrix $V(t)^{-1}$ is then a non-negative matrix with row
sum $1$. We can now see $-U(t)\,V(t)^{-1}$ has the row sum
$-2(A+\int_{t}^{T}\phi_{s}\,ds)$. ∎

######  Remark 6.4.

Case (i) in Theorem 6.3 pertains to the differential game of $N$ risk-neutral
players, while case (ii) is concerned with the two-player differential game.
The following Theorem 6.5 addresses the two-player stochastic differential
game.

One crucial characteristic of the solution, shown in Theorem 6.3, is the
uniform row sum of the solution $\mathcal{X}_{t}$ across all $t$. This
attribute helps solve the stochastic Riccati equation and further enables the
derivation of additional properties of the solution.

######  Theorem 6.5.

If $N=2$, then the BSRE (43) accepts a unique bounded solution $\mathcal{X}$.
Furthermore, for all $t$ the matrix $\mathcal{X}_{t}$ is non-positive with row
sum $-2\,\mathbb{E}_{t}[A+\int_{t}^{T}\phi_{s}\,ds]$, and is strictly
diagonally dominant of column entries.

###### PROOF.

Since the solution exhibits a uniform row sum of
$-2A-2\int_{t}^{T}\phi_{s}\,ds$ for all $t$ in the deterministic case, we
consider the ansatz that the row sum becomes
$\mathcal{S}_{t}:=\mathbb{E}_{t}[-2A-2\int_{t}^{T}\phi_{s}\,ds]$ in the
stochastic context. Specifically, we propose the form

$\mathcal{X}_{t}=\begin{bmatrix}\mathcal{X}_{11}(t)&\mathcal{X}_{12}(t)\\\
\mathcal{X}_{21}(t)&\mathcal{X}_{22}(t)\end{bmatrix}=\begin{bmatrix}\chi_{1}(t)&\mathcal{S}_{t}-\chi_{1}(t)\\\
\mathcal{S}_{t}-\chi_{2}(t)&\chi_{2}(t)\end{bmatrix}$ (54)

where $\chi_{1},\chi_{2}$ are some processes to be determined. Noting that
$\mathcal{B}$ can be represented as

$\mathcal{B}_{t}=\begin{bmatrix}\mathfrak{b}_{1}(t)&-\mathfrak{b}_{1}(t)\\\
-\mathfrak{b}_{2}(t)&\mathfrak{b}_{2}(t)\end{bmatrix}$

for some bounded non-negative processes $\mathfrak{b}_{1}$ and
$\mathfrak{b_{2}}$, we compute

$\displaystyle\mathcal{X}_{t}\,\mathcal{B}_{t}\,\mathcal{X}_{t}$
$\displaystyle\hskip
17.07182pt=\begin{bmatrix}(\mathfrak{b}_{1}\chi_{1}-\mathfrak{b}_{2}\mathcal{S}+\mathfrak{b}_{2}\chi_{1})\cdot(\chi_{1}-\mathcal{S}+\chi_{2})&-(\mathfrak{b}_{1}\chi_{1}-\mathfrak{b}_{2}\mathcal{S}+\mathfrak{b}_{2}\chi_{1})\cdot(\chi_{1}-\mathcal{S}+\chi_{2})\\\
-(\mathfrak{b}_{1}\chi_{2}-\mathfrak{b}_{1}\mathcal{S}+\mathfrak{b}_{2}\chi_{2})\cdot(\chi_{1}-\mathcal{S}+\chi_{2})&(\mathfrak{b}_{1}\chi_{2}-\mathfrak{b}_{1}\mathcal{S}+\mathfrak{b}_{2}\chi_{2})\cdot(\chi_{1}-\mathcal{S}+\chi_{2})\end{bmatrix},$

where the $t$ variable is omitted for conciseness whenever appropriate. Since
$\mathbb{E}_{t}[-2A-2\int_{0}^{T}\phi_{s}\,ds]$ is a martingale, we define
$\mathcal{H}\in\mathbb{H}^{2}$ such that

$\mathbb{E}_{t}\big{[}-2A-2\int_{0}^{T}\phi_{s}\,ds\big{]}=\int_{0}^{t}\mathcal{H}_{s}\,dW_{s}.$

Plugging (54) back to (43), suppose the BSDE

$d\chi_{1}(t)=\Big{[}-\big{(}\mathfrak{b}_{1}\chi_{1}-\mathfrak{b}_{2}\mathcal{S}+\mathfrak{b}_{2}\chi_{1}\big{)}\big{(}\chi_{1}-\mathcal{S}+\chi_{2}\big{)}+2\phi_{t}\Big{]}\,dt+\mathcal{Z}_{11}(t)\,dW_{t}$

with $\chi_{1}(T)=-2A$ is well-posed, we can have

$\displaystyle d\mathcal{X}_{12}(t)$
$\displaystyle=d\mathcal{S}_{t}-d\chi_{1}(t)$
$\displaystyle=2\phi_{t}\,dt+\mathcal{H}_{t}\,dW_{t}-\Big{[}-\big{(}\mathfrak{b}_{1}\chi_{1}-\mathfrak{b}_{2}\mathcal{S}+\mathfrak{b}_{2}\chi_{1}\big{)}\big{(}\chi_{1}-\mathcal{S}+\chi_{2}\big{)}+2\phi_{t}\Big{]}\,dt-\mathcal{Z}_{11}(t)\,dW_{t}$
$\displaystyle=\big{(}\mathfrak{b}_{1}\chi_{1}-\mathfrak{b}_{2}\mathcal{S}+\mathfrak{b}_{2}\chi_{1}\big{)}\big{(}\chi_{1}-\mathcal{S}+\chi_{2}\big{)}\,dt+\big{(}\mathcal{H}_{t}-\mathcal{Z}_{11}(t)\big{)}\,dW_{t}$

and $\mathcal{X}_{12}(T)=-2A+2A=0$. Since it suffices to let
$\mathcal{Z}_{12}(t)=\mathcal{H}_{t}-\mathcal{Z}_{11}(t)$, the hypothesis on
the position of $\mathcal{X}_{12}(t)$ has been checked; the one of
$\mathcal{X}_{21}(t)$ can be justified via similar calculations. We can now
shift our focus to the BSDEs for $\chi_{1}$ and $\chi_{2}$:

$\displaystyle d\chi_{1}(t)$
$\displaystyle=\Big{[}-\big{(}\mathfrak{b}_{1}\chi_{1}-\mathfrak{b}_{2}\mathcal{S}+\mathfrak{b}_{2}\chi_{1}\big{)}\cdot\big{(}\chi_{1}-\mathcal{S}+\chi_{2}\big{)}+2\phi_{t}\Big{]}\,dt+\mathcal{Z}_{11}(t)\,dW_{t};\quad\chi_{1}(T)=-2A,$
(55) $\displaystyle d\chi_{2}(t)$
$\displaystyle=\Big{[}-\big{(}\mathfrak{b}_{1}\chi_{2}-\mathfrak{b}_{1}\mathcal{S}+\mathfrak{b}_{2}\chi_{2}\big{)}\cdot\big{(}\chi_{1}-\mathcal{S}+\chi_{2}\big{)}+2\phi_{t}\Big{]}\,dt+\mathcal{Z}_{22}(t)\,dW_{t};\quad\chi_{2}(T)=-2A.$

Setting $\vartheta_{1}:=\chi_{1}+\chi_{2}-\mathcal{S}$ and
$\vartheta_{2}:=\chi_{1}-\mathcal{S}$, the system (55) can be equivalently
transformed to

$\displaystyle d\vartheta_{1}(t)$
$\displaystyle=\Big{[}-\big{(}\mathfrak{b}_{1}(t)+\mathfrak{b}_{2}(t)\big{)}\,\vartheta_{1}(t)^{2}+2\phi_{t}\Big{]}\,dt+\tilde{\mathcal{Z}}_{1}(t)\,dW_{t};\quad\vartheta_{1}(T)=-2A,$
(56) $\displaystyle d\vartheta_{2}(t)$
$\displaystyle=-\vartheta_{1}(t)\,\Big{[}\big{(}\mathfrak{b}_{1}(t)+\mathfrak{b}_{2}(t)\big{)}\,\vartheta_{2}(t)+\mathfrak{b}_{1}(t)\,\mathcal{S}_{t}\Big{]}\,dt+\tilde{\mathcal{Z}}_{2}(t)\,dW_{t};\quad\vartheta_{2}(T)=0.$

Notice the first equation of (56) is decoupled and its well-posedness is
studied in [26]. We thus know there exists a unique bounded solution
$\vartheta_{1}$, satisfying $\vartheta_{1}(t)\in[-2(\bar{A}+\bar{\phi}\,T),0]$
for all $t$. Whereas the second equation in (56) is linear, we can write the
solution explicitly as

$\vartheta_{2}(t)=\mathbb{E}_{t}\Big{[}\int_{t}^{T}\vartheta_{1}(u)\,\mathfrak{b}_{1}(u)\,\mathcal{S}_{u}\,e^{\int_{t}^{u}\vartheta_{1}(s)\,(\mathfrak{b}_{1}(s)+\mathfrak{b}_{2}(s))\,ds}\,du\Big{]}\in[0,4\bar{\mathfrak{b}}\,(\bar{A}+\bar{\phi}\,T)^{2}\,T],$

where $\bar{\mathfrak{b}}$ is the upper bound for the process
$\mathfrak{b}_{1}$. Bounded solutions $\chi_{1}$ and $\chi_{2}$ can then be
constructed based on $\vartheta_{1}$ and $\vartheta_{2}$. On the other hand,
given the process $\vartheta_{1}$, the process $\chi_{1}$ also satisfies the
linear BSDE

$d\chi_{1}(t)=\Big{[}-\vartheta_{1}\,(\mathfrak{b}_{1}+\mathfrak{b}_{2})\,\chi_{1}+\big{(}2\phi+\mathfrak{b}_{2}\,\mathcal{S}\,\vartheta_{1}\big{)}\Big{]}\,dt+\mathcal{Z}_{11}(t)\,dW_{t};\quad\chi_{1}(T)=-2A,$

the unique solution of which reads

$\displaystyle\chi_{1}(t)$
$\displaystyle=-\,\mathbb{E}\Big{[}2A\,e^{\int_{t}^{T}\vartheta_{1}(u)\,(\mathfrak{b}_{1}(u)+\mathfrak{b}_{2}(u))\,du}+\int_{t}^{T}\big{(}2\phi_{u}+\mathfrak{b}_{2}(u)\,\mathcal{S}(u)\,\vartheta_{1}(u)\big{)}\,e^{\int_{t}^{u}\vartheta_{1}(s)\,(\mathfrak{b}_{1}(s)+\mathfrak{b}_{2}(s))\,ds}\,du\Big{]}$
$\displaystyle\leq 0.$

Combining the results on $\chi_{1}$, $\vartheta_{1}$, and $\vartheta_{2}$, in
the matrix form (54) we can conclude that: (i) entries of the first row are
non-positive; (2) the difference $\mathcal{X}_{11}-\mathcal{X}_{21}\leq 0$.
With a symmetric discussion for the $\chi_{2}$, the following can be similarly
achieved: (iii) entries of the second row are non-positive; (iv) the
difference $\mathcal{X}_{22}-\mathcal{X}_{12}\leq 0$. For the uniqueness, it
suffices to see that any two bounded solutions solve the same Lipschitz BSDE,
which only admits a unique solution. ∎

Since we have solved the two-dimensional case completely, the infinite-player
game studied in the beginning of this section then has a Nash equilibrium.

######  Proposition 6.6.

The infinite-player game described in Proposition 6.1 has a Nash equilibrium.

The remaining of this section is devoted to the economic interpretation. We
regard the price impact as the price movement triggered by order imbalances.
Given some $q_{0}\in\mathbb{R}$, consider the following condition:

$q_{0}^{i}=q_{0}\text{\, for all \,}i.$ (57)

Since (57) indicates that all agents share the same initial inventory, the
ordering property 5.2 tells us that they also share the same strategy in the
equilibrium. Hence, we treat the price competition under (57) as weak, and it
holds for all $i$ that

$dQ_{t}^{i}=b_{t}\,\Lambda(0)\,dt-a_{t}\,\Lambda(0)\,dt$

and $Y^{i}$ can be subsequently computed explicitly. For simplicity, here we
let $q_{0}=0$ and $\phi_{t}=\phi$, with $\phi$ and $A$ being non-negative
constants. Also, we let $\Lambda$ be the exponential function as in Example
3.12. The ask price at time $t$ then reads

$\displaystyle S_{t}+\frac{1}{\gamma}$
$\displaystyle+2\big{[}A+\phi(T-t)\big{]}\int_{0}^{t}(a_{u}-b_{u})\,du$ (58)
$\displaystyle+2A\,\mathbb{E}_{t}\Big{[}\int_{t}^{T}(a_{u}-b_{u})\,du\Big{]}+2\phi\,\mathbb{E}_{t}\Big{[}\int_{t}^{T}(T-u)(a_{u}-b_{u})\,du\Big{]},$

where integral terms are considered as the price impact. The impact is further
linear if $\phi=0$. We regard the integral in the first line of (58) as the ex
post impact since it encapsulates the influence of previous orders.
Additionally, integrals in the second line of (58)—indicating the adjustment
in correspondence to the expected future imbalances—are considered as the ex
ante impact.

We then look at the two-player game for convenience. First, note that the
decoupling field can be seen as a generalization of the derivative of the
terminal payoff. More specifically, fixing any $s\in[0,T]$, we can truncate
the solution
$(\boldsymbol{Q}_{t},\boldsymbol{Y}_{t},\boldsymbol{M}_{t})_{t\in[0,T]}$ with
respect to time, resulting in
$(\boldsymbol{Q}_{t},\boldsymbol{Y}_{t},\boldsymbol{M}_{t})_{t\in[0,s]}$. In
reference to agent $1$, the corresponding strategy
$(\psi^{1,a}(\boldsymbol{Y}_{t}),\psi^{1,b}(\boldsymbol{Y}_{t}))_{t\in[0,s]}$
will solve an optimal control problem if $u^{1}(s,(q^{1},Q^{2}_{s}))$, as a
function of $q^{1}$, can be regarded as the derivative of a new terminal
penalty at time $s$. The stochastic maximum principle implies that a
sufficient condition for this is: the terminal penalty is concave with respect
to $q^{1}$. To verify the concavity, the definition of the matrix $\chi$ and
Theorem 6.5 yield

$u^{1}(s,(\tilde{q}^{1},q^{2}))-u^{1}(s,(q^{1},q^{2}))=-C\,(\tilde{q}^{1}-q^{1})$
(59)

for any $\tilde{q}^{1},q^{1}\in\mathbb{R}$, where $C$ is some non-negative
random variable. If $P(s,q^{1},q^{2})$ is some function such that $\partial
P/\partial q^{1}=u^{1}$, it can be inferred from (59) that $P$ is concave in
$q^{1}$. Hence, the process
$(\psi^{1,a}(\boldsymbol{Y}_{t}),\psi^{1,b}(\boldsymbol{Y}_{t}))_{t\in[0,s]}$
solves the stochastic control problem with $P(s,\cdot,Q_{s}^{2})$ as the
terminal penalty. On the other hand, due the symmetry
$u^{1}(s,(q^{1},q^{2}))=u^{2}(s,(q^{2},q^{1}))$, the new terminal penalty for
agent $2$ can be $P(s,\cdot,Q_{s}^{1})$. It suffices to look at one player.

We introduce the $\mathcal{F}_{s}$-measurable random variable

$\displaystyle\mathfrak{q}:=\mathbb{E}_{s}\Big{[}A\,\int_{s}^{T}\Lambda(0)\,(a_{u}-b_{u})\,du+\int_{s}^{T}$
$\displaystyle\phi_{r}\,\big{(}\int_{s}^{r}\Lambda(0)\,(a_{u}-b_{u})\,du\big{)}\,dr\Big{]}$
(60)
$\displaystyle\Big{/}\;\mathbb{E}_{s}\Big{[}A+\int_{s}^{T}\phi_{r}\,dr\Big{]},$

where the penalty coefficients are forced to be non-zero. Note that
$\mathfrak{q}$ represents a weighted average of future order imbalances,
determined by the penalty coefficients, and $\mathfrak{q}=0$ if $s=T$. By the
definition of decoupling fields, the value of
$u(s,(\mathfrak{q},\mathfrak{q}))$ is determined by $(Y_{s}^{1},Y_{s}^{2})$
that solves the FBSDE

$\left\\{\begin{aligned}
dQ_{t}^{i}&=b_{t}\,\Lambda\big{(}\psi^{i,b}(\boldsymbol{Y}_{t})-\bar{\psi}^{i,b}(\boldsymbol{Y}_{t})\big{)}\,dt-
a_{t}\,\Lambda\big{(}\psi^{i,a}(\boldsymbol{Y}_{t})-\bar{\psi}^{i,a}(\boldsymbol{Y}_{t})\big{)}\,dt,\\\
\,dY_{t}^{i}&=2\phi_{t}\,Q_{t}^{i}\,dt+dM_{t}^{i},\\\
Q_{s}^{i}&=\mathfrak{q},\quad Y_{T}^{i}=-2A\,Q_{T}^{i},\end{aligned}\right.$
(61)

for all $i$. Because the initial inventories are identical, it follows

$Y_{s}^{i}=-\mathbb{E}_{s}\Big{[}A\,\big{(}\mathfrak{q}+\int_{s}^{T}\Lambda(0)\,(b_{r}-a_{r})\,dr\big{)}+\int_{s}^{T}\phi_{r}\,\big{(}\mathfrak{q}+\int_{s}^{r}\Lambda(0)\,(b_{u}-a_{u})\,du\big{)}\Big{]}.$

Due to the choice (60), we can now see $Y_{s}^{i}=0$ and thus
$u(s,(\mathfrak{q},\mathfrak{q}))$ is a vector with zeros. The concavity
implies that the agent is still penalizing the excessive inventories at time
$s$. Moreover, conditional on $Q_{s}^{2}=\mathfrak{q}$, the penalty is applied
to the deviation from $\mathfrak{q}$ as a generalization of $0$ at time $T$.
For the case that $q^{2}>q$, Theorem 6.5 again yields

$u^{1}(s,(\mathfrak{q},q^{2}))=u^{1}(s,(\mathfrak{q},\mathfrak{q}))-C\,(q^{2}-\mathfrak{q})\leq
u^{1}(s,(\mathfrak{q},\mathfrak{q}))=0,$

where $C$ is some non-negative random variable. Notice the (generalized)
derivative of $u^{1}$ with respect to $q^{1}$ is bounded away from $0$ if the
same condition holds for coefficients $\phi$ and $A$. Then, there exists a
random variable $q^{*}\leq\mathfrak{q}$ such that $u^{1}(s,(q^{*},q^{2}))=0$.
In the economic context, it signifies that if the inventory level of the other
player is higher than the benchmark value $\mathfrak{q}$ at time $s$, the
target of her own inventory will be smaller than $\mathfrak{q}$. The case
$q^{2}<\mathfrak{q}$ can be similarly discussed.

##  7\. Heterogeneous Risk Coefficients

While the preceding sections mainly focus on homogeneous agents, this final
section delves into heterogeneous risk coefficients under the linear setting.
In addition to the well-posedness result, we will present an example to show
that the desirable ordering property in the homogeneous case breaks down.

We consider a two-player game in the linear setting of Section 2. Recall that
the inventory and cash of the agent $1$ are modelled by $X_{t}^{1}$ and
$Q_{t}^{1}$ accordingly:

$\displaystyle
X_{t}^{1}=\int_{0}^{t}(S_{u}+\delta_{u}^{1,a})\,a_{u}\,(\,\zeta-\gamma\,(\delta_{u}^{1,a}-\delta_{u}^{2,a})\,)\,du-\int_{0}^{t}(S_{u}-\delta_{u}^{1,b})\,b_{u}\,(\,\zeta-\gamma\,(\delta_{u}^{1,b}-\delta_{u}^{2,b})\,)\,du,$
$\displaystyle
Q_{t}^{1}=q_{0}^{1}-\int_{0}^{t}a_{u}\,(\,\zeta-\gamma\,(\delta_{u}^{1,a}-\delta_{u}^{2,a})\,)\,du+\int_{0}^{t}b_{u}\,(\,\zeta-\gamma\,(\delta_{u}^{1,b}-\delta_{u}^{2,b})\,)\,du.$

Note that $\bar{\delta}^{1,a}$ can be replaced by $\delta^{2,a}$. The ones of
agent $2$ are symmetric. Controlling
$\boldsymbol{\delta}^{i}\in\mathbb{H}^{2}\times\mathbb{H}^{2}$, the agent
$i\in\\{1,2\\}$ aims at maximizing the objective functional

$\displaystyle J(\boldsymbol{\delta}^{i};\boldsymbol{\delta}^{-i}):$
$\displaystyle=\mathbb{E}\Big{[}X_{T}^{i}+S_{T}\,Q_{T}^{i}-\int_{0}^{T}\phi_{t}^{i}\,\big{(}Q_{t}^{i}\big{)}^{2}\,dt-A^{i}\,\big{(}Q_{T}^{i}\big{)}^{2}\Big{]}$
$\displaystyle=\mathbb{E}\Big{[}\int_{0}^{T}\delta_{t}^{i,a}\,a_{t}\,(\,\zeta-\gamma\,(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})\,)\,dt+\int_{0}^{T}\delta_{t}^{i,b}\,b_{t}\,(\,\zeta-\gamma\,(\delta_{t}^{i,b}-\bar{\delta}_{t}^{i,b})\,)\,dt$
$\displaystyle\hskip
85.35826pt-\int_{0}^{T}\phi_{t}^{i}\,\big{(}Q_{t}^{i}\big{)}^{2}\,dt-A^{i}\,\big{(}Q_{T}^{i}\big{)}^{2}\Big{]},$

where penalties $\phi^{i}:=(\phi_{t}^{i})_{t\in[0,T]}\in\mathbb{H}^{2}$ and
$A^{i}\in L^{2}(\Omega,\mathcal{F}_{T})$ are non-negative, satisfying
$\phi_{t}^{i}\leq\bar{\phi}$ and $A^{i}\leq\bar{A}$ for some constants
$\bar{\phi},\bar{A}>0$. The goal is to find a Nash equilibrium in the same
sense as (3).

######  Theorem 7.1.

A strategy profile $(\boldsymbol{\delta}^{j})_{1\leq j\leq
2}\in(\mathbb{H}^{2}\times\mathbb{H}^{2})^{2}$ forms a Nash equilibrium if and
only if

$\displaystyle\delta_{t}^{1,a}$
$\displaystyle=\frac{2}{3}\,Y_{t}^{1}+\frac{1}{3}\,Y_{t}^{2}+\frac{\zeta}{\gamma},\quad\delta_{t}^{1,b}=-\frac{2}{3}\,Y_{t}^{1}-\frac{1}{3}\,Y_{t}^{2}+\frac{\zeta}{\gamma},$
(62) $\displaystyle\delta_{t}^{2,a}$
$\displaystyle=\frac{1}{3}\,Y_{t}^{1}+\frac{2}{3}\,Y_{t}^{2}+\frac{\zeta}{\gamma},\quad\delta_{t}^{2,b}=-\frac{1}{3}\,Y_{t}^{1}-\frac{2}{3}\,Y_{t}^{2}+\frac{\zeta}{\gamma},$

where $(Y^{1},Y^{2})$ solves the following system of FBSDEs:

$\left\\{\begin{aligned}
\;&dQ_{t}^{1}=\zeta\,(b_{t}-a_{t})\,dt+\gamma\,(a_{t}+b_{t})\,(Y_{t}^{1}-Y_{t}^{2})/3\,dt,\\\
&dY_{t}^{1}=2\phi_{t}^{1}\,Q_{t}^{1}\,dt+dM^{1}_{t},\\\
&Q_{0}^{1}=q_{0}^{1},\quad Y_{T}^{1}=-A^{1}\,Q_{T}^{1};\end{aligned}\right.$
(63) $\left\\{\begin{aligned}
\;&dQ_{t}^{2}=\zeta\,(b_{t}-a_{t})\,dt+\gamma\,(a_{t}+b_{t})\,(Y_{t}^{2}-Y_{t}^{1})/3\,dt,\\\
&dY_{t}^{2}=2\phi_{t}^{2}\,Q_{t}^{2}\,dt+dM^{2}_{t},\\\
&Q_{0}^{2}=q_{0}^{2},\quad Y_{T}^{2}=-A^{2}\,Q_{T}^{2},\end{aligned}\right.$
(64)

with $M^{i}_{t}\in\mathbb{M}$ for all $i\in\\{1,2\\}$.

###### PROOF.

With reference to the homogeneous case, the heterogeneity of the penalty
parameters does not affect the concavity of the objective functional. By
Theorem 2.4, strategy profile $(\boldsymbol{\delta}^{j})_{1\leq j\leq
2}\in(\mathbb{H}^{2}\times\mathbb{H}^{2})^{2}$ forms a Nash equilibrium if and
only if it solves

$\left\\{\begin{aligned}
\;&dQ_{t}^{i}=-a_{t}\,\big{(}\zeta+\gamma\bar{\delta}^{i,a}_{t}-\gamma\delta_{t}^{i,a}\big{)}dt+b_{t}\,\big{(}\zeta+\gamma\bar{\delta}^{i,b}_{t}-\gamma\delta_{t}^{i,b}\big{)}dt,\\\
&d\delta_{t}^{i,a}=d\bar{\delta}^{i,a}_{t}/2+\phi_{t}\,Q_{t}^{i}\,dt-d\tilde{M}_{t}^{i},\\\
&d\delta_{t}^{i,b}=d\bar{\delta}^{i,b}_{t}/2-\phi_{t}\,Q_{t}^{i}\,dt+d\tilde{M}_{t}^{i},\\\
&Q_{0}^{i}=q_{0}^{i},\quad\delta_{T}^{i,a}=\zeta/(2\gamma)+\bar{\delta}^{i,a}_{T}/2-A\,Q_{T}^{i},\quad\delta_{T}^{i,b}=\zeta/(2\gamma)+\bar{\delta}^{i,b}_{T}/2+A\,Q_{T}^{i},\end{aligned}\right.$
(65)

with $\tilde{M}^{i}\in\mathbb{M}$ for all $i\in\\{1,2\\}$. System (63)-(64)
can be obtain from (65) via the transform (62). ∎

We introduce the following processes and random variables:

$\displaystyle\boldsymbol{H}_{t}=\frac{\gamma}{3}\,(a_{t}+b_{t})\,\begin{bmatrix}1&-1\\\
-1&1\end{bmatrix},\hskip
11.38092pt\boldsymbol{G}_{t}=\begin{pmatrix}\,\zeta\,(b_{t}-a_{t})\,\\\
\,\zeta\,(b_{t}-a_{t})\,\end{pmatrix},$
$\displaystyle\boldsymbol{D}_{t}=\begin{bmatrix}2\,\phi_{t}^{1}&0\\\
0&2\,\phi_{t}^{2}\end{bmatrix},\hskip
11.38092pt\boldsymbol{A}=\begin{bmatrix}2A^{1}&0\\\ 0&2A^{2}\end{bmatrix}.$

System (63)-(64) can then be written neatly in the vector form

$\left\\{\begin{aligned}
\;&d\boldsymbol{Q}_{t}=\boldsymbol{G}_{t}\,dt+\boldsymbol{H}_{t}\,\boldsymbol{Y}_{t}\,dt,\\\
&d\boldsymbol{Y}_{t}=\boldsymbol{D}_{t}\,\boldsymbol{Q}_{t}\,dt+d\boldsymbol{M}_{t},\\\
&\boldsymbol{Q}_{0}=\boldsymbol{q}_{0},\quad\boldsymbol{Y}_{T}=-\boldsymbol{A}\,\boldsymbol{Q}_{T}.\end{aligned}\right.$
(66)

######  Theorem 7.2.

The FBSDE (66) accepts a unique solution
$(\boldsymbol{Q},\boldsymbol{Y},\boldsymbol{M})\in(\mathbb{S}^{2}\times\mathbb{S}^{2}\times\mathbb{M})^{2}$.
Moreover, the solution has the relation

$\boldsymbol{Y}_{t}=\boldsymbol{R}_{t}\,\boldsymbol{Q}_{t}+\boldsymbol{P}_{t},$
(67)

where the matrix-valued process $\boldsymbol{R}$ solves the BSRE

$d\boldsymbol{R}_{t}=\big{(}\boldsymbol{D}_{t}-\boldsymbol{R}_{t}\,\boldsymbol{H}_{t}\,\boldsymbol{R}_{t}\big{)}\,dt+\boldsymbol{Z}_{t}\,dW_{t},\quad\boldsymbol{R}_{T}=-\boldsymbol{A}.$
(68)

###### PROOF.

Given the linearity of system (66), we think of the linear ansatz (67). By
plugging (67) back to (66) and matching the coefficients, we can see
$\boldsymbol{R}$ solves the BSRE (68) and $\boldsymbol{P}$ satisfies the
linear BSDE

$d\boldsymbol{P}_{t}=-\big{(}\boldsymbol{R}_{t}\,\boldsymbol{H}_{t}\,\boldsymbol{P}_{t}+\boldsymbol{R}_{t}\,\boldsymbol{G}_{t}\big{)}+\tilde{\boldsymbol{Z}}_{t}\,dW_{t},\quad\boldsymbol{P}_{T}=0.$
(69)

Let us look at (68) and first remark that both $\boldsymbol{D}_{t}$ and
$\boldsymbol{A}$ are symmetric positive semi-definite matrices for all $t$.
The same is true for $\boldsymbol{H}_{t}$ because it is a symmetric
$M_{0}$-matrix. Equation (68) is similar to the type of BSRE studied in [34],
from which we learn there exists a unique bounded solution $\boldsymbol{R}$.
Consequently, the linear BSDE (69) has bounded coefficients. The resulting
Lipschitz property further guarantees the well-posedness of $\boldsymbol{P}$.
A solution $(\boldsymbol{Q},\boldsymbol{Y},\boldsymbol{M})$ can be obtained by
inserting (67) back to (66). Realizing that (67) is also a regular decoupling
field, the solution is also unique according to Theorem 5.9. ∎

######  Example 7.3.

In this example, we will see the ordering property may break down under
heterogeneity risk preferences. The following conditions are examined:

  * •

order flows satisfy $a_{t}=b_{t}=3/2$ for all $t$;

  * •

random variables $A^{1}$ and $A^{2}$ are deterministic, with $A^{1}>A^{2}>0$;

  * •

running penalties satisfy $\phi_{t}^{1}=\phi_{t}^{2}=1$;

  * •

it holds $0<q_{0}^{1}<q_{0}^{2}<\frac{A^{1}}{A^{2}}\,q_{0}^{1}$ and $T$ is
large enough.

Due to the lack of noises, BSRE (68) degenerates to the Riccati equation

$\displaystyle d\boldsymbol{R}_{t}$
$\displaystyle=\big{(}\boldsymbol{D}_{t}-\boldsymbol{R}_{t}\,\boldsymbol{H}_{t}\,\boldsymbol{R}_{t}\big{)}\,dt$
(70) $\displaystyle=\bigg{(}\begin{bmatrix}1&0\\\
0&1\end{bmatrix}-\boldsymbol{R}_{t}\,\gamma\,\begin{bmatrix}1&0\\\
-1&0\end{bmatrix}\,\begin{bmatrix}1&0\\\
-1&0\end{bmatrix}^{*}\,\boldsymbol{R}_{t}\bigg{)}\,dt$

such that $\boldsymbol{R}_{T}=-\boldsymbol{A}$. We apply the asymptotic result
studied in Theorem 2.1 by [38], after checking that:

  * •

matrix $\begin{bmatrix}1&0\\\ -1&0\end{bmatrix}$ is stable because it is an
$M_{0}$-matrix;

  * •

$\bigg{(}\begin{bmatrix}1&0\\\ 0&1\end{bmatrix},\boldsymbol{0}\bigg{)}$ is
observable since matrix $\begin{bmatrix}1&0&0&0\\\ 0&1&0&0\end{bmatrix}$ is of
rank $2$.

Consequently, we know there exists a matrix
$\boldsymbol{R}_{\infty}\in\mathbb{R}^{2\times 2}$ that is the unique negative
definite solution of algebraic Riccati equation

$\displaystyle 0$ $\displaystyle=\begin{bmatrix}1&0\\\
0&1\end{bmatrix}-\boldsymbol{R}\,\gamma\,\begin{bmatrix}1&0\\\
-1&0\end{bmatrix}\,\begin{bmatrix}1&0\\\
-1&0\end{bmatrix}^{*}\,\boldsymbol{R}$ (71)
$\displaystyle=\begin{bmatrix}1&0\\\
0&1\end{bmatrix}-\boldsymbol{R}\,\gamma\,\begin{bmatrix}1&-1\\\
-1&1\end{bmatrix}\,\boldsymbol{R}.$

Moreover, it holds that

$\boldsymbol{R}_{\infty}=\lim_{t\to-\infty}\boldsymbol{R}_{t}.$ (72)

We claim that $\boldsymbol{R}_{\infty}$ is both symmetric and persymmetric.
Indeed, given $\boldsymbol{R}_{\infty}$ as the solution of (72), consider the
permutation transform

$\tilde{\boldsymbol{R}}_{\infty}:=\begin{bmatrix}0&1\\\
1&0\end{bmatrix}\,\boldsymbol{R}_{\infty}\,\begin{bmatrix}0&1\\\
1&0\end{bmatrix}.$

Because direct calculation yields

$\displaystyle 0$ $\displaystyle=\begin{bmatrix}0&1\\\
1&0\end{bmatrix}\,\begin{bmatrix}1&0\\\
0&1\end{bmatrix}\,\begin{bmatrix}0&1\\\ 1&0\end{bmatrix}-\begin{bmatrix}0&1\\\
1&0\end{bmatrix}\,\boldsymbol{R}_{\infty}\,\gamma\,\begin{bmatrix}1&-1\\\
-1&1\end{bmatrix}\,\boldsymbol{R}_{\infty}\,\begin{bmatrix}0&1\\\
1&0\end{bmatrix}$ $\displaystyle=\begin{bmatrix}1&0\\\
0&1\end{bmatrix}-\tilde{\boldsymbol{R}}_{\infty}\,\gamma\,\begin{bmatrix}1&-1\\\
-1&1\end{bmatrix}\,\tilde{\boldsymbol{R}}_{\infty},$

matrix $\tilde{\boldsymbol{R}}_{\infty}$ solves the same equation (72) as
$\boldsymbol{R}_{\infty}$. The symmetric and persymmetric property thus can be
deduced from the uniqueness of the solution. We then write

$\boldsymbol{R}_{\infty}=\begin{bmatrix}\mathfrak{r}_{1}&\mathfrak{r}_{2}\\\
\mathfrak{r}_{2}&\mathfrak{r}_{1}\end{bmatrix}.$

Note that $\mathfrak{r}_{1}-\mathfrak{r}_{2}<0$ due to the negative
definiteness. Defining

$\begin{bmatrix}\,Y_{\infty}^{1}\,\\\
\,Y_{\infty}^{2}\,\end{bmatrix}:=\boldsymbol{R}_{\infty}\,\begin{bmatrix}\,q_{0}^{1}\,\\\
\,q_{0}^{2}\,\end{bmatrix},$

we can derive

$Y_{\infty}^{1}-Y_{\infty}^{2}=(\mathfrak{r}_{1}-\mathfrak{r}_{2})\,(q_{0}^{1}-q_{0}^{2})>0.$

Since the condition $a_{t}=b_{t}$ yields
$\boldsymbol{G}_{t}=\boldsymbol{P}_{t}=0$, the decoupling field (67) tells us

$\begin{bmatrix}\,Y_{0}^{1}\,\\\
\,Y_{0}^{2}\,\end{bmatrix}=\boldsymbol{R}_{0}\,\begin{bmatrix}\,q_{0}^{1}\,\\\
\,q_{0}^{2}\,\end{bmatrix}.$

Combined with the asymptotic behavior (72), we can find $T$ large enough such
that $Y_{0}^{1}-Y_{0}^{2}>0$.

To show it can not be true that $Y_{t}^{1}-Y_{t}^{2}\geq 0$ for all $t$, let
us assume $\Delta Y_{t}:=Y_{t}^{1}-Y_{t}^{2}\geq 0$ and it follows by
definition

$Q_{T}^{1}=q_{0}^{1}+\int_{0}^{T}\Delta Y_{t}\,dt.$

If we similarly define $\Delta Q_{t}:=Q_{t}^{1}-Q_{t}^{2}$, it turns out
$(\Delta Q,\Delta Y)$ solves

$\left\\{\begin{aligned} \;&d\Delta Q_{t}=2\,\Delta Y_{t}\,dt,\\\ &d\Delta
Y_{t}=\Delta Q_{t}\,dt,\\\ &\Delta Q_{0}=q_{0}^{1}-q_{0}^{2},\quad\Delta
Y_{T}=-2A^{2}\,\Delta Q_{T}-2(A^{1}-A^{2})\,Q_{T}^{1}.\end{aligned}\right.$

We then find

$\displaystyle\Delta Y_{T}$
$\displaystyle=-2A^{2}\,\Big{(}q_{0}^{1}-q_{0}^{2}+2\int_{0}^{T}\Delta
Y_{t}\,dt\Big{)}-2(A^{1}-A^{2})\,\Big{(}q_{0}^{1}+\int_{0}^{T}\Delta
Y_{t}\,dt\Big{)}$
$\displaystyle=-2A^{2}\,\Big{(}\frac{A^{1}}{A^{2}}\,q_{0}^{1}-q_{0}^{2}\Big{)}-2(A^{1}+A^{2})\,\int_{0}^{T}\Delta
Y_{t}\,dt$ $\displaystyle<0.$

Therefore, although $\Delta Y_{0}>0$, it can not hold for all $t$ and there
exist some time $s<T$ such that $\Delta Y_{s}<0$. Finally, utilizing (62), the
sign of $\Delta Y$ determines which agent quotes a better price. For example,
the positive value of $\Delta Y_{0}$ indicates that agent $1$ provides a worse
sell price (but a better buy price) than agent $2$. For an economical
interpretation, we first remark that heterogeneity lies in the terminal
penalty coefficient. As $T$ becomes sufficiently large, the solution
$\boldsymbol{R}_{0}$ converges to $\boldsymbol{R}_{\infty}$, solving the
algebraic equation (71) independent of the terminal condition. Essentially,
the influence exerted by the terminal penalty diminishes as we move backward
in time. Consequently, agents are approximately homogeneous at the initial
time. As discussed earlier, the better buy price offered by agent $1$ results
from the lower inventory level. However, as $T$ approaches, this approximation
loses validity, leading to the ‘crossing’ of quotes.

## References

  * [1] Robert Almgren and Neil Chriss “Optimal execution of portfolio transactions” In _Journal of Risk_ 3, 2001, pp. 5–40
  * [2] Robert F Almgren “Optimal execution with nonlinear impact functions and trading-enhanced risk” In _Applied mathematical finance_ 10.1 Taylor & Francis, 2003, pp. 1–18
  * [3] Stefan Ankirchner, Alexander Fromm, Thomas Kruse and Alexandre Popier “Optimal position targeting via decoupling fields” In _The Annals of Applied Probability_ 30.2 Institute of Mathematical Statistics, 2020, pp. 644–672
  * [4] Marco Avellaneda and Sasha Stoikov “High-frequency trading in a limit order book” In _Quantitative Finance_ 8.3 Taylor & Francis, 2008, pp. 217–224
  * [5] Peter Bank, H Mete Soner and Moritz Voß “Hedging with temporary price impact” In _Mathematics and financial economics_ 11.2 Springer, 2017, pp. 215–239
  * [6] Erhan Bayraktar and Michael Ludkovski “Liquidation in limit order books with controlled intensity” In _Mathematical Finance_ 24.4 Wiley Online Library, 2014, pp. 627–650
  * [7] Nataliya Bershova and Dmitry Rakhlin “The non-linear market impact of large trades: Evidence from buy-side order flow” In _Quantitative finance_ 13.11 Taylor & Francis, 2013, pp. 1759–1778
  * [8] Dimitris Bertsimas and Andrew W Lo “Optimal control of execution costs” In _Journal of financial markets_ 1.1 Elsevier, 1998, pp. 1–50
  * [9] Frédéric Bucci, Michael Benzaquen, Fabrizio Lillo and Jean-Philippe Bouchaud “Crossover from linear to square-root market impact” In _Physical review letters_ 122.10 APS, 2019, pp. 108302
  * [10] Luciano Campi and Diego Zabaljauregui “Optimal market making under partial information with general intensities” In _Applied Mathematical Finance_ 27.1-2 Taylor & Francis, 2020, pp. 1–45
  * [11] René Carmona “Lectures on BSDEs, Stochastic Control, and Stochastic Differential Games with Financial Applications” SIAM, 2016
  * [12] René Carmona and François Delarue “Probabilistic theory of mean field games with applications I-II” Springer, 2018
  * [13] Álvaro Cartea, Sebastian Jaimungal and José Penalva “Algorithmic and high-frequency trading” Cambridge University Press, 2015
  * [14] Francis Clarke “Functional analysis, calculus of variations and optimal control” Springer, 2013
  * [15] Frank H Clarke “Optimization and nonsmooth analysis” SIAM, 1990
  * [16] Christian Clason “Nonsmooth analysis and optimization” In _arXiv preprint arXiv:1708.04180_ , 2017
  * [17] Rama Cont and Wei Xiong “Dynamics of market making algorithms in dealer markets: Learning and tacit collusion” In _Mathematical Finance_ Wiley Online Library, 2022
  * [18] Gerhard Freiling “A survey of nonsymmetric Riccati equations” In _Linear algebra and its applications_ 351 Elsevier, 2002, pp. 243–270
  * [19] Alexander Fromm “Theory and applications of decoupling fields for forward-backward stochastic differential equations” Humboldt-Universität zu Berlin, Mathematisch-Naturwissenschaftliche Fakultät II, 2015
  * [20] Alexander Fromm and Peter Imkeller “Existence, uniqueness and regularity of decoupling fields to multidimensional fully coupled FBSDEs” In _arXiv preprint arXiv:1310.0499_ , 2013
  * [21] Marek Galewski and Marius Rădulescu “On a global implicit function theorem for locally Lipschitz maps via non-smooth critical point theory” In _Quaestiones Mathematicae_ 41.4 Taylor & Francis, 2018, pp. 515–528
  * [22] Olivier Guéant “Optimal market making” In _Applied Mathematical Finance_ 24.2 Taylor & Francis, 2017, pp. 112–154
  * [23] Olivier Guéant “The Financial Mathematics of Market Liquidity” CRC Press, 2016
  * [24] Olivier Guéant, Charles-Albert Lehalle and Joaquin Fernandez-Tapia “Dealing with the inventory risk: a solution to the market making problem” In _Mathematics and financial economics_ 7 Springer, 2013, pp. 477–507
  * [25] Olivier Guéant, Charles-Albert Lehalle and Joaquin Fernandez-Tapia “Optimal portfolio liquidation with limit orders” In _SIAM Journal on Financial Mathematics_ 3.1 SIAM, 2012, pp. 740–764
  * [26] Ivan Guo, Shijia Jin and Kihun Nam “Macroscopic Market Making” In _arXiv preprint arXiv:2307.14129_ , 2023
  * [27] Thomas Ho and Hans R Stoll “Optimal dealer pricing under transactions and return uncertainty” In _Journal of Financial economics_ 9.1 Elsevier, 1981, pp. 47–73
  * [28] Roger A Horn and Charles R Johnson “Topics in matrix analysis” Cambridge university press, 1994
  * [29] Kaitong Hu, Zhenjie Ren and Nizar Touzi “On path-dependent multidimensional forward-backward SDEs” In _Numerical Algebra, Control and Optimization_ American Institute of Mathematical Sciences, 2022
  * [30] Dariusz Idczak “A global implicit function theorem and its applications to functional equations” In _arXiv preprint arXiv:1401.4049_ Citeseer, 2014
  * [31] Jialiang Luo and Harry Zheng “Dynamic Equilibrium of Market Making with Price Competition” In _Dynamic Games and Applications_ 11.3 Springer, 2021, pp. 556–579
  * [32] Jin Ma, Zhen Wu, Detao Zhang and Jianfeng Zhang “On well-posedness of forward–backward SDEs—A unified approach” In _The Annals of Applied Probability_ 25.4 Institute of Mathematical Statistics, 2015, pp. 2168–2214
  * [33] Shige Peng “A general stochastic maximum principle for optimal control problems” In _SIAM Journal on control and optimization_ 28.4 SIAM, 1990, pp. 966–979
  * [34] Shige Peng “Stochastic hamilton–jacobi–bellman equations” In _SIAM Journal on Control and Optimization_ 30.2 SIAM, 1992, pp. 284–304
  * [35] Shige Peng and Zhen Wu “Fully coupled forward-backward stochastic differential equations and applications to optimal control” In _SIAM Journal on Control and Optimization_ 37.3 SIAM, 1999, pp. 825–843
  * [36] Max O Souza and Yuri Thamsten “On regularized optimal execution problems and their singular limits” In _Applied Mathematical Finance_ Taylor & Francis, 2022, pp. 1–31
  * [37] James M Varah “A lower bound for the smallest singular value of a matrix” In _Linear Algebra and its applications_ 11.1 Elsevier, 1975, pp. 3–5
  * [38] William M Wonham “On a matrix Riccati equation of stochastic control” In _SIAM Journal on Control_ 6.4 SIAM, 1968, pp. 681–697

##  8\. Proofs: Stochastic Maximum Principle

This section is devoted to the game version of the stochastic maximum
principle. We follow the argument in [11] that is devoted to the control
problem, and start with the necessary condition. Fix an admissible strategy
profile $(\boldsymbol{\delta}^{j})_{1\leq j\leq
N}\in(\mathbb{A}\times\mathbb{A})^{N}$ and denote by
$(Q^{1},\dots,Q^{N})\in(\mathbb{S}^{2})^{N}$ the corresponding controlled
inventory. Next, in view of player $i$, let us consider
$\boldsymbol{\beta}\in\mathbb{A}\times\mathbb{A}$, which is uniformly bounded
and satisfies
$\boldsymbol{\delta}^{i}+\boldsymbol{\beta}\in\mathbb{A}\times\mathbb{A}$. The
direction $\boldsymbol{\beta}$ will be used to compute the Gâteaux derivative
of $J^{i}$. For each $\epsilon>0$ small enough, define an admissible control
$\boldsymbol{\delta}^{i,\epsilon}\in\mathbb{A}\times\mathbb{A}$ as
$\boldsymbol{\delta}^{i,\epsilon}_{t}=\boldsymbol{\delta}^{i}_{t}+\epsilon\,\boldsymbol{\beta}_{t}$,
and the corresponding controlled state is denoted by
$Q^{i,\epsilon}\in\mathbb{S}^{2}$. Let $V$ be the solution of the equation

$dV_{t}=\partial_{\boldsymbol{\delta}}\theta(t,\boldsymbol{\delta}^{i}_{t};\boldsymbol{\delta}^{-i}_{t})\cdot\boldsymbol{\beta}_{t}\,dt,$

with the initial condition $V_{0}=0$, where

$\theta(t,\boldsymbol{\delta}^{i};\boldsymbol{\delta}^{-i}):=-a_{t}\,\Lambda(\delta^{i,a}-\bar{\delta}^{i,a})+b_{t}\,\Lambda(\delta^{i,b}-\bar{\delta}^{i,b}).$

It is clear that $V$ is uniformly bounded in any finite time horizon.

######  Lemma 8.1.

The functional $\boldsymbol{\delta}\hookrightarrow
J^{i}(\boldsymbol{\delta};\boldsymbol{\delta}^{-i})$ is Gâteaux differentiable
and the derivative reads

$\displaystyle\frac{d}{d\epsilon}J^{i}(\boldsymbol{\delta}^{i}+\epsilon\boldsymbol{\beta};\boldsymbol{\delta}^{-i})\big{|}_{\epsilon=0}:$
$\displaystyle=\lim_{\epsilon\searrow
0}\frac{1}{\epsilon}\,\big{[}J^{i}(\boldsymbol{\delta}^{i}+\epsilon\boldsymbol{\beta};\boldsymbol{\delta}^{-i})-J^{i}(\boldsymbol{\delta}^{i};\boldsymbol{\delta}^{-i})\big{]}$
(73)
$\displaystyle=\mathbb{E}\Big{[}\int_{0}^{T}\big{(}\partial_{q}f(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})\,V_{t}+\partial_{\boldsymbol{\delta}}f(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})\cdot\boldsymbol{\beta}_{t}\big{)}\,dt-2\,A^{i}\,V_{T}\,Q_{T}^{i}\Big{]},$

where $f$ is the running payoff given as

$f(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})=b_{t}\,\delta_{t}^{i,b}\,\Lambda(\delta_{t}^{i,b}-\bar{\delta}_{t}^{i,b})+a_{t}\,\delta_{t}^{i,a}\,\Lambda(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})-\phi_{t}^{i}\,\big{(}Q_{t}^{i}\big{)}^{2}.$

###### PROOF.

Regarding the ask side of the running payoff, we compute that

$\displaystyle\frac{1}{\epsilon}$
$\displaystyle\mathbb{E}\Big{[}\int_{0}^{T}a_{t}\,(\delta_{t}^{i,a}+\epsilon\,\beta_{t}^{a})\,\Lambda(\delta_{t}^{i,a}+\epsilon\,\beta_{t}^{a}-\bar{\delta}_{t}^{i,a})\,dt-\int_{0}^{T}a_{t}\,\delta_{t}^{i,a}\,\Lambda(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})\,dt\Big{]}$
$\displaystyle=\mathbb{E}\Big{[}\int_{0}^{T}a_{t}\,\delta_{t}^{i,a}\,\frac{1}{\epsilon}\,\big{(}\Lambda(\delta_{t}^{i,a}+\epsilon\,\beta_{t}^{a}-\bar{\delta}_{t}^{i,a})-\Lambda^{a}(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})\big{)}\,dt\Big{]}+\mathbb{E}\Big{[}\int_{0}^{T}a_{t}\,\beta_{t}^{a}\,\Lambda(\delta_{t}^{i,a}+\epsilon\,\beta_{t}^{a}-\bar{\delta}_{t}^{i,a})\,dt\Big{]}.$

To perform the limiting procedure, it suffices to notice

$\displaystyle\int_{0}^{T}a_{t}\,\delta_{t}^{i,a}\,\frac{1}{\epsilon}\,\big{(}\Lambda(\delta_{t}^{i,a}+\epsilon\,\beta_{t}^{a}-\bar{\delta}_{t}^{i,a})-\Lambda^{a}(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})\big{)}\,dt$
$\displaystyle\leq C\,\int_{0}^{T}a_{t}\,|\delta_{t}^{i,a}|\,dt,$
$\displaystyle\int_{0}^{T}a_{t}\,\beta_{t}^{a}\,\Lambda(\delta_{t}^{i,a}+\epsilon\,\beta_{t}^{a}-\bar{\delta}_{t}^{i,a})\,dt$
$\displaystyle\leq C\,\int_{0}^{T}a_{t}\,dt,$ $\displaystyle
a_{t}\,\delta_{t}^{i,a}\,\frac{1}{\epsilon}\,\big{(}\Lambda(\delta_{t}^{i,a}+\epsilon\,\beta_{t}^{a}-\bar{\delta}_{t}^{i,a})-\Lambda^{a}(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})\big{)}$
$\displaystyle\leq C\,a_{t}\,|\delta_{t}^{i,a}|,$ $\displaystyle
a_{t}\,\beta_{t}^{a}\,\Lambda(\delta_{t}^{i,a}+\epsilon\,\beta_{t}^{a}-\bar{\delta}_{t}^{i,a})$
$\displaystyle\leq C\,a_{t},$

and the integrability on the right hand sides. Hence, the dominated
convergence theorem yields

$\displaystyle\lim_{\epsilon\to
0}\frac{1}{\epsilon}\,\mathbb{E}\Big{[}\int_{0}^{T}a_{t}\,(\delta_{t}^{i,a}$
$\displaystyle+\epsilon\,\beta_{t}^{a})\,\Lambda(\delta_{t}^{i,a}+\epsilon\,\beta_{t}^{a}-\bar{\delta}_{t}^{i,a})\,dt-\int_{0}^{T}a_{t}\,\delta_{t}^{i,a}\,\Lambda(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})\,dt\Big{]}$
$\displaystyle=\mathbb{E}\Big{[}\int_{0}^{T}a_{t}\,\delta_{t}^{i,a}\,\big{(}\Lambda(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})\big{)}^{\prime}\,\beta_{t}^{a}\,dt\Big{]}+\mathbb{E}\Big{[}\int_{0}^{T}a_{t}\,\beta_{t}^{a}\,\Lambda(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})\,dt\Big{]}$
$\displaystyle=\mathbb{E}\Big{[}\int_{0}^{T}\partial_{\delta^{a}}f(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})\,\beta_{t}^{a}\,dt\Big{]}.$

Since the bid side can be computed in the same way, the term with
$\partial_{\boldsymbol{\delta}}f$ in (73) is thus verified. With respect to
the running inventory penalty, we further calculate that

$\displaystyle\frac{1}{\epsilon}\,\mathbb{E}\Big{[}\int_{0}^{T}\phi_{t}^{i}$
$\displaystyle\big{(}Q_{t}^{i,\epsilon}\big{)}^{2}\,dt-\int_{0}^{T}\phi^{i}_{t}\,\big{(}Q_{t}^{i}\big{)}^{2}\,dt\Big{]}$
$\displaystyle=\mathbb{E}\Big{[}\int_{0}^{T}\phi_{t}^{i}\,\big{(}Q_{t}^{i,\epsilon}+Q_{t}^{i}\big{)}\,\Big{(}-\int_{0}^{t}a_{u}\,\frac{1}{\epsilon}\,\big{[}\Lambda(\delta_{u}^{i,a}+\epsilon\,\beta_{u}^{a}-\bar{\delta}_{u}^{i,a})-\Lambda(\delta_{u}^{i,a}-\bar{\delta}_{u}^{i,a})\big{]}\,du$
$\displaystyle\hskip
116.65646pt+\int_{0}^{t}b_{u}\,\frac{1}{\epsilon}\,\big{[}\Lambda(\delta_{u}^{i,b}+\epsilon\,\beta_{u}^{b}-\bar{\delta}_{u}^{i,b})-\Lambda(\delta_{u}^{i,b}-\bar{\delta}_{u}^{i,b})\big{]}\,du\Big{)}\,dt\Big{]}.$

Similarly, being aware of the boundedness of $\phi^{i}$, $Q^{i}$,
$Q^{i,\epsilon}$, $a$, $b$, as well as

$\frac{1}{\epsilon}\,\big{[}\Lambda(\delta_{u}^{i,a}+\epsilon\,\beta_{u}^{a}-\bar{\delta}_{u}^{i,a})-\Lambda(\delta_{u}^{i,a}-\bar{\delta}_{u}^{i,a})\big{]}\text{\quad
and
\quad}\frac{1}{\epsilon}\,\big{[}\Lambda(\delta_{u}^{i,b}+\epsilon\,\beta_{u}^{b}-\bar{\delta}_{u}^{i,b})-\Lambda(\delta_{u}^{i,b}-\bar{\delta}_{u}^{i,b})\big{]},$

we can conclude by the dominated convergence theorem that

$\lim_{\epsilon\to
0}\,\frac{1}{\epsilon}\,\mathbb{E}\Big{[}\int_{0}^{T}\phi_{t}^{i}\,\big{(}Q_{t}^{i,\epsilon}\big{)}^{2}\,dt-\int_{0}^{T}\phi^{i}_{t}\,\big{(}Q_{t}^{i}\big{)}^{2}\,dt\Big{]}=2\,\mathbb{E}\Big{[}\int_{0}^{T}\phi_{t}^{i}\,Q_{t}^{i}\,V_{t}\,dt\Big{]}.$

While the terminal penalty part can be justified in the same way, the proof is
complete. ∎

Let $(Y^{i},M^{i})$ be the adjoint processes associated with
$(\boldsymbol{\delta}^{j})_{1\leq j\leq N}$, i.e., processes $(Y,Z)$ solve the
BSDE

$dY_{t}^{i}=2\phi_{t}^{i}\,Q_{t}^{i}\,dt+dM_{t}^{i},\quad
Y_{T}^{i}=-2A^{i}\,Q_{T}^{i}.$

The following duality relation provides an expression of the Gâteaux
derivative of the cost functional in terms of the Hamiltonian of the system.

######  Lemma 8.2.

The duality relation is given by

$\mathbb{E}[Y_{T}^{i}\,V_{T}]=\mathbb{E}\Big{[}\int_{0}^{T}\big{(}\partial_{\boldsymbol{\delta}}\theta(t,\boldsymbol{\delta}^{i}_{t};\boldsymbol{\delta}^{-i}_{t})\cdot\boldsymbol{\beta}_{t}\,Y_{t}^{i}-\partial_{q}f(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})\,V_{t}\big{)}\,dt\Big{]}.$

Such duality further implies

$\frac{d}{d\epsilon}\,J^{i}(\boldsymbol{\delta}^{i}+\epsilon\boldsymbol{\beta};\boldsymbol{\delta}^{-i})\big{|}_{\epsilon=0}=\mathbb{E}\Big{[}\int_{0}^{T}\partial_{\boldsymbol{\delta}}H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}_{t}^{-i})\cdot\boldsymbol{\beta}_{t}\,dt\Big{]}.$

###### PROOF.

The integration by parts yields

$\displaystyle Y_{T}^{i}\,V_{T}$
$\displaystyle=Y_{0}^{i}\,V_{0}+\int_{0}^{T}Y_{t}^{i}\,dV_{t}+\int_{0}^{T}V_{t}\,dY_{t}^{i}$
$\displaystyle=\int_{0}^{T}Y_{t}\,\partial_{\boldsymbol{\delta}}\theta(t,\boldsymbol{\delta}^{i}_{t};\boldsymbol{\delta}^{-i}_{t})\cdot\boldsymbol{\beta}_{t}\,dt+\int_{0}^{T}2\,\phi_{t}^{i}\,V_{t}\,Q_{t}^{i}\,dt+\int_{0}^{T}V_{t}\,dM_{t}^{i}.$

The duality is obtain through taking the expectation. Using this relation, we
can further compute

$\displaystyle\frac{d}{d\epsilon}\,J^{i}(\boldsymbol{\delta}^{i}+\epsilon\boldsymbol{\beta}$
$\displaystyle;\boldsymbol{\delta}^{-i})\big{|}_{\epsilon=0}$
$\displaystyle=\mathbb{E}\Big{[}\int_{0}^{T}\big{(}\partial_{q}f(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})\,V_{t}+\partial_{\boldsymbol{\delta}}f(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})\cdot\boldsymbol{\beta}_{t}\big{)}\,dt+V_{T}\,Y_{T}^{i}\Big{]}$
$\displaystyle=\mathbb{E}\Big{[}\int_{0}^{T}Y_{t}^{i}\,\partial_{\boldsymbol{\delta}}\theta(t,\boldsymbol{\delta}^{i}_{t};\boldsymbol{\delta}^{-i}_{t})\cdot\boldsymbol{\beta}_{t}\,dt+\int_{0}^{T}\partial_{\boldsymbol{\delta}}f(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})\cdot\boldsymbol{\beta}_{t}\big{)}\,dt\Big{]}$
$\displaystyle=\mathbb{E}\Big{[}\int_{0}^{T}\partial_{\boldsymbol{\delta}}H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}_{t}^{-i})\cdot\boldsymbol{\beta}_{t}\,dt\Big{]}.$

∎

With these preliminary steps, we have now arrived at the necessary condition
for optimality.

######  Theorem 8.3 (Necessary condition).

If an admissible strategy profile $(\boldsymbol{\delta}^{j})_{1\leq j\leq N}$
is a Nash equilibrium, $(Q^{j})_{1\leq j\leq N}$ are the corresponding
controlled inventories, and $(Y^{j},M^{j})_{1\leq j\leq N}$ are the associated
adjoint processes, then it holds for any player $i$ that

$H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}_{t}^{-i})\geq
H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\beta};\boldsymbol{\delta}_{t}^{-i}),$

$dt\times d\mathbb{P}-a.s.$ for any
$\boldsymbol{\beta}\in[-\xi,\xi]\times[-\xi,\xi]$.

###### PROOF.

Fix any $i\in\\{1,\dots,N\\}$. Due to the convexity of the admissible space,
given any admissible and bounded
$\boldsymbol{\beta}\in\mathbb{A}\times\mathbb{A}$, we can choose the
perturbation
$\boldsymbol{\delta}^{i,\epsilon}=\boldsymbol{\delta}^{i}+\epsilon\,(\boldsymbol{\beta}-\boldsymbol{\delta}^{i})$
which is still admissible. Since $(\boldsymbol{\delta}^{j})_{1\leq j\leq N}$
is a Nash equilibrium, the following inequality should hold:

$\frac{d}{d\epsilon}\,J^{i}(\boldsymbol{\delta}^{i}+\epsilon\,(\boldsymbol{\beta}-\boldsymbol{\delta}^{i});\boldsymbol{\delta}^{-i})\big{|}_{\epsilon=0}=\mathbb{E}\Big{[}\int_{0}^{T}\partial_{\boldsymbol{\delta}}H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}_{t}^{-i})\cdot(\boldsymbol{\beta}_{t}-\boldsymbol{\delta}_{t}^{i})\,dt\Big{]}\leq
0.$

From this we can see

$\partial_{\boldsymbol{\delta}}H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}_{t}^{-i})\cdot(\boldsymbol{\beta}-\boldsymbol{\delta}_{t}^{i})\leq
0,$ (74)

$dt\times d\mathbb{P}-$a.s. for all
$\boldsymbol{\beta}\in[-\xi,\xi]\times[-\xi,\xi]$. Look at the ask side for
example. Regarding the condition (74), it can only happen at time $t$ if one
of the following three cases holds:

$\displaystyle\delta_{t}^{i,a}\in(-\xi,+\xi)\text{\quad with
\quad}\partial_{\boldsymbol{\delta}^{a}}H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}_{t}^{-i})=0,$
$\displaystyle\delta_{t}^{i,a}=-\xi\leq\beta^{a}\text{\quad with
\quad}\partial_{\boldsymbol{\delta}^{a}}H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}_{t}^{-i})\leq
0,$ $\displaystyle\delta_{t}^{i,a}=+\xi\geq\beta^{a}\text{\quad with
\quad}\partial_{\boldsymbol{\delta}^{a}}H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}_{t}^{-i})\geq
0.$

Suppose that

$\partial_{\boldsymbol{\delta}^{a}}H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}_{t}^{-i})=-a_{t}\,\big{(}\Lambda^{a}(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})\big{)}^{\prime}\,\Big{[}Y_{t}-\delta_{t}^{i,a}-\frac{\Lambda(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})}{\big{(}\Lambda(\delta_{t}^{i,a}-\bar{\delta}_{t}^{i,a})\big{)}^{\prime}}\Big{]}=0,$

the monotonicity of
$\delta+\frac{\Lambda(\delta)}{(\Lambda(\delta))^{\prime}}$ (see Assumption
3.1) implies that $\delta_{t}^{i,a}$ maximizes the corresponding part of the
Hamiltonian. In another case, think of $\delta_{t}^{i,a}=-\xi$ with
$\partial_{\boldsymbol{\delta}^{a}}H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}_{t}^{-i})\leq
0$. Since $H^{i}$ is first increasing and then decreasing, this property
infers that $\delta_{t}^{i,a}=-\xi$ is the maximizer of the Hamiltonian on the
interval $[-\xi,\xi]$. While the remaining case and the bid side can be
discussed in a similar fashion, all three cases imply that the strategy
$\boldsymbol{\delta}^{i}$ maximizes the $H^{i}$ along the optimal paths. The
proof is complete since the previous discussion holds for any player. ∎

The sufficient condition is more straightforward to derive, stated as below.

######  Theorem 8.4 (Sufficient condition).

Let $(\boldsymbol{\delta}^{j})_{1\leq j\leq N}$ be an admissible strategy
profile, $(Q^{1}\,\dots,Q^{N})$ be the corresponding controlled inventories,
and $(Y^{j},M^{j})_{1\leq j\leq N}$ be the associated adjoint processes. If it
holds $dt\times d\mathbb{P}$-a.s. that

$H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}_{t}^{-i})=\sup_{\boldsymbol{\beta}\in[-\xi,\xi]^{2}}H^{i}(t,Q_{t}^{i},Y_{t}^{i},\boldsymbol{\beta};\boldsymbol{\delta}_{t}^{-i}),$

for all $i\in\\{1,\dots,N\\}$, then $(\boldsymbol{\delta}^{j})_{1\leq j\leq
N}$ is a Nash equilibrium.

###### PROOF.

Fix any $i\in\\{1,\dots,N\\}$. Let
$\boldsymbol{\beta}\in\mathbb{A}\times\mathbb{A}$ be a generic admissible
strategy and $Q^{i^{\prime}}$ be the state process associated with the profile
$(\boldsymbol{\beta},\boldsymbol{\delta}^{-i})$. Due to the concaveness of the
terminal penalty, we have

$\displaystyle\mathbb{E}\big{[}-A^{i}\,\big{(}Q_{T}^{i}\big{)}^{2}+A^{i}\,\big{(}Q^{i^{\prime}}_{T}\big{)}^{2}\big{]}$
$\displaystyle\geq\mathbb{E}\big{[}-2A^{i}\,Q^{i}_{T}\,(Q^{i}_{T}-Q^{i^{\prime}}_{T})\big{]}$
$\displaystyle=\mathbb{E}\big{[}Y_{T}^{i}\,(Q_{T}^{i}-Q^{i^{\prime}}_{T})\big{]}$
$\displaystyle=\mathbb{E}\Big{[}\int_{0}^{T}2\phi_{t}^{i}\,Q_{t}^{i}\,(Q_{t}^{i}-Q^{i^{\prime}}_{t})\,dt+\int_{0}^{T}Y_{t}^{i}\,d(Q_{t}^{i}-Q^{i^{\prime}}_{t})\Big{]}.$

By the definition of the Hamiltonian, one can also obtain

$\displaystyle\mathbb{E}\int_{0}^{T}\big{[}f(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})$
$\displaystyle-f(t,Q_{t}^{i^{\prime}},\boldsymbol{\beta}_{t};\boldsymbol{\delta}^{-i}_{t})\big{]}\,dt$
$\displaystyle=\mathbb{E}\int_{0}^{T}\big{[}H(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})-H(t,Q_{t}^{i^{\prime}},\boldsymbol{\beta}_{t};\boldsymbol{\delta}^{-i}_{t})\big{]}\,dt-\mathbb{E}\int_{0}^{T}Y_{t}^{i}\,d(Q_{t}^{i}-Q^{i^{\prime}}_{t}),$

where we recall that $f$ is the running payoff of agent $i$. Combining two
results above, we can get

$\displaystyle
J^{i}(\boldsymbol{\delta}^{i};\boldsymbol{\delta}^{-i})-J^{i}(\boldsymbol{\beta};\boldsymbol{\delta}^{-i})$
$\displaystyle=\mathbb{E}\int_{0}^{T}\big{[}f(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})-f(t,Q_{t}^{i^{\prime}},\boldsymbol{\beta}_{t};\boldsymbol{\delta}^{-i}_{t})\big{]}\,dt+\mathbb{E}\big{[}-2A^{i}\,(Q_{T}^{i})^{2}+2A^{i}\,(Q^{i^{\prime}}_{T})^{2}\big{]}$
$\displaystyle\geq\mathbb{E}\int_{0}^{T}\big{[}H(t,Q_{t}^{i},\boldsymbol{\delta}_{t}^{i};\boldsymbol{\delta}^{-i}_{t})-H(t,Q_{t}^{i^{\prime}},\boldsymbol{\beta}_{t};\boldsymbol{\delta}^{-i}_{t})+2\phi_{t}^{i}\,Q_{t}^{i}\,(Q_{t}^{i}-Q^{i^{\prime}}_{t})\big{]}\,dt$
$\displaystyle\geq\mathbb{E}\int_{0}^{T}\big{[}H(t,Q_{t}^{i},\boldsymbol{\beta}_{t};\boldsymbol{\delta}^{-i}_{t})-H(t,Q_{t}^{i^{\prime}},\boldsymbol{\beta}_{t};\boldsymbol{\delta}^{-i}_{t})+2\phi_{t}^{i}\,Q_{t}^{i}\,(Q_{t}^{i}-Q^{i^{\prime}}_{t})\big{]}\,dt$
$\displaystyle=\mathbb{E}\int_{0}^{T}\phi_{t}^{i}\,\big{(}Q_{t}^{i}-Q^{i^{\prime}}_{t}\big{)}^{2}\,dt$
$\displaystyle\geq 0.$

The proof is complete since the previous discussion holds for any player. ∎

##  9\. Proofs: Non-smooth Analysis and Implicit Function

Let $X$ be a separable Banach space. Given function $F:X\to\mathbb{R}$, the
generalized directional derivative in $x\in X$ with respect to the direction
$h\in X$ is given by

$F^{\circ}(x;h):=\limsup_{\begin{subarray}{c}y\to x\\\ t\to
0\end{subarray}}\frac{1}{t}\,\big{[}F(y+t\,h)-F(y)\big{]}.$

In addition, we review the Rademacher’s theorem.

######  Theorem 9.1 (Rademacher).

For some $m,l\in\mathbb{N}$, let $U\subseteq\mathbb{R}^{m}$ be open and
$F:U\to\mathbb{R}^{l}$ be locally Lipschitz continuous. Then $F$ is Fréchet
differentiable in almost every $x\in U$.

In view of Rademacher’s theorem, we introduce the generalized derivative and
the generalized Jacobian matrix as follows.

######  Definition 9.2.

(1) For a locally Lipschitz function $F:X\to\mathbb{R}$, its (Clarke)
generalized derivative in $x\in X$ is defined by

$\partial F(x):=\big{\\{}x^{*}\in X^{*}\,:\,\langle x^{*},h\rangle_{X}\leq
F^{\circ}(x;h),\;\forall h\in X\big{\\}}.$

Here, $\langle\cdot,\cdot\rangle_{X}$ denotes the duality product between $X$
and $X^{*}$. Specifically, if $X=\mathbb{R}^{m}$, then $F$ is differentiable
on $\mathbb{R}^{m}\,\backslash\,E_{F}$ for a set
$E_{F}\subseteq\mathbb{R}^{m}$ of Lebesgue measure $0$. Then, the generalized
derivative can be characterized by

$\partial F(x)=\text{co }\big{\\{}\lim_{n\to\infty}\nabla
F(x_{n})\,:\,x_{n}\to x,\,x_{n}\notin E_{F}\,\big{\\}},$

where co $A$ denotes the convex hull of $A\subseteq\mathbb{R}^{m}$.

(2) For a locally Lipschitz function $F:\mathbb{R}^{m}\to\mathbb{R}^{l}$, its
(Clarke) generalized Jacobian in $x\in\mathbb{R}^{m}$ is defined as

$\partial F(x):=\text{co }\big{\\{}\lim_{n\to\infty}\nabla
F(x_{n})\,:\,x_{n}\to x,\,x_{n}\notin E_{F}\,\big{\\}}.$

Here, we let $F$ be differentiable on $\mathbb{R}^{m}\,\backslash\,E_{F}$ for
a set $E_{F}\subseteq\mathbb{R}^{m}$ of Lebesgue measure $0$.

######  Remark 9.3.

Unlike some literature, we will not use $\partial_{C}$ to denote the
generalized derivative. Instead, we will rely on the subscript notation to
represent the partial derivative. For illustration,

$\displaystyle\partial_{x}F(x,y)$ denotes the generalized derivative of $F$
with respect to its first variable; $\displaystyle\partial F(x,y)$
$\displaystyle\text{ denotes the generalized derivative of $F$ with respect to
the variable }(x,y).$

It is worth noting that in the smooth case, the generalized derivatives and
classical derivatives coincide, so there should be no ambiguity in our
notation.

With these prerequisites in place, we now proceed to present the proof of the
implicit function theorem.

###### PROOF.

(Propostion 3.9) Given the first two conditions, one already knows there
exists a unique locally Lipschitz function $f:\mathbb{R}^{m}\to\mathbb{R}^{n}$
such that $F(x,y)=0$ and $x=f(y)$ are equivalent in the set
$\mathbb{R}^{n}\times\mathbb{R}^{m}$, according to Theorem 4 in [21]. It
suffices to verify the global Lipschitz property.

First, the proof of Theorem 4 is accomplished by concatenating the local
result given by Theorem 2 in [21]. For every $y\in\mathbb{R}^{m}$, the
resulting function $f$ is Lipschitz in a neighborhood containing $y$, and thus
is locally Lipschitz overall. Subsequently, the proof of Theorem 2 is built
upon the inverse function theorem, referring to Theorem 7.1.1 in [15].
Illustrated by the Corollary on page 256 in [15], we can define the function
$\tilde{F}:\mathbb{R}^{m+n}\to\mathbb{R}^{m+n}$ by

$\tilde{F}(x,y)=\big{[}y,\;F(x,y)\big{]},$

and apply the inverse function theorem to $\tilde{F}$ to justify the local
version of the implicit function theorem (i.e., Theorem 2 in [21]). The
function $f$ is then Lipschitz (in a neighborhood) because of the Lipschitz
inverse according to Theorem 7.1.1. Finally, the Lipschitz coefficient of the
Lipschitz inverse is $1/\delta$, provided that the distance between the
distance between $\partial\tilde{F}(\hat{x},\hat{y})\,\tilde{S}$—letting
$\tilde{S}$ signify the unit sphere in $\mathbb{R}^{m+n}$—and $0$ is
$2\delta$, for a $(\hat{x},\hat{y})\in\mathbb{R}^{n}\times\mathbb{R}^{m}$
fixed at the beginning.

Therefore, given that there exists a constant $\upsilon>0$ such that the
distance between $\partial\tilde{F}(x,y)\,\tilde{S}$ and $0$ is at least
$\upsilon$ for all $(x,y)\in\mathbb{R}^{n}\times\mathbb{R}^{m}$, the above
procedure tells us that the function $f$ is Lipschitz in every small
neighborhood, with the Lipschitz coefficient being bounded uniformly by
$2/\upsilon$. It follows $f$ is (globally) Lipschitz. ∎

We finish with the calculus rules of the generalized derivative.

######  Theorem 9.4 ([16], [15]).

Let $F:X\to\mathbb{R}$ be locally Lipschitz continuous near $x\in X$ and
$\kappa\in\mathbb{R}$, then

$\partial(\kappa\,F)(x)=\kappa\partial F(x).$

Let $G:X\to\mathbb{R}$ also be locally Lipschitz continuous near $x\in X$,
then

$\partial(F+G)(x)\subseteq\partial F(x)+\partial
G(x):=\big{\\{}f+g\,:\,f\in\partial F(x),\,g\in\partial G(x)\big{\\}}.$

The above are also true for the case $F,G:\mathbb{R}^{m}\to\mathbb{R}^{l}$.