between the current case and Case 3 in Section 5.2.3 is that the jumps on the
contours $\partial D_{\lambda_{j,k}}^{\epsilon}$ and $\partial
D_{\overline{\lambda_{j,k}}}^{\epsilon}$ for all $1\leq k\leq N_{j}$ do not
decay to the identity matrix uniformly. Other than this, the two cases are
quite alike, in the sense that, all jumps enclosing eigenvalues
$\lambda_{l,k}$ and $\overline{\lambda_{l,k}}$ with $1\leq l\leq j-1$ decay to
the identity matrix uniformly, and the ones enclosing $\lambda_{l,k}$ and
$\overline{\lambda_{l,k}}$ with $l>j$ grow uniformly, as $t\to+\infty$. As a
result, both deformations (5.9) and (5.13) apply to the current case. Thus, we
start the discussion of Case 4 with RHP 9. However, unlike the previous case,
one cannot take the limit $t\to+\infty$ of RHP 9 and arrive at Equation
(5.18), because the jumps surrounding eigenvalues $\lambda_{j,k}$ are
oscillatory as $t$ changes. In fact, RHP 9 do not have a limit in Case 4.

To continue our calculation, let us first build a parametrix for RHP 9 by
ignoring all decaying jumps.

######  Riemann-Hilbert Problem 10 (Jump form).

Let $\mathbf{P}_{\mathrm{s}+,j}(\lambda;t,z)$ be a $2\times 2$ matrix function
on $\lambda\in\mathbb{C}$. Suppose that
$\mathbf{P}_{\mathrm{s}+,j}(\lambda;t,z)$ has the asymptotics
$\mathbf{P}_{\mathrm{s}+,j}(\lambda;t,z)\to\mathbb{I}$ as $\lambda\to\infty$,
is sectional analytic and satisfies the following jumps,

(5.21) $\displaystyle\mathbf{P}_{\mathrm{s}+,j}^{+}(\lambda;t,z)$
$\displaystyle=\mathbf{P}_{\mathrm{s}+,j}^{-}(\lambda;t,z)\delta_{j+1}(\lambda)^{\sigma_{3}}\mathbf{V}_{j,k}(\lambda;t,z)^{-1}\delta_{j+1}(\lambda)^{-\sigma_{3}}\,,\quad$
$\displaystyle\lambda\in\partial D_{\lambda_{j,k}}^{\epsilon}\,,\quad 1\leq
k\leq N_{j}\,,$ $\displaystyle\mathbf{P}_{\mathrm{s}+,j}^{+}(\lambda;t,z)$
$\displaystyle=\mathbf{P}_{\mathrm{s}+,j}^{-}(\lambda;t,z)\delta_{j+1}(\lambda)^{\sigma_{3}}\mathbf{V}_{j,k}^{\dagger}(\lambda;t,z)\delta_{j+1}(\lambda)^{-\sigma_{3}}\,,\quad$
$\displaystyle\lambda\in\partial
D_{\overline{\lambda_{j,k}}}^{\epsilon}\,,\quad 1\leq k\leq N_{j}\,,$

where $\mathbf{V}_{j,k}(\lambda;t,z)$ is still given in RHP 2.

The subscript “$\mathrm{s}+$” denotes the case in a stable medium with
asymptotics $t\to+\infty$. For convenience, we explicitly express the jump
matrices of RHP 10 below

(5.22)
$\displaystyle\delta_{j+1}(\lambda)^{\sigma_{3}}\mathbf{V}_{j,k}(\lambda;t,z)^{-1}\delta_{j+1}(\lambda)^{-\sigma_{3}}$
$\displaystyle=\begin{pmatrix}1&0\\\
-\frac{\omega_{j,k}\delta_{j+1}(\lambda)^{-2}\mathrm{e}^{-2\mathrm{i}\theta}}{\lambda-\lambda_{j,k}}&1\end{pmatrix}\,,$
$\displaystyle\delta_{j+1}(\lambda)^{\sigma_{3}}\mathbf{V}_{j,k}^{\dagger}(\lambda;t,z)\delta_{j+1}(\lambda)^{-\sigma_{3}}$
$\displaystyle=\begin{pmatrix}1&\frac{\overline{\omega_{j,k}}\delta_{j+1}(\lambda)^{2}\mathrm{e}^{2\mathrm{i}\theta}}{\lambda-\overline{\lambda_{j,k}}}\\\
0&1\end{pmatrix}\,.$

Note that the $(2,1)$ entry of
$\delta_{j+1}(\lambda)^{\sigma_{3}}\mathbf{V}_{j,k}(\lambda;t,z)^{-1}\delta_{j+1}(\lambda)^{-\sigma_{3}}$
and the $(1,2)$ entry of
$\delta_{j+1}(\lambda)^{\sigma_{3}}\mathbf{V}_{j,k}^{\dagger}(\lambda;t,z)\delta_{j+1}(\lambda)^{-\sigma_{3}}$
have simple poles at $\lambda=\lambda_{j,k}$ and
$\lambda=\overline{\lambda_{j,k}}$, respectively. Hence, one can convert the
above RHP with jumps into an equivalent RHP with residue conditions, just like
the equivalence between RHPs 1 and 2. Comparing the residue form of RHP 10
with RHP 1, one immediately realizes that RHP 10 describes a DSG containing
$N_{j}$ solitons. The eigenvalues are still the same $\Lambda_{j}$ in the
upper half plane, but the corresponding norming constants are modified as the
set $\Omega_{j}^{(+)}$ containing $\omega_{j,k}^{(+)}$. Therefore, the
parametrix $\mathbf{P}_{\mathrm{s}+,j}(\lambda;t,z)$ is solvable, and is
bounded as $t\to+\infty$ via Theorem 2 with $z=V_{j}t$.

One looks at the difference between the original matrix function
$\mathbf{M}^{(3)}(\lambda;t,z)$ from RHP 9 and its parametrix
$\mathbf{P}_{\mathrm{s}+,j}(\lambda;t,z)$ from RHP 10, by defining an error
function

(5.23)
$\mathbf{E}_{\mathrm{s}+,j}(\lambda;t,z)\coloneqq\mathbf{M}_{j}^{(3)}(\lambda;t,z)\mathbf{P}_{\mathrm{s}+,j}(\lambda;t,z)^{-1}\,,\qquad\mbox{for
a fixed }j=1,2,\dots,J\,.$

Then, RHPs 9 and 10 yield a RHP for the error function
$\mathbf{E}_{\mathrm{s}+,j}(\lambda;t,z)$.

######  Riemann-Hilbert Problem 11.

Let $\mathbf{E}_{\mathrm{s}+,j}(\lambda;t,z)$ be a $2\times 2$ matrix function
on $\lambda\in\mathbb{C}$. It has asymptotics
$\mathbf{E}_{\mathrm{s}+,j}(\lambda;t,z)\to\mathbb{I}$ as $\lambda\to\infty$,
is sectional analytic and satisfies jumps, for $1\leq l<j$ and $1\leq k\leq
N_{l}$

(5.24) $\displaystyle\mathbf{E}_{\mathrm{s}+,j}^{+}(\lambda;t,z)$
$\displaystyle=\mathbf{E}_{\mathrm{s}+,j}^{-}(\lambda;t,z)\mathbf{P}_{\mathrm{s}+,j}\delta_{j+1}(\lambda)^{\sigma_{3}}\mathbf{V}_{l,k}(\lambda;t,z)^{-1}\delta_{j+1}(\lambda)^{-\sigma_{3}}\mathbf{P}_{\mathrm{s}+,j}^{-1}\,,\quad$
$\displaystyle\lambda\in\partial D_{\lambda_{l,k}}^{\epsilon}\,,$
$\displaystyle\mathbf{E}_{\mathrm{s}+,j}^{+}(\lambda;t,z)$
$\displaystyle=\mathbf{E}_{\mathrm{s}+,j}^{-}(\lambda;t,z)\mathbf{P}_{\mathrm{s}+,j}\delta_{j+1}(\lambda)^{\sigma_{3}}\mathbf{V}_{l,k}^{\dagger}(\lambda;t,z)\delta_{j+1}(\lambda)^{-\sigma_{3}}\mathbf{P}_{\mathrm{s}+,j}^{-1}\,,\quad$
$\displaystyle\lambda\in\partial D_{\overline{\lambda_{l,k}}}^{\epsilon}\,,$

and for $j+1\leq l\leq J$ and $1\leq k\leq N_{l}$

(5.25) $\displaystyle\mathbf{E}_{\mathrm{s}+,j}^{+}(\lambda;t,z)$
$\displaystyle=\mathbf{E}_{\mathrm{s}+,j}^{-}(\lambda;t,z)\mathbf{P}_{\mathrm{s}+,j}\mathbf{W}_{j,l,k}^{(1)}(\lambda;t,z)\mathbf{P}_{\mathrm{s}+,j}^{-1}\,,\qquad$
$\displaystyle\lambda\in\partial D_{\lambda_{l,k}}^{\epsilon}\,,$
$\displaystyle\mathbf{E}_{\mathrm{s}+,j}^{+}(\lambda;t,z)$
$\displaystyle=\mathbf{E}_{\mathrm{s}+,j}^{-}(\lambda;t,z)\mathbf{P}_{\mathrm{s}+,j}[\mathbf{W}_{j,l,k}^{(1)}]^{-\dagger}(\lambda;t,z)\mathbf{P}_{\mathrm{s}+,j}^{-1}\,,\qquad$
$\displaystyle\lambda\in\partial D_{\overline{\lambda_{l,k}}}^{\epsilon}\,,$

Note that $\mathbf{P}_{\mathrm{s}+,j}(\lambda;t,z)$ is bounded as
$t\to+\infty$. Thus, all jump matrices for
$\mathbf{E}_{\mathrm{s}+,j}(\lambda;t,z)$ in RHP 11 decay to the identity
matrix exponentially and uniformly as $t\to+\infty$. Consequently, one writes

(5.26) $\mathbf{E}_{\mathrm{s}+,j}(\lambda;t,z)=\mathbb{I}+o(1)\,,\qquad
t\to+\infty\,,$

for every fixed $j=1,2,\dots,J$. Because of the boundedness of
$\mathbf{P}_{\mathrm{s}+,j}(\lambda;t,z)$, the definition (5.23) yields

(5.27)
$\mathbf{M}_{j}^{(3)}(\lambda;t,z)=(\mathbb{I}+o(1))\mathbf{P}_{\mathrm{s}+,j}(\lambda;t,z)\,,\qquad
t\to+\infty\,.$

Therefore, the solution to RHP 2 can be reconstructed similarly to Equation
(5.19), as follows,

(5.28)
$\mathbf{M}(\lambda;t,z)=(\mathbb{I}+o(1))\mathbf{P}_{\mathrm{s}+,j}(\lambda;t,z)\delta_{j+1}(\lambda)^{\sigma_{3}}\,,\qquad|\lambda|\gg
1\mbox{ or }|\lambda|\ll 1\,.$

Recall that $\delta_{j+1}(\lambda)\to 1$ as $\lambda\to\infty$, and
$\delta_{j+1}(\lambda)^{\sigma_{3}}$ commutes with $\sigma_{3}$. So, by Lemma
1 and Theorem 2 one concludes

(5.29) $\displaystyle q(t,z;\Lambda,\Omega)$
$\displaystyle=q(t,z;\Lambda_{j},\Omega_{j}^{(+)})+o(1)\,,$ $\displaystyle
D(t,z;\Lambda,\Omega)$
$\displaystyle=D(t,z;\Lambda_{j},\Omega_{j}^{(+)})+o(1)\,,$ $\displaystyle
P(t,z;\Lambda,\Omega)$
$\displaystyle=P(t,z;\Lambda_{j},\Omega_{j}^{(+)})+o(1)\,,$

where the solution
$\\{q_{\mathrm{s}+,j},D_{\mathrm{s}+,j},P_{\mathrm{s}+,j}\\}$ is an
$N_{j}$-DSG derived from RHP 10 with eigenvalues $\Lambda_{j}$ and modified
norming constants $\Omega_{j}^{(+)}$.

###  5.3. Soliton asymptotics as $t\to-\infty$ in a stable medium

Figure 12. Similarly to Figures 10 and 11, but with $t<0$ representing the
four cases in Table 1 . Top row: situations $\xi\leq 0<V_{1}$ (left) and
$0<\xi<V_{1}$ (right) corresponding to Case 1. Bottom left: $\xi>V_{J}$
corresponding to Case 2. Bottom center: $V_{j-1}<\xi<V_{j}$ slightly different
but equivalent to Case 3. Bottom right: $\xi=V_{j}$ corresponding to Case 4.

One needs to consider the four cases in Table 1, but with the limit
$t\to-\infty$. Because $t<0$ instead of $t>0$ from Section 5.2, and because of
$\theta(\lambda;t,z)=\lambda t+\xi t/(2\lambda)$, the quantity
$\Re(\mathrm{i}\theta(t,z))$ flips its sign comparing to the previous
$t\to+\infty$ discussion.

The illustrative plots are given in Figure 12, which should be compared with
Figures 10 and 11. It is then necessary to notice the similarities and
differences between the situations $t\to-\infty$ and $t\to+\infty$. Built on
the detailed discussion on the case of $t\to+\infty$ in Section 5.2, the
results of the $t\to-\infty$ case can be easily obtained, in similar
treatments. As such, most details will be omitted in order to avoid
repetitions.

It is important to point out that we have chosen $V_{j-1}<\xi<V_{j}$ in Figure
12(bottom center), instead of $V_{j}<\xi<V_{j+1}$ presented in Table 1. We
shift the index in order to simplify calculations in Section 5.3.3, as will be
discussed in detail there. The shifted case is equivalent to Case 3 in Table
1. One just needs to write $j=2,3,\dots,J$ instead of $j=1,2,\dots,J-1$ in
Section 5.2.3.

####  5.3.1. The case of $\xi<V_{1}$

There are two subcases $\xi\leq 0<V_{1}$ and $0<\xi<V_{1}$ as shown in Figure
12(top left and top right). Physically, they describe solutions inside and
outside of the light cone, but mathematically, they can be treated in
identical ways. For small enough radii $\epsilon>0$, all jumps in RHP 2 are
decaying exponentially and uniformly to the identity matrix as $t\to-\infty$.
Therefore, the solution can be written as
$\mathbf{M}(\lambda;t,z)=\mathbb{I}+o(1)$ as $t\to-\infty$. Consequently, one
has

(5.30)
$q(t,z)=o(1)\,,\qquad\boldsymbol{\rho}(t,z)=D_{-}\sigma_{3}+o(1)\,,\qquad
t\to-\infty\,.$

####  5.3.2. The case of $\xi>V_{J}$

As can be calculated and also observed from Figure 12(bottom left), for small
enough $\epsilon>0$, all jumps in RHP 2 grow exponentially as $t\to-\infty$.
Thus, one needs to reuse the technique from Section 5.2.1 to order to turn the
growing jumps into decaying ones. We define

(5.31)
$\mathbf{M}^{(5)}(\lambda;t,z)\coloneqq\begin{cases}\mathbf{M}(\lambda;t,z)\delta_{1}^{\sigma_{3}}\,,&\lambda\in\mathbb{C}\setminus\bigg{\\{}\bigcup\limits_{j=1}^{J}\bigcup\limits_{k=1}^{N_{j}}\big{(}D_{\lambda_{j,k}}^{\epsilon}\bigcup
D_{\overline{\lambda_{j,k}}}^{\epsilon}\big{)}\bigg{\\}}\,,\\\
\mathbf{M}(\lambda;t,z)\big{(}\lambda-\overline{\lambda_{j,k}}\big{)}^{\sigma_{3}}(\mathbf{Y}_{j,k}^{(1)}(\lambda;t,z))^{-1}\,,&\lambda\in
D_{\lambda_{j,k}}^{\epsilon}\,,\\\
\mathbf{M}(\lambda;t,z)\big{(}\lambda-\lambda_{j,k}\big{)}^{-\sigma_{3}}(\mathbf{Y}_{j,k}^{(1)})^{\dagger}(\lambda;t,z)\,,&\lambda\in
D_{\overline{\lambda_{j,k}}}^{\epsilon}\,,\end{cases}$

where for $j=1,2,\dots,J$ and $k=1,2,\dots,N_{j}$, and
$\mathbf{Y}_{j,k}^{(1)}(\lambda;t,z)$ is defined in Equation (5.4). Then, the
newly defined matrix function $\mathbf{M}^{(5)}(\lambda;t,z)$ satisfies the
following RHP.

######  Riemann-Hilbert Problem 12.

Let $\mathbf{M}^{(5)}(\lambda;t,z)$ be a $2\times 2$ sectional analytic matrix
function on $\lambda\in\mathbb{C}$. It has asymptotics
$\mathbf{M}^{(5)}(\lambda;t,z)\to\mathbb{I}$ as $\lambda\to\infty$. For
$j=1,2,\dots,J$ and $k=1,2,\dots,N_{j}$, it has jumps

(5.32) $\displaystyle\mathbf{M}^{(5)+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}^{(5)-}(\lambda;t,z)\begin{pmatrix}1&-\omega_{j,k}^{-1}(\lambda-\lambda_{j,k})\frac{\overline{p_{1}}(\lambda_{j,k})^{2}}{p_{1}(\lambda)^{2}}\mathrm{e}^{2\mathrm{i}\theta}\\\
0&1\end{pmatrix}\,,\qquad\lambda\in\partial D_{\lambda_{j,k}}^{\epsilon}\,,$
$\displaystyle\mathbf{M}^{(5)+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}^{(5)-}(\lambda;t,z)\begin{pmatrix}1&0\\\
\overline{\omega_{j,k}}^{-1}(\lambda-\overline{\lambda_{j,k}})\frac{p_{1}(\overline{\lambda_{j,k}})^{2}}{\overline{p_{1}}(\lambda)}\mathrm{e}^{-2\mathrm{i}\theta}&1\end{pmatrix}\,,\qquad\lambda\in\partial
D_{\overline{\lambda_{j,k}}}^{\epsilon}\,.$

Clearly, all jumps are decaying exponentially to the identity matrix as
$t\to-\infty$. Therefore, one concludes

(5.33) $q(t,z)=o(1)\,,\qquad D(t,z)=D_{-}(1+o(1))\,,\qquad
P(t,z)=o(1)\,,\qquad t\to-\infty\,.$

####  5.3.3. The case of $V_{j-1}<\xi<V_{j}$

First of all, let us fix the value of $j=2,3,\dots,J$. One observes from
Figure 12(bottom center) that all jumps inside the red circles are growing
exponentially while the ones outsides are decaying as $t\to-\infty$. Hence,
one needs to modify the jumps with $\mathbf{V}_{l,k}$ with $l=1,2,\dots,j-1$
and $k=1,2,\dots,N_{l}$.

The reason why we have shifted the index from Case 3 in Table 1 is to simplify
expressions in this subsection. In particular, with the original expression
$V_{j}<\xi<V_{j+1}$, one would need to employ quantities
$\delta_{j+1}(\lambda)$ and $p_{j+1}(\lambda)$ in calculations, but with the
shifted one, one can just use $\delta_{j}(\lambda)$ and $p_{j}(\lambda)$ as
shown below.

One defines

(5.34)
$\mathbf{M}_{j}^{(6)}(\lambda;t,z)\coloneqq\begin{cases}\mathbf{M}(\lambda;t,z)\delta_{1}(\lambda)^{-\sigma_{3}}\delta_{j}(\lambda)^{\sigma_{3}}\,,&\lambda\in\mathbb{C}\setminus\bigcup\limits_{l=1}^{j-1}\bigcup\limits_{k=1}^{N_{l}}\big{(}D_{\lambda_{l,k}}^{\epsilon}\bigcup
D_{\overline{\lambda_{l,k}}}^{\epsilon}\big{)}\,,\\\
\mathbf{M}(\lambda;t,z)(\lambda-\overline{\lambda_{l,k}})^{\sigma_{3}}\mathbf{Z}_{j,l,k}^{(1)}(\lambda;t,z)^{-1}\,,&\lambda\in
D_{\lambda_{l,k}}^{\epsilon}\,,\\\
\mathbf{M}(\lambda;t,z)(\lambda-\lambda_{l,k})^{-\sigma_{3}}(\mathbf{Z}_{j,l,k}^{(1)})^{\dagger}(\lambda;t,z)\,,&\lambda\in
D_{\overline{\lambda_{l,k}}}^{\epsilon}\,,\end{cases}$

where for $l=1,2,\dots,j-1$ and $k=1,2,\dots,N_{l}$,

(5.35)
$\mathbf{Z}_{j,l,k}^{(1)}(\lambda;t,z)\coloneqq\begin{pmatrix}(\lambda-\overline{\lambda_{l,k}})\frac{p_{j}\overline{p_{j}}}{p_{1}\overline{p_{1}}}\bigg{[}\frac{\overline{p_{1}}(\lambda)^{2}}{\overline{p_{j}}(\lambda)^{2}}-\frac{\overline{p_{1}}(\lambda_{l,k})^{2}}{\overline{p_{j}}(\lambda_{l,k})^{2}}\bigg{]}&\omega_{l,k}^{-1}\frac{\lambda-\lambda_{l,k}}{\lambda-\overline{\lambda_{l,k}}}\frac{\overline{p_{1}}(\lambda_{l,k})^{2}}{\overline{p_{j}}(\lambda_{l,k})^{2}}\frac{p_{j}\overline{p_{j}}}{p_{1}\overline{p_{1}}}\mathrm{e}^{2\mathrm{i}\theta(\lambda;t,z)}\\\
-\omega_{l,k}\frac{p_{1}\overline{p_{j}}}{\overline{p_{1}}p_{j}}\frac{\lambda-\overline{\lambda_{l,k}}}{\lambda-\lambda_{l,k}}\mathrm{e}^{-2\mathrm{i}\theta(\lambda;t,z)}&\frac{1}{\lambda-\overline{\lambda_{l,k}}}\frac{p_{1}\overline{p_{j}}}{\overline{p_{1}}p_{j}}\end{pmatrix}\,.$

Then, the new function $\mathbf{M}_{j}^{(6)}(\lambda;t,z)$ satisfies the
following RHP.

######  Riemann-Hilbert Problem 13.

Let $\mathbf{M}_{j}^{(6)}(\lambda;t,z)$ be a $2\times 2$ sectional analytic
matrix function on $\lambda\in\mathbb{C}$, with asymptotics
$\mathbf{M}_{j}^{(6)}(\lambda;t,z)\to\mathbb{I}$ as $\lambda\to\infty$, and
satisfying the following jumps, for $l=1,2,\dots,j-1$ and $k=1,2,\dots,N_{l}$,

(5.36) $\displaystyle\mathbf{M}_{j}^{(6)+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}_{j}^{(6)-}(\lambda;t,z)\delta_{1}^{\sigma_{3}}\delta_{j}^{-\sigma_{3}}\mathbf{V}_{l,k}(\lambda;t,z)^{-1}(\lambda-\overline{\lambda_{l,k}})^{\sigma_{3}}\mathbf{Z}_{j,l,k}^{(1)}(\lambda;t,z)^{-1}\,,\qquad\lambda\in\partial
D_{\lambda_{l,k}}^{\epsilon}\,,$
$\displaystyle\mathbf{M}_{j}^{(6)+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}_{j}^{(6)-}(\lambda;t,z)\delta_{1}^{\sigma_{3}}\delta_{j}^{-\sigma_{3}}\mathbf{V}_{l,k}^{\dagger}(\lambda;t,z)(\lambda-\lambda_{l,k})^{-\sigma_{3}}(\mathbf{Z}_{j,l,k}^{(1)})^{\dagger}(\lambda;t,z)\,,\qquad\lambda\in\partial
D_{\overline{\lambda_{l,k}}}^{\epsilon}\,,\,,$

and for $l=j,j+1,\dots,J$ and $k=1,2,\dots,N_{l}$,

(5.37) $\displaystyle\mathbf{M}_{j}^{(6)+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}_{j}^{(6)-}(\lambda;t,z)\delta_{1}^{\sigma_{3}}\delta_{j}^{-\sigma_{3}}\mathbf{V}_{l,k}(\lambda;t,z)^{-1}\delta_{j}^{\sigma_{3}}\delta_{1}^{-\sigma_{3}}\,,\qquad\lambda\in\partial
D_{\lambda_{l,k}}^{\epsilon}\,,$
$\displaystyle\mathbf{M}_{j}^{(6)+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}_{j}^{(6)-}(\lambda;t,z)\delta_{1}^{\sigma_{3}}\delta_{j}^{-\sigma_{3}}\mathbf{V}_{l,k}^{\dagger}(\lambda;t,z)\delta_{j}^{\sigma_{3}}\delta_{1}^{-\sigma_{3}}\,,\qquad\lambda\in\partial
D_{\overline{\lambda_{l,k}}}^{\epsilon}\,,$

where $\mathbf{V}_{l,k}(\lambda;t,z)$ is given in RHP 2.

In particular, some useful jumps are given explicitly below

(5.38)
$\displaystyle\delta_{1}^{\sigma_{3}}\delta_{j}^{-\sigma_{3}}\mathbf{V}_{l,k}(\lambda;t,z)^{-1}(\lambda-\overline{\lambda_{l,k}})^{\sigma_{3}}\mathbf{Z}_{j,l,k}^{(1)}(\lambda;t,z)^{-1}$
$\displaystyle=\begin{pmatrix}1&-\frac{\lambda-\lambda_{l,k}}{\omega_{l,k}}\frac{\overline{p_{1}}(\lambda_{l,k})^{2}}{\overline{p_{j}}(\lambda_{l,k})^{2}}\frac{p_{j}(\lambda)^{2}}{p_{1}(\lambda)^{2}}\mathrm{e}^{2\mathrm{i}\theta}\\\
0&1\end{pmatrix}\,,$
$\displaystyle\delta_{1}^{\sigma_{3}}\delta_{j}^{-\sigma_{3}}\mathbf{V}_{l,k}^{\dagger}(\lambda;t,z)(\lambda-\lambda_{l,k})^{-\sigma_{3}}(\mathbf{Z}_{j,l,k}^{(1)})^{\dagger}(\lambda;t,z)$
$\displaystyle=\begin{pmatrix}1&0\\\
\frac{\lambda-\overline{\lambda_{l,k}}}{\overline{\omega_{l,k}}}\frac{p_{1}(\overline{\lambda_{l,k}})^{2}}{p_{j}(\overline{\lambda_{l,k}})^{2}}\frac{\overline{p_{j}}(\lambda)^{2}}{\overline{p_{1}}(\lambda)^{2}}\mathrm{e}^{-2\mathrm{i}\theta}&1\end{pmatrix}\,.$

It is easy to verify that all jumps in RHP 13 decay exponentially to the
identity matrix as $t\to-\infty$, so one can write
$\mathbf{M}_{j}^{(6)}(\lambda;t,z)=\mathbb{I}(1+o(1))$. Therefore,

(5.39) $q(t,z)=o(1)\,,\qquad P(t,z)=o(1)\,,\qquad D(t,z)=D_{-}+o(1)\,,\qquad
t\to-\infty\,.$

The calculation is valid for every $j=2,3,\dots,J$, so we have proved that
when $\xi$ is between adjacent velocities, the solution decays to the ZBG as
$t\to-\infty$.

####  5.3.4. The case of $\xi=V_{j}$

First, let us fix the value of $j=1,2,\dots,J$. According to Figure 12(bottom
right), this case is analogous to the one in Section 5.2.4, in the sense that
all the jumps of RHP 2 with indices $l<j$ are growing as $t\to-\infty$ and the
ones with $l>j$ are decaying to the identity matrix in the asymptotics.
Different from the previous case, here the jumps with $l=j$ stay bounded. The
steps to treat the growing jumps is identical to the previous case, so we omit
the calculations and start with RHP 13, which is similar to what happens in
Section 5.2.4. As a result, all jumps with indices $l\neq j$ decay
exponentially to the identity matrix as $t\to-\infty$. By ignoring all
decaying jumps, one arrives at the following leading-order problem.

######  Riemann-Hilbert Problem 14 (Jump form).

Let $j=1,2,\dots,J$ be fixed. Let $\mathbf{P}_{\mathrm{s}-,j}(\lambda;t,z)$ be
a $2\times 2$ sectional analytical matrix function on $\mathbb{C}$.
$\mathbf{P}_{\mathrm{s}-,j}(\lambda;t,z)\to\mathbb{I}$ as $\lambda\to\infty$.
This matrix function also satisfies the following jumps

(5.40) $\displaystyle\mathbf{P}_{\mathrm{s}-,j}^{+}$
$\displaystyle=\mathbf{P}_{\mathrm{s}-,j}^{-}\delta_{1}^{\sigma_{3}}\delta_{j}^{-\sigma_{3}}\mathbf{V}_{j,k}^{-1}\delta_{j}^{\sigma_{3}}\delta_{1}^{-\sigma_{3}}=\mathbf{P}_{\mathrm{s}-,j}^{-}\begin{pmatrix}1&0\\\
-\frac{\omega_{j,k}}{\lambda-\lambda_{j,k}}\frac{\delta_{j}(\lambda)^{2}}{\delta_{1}(\lambda)^{2}}\mathrm{e}^{-2\mathrm{i}\theta}&1\end{pmatrix}\,,\qquad\lambda\in\partial
D_{\lambda_{j,k}}^{\epsilon}\,,$ $\displaystyle\mathbf{P}_{\mathrm{s}-,j}^{+}$
$\displaystyle=\mathbf{P}_{\mathrm{s}-,j}^{-}\delta_{1}^{\sigma_{3}}\delta_{j}^{-\sigma_{3}}\mathbf{V}_{j,k}^{\dagger}\delta_{j}^{\sigma_{3}}\delta_{1}^{-\sigma_{3}}=\mathbf{P}_{\mathrm{s}-,j}^{-}\begin{pmatrix}1&\frac{\overline{\omega_{j,k}}}{\lambda-\overline{\lambda_{j,k}}}\frac{\delta_{1}(\lambda)^{2}}{\delta_{j}(\lambda)^{2}}\mathrm{e}^{2\mathrm{i}\theta}\\\
0&1\end{pmatrix}\,,\qquad\lambda\in\partial
D_{\overline{\lambda_{j,k}}}^{\epsilon}\,,$

The subscript “$\mathrm{s}-$” denote the asymptotics in a stable medium as
$t\to-\infty$. If needed, the above RHP can be turned in to a residue form.
According to Theorem 2, solution $\mathbf{P}_{\mathrm{s}-,j}(\lambda;t,z)$ is
bounded as $t\to-\infty$. One then considers the following error function,

(5.41)
$\mathbf{E}_{\mathrm{s}-,j}(\lambda;t,z)\coloneqq\mathbf{M}_{j}^{(6)}(\lambda;t,z)\mathbf{P}_{\mathrm{s}-,j}(\lambda;t,z)^{-1}\,,\qquad\mbox{
for a fixed }j=1,2,\dots,J\,,$

Similarly to the matrix $\mathbf{E}_{\mathrm{s}+,j}(\lambda;t,z)$ from
Equation (5.23), it is straightforward to verify that
$\mathbf{E}_{\mathrm{s}-,j}(\lambda;t,z)$ satisfies a RHP whose jump matrices
decay uniformly to the identity matrix as $t\to-\infty$. Hence, one writes
$\mathbf{E}_{\mathrm{s}-,j}(\lambda;t,z)=\mathbb{I}+o(1)$, as $t\to-\infty$.
Consequently, it is
$\mathbf{M}_{j}^{(6)}(\lambda;t,z)=(\mathbb{I}+o(1))\mathbf{P}_{\mathrm{s}-,j}(\lambda;t,z)$,
as $t\to-\infty$. Tracing back all the deformations in Section 5.3.3, one
finally obtains as $t\to-\infty$,
$\mathbf{M}(\lambda;t,z)=(\mathbb{I}+o(1))\mathbf{P}_{\mathrm{s}-,j}(\lambda;t,z)\delta_{1}(\lambda)^{\sigma_{3}}\delta_{j}(\lambda)^{-\sigma_{3}}$,
with $|\lambda|\gg 1\mbox{ or }|\lambda|\ll 1$. Recall $\delta_{1}(\lambda)\to
1$ and $\delta_{j}(\lambda)\to 1$ as $\lambda\to\infty$, and
$\delta_{1}(\lambda)^{\sigma_{3}}\delta_{j}(\lambda)^{-\sigma_{3}}$ commutes
with $\sigma_{3}$, Lemma 1 and Theorem 2 yields

(5.42) $\displaystyle q(t,z;\Lambda,\Omega)$
$\displaystyle=q(t,z;\Lambda_{j},\Omega_{j}^{(-)})+o(1)\,,$ $\displaystyle
D(t,z;\Lambda,\Omega)$
$\displaystyle=D(t,z;\Lambda_{j},\Omega_{j}^{(-)})+o(1)\,,$ $\displaystyle
P(t,z;\Lambda,\Omega)$
$\displaystyle=P(t,z;\Lambda_{j},\Omega_{j}^{(-)})+o(1)\,,$

where the solution
$\\{q_{\mathrm{s}-,j},D_{\mathrm{s}-,j},P_{\mathrm{s}-,j}\\}$ is an
$N_{j}$-DSG derived from RHP 14 with eigenvalues $\Lambda_{j}$ and modified
norming constants $\Omega_{j}^{(-)}\coloneqq\\{\omega_{j,k}^{(-)}\\}$ defined
in Theorem 3.

###  5.4. Soliton asymptotics in an unstable medium

It is time to discuss the soliton asymptotics in an unstable medium with
$J\geq 2$.

Table 2. Four cases in the long-time asymptotics in an unstable medium. Unstable Cases |  Relations between $\xi$ and velocities |  Description
---|---|---
Case 1 | $\vphantom{\Big{|}}\xi<V_{\mathrm{u},J}<\dots<V_{\mathrm{u},1}$ |  $\xi$ is smaller than all soliton velocities
Case 2 | $\vphantom{\Big{|}}V_{\mathrm{u},J}<\dots<V_{\mathrm{u},1}<\xi$ |  $\xi$ is larger than all soliton velocities
Case 3 |  $\vphantom{\Big{|}}V_{\mathrm{u},j+1}<\xi<V_{\mathrm{u},j}$ with $j=1,\dots,J-1$ |  $\xi$ is in-between adjacent velocities
Case 4 |  $\vphantom{\Big{|}}V_{\mathrm{u},j+1}<\xi=V_{\mathrm{u},j}<V_{\mathrm{u},j-1}$ with $j=1,\dots,J$ |  $\xi$ coincides with velocity $V_{\mathrm{u},j}$

According to Theorem 2, let us denote $V_{\mathrm{s},j}=2r_{j}^{2}$ as the
velocity in a stable medium, while
$V_{\mathrm{u},j}=-2r_{j}^{2}=-V_{\mathrm{s},j}$ as the one in an unstable
medium. The crucial component in the Deift-Zhou’s nonlinear steepest descent
method is the sign of $\Re(\mathrm{i}\theta(\lambda_{j,k};t,z))$, which are
given in the stable and unstable case, respectively, as follows,

(5.43)
$\displaystyle\Re(\mathrm{i}\theta_{\mathrm{s},j,k}(\xi,t))=\frac{y_{j,k}}{2r_{j}^{2}}(\xi-
V_{\mathrm{s},j})t\,,\qquad\Re(\mathrm{i}\theta_{\mathrm{u},j,k}(\xi,t))=-\frac{y_{j,k}}{2r_{j}^{2}}(\xi-
V_{\mathrm{u},j})t\,,$

where $y_{j,k}=\Im(\lambda_{j,k})>0$ and $z=\xi t$. As a result, it can be
seen

(5.44)
$\Re(\mathrm{i}\theta_{\mathrm{u},j,k}(\xi,t))=\Re(\mathrm{i}\theta_{\mathrm{s},j,k}(-\xi,t))\,.$

Therefore, the asymptotics in an unstable medium along the line $z=\xi t$ can
be regarded as asymptotics in a stable medium along a different line $z=-\xi
t$. This relation allows us to compute the asymptotics in the unstable case
easily from known results of the stable case.

Clearly, the soliton velocities $V_{\mathrm{u},j}=-2r_{j}^{2}$ defined in
Theorem 2 are ordered below,

(5.45) $V_{J}<V_{J-1}<\dots<V_{1}<0\,.$

Similarly to what happens in a stable, there are four fundamental cases given
in Table 2.

In the asymptotics $t\to+\infty$, one can relate the cases between stable
medium (cf. Table 1) and unstable medium (cf. Table 2) below.

  1. Unstable Case 1

Equation (5.44) implies that the unstable Case 1
($\xi<V_{\mathrm{u},J}<\dots<V_{\mathrm{u},1}$) can be calculated equivalently
to the stable Case 2 ($V_{\mathrm{s},1}<\dots<V_{\mathrm{s},J}<-\xi$) in the
same asymptotics $t\to+\infty$, where $z=\xi t$ with $\xi<0$. Hence, one
arrives at $q(t,z)=o(1)$ and
$\boldsymbol{\rho}(t,z)=D_{-}\sigma_{3}(\mathbb{I}+o(1))$ following Section
5.2.2.

  2. Unstable Case 2

Equation (5.44) implies that the unstable Case 2
($V_{\mathrm{u},J}<\dots<V_{\mathrm{u},1}<\xi$) is equivalent to the stable
Case 1 ($-\xi<V_{\mathrm{s},1}<\dots<V_{\mathrm{s},J}$). Therefore, one
arrives at the same asymptotic result $q(t,z)=o(1)$ and
$\boldsymbol{\rho}(t,z)=D_{-}\sigma_{3}(\mathbb{I}+o(1))$ as $t\to+\infty$
following Section 5.2.1.

  3. Unstable Case 3

Equation (5.44) implies that the unstable Case 3
($V_{\mathrm{u},j+1}<\xi<V_{\mathrm{u},j}$) is equivalent to the stable Case 3
($V_{\mathrm{s},j}<-\xi<V_{\mathrm{s},j+1}$). Therefore, one arrives at the
same asymptotic result $q(t,z)=o(1)$ and
$\boldsymbol{\rho}(t,z)=D_{-}\sigma_{3}(\mathbb{I}+o(1))$ as $t\to+\infty$
following Section 5.2.3.

  4. Unstable Case 4

Equation (5.44) implies that the unstable Case 4 ($\xi=V_{\mathrm{u},j}$) is
equivalent to the stable Case 4 ($-\xi=V_{\mathrm{s},j}$). Therefore, by
repeating similar calculations as in Section 5.2.4, one arrives at the formula
for the asymptotic solution with $\Lambda_{j}$ and $\Omega_{j}^{(+)}$ in
Theorem 3.

The other asymptotics $t\to-\infty$ between the stable and unstable cases can
also be related. The details are omitted for brevity. The asymptotic DSG is
given in Theorem 3. Finally, this completes the proof of Theorem 3 for all
cases.

##  6\. High-order solitons

In this section, we discuss various aspects about the high-order solitons.
First, we show how to fuse $N\geq 2$ simple poles of RHP 1 in order to get an
$N$th order pole. As a result, the corresponding $N$-soliton solution merge
into a single $N$th-order soliton.

###  6.1. Derivation of an $N$th order soliton from fusion

Suppose RHP 3 with the eigenvalue $\lambda_{\circ}$ and norming constants
$\\{\omega_{\circ,k}\\}_{k=0}^{N-1}$ given. We would like to find an explicit
limiting procedure, such that when merging all eigenvalues $\Lambda$ from RHP
2, the $N$ simple poles $\omega_{j,k}/(\lambda-\lambda_{j,k})$ fuse into the
$N$th-order pole $p_{\circ}(\lambda)/(\lambda-\lambda_{\circ})^{N}$ appearing
in RHP 3.

Therefore, we start with RHP 2 for the $N$-soliton solutions. To better take
limits, we rewrite this RHP in an equivalent form, in which the jump contours
$\partial D_{\lambda_{j,k}}^{\epsilon}$ are modified into $\partial
D_{\Lambda}^{\epsilon}$, which is a single contour surrounding all eigenvalues
in the upper half plane. Moreover, we would like to relabel all eigenvalues as
$\lambda_{j,k}\mapsto\lambda_{s}$ for $s=1,2,\dots,N$, and same for the
norming constants $\omega_{j,k}\mapsto\omega_{s}$. As a result, the new RHP is
given by below.

######  Riemann-Hilbert Problem 15 (Equivalent jump-form of RHP 2).

Let $\mathbf{M}(\lambda;t,z)$ be a sectional analytic $2\times 2$ matrix
function, with asymptotics $\mathbf{M}(\lambda;t,z)\to\mathbb{I}$ as
$\lambda\to\infty$ and jumps

(6.1) $\displaystyle\mathbf{M}^{+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}^{-}(\lambda;t,z)\mathbf{V}(\lambda;t,z)^{-1}\,,\qquad$
$\displaystyle\partial D_{\Lambda}^{\epsilon}\,,$
$\displaystyle\mathbf{M}^{+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}^{-}(\lambda;t,z)\mathbf{V}^{\dagger}(\lambda;t,z)\,,\qquad$
$\displaystyle\partial D_{\overline{\Lambda}}^{\epsilon}\,,$

where

(6.2) $\mathbf{V}\coloneqq\begin{pmatrix}1&0\\\
\sum\limits_{s=1}^{N}\frac{\omega_{j}}{\lambda-\lambda_{j}}\mathrm{e}^{-2\mathrm{i}\theta}&1\end{pmatrix}\,.$

Recall that we would like to take limits $\lambda_{s}\to\lambda_{\circ}$, so
it is natural to write

(6.3) $\lambda_{s}=\lambda_{\circ}+\epsilon_{s}\,,$

where the target eigenvalue $\lambda_{\circ}$ is a fixed complex number in the
upper half plane. As a result, the relevant quantity in the above RHP becomes

(6.4)
$\sum_{s=1}^{N}\frac{\omega_{s}}{\lambda-\lambda_{s}}=\sum_{s=1}^{N}\frac{\omega_{s}}{\lambda-\lambda_{\circ}-\epsilon_{s}}\,.$

We next show that the fusion of eigenvalues must be done in a special way.
Otherwise, the result is quite trivial.

######  Remark 10 (A naive limit).

We demonstrate here that taking the limit $\epsilon_{s}\to 0$ for all
$s=1,2,\dots,N$ directly will yield a trivial result. Note that this naive
limit does not rescale the norming constants $\omega_{s}$. Then,

(6.5)
$\sum_{s}\frac{\omega_{s}}{\lambda-\lambda_{s}}=\sum_{s}\frac{\omega_{s}}{\lambda-\lambda_{\circ}-\epsilon_{s}}\to\sum_{s}\frac{\omega_{s}}{\lambda-\lambda_{\circ}}=\frac{\sum_{s}\omega_{s}}{\lambda-\lambda_{\circ}}\,,$

as $\epsilon_{s}\to 0$ for all $s$. So, the $N$ simple poles
$\frac{\omega_{s}}{\lambda-\lambda_{s}}$ fuse into one simple pole at
$\lambda=\lambda_{\circ}$. Correspondingly, the $N$-soliton solution becomes
an one-soiton solution, with its eigenvalue $\lambda=\lambda_{\circ}$ and
norming constant $\sum_{s}\omega_{s}$. Thus, in order to get the $N$th-order
soliton, one has to rescale the norming constants $\omega_{s}$, i.e., assuming
$\omega_{s}=\omega_{s}(\epsilon_{1},\dots\epsilon_{N})$. It remains to show
that such rescaling exists.

Let us calculate the sum from Equation (6.4) more carefully

(6.6)
$\sum_{s=1}^{N}\frac{\omega_{s}}{\lambda-\lambda_{s}}=\frac{\sum_{s=1}^{N}\omega_{s}p_{s}(\lambda)}{\prod_{s=1}^{N}(\lambda-\lambda_{s})}\,,\quad\mbox{where
we define}\quad
p_{s}(\lambda)\coloneqq\frac{\prod_{l=1}^{N}(\lambda-\lambda_{l})}{\lambda-\lambda_{s}}\,.$

Note that the quantity $p_{s}(\lambda)$ is a polynomial of $\lambda$ of degree
$N-1$. Thus, $\sum_{s}\omega_{s}p_{s}(\lambda)$ is also a polynomial of degree
at most $N-1$. The leading term of $p_{s}(\lambda)$ is simply $\lambda^{N-1}$,
and the leading term of $\sum_{s}\omega_{s}p_{s}(\lambda)$ is
$\lambda^{N-1}\sum_{s}\omega_{s}$. For now, let us assume
$\sum_{s}\omega_{s}\in\mathbb{C}\setminus\\{0\\}$, so that the degree of
$\sum_{s}\omega_{s}p_{s}(\lambda)$ is exactly $N-1$.

######  Lemma 2.

For the given polynomial $p_{\circ}(\lambda)$ from RHP 3, the following
equation for unknowns $\\{\omega_{s}\\}_{s=1}^{N}$

(6.7) $\sum_{s=1}^{N}\omega_{s}p_{s}(\lambda)=p_{\circ}(\lambda)$

with $\lambda\in\mathbb{C}$ has a unique solution.

###### Proof.

Recall that the given polynomial $p_{\circ}(\lambda)$ from RHP 3 has the
coefficients $\\{\omega_{\circ,k}\\}_{k=0}^{N-1}$ with
$\omega_{\circ,k}\in\mathbb{C}$ and $\omega_{\circ,0}\neq 0$. We show that the
unknowns $\\{\omega_{l}\\}_{l=1}^{N}$ can be solved from Equation (6.7)
explicitly. It is important to notice that the polynomial $p_{s}(\lambda)$ of
degree $N-1$ has roots
$\\{\lambda_{1},\dots,\lambda_{s-1},\lambda_{s+1},\dots,\lambda_{N}\\}$. Thus,
substituting $\lambda=\lambda_{s}$ into Equation (6.7) yields
$\omega_{s}p_{s}(\lambda_{s})=p_{\circ}(\lambda_{s})$, which implies

(6.8)
$\omega_{s}=\frac{p_{\circ}(\lambda_{s})}{p_{s}(\lambda_{s})}=\frac{\sum_{k=0}^{N-1}\omega_{\circ,k}\epsilon_{s}^{k}}{\prod_{k\neq
s}(\epsilon_{s}-\epsilon_{k})}\,,\qquad s=1,2,\dots,N-1\,.$

Hence, the unknowns $\\{\omega_{1},\dots,\omega_{N}\\}$ are solved explicitly
and uniquely. ∎

Lemma 2 implies that the rescaling
$\omega_{s}=\omega_{s}(\epsilon_{1},\dots,\epsilon_{N})$ exists. Using this,
we next show that the $N$-soliton solution can be transformed into an $N$th-
oder soliton. Lemma 2 implies that Equation (6.6) becomes
$\sum_{s=1}^{N}\frac{\omega_{s}}{\lambda-\lambda_{s}}=p_{\circ}(\lambda)/\prod_{s=1}^{N}(\lambda-\lambda_{s})$.
Then, taking limits $\lambda_{s}\to\lambda_{\circ}$ for all $s=1,2,\dots,N$
with the particular $\omega_{s}$ from Equation (6.8), the sum becomes
$\sum_{s=1}^{N}\frac{\omega_{s}}{\lambda-\lambda_{s}}\to
p_{\circ}(\lambda)/(\lambda-\lambda_{\circ})^{N}$. One can treat the poles in
the lower half plane in a similar way, what we omit here. Hence, RHP 15
immediately becomes RHP 3, where the contour $\partial D_{\Lambda}^{\epsilon}$
becomes $\partial D_{\lambda_{\circ}}^{\epsilon}$ in the upper half plane.

Finally, we conclude that by taking the limit $\lambda_{s}\to\lambda_{\circ}$
for all $s=1,2,\dots,N$ of the $N$-soliton solution of MBEs, with particular
rescaling of the norming constants $\\{\omega_{s}\\}_{s=1}^{N}$ described in
Equation (6.8), the $N$-soliton solution becomes an $N$th-order soliton,
described by RHP 3.

###  6.2. Derivation of the pole form of RHP 3

In this section, we rewrite RHP 3 for the $N$th-order soliton in an equivalent
form, where the jump conditions become appropriate pole conditions, mimicking
the residue conditions in RHP 1 for the $N$-soliton solutions. To start, let
use define a new matrix function

(6.9)
$\mathbf{N}_{\circ}(\lambda;t,z)\coloneqq\begin{cases}\mathbf{M}_{\circ}(\lambda;t,z)\,,&\lambda\in\mathbb{C}\setminus\\{D_{\lambda_{\circ}}^{\epsilon},D_{\overline{\lambda_{\circ}}}^{\epsilon}\\}\,,\\\
\mathbf{M}_{\circ}(\lambda;t,z)\mathbf{V}_{\circ}(\lambda;t,z)\,,&\lambda\in
D_{\lambda_{\circ}}^{\epsilon}\,,\\\
\mathbf{M}_{\circ}(\lambda;t,z)\mathbf{V}_{\circ}^{-\dagger}(\lambda;t,z)\,,&\lambda\in
D_{\overline{\lambda_{\circ}}}^{\epsilon}\,,\end{cases}$

where $\mathbf{M}_{\circ}(\lambda;t,z)$ and $\mathbf{V}_{\circ}(\lambda;t,z)$
are defined in RHP 3. It is easy to verify that:
$\mathbf{N}_{\circ}(\lambda;t,z)\to\infty$ as $\lambda\to\infty$;
$\mathbf{N}_{\circ}(\lambda;t,z)$ does not admit any jumps on $\partial
D_{\lambda_{\circ}}^{\epsilon}$ and $\partial
D_{\overline{\lambda_{\circ}}}^{\epsilon}$; and
$\mathbf{N}_{\circ}(\lambda;t,z)$ is meromorphic on $\mathbb{C}$, with poles
at $\lambda=\lambda_{\circ}$ and $\lambda=\overline{\lambda_{\circ}}$. The
only thing left is to derive the pole conditions for
$\mathbf{N}_{\circ}(\lambda;t,z)$. Recall the pole operator
$\mathop{\mathcal{P}_{n}}$ from Definition 2, which can be used to compute the
Laurent coefficient of $\mathbf{N}_{\circ}(\lambda;t,z)$ at
$\lambda=\lambda_{\circ}$ and $\lambda=\overline{\lambda_{\circ}}$. Here, we
show detailed steps for the case of $\lambda=\lambda_{\circ}$. The other point
can be addressed similarly.

(6.10)
$\displaystyle\mathop{\mathcal{P}_{-n}}_{\lambda=\lambda_{\circ}}\mathbf{N}_{\circ}(\lambda;t,z)=\mathop{\mathcal{P}_{-n}}_{\lambda=\lambda_{\circ}}\begin{pmatrix}\mathbf{M}_{\circ,1}+\frac{p_{\circ}(\lambda)}{(\lambda-\lambda_{\circ})^{N}}\mathrm{e}^{-2\mathrm{i}\theta}\mathbf{M}_{\circ,2}&\mathbf{M}_{\circ,2}\end{pmatrix}\,,$

where the subscripts $1$ and $2$ denote the first and second columns of the
matrix, respectively, and $n=1,2,\dots,N$. Also, we suppress the variable
dependence of $\mathbf{M}_{\circ}$ in the above and below equations for
simplicity. Recall that $\mathbf{M}_{\circ}$ is analytic at
$\lambda=\lambda_{\circ}$, so by distributing $\mathop{\mathcal{P}_{-n}}$ in
each column, one gets

(6.11)
$\displaystyle\mathop{\mathcal{P}_{-n}}_{\lambda=\lambda_{\circ}}\mathbf{N}_{\circ}(\lambda;t,z)$
$\displaystyle=\mathop{\mathcal{P}_{-n}}_{\lambda=\lambda_{\circ}}\mathbf{M}_{\circ}\begin{pmatrix}0&0\\\
\frac{p_{0}(\lambda)}{(\lambda-\lambda_{\circ})^{N}}\mathrm{e}^{-2\mathrm{i}\theta}&0\end{pmatrix}\,,\qquad
n=1,2,\dots,N\,.$

It is necessary to calculate
$\mathop{\mathcal{P}_{-n}}\limits_{\lambda=\lambda_{\circ}}\mathbf{M}_{\circ}\frac{p_{\circ}(\lambda)}{(\lambda-\lambda_{\circ})^{N}}\mathrm{e}^{-2\mathrm{i}\theta}$.
The analyticity of
$p_{0}(\lambda)\mathrm{e}^{-2\mathrm{i}\theta(\lambda;t,z)}\mathbf{M}_{\circ}(\lambda;t,z)$
at $\lambda=\lambda_{\circ}$ yields
$\mathop{\mathcal{P}_{-n}}\limits_{\lambda=\lambda_{\circ}}\frac{p_{\circ}(\lambda)\mathrm{e}^{-2\mathrm{i}\theta}}{(\lambda-\lambda_{\circ})^{N}}\mathbf{M}_{\circ}=\mathop{\mathcal{P}_{N-n}}\limits_{\lambda=\lambda_{\circ}}p_{\circ}(\lambda)\mathrm{e}^{-2\mathrm{i}\theta}\mathbf{M}_{\circ}$,
for $n=1,2,\dots,N$. Using the expression for $p_{\circ}(\lambda)$ in Equation
(2.27), the computation continues

(6.12)
$\displaystyle\mathop{\mathcal{P}_{-n}}_{\lambda=\lambda_{\circ}}\frac{p_{\circ}(\lambda)\mathrm{e}^{-2\mathrm{i}\theta}}{(\lambda-\lambda_{\circ})^{N}}\mathbf{M}_{\circ}$
$\displaystyle=\sum_{i=0}^{N-n}\mathop{\mathcal{P}_{N-n}}_{\lambda=\lambda_{\circ}}\omega_{\circ,i}(\lambda-\lambda_{\circ})^{i}\mathrm{e}^{-2\mathrm{i}\theta}\mathbf{M}_{\circ}+\sum_{i=N-n+1}^{N-1}\mathop{\mathcal{P}_{N-n}}_{\lambda=\lambda_{\circ}}\omega_{\circ,i}(\lambda-\lambda_{\circ})^{i}\mathrm{e}^{-2\mathrm{i}\theta}\mathbf{M}_{\circ}$

Clearly, terms in the second summation are all zeros. Thus, the equation
continues

(6.13)
$\displaystyle\mathop{\mathcal{P}_{-n}}_{\lambda=\lambda_{\circ}}\frac{p_{\circ}(\lambda)\mathrm{e}^{-2\mathrm{i}\theta}}{(\lambda-\lambda_{\circ})^{N}}\mathbf{M}_{\circ}=\sum_{i=0}^{N-n}\mathop{\mathcal{P}_{N-n}}_{\lambda=\lambda_{\circ}}\omega_{\circ,i}(\lambda-\lambda_{\circ})^{i}\mathrm{e}^{-2\mathrm{i}\theta}\mathbf{M}_{\circ}=\sum_{i=0}^{N-n}\omega_{\circ,i}\mathop{\mathcal{P}_{N-n-i}}_{\lambda=\lambda_{\circ}}\mathrm{e}^{-2\mathrm{i}\theta}\mathbf{M}_{\circ}\,,$

for $n=1,2,\dots,N$. Due to analyticity, the Taylor expansion of the product
$\mathrm{e}^{-2\mathrm{i}\theta}\mathbf{M}_{\circ}$ at
$\lambda=\lambda_{\circ}$ is

(6.14)
$\mathrm{e}^{-2\mathrm{i}\theta}\mathbf{M}_{\circ}=\sum_{k=0}^{\infty}\frac{(\lambda-\lambda_{\circ})^{k}}{k!}\frac{\partial^{k}}{\partial\lambda^{k}}(\mathrm{e}^{-2\mathrm{i}\theta}\mathbf{M}_{\circ})_{\lambda_{\circ}}\,.$

Thus, each Laurent coefficient in Equation (6.13) can be calculated as (with
help of the product rule)

(6.15)
$\mathop{\mathcal{P}_{N-n-i}}_{\lambda=\lambda_{\circ}}\mathrm{e}^{-2\mathrm{i}\theta}\mathbf{M}_{\circ}=\frac{1}{(N-n-i)!}\sum_{k=0}^{N-n-i}{{N-n-i}\choose
k}\frac{\partial^{N-n-i-k}}{\partial\lambda^{N-n-
i-k}}(\mathrm{e}^{-2\mathrm{i}\theta})_{\lambda_{\circ}}\frac{\partial^{k}}{\partial\lambda^{k}}(\mathbf{M}_{\circ})_{\lambda_{\circ}}\,.$

Therefore, one can continue the calculation of the pole condition from
Equation (6.11)

(6.16)
$\displaystyle\mathop{\mathcal{P}_{-n}}_{\lambda=\lambda_{\circ}}\mathbf{N}_{\circ}$
$\displaystyle=\mathop{\mathcal{P}_{-n}}_{\lambda=\lambda_{\circ}}\mathbf{M}_{\circ}\frac{p_{\circ}(\lambda)}{(\lambda-\lambda_{\circ})^{N}}\mathrm{e}^{-2\mathrm{i}\theta}\begin{pmatrix}0&0\\\
1&0\end{pmatrix}$
$\displaystyle=\sum_{i=0}^{N-n}\omega_{\circ,i}\mathop{\mathcal{P}_{N-n-i}}_{\lambda=\lambda_{\circ}}\mathrm{e}^{-2\mathrm{i}\theta}\mathbf{M}_{\circ}\begin{pmatrix}0&0\\\
1&0\end{pmatrix}$
$\displaystyle=\sum_{i=0}^{N-n}\frac{\omega_{\circ,i}}{(N-n-i)!}\sum_{k=0}^{N-n-i}{{N-n-i}\choose
k}\frac{\partial^{N-n-i-k}}{\partial\lambda^{N-n-
i-k}}(\mathrm{e}^{-2\mathrm{i}\theta})_{\lambda_{\circ}}\frac{\partial^{k}}{\partial\lambda^{k}}(\mathbf{M}_{\circ})_{\lambda_{\circ}}\begin{pmatrix}0&0\\\
1&0\end{pmatrix}\,.$

Next, we need to convert the derivatives of $\mathbf{M}_{\circ}$ back to the
ones of $\mathbf{N}_{\circ}$ on the right-hand side of the above equation.
Recall that near $\lambda=\lambda_{\circ}$, it is
$\mathbf{M}_{\circ}=\mathbf{N}_{\circ}\mathbf{V}_{\circ}^{-1}$, so

(6.17)
$\frac{\partial^{k}}{\partial\lambda^{k}}(\mathbf{M}_{\circ})=\frac{\partial^{k}}{\partial\lambda^{k}}(\mathbf{N}_{\circ}\mathbf{V}_{\circ}^{-1})=\sum_{j=0}^{k}{k\choose
j}\frac{\partial^{j}}{\partial\lambda^{j}}(\mathbf{N}_{\circ})\frac{\partial^{k-j}}{\partial\lambda^{k-j}}(\mathbf{V}_{\circ}^{-1})\,,\qquad\lambda\in
D_{\lambda_{\circ}}^{\epsilon}\setminus\\{\lambda_{\circ}\\}\,.$

Recall also that $\mathbf{V}_{\circ}^{-1}=\begin{pmatrix}1&0\\\
-\frac{p_{\circ}\mathrm{e}^{-2\mathrm{i}\theta}}{(\lambda-\lambda_{\circ})^{N}}&1\end{pmatrix}$,
so that
$\frac{\partial^{k-j}}{\partial\lambda^{k-j}}(\mathbf{V}_{\circ}^{-1})=\begin{pmatrix}0&0\\\
-\frac{\partial^{k-j}}{\partial\lambda^{k-j}}\big{(}\frac{p_{\circ}\mathrm{e}^{-2\mathrm{i}\theta}}{(\lambda-\lambda_{\circ})^{N}}\big{)}&0\end{pmatrix}$
when $k-j>0$. Such lower triangular matrices with zero diagonals when
multiplying by $\begin{pmatrix}0&0\\\ 1&0\end{pmatrix}$ yields the zero
matrix. As a result, the product in Equation (6.16) becomes

(6.18)
$\displaystyle\frac{\partial^{k}\mathbf{M}_{\circ}}{\partial\lambda^{k}}\begin{pmatrix}0&0\\\
1&0\end{pmatrix}$
$\displaystyle=\frac{\partial^{k}\mathbf{N}_{\circ}}{\partial\lambda^{k}}\begin{pmatrix}0&0\\\
1&0\end{pmatrix}\,,\qquad\lambda\in
D_{\lambda_{\circ}}^{\epsilon}\setminus\\{\lambda_{\circ}\\}\,.$

Then, the combination of Equations (6.16) and (6.18) yields

(6.19)
$\displaystyle\mathop{\mathcal{P}_{-n}}_{\lambda=\lambda_{\circ}}\mathbf{N}_{\circ}$
$\displaystyle=\lim_{\lambda\to\lambda_{\circ}}\sum_{k=0}^{N-n}\sum_{i=0}^{N-n-k}\frac{\omega_{\circ,i}}{(N-n-
i-k)!\,k!}\frac{\partial^{N-n-
i-k}\mathrm{e}^{-2\mathrm{i}\theta}}{\partial\lambda^{N-n-
i-k}}\frac{\partial^{k}\mathbf{N}_{\circ}}{\partial\lambda^{k}}\begin{pmatrix}0&0\\\
1&0\end{pmatrix}$

Similarly, one can derive the pole conditions for $N_{\circ}(\lambda;t,z)$ at
$\lambda=\overline{\lambda_{\circ}}$

(6.20)
$\mathop{\mathcal{P}_{-n}}_{\lambda=\overline{\lambda_{\circ}}}\mathbf{N}_{\circ}=-\lim_{\lambda\to\lambda_{\circ}}\sum_{k=0}^{N-n}\sum_{i=0}^{N-n-k}\frac{\overline{\omega_{\circ,i}}}{(N-n-
i-k)!\,k!}\frac{\partial^{N-n-
i-k}\mathrm{e}^{2\mathrm{i}\theta}}{\partial\lambda^{N-n-
i-k}}\frac{\partial^{k}\mathbf{N}_{\circ}}{\partial\lambda^{k}}\begin{pmatrix}0&1\\\
0&0\end{pmatrix}\,.$

Finally, let us rename the matrix $\mathbf{N}_{\circ}(\lambda;t,z)$ back to
$\mathbf{M}_{\circ}(\lambda;t,z)$. Then, RHP 4 is obtained.

###  6.3. $N$th-order soliton solution formula

In this section, we show how to derive the $N$th-order soliton solution
formula from RHP 4. For brevity, we suppress the variable dependence
$(\lambda;t,z)$ in this section. Due to the particular form of the RHP, we
assume the following ansatz for the matrix $\mathbf{M}_{\circ}(\lambda;t,z)$

(6.21)
$\mathbf{M}_{\circ}=\mathbb{I}+\sum_{n=1}^{N}\frac{1}{(\lambda-\lambda_{\circ})^{n}}\mathop{\mathcal{P}_{-n}}\limits_{\lambda=\lambda_{\circ}}\mathbf{M}_{\circ}+\sum_{n=1}^{N}\frac{1}{(\lambda-\overline{\lambda_{\circ}})^{n}}\mathop{\mathcal{P}_{-n}}_{\lambda=\overline{\lambda_{\circ}}}\mathbf{M}_{\circ}\,.$

As discussed in Remark 2, one only needs to reconstruct the top row, so let us
denote

(6.22) $X\coloneqq M_{\circ,1,1}\,,\qquad Y\coloneqq M_{\circ,1,2}\,.$

Thus, the top row of Equation (6.21) can be expressed explicitly as

(6.23)
$\begin{pmatrix}X&Y\end{pmatrix}=\begin{pmatrix}1&0\end{pmatrix}+\sum_{n=1}^{N}\frac{1}{(\lambda-\lambda_{\circ})^{n}}\mathop{\mathcal{P}_{-n}}_{\lambda=\lambda_{\circ}}\begin{pmatrix}X&Y\end{pmatrix}+\sum_{n=1}^{N}\frac{1}{(\lambda-\overline{\lambda_{\circ}})^{n}}\mathop{\mathcal{P}_{-n}}_{\lambda=\overline{\lambda_{\circ}}}\begin{pmatrix}X&Y\end{pmatrix}\,.$

Substituting both pole conditions into the above equation, one gets

(6.24) $\displaystyle X$
$\displaystyle=1+\sum_{n=1}^{N}\frac{1}{(\lambda-\lambda_{\circ})^{n}}\lim_{\lambda\to\lambda_{\circ}}\sum_{k=0}^{N-n}\sum_{i=0}^{N-n-k}\frac{\omega_{\circ,i}}{(N-n-
i-k)!\,k!}\frac{\partial^{N-n-
i-k}\mathrm{e}^{-2\mathrm{i}\theta}}{\partial\lambda^{N-n-
i-k}}\frac{\partial^{k}Y}{\partial\lambda^{k}}\,,$ $\displaystyle Y$
$\displaystyle=-\sum_{n=1}^{N}\frac{1}{(\lambda-\overline{\lambda_{\circ}})^{n}}\lim_{\lambda\to\overline{\lambda_{\circ}}}\sum_{k=0}^{N-n}\sum_{i=0}^{N-n-k}\frac{\overline{\omega_{\circ,i}}}{(N-n-
i-k)!\,k!}\frac{\partial^{N-n-
i-k}\mathrm{e}^{2\mathrm{i}\theta}}{\partial\lambda^{N-n-
i-k}}\frac{\partial^{k}X}{\partial\lambda^{k}}\,.$

Clearly, $X$ is analytic at $\lambda=\overline{\lambda_{\circ}}$ and $Y$ is
analytic at $\lambda=\lambda_{\circ}$, so the limits can be computed as

(6.25) $\displaystyle X$
$\displaystyle=1+\sum_{n=1}^{N}\frac{1}{(\lambda-\lambda_{\circ})^{n}}\sum_{k=0}^{N-n}\sum_{i=0}^{N-n-k}\frac{\omega_{\circ,i}}{(N-n-
i-k)!\,k!}\frac{\partial^{N-n-
i-k}\mathrm{e}^{-2\mathrm{i}\theta}}{\partial\lambda^{N-n-
i-k}}(\lambda_{\circ})\frac{\partial^{k}Y}{\partial\lambda^{k}}(\lambda_{\circ})\,,$
$\displaystyle Y$
$\displaystyle=-\sum_{n=1}^{N}\frac{1}{(\lambda-\overline{\lambda_{\circ}})^{n}}\sum_{k=0}^{N-n}\sum_{i=0}^{N-n-k}\frac{\overline{\omega_{\circ,i}}}{(N-n-
i-k)!\,k!}\frac{\partial^{N-n-
i-k}\mathrm{e}^{2\mathrm{i}\theta}}{\partial\lambda^{k}}(\overline{\lambda_{\circ}})\frac{\partial^{k}X}{\partial\lambda^{k}}(\overline{\lambda_{\circ}})\,.$

We next change the order of sums as
$\sum_{n=1}^{N}\sum_{k=0}^{N-n}=\sum_{k=0}^{N-1}\sum_{n=1}^{N-k}$ and use the
quantity $\gamma_{k}(\lambda)$ defined in Theorem 5 to simplify the equations,

(6.26)
$X(\lambda)=1+\sum_{k=0}^{N-1}\gamma_{k}(\lambda)\frac{\partial^{k}Y}{\partial\lambda^{k}}(\lambda_{\circ})\,,\qquad
Y(\lambda)=-\sum_{k=0}^{N-1}\overline{\gamma_{k}}(\lambda)\frac{\partial^{k}X}{\partial\lambda^{k}}(\overline{\lambda_{\circ}})\,.$

Let us define more quantities in order to solve the above system

(6.27) $X_{s}\coloneqq X^{(s)}(\overline{\lambda_{\circ}})\,,\qquad
Y_{s}\coloneqq Y^{(s)}(\lambda_{\circ})\,,\qquad s=0,1,\dots,N-1\,.$

Substituting $\lambda=\lambda_{\circ}$ and
$\lambda=\overline{\lambda_{\circ}}$ in Equation (6.26), one gets

(6.28)
$X_{0}=1+\sum_{k=0}^{N-1}\gamma_{k}(\overline{\lambda_{\circ}})Y_{k}\,,\qquad
Y_{0}=-\sum_{k=0}^{N-1}\overline{\gamma_{k}}(\lambda_{\circ})X_{k}\,.$

Moreover, by taking derivatives of the system (6.26) $s$ times and
substituting $\lambda=\lambda_{\circ}$ and
$\lambda=\overline{\lambda_{\circ}}$, one gets more equations

(6.29)
$X_{s}=\sum_{k=0}^{N-1}\gamma_{k}^{(s)}(\overline{\lambda_{\circ}})Y_{k}\,,\qquad
Y_{s}=-\sum_{k=0}^{N-1}\overline{\gamma_{k}^{(s)}}(\lambda_{\circ})X_{k}\,,\qquad
s=1,2,\dots,N-1\,.$

Equations (6.28) and (6.29) form a linear system for
$\\{X_{s},Y_{s}\\}_{s=0}^{N-1}$, which can be solved in the same manner as in
Section 3.2. The final solution formula is given in Theorem 5.

##  7\. Soliton gas

In this section, we explain how one can take the formal limit $N\to+\infty$ of
an $N$-DSG and derive an $\infty$-DSG, i.e., soliton gas. The derivation is
similar to [27].

Recall that $N$-DSG is described by RHP 2 with $J=1$. It is easy to verify
that this RHP can be equivalently turned into the following one with an
enlarged jump contour.

######  Riemann-Hilbert Problem 16.

Suppose $J=1$ and $N\geq 1$, and the sets $\Lambda$, $\Omega$ are given
according to Definition 1. Let $\lambda\mapsto\mathbf{M}^{(7)}(\lambda;t,z)$
be a $2\times 2$ matrix function on $\mathbb{C}$, with asymptotics
$\mathbf{M}^{(7)}(\lambda;t,z)\to\mathbb{I}$ as $\lambda\to\infty$. It is also
sectional analytic with jumps

(7.1) $\displaystyle\mathbf{M}^{(7)+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}^{(7)-}(\lambda;t,z)\mathbf{V}^{(7)}(\lambda;t,z)^{-1}\,,\qquad$
$\displaystyle\lambda\in\partial D_{\Lambda}\,,$
$\displaystyle\mathbf{M}^{(7)+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}^{(7)-}(\lambda;t,z)[\mathbf{V}^{(7)}]^{\dagger}(\lambda;t,z)\,,\qquad$
$\displaystyle\lambda\in\partial D_{\overline{\Lambda}}\,,$

where the jump matrix $\mathbf{V}^{(7)}(\lambda;t,z)$ is given by

(7.2) $\mathbf{V}^{(7)}(\lambda;t,z)\coloneqq\begin{pmatrix}1&0\\\
\mathrm{e}^{-2\mathrm{i}\theta(\lambda;t,z)}\sum_{k=1}^{N}\frac{\omega_{1,k}}{\lambda-\lambda_{1,k}}&1\end{pmatrix}\,,$

and the jump contours $\partial D_{\Lambda}$ and $\partial
D_{\overline{\Lambda}}$ are illustrated in Figure 13(center) with
counterclockwise orientation.

We point out that the region $D_{\Lambda}\subset\mathbb{C}^{+}$ encloses all
the discrete eigenvalues $\Lambda$ in the upper half plane, and similarly,
$D_{\overline{\Lambda}}\subset\mathbb{C}^{-}$ encloses $\overline{\Lambda}$ in
the lower half plane. Here, $\overline{\Lambda}$ denote the complex conjugate
of the set $\Lambda$, not its closure. Therefore, the two jump contours
$\partial D_{\Lambda}$ and $\partial D_{\overline{\Lambda}}$ never intersect.
Furthermore, in the rest of this work, we assume implicitly that $\partial
D_{\Lambda}$ is sufficiently close to $\Lambda$ and contains an arc of the
circle with radius $r_{1}$, as illustrated in Figure 13(left and center).

Figure 13.  Left and center: The transformation from the $N$-DSG RHP 2 with
function $\mathbf{M}(\lambda;t,z)$ and $J=1$ to an equivalent RHP 16 with
combined jump contours. Right: An illustration of: (i) the contour
$\mathrm{A}$ and $\overline{\mathrm{A}}$ from Equation (2.36) shown in red;
(ii) the regions $D_{\mathrm{A}}$ and $D_{\overline{\mathrm{A}}}$ with
boundaries shown in blue, oriented counterclockwise. The blue dot are
representations of discrete eigenvalues $\Lambda$ and $\overline{\Lambda}$.

With the new setup, we are now ready to take the limit of RHP 16 as
$N\to+\infty$ while keeping $J=1$. We would like to take this limit in a
controllable way, in the sense that all discrete eigenvalues in
$\Lambda_{1}=\Lambda$ are accumulating on an arc centered at the origin.
Recall the arc $\mathrm{A}$ defined in Definition 3, so we know that
$\Lambda\subset\mathrm{A}$. Correspondingly, one defines a new region
$D_{\mathrm{A}}\subset\mathbb{C}^{+}$ which contains the whole arc
$\mathrm{A}$ and does not intersects with the real axis. Please see an
illustration of the contour $\mathrm{A}$ and corresponding region
$D_{\mathrm{A}}$ shown in Figure 13(right).

It is clear that in RHP 16 the poles are encoded in the jump matrix
$\mathbf{V}^{(7)}(\lambda;t,z)$ via the sum in the $(1,2)$ component. In order
to take the limit $N\to+\infty$, we rewrite the relevant function as follows,

(7.3)
$\sum_{k=1}^{N}\frac{\omega_{1,k}}{\lambda-\lambda_{1,k}}=\int_{\mathrm{A}}\omega_{1}(s;N)\frac{\,\mathrm{d}s}{s-\lambda}\,,\qquad\omega_{1}(\lambda;N)\coloneqq-\sum_{k=1}^{N}\omega_{1,k}\delta(\lambda-\lambda_{1,k})\,,\qquad\lambda\not\in\mathrm{A}\,,$

where $\delta(\cdot)$ is the Dirac delta distribution defined on the arc
$\mathrm{A}$, oriented left-to-right i.e.,
$\int_{\mathrm{A}}f(s)\delta(s-s_{0})\,\mathrm{d}s=f(s_{0})$, provided that
$f(z)$ is continuous at every interior point $z=s_{0}\in\mathrm{A}$. In
Equation (7.3), the norming constants $\omega_{1,k}$ can be recognized as
weights of each Dirac delta, and the sum is not normalized, i.e,
$\int_{\mathrm{A}}\omega_{1}(s;N)\,\mathrm{d}s\neq 1$. Hence, the function
$\omega_{1}(s;N)$ can be considered as an (unnormalized) complex distribution
on the arc $\mathrm{A}$. Therefore, taking the limit $N\to+\infty$ of the
left-hand side of Equation (7.3) is equivalent to calculating
$\lim_{N\to+\infty}\omega_{1}(\lambda;N)$. Formally speaking, the latter limit
describes infinitely many Dirac deltas accumulating on the arc $\mathrm{A}$
with their own complex-valued weights, which can be considered as an
interpolation of an (unnormalized) continuous distribution. So, the formal
result can be written as

(7.4) $\omega(\lambda)\coloneqq
2\pi\mathrm{i}\lim_{N\to+\infty}\omega_{1}(\lambda;N)\,,\qquad\lambda\in\mathrm{A}_{1}\,,$

where the function $\omega(\lambda)$ is a complex-valued function. The factor
$2\pi\mathrm{i}$ in the above equation is for convenience for later usages.
Also, we assume $\omega\in L_{2}(\mathrm{A})$ for convenience. It is important
to point out that the support of $\omega(\lambda)$ is an arbitrary subset of
$\mathrm{A}$. Please see Figure 14 for demonstration.

Figure 14.  An illustration of the eigenvalue distribution in the limit
$N\to+\infty$. In the right panel, the thick blue arcs represent
$\operatorname{supp}\omega\subset\mathrm{A}$ and
$\operatorname{supp}\overline{\omega}\subset\overline{\mathrm{A}}$.

Hence, the limit produces

(7.5)
$\displaystyle\lim_{N\to+\infty}\sum_{k=1}^{N}\frac{\omega_{1,k}}{\lambda-\lambda_{1,k}}$
$\displaystyle=\frac{1}{2\pi\mathrm{i}}\int_{\mathrm{A}}\frac{\omega(s)\,\mathrm{d}s}{s-\lambda}\,,\qquad\lambda\not\in\mathrm{A}\,,$
$\displaystyle\lim_{N\to+\infty}\sum_{k=1}^{N}\frac{\overline{\omega_{1,k}}}{\lambda-\overline{\lambda_{1,k}}}$
$\displaystyle=\frac{1}{2\pi\mathrm{i}}\int_{\overline{\mathrm{A}}}\frac{\overline{\omega}(s)\,\mathrm{d}s}{s-\lambda}\,,\qquad\lambda\not\in\overline{\mathrm{A}}\,.$

where the integration contour is oriented left-to-right. To simplify the above
integrals, one defines the following Cauchy transforms

(7.6)
$\mathscr{C}_{\omega}(\lambda)\coloneqq\frac{1}{2\pi\mathrm{i}}\int_{\mathrm{A}}\frac{\omega(s)\,\mathrm{d}s}{s-\lambda}.$

Then, the above limits can be rewritten as
$\lim\limits_{N\to+\infty}\sum_{k=1}^{N}\frac{\omega_{1,k}}{\lambda-\lambda_{1,k}}=\mathscr{C}_{\omega}(\lambda)$,
and
$\lim\limits_{N\to+\infty}\sum_{k=1}^{N}\frac{\overline{\omega_{1,k}}}{\lambda-\overline{\lambda_{1,k}}}=\mathscr{C}_{\overline{\omega}}(\lambda)$.
It is easy to verify that
$\overline{\mathscr{C}_{\omega}}(\lambda)=-\mathscr{C}_{\overline{\omega}}(\lambda)$.
Consequently, the jumps in RHP 16 in the limit becomes

(7.7) $\displaystyle\lim_{N\to+\infty}\mathbf{V}^{(7)}(\lambda;t,z)^{-1}$
$\displaystyle=\lim_{N\to+\infty}\begin{pmatrix}1&0\\\
-\mathrm{e}^{-2\mathrm{i}\theta(\lambda;t,z)}\sum_{k=1}^{N}\frac{\omega_{1,k}}{\lambda-\lambda_{1,k}}&1\end{pmatrix}=\begin{pmatrix}1&0\\\
-\mathrm{e}^{-2\mathrm{i}\theta(\lambda;t,z)}\mathscr{C}_{\omega}(\lambda)&1\end{pmatrix}\,,$
$\displaystyle\lim_{N\to+\infty}[\mathbf{V}^{(7)}]^{\dagger}(\lambda;t,z)$
$\displaystyle=\lim_{N\to+\infty}\begin{pmatrix}1&\mathrm{e}^{2\mathrm{i}\theta(\lambda;t,z)}\sum_{k=1}^{N}\frac{\overline{\omega_{1,k}}}{\lambda-\overline{\lambda_{1,k}}}\\\
0&1\end{pmatrix}=\begin{pmatrix}1&\mathrm{e}^{2\mathrm{i}\theta(\lambda;t,z)}\mathscr{C}_{\overline{\omega}}(\lambda)\\\
0&1\end{pmatrix}\,.$

Following these calculations, RHP 16 as $N\to+\infty$ becomes the following
one.

######  Riemann-Hilbert Problem 17 (Cauchy-form of soliton gas).

Let an arc $\mathrm{A}$ and a function $\omega(\lambda)$ be given according to
Definition 3. Let $\lambda\mapsto\mathbf{M}^{(8)}(\lambda;t,z)$ be a $2\times
2$ matrix function on $\mathbb{C}$, with asymptotics
$\mathbf{M}^{(8)}(\lambda;t,z)\to\mathbb{I}$ as $\lambda\to\infty$. It is also
sectional analytic with jumps

(7.8) $\displaystyle\mathbf{M}^{(8)+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}^{(8)-}(\lambda;t,z)\mathbf{V}^{(8)}(\lambda;t,z)^{-1}\,,\qquad$
$\displaystyle\lambda\in\partial D_{\mathrm{A}}\,,$
$\displaystyle\mathbf{M}^{(8)+}(\lambda;t,z)$
$\displaystyle=\mathbf{M}^{(8)-}(\lambda;t,z)[\mathbf{V}^{(8)}]^{\dagger}(\lambda;t,z)\,,\qquad$
$\displaystyle\lambda\in\partial D_{\overline{\mathrm{A}}}\,,$

where the jump matrix $\mathbf{V}^{(8)}(\lambda;t,z)$ is given by

(7.9) $\mathbf{V}^{(8)}(\lambda;t,z)\coloneqq\begin{pmatrix}1&0\\\
\mathrm{e}^{-2\mathrm{i}\theta(\lambda;t,z)}\mathscr{C}_{\omega}(\lambda)&1\end{pmatrix}\,,$

and the jump contours $\partial D_{\mathrm{A}}$ and $\partial
D_{\overline{\mathrm{A}}}$ are illustrated in Figure 14(right) with left-to-
right orientation.

Finally, one can convert the above RHP into a simpler form, by closing the
jump contour $\partial D_{\mathrm{A}}$ onto the arc $\mathrm{A}$ and $\partial
D_{\overline{\mathrm{A}}}$ onto $\overline{\mathrm{A}}$, respectively. Then,
RHP 17 is converted to RHP 5.

## Acknowledgments

This project was partially supported by the NSF of China (No. 12201526), the
Fundamental Research Funds for the Central Universities (No. 20720220040), and
the Natural Science Foundation of Fujian Province of China (No. 2022J01032).

## References

  * [1] O
  * [2] A. Abeya, B. Prinari, G. Biondini, and P. G. Kevrekidis, Solitons and soliton interactions in repulsive spinor Bose–Einstein condensates with nonzero background, Eur. Phys. J. Plus 136, 1126 (2021).
  * [3] M. J. Ablowitz, D. J. Kaup, A. C. Newell, and H. Segur, Method for Solving the Sine-Gordon Equation, Phys. Rev. Lett. 30, 1262 (l973).
  * [4] M. J. Ablowitz, D. J. Kaup, A. C. Newell, and H. Segur. The inverse scattering transform-Fourier analysis for nonlinear problems, Stud. Appl. Math. 53, 249–315 (1974).
  * [5] M. J. Ablowitz, D. J. Kaup, and A. C. Newell, Coherent pulse propagation, a dispersive irreversible phenomenon, J. Math. Phys. 15, 1852–1858 (1974).
  * [6] M.J. Ablowitz, B. Prinari, A.D. Trubatch, Discrete and Continuous Nonlinear Schrödinger Systems, Cambridge University Press, 2004.
  * [7] M. J. Ablowitz and H. Segur, Solitons and the inverse scattering transform, SIAM, Philadelphia, PA, 1981.
  * [8] L. Allen and J. H. Eberly, Optical resonance and two-level atoms, Wiley, New York, NY, 1975.
  * [9] R. Beals and R. R. Coifman, Scattering and inverse scattering for first order systems, Comm. Pure Appl. Math. 37, 39–90 (1984).
  * [10] R. Beals and R. R. Coifman, Inverse scattering and evolution equations, Comm. Pure Appl. Math. 38, 29–42 (1985).
  * [11] D. Bilman and P. D. Miller, Broader universality of rogue waves of infinite order, Physica D 435, 133289 (2022)
  * [12] D. Bilman, L. Ling, and P. D. Miller, Extreme superposition: Rogue waves of infinite order and the Painlevé-III hierarchy, Duke Math. J. 169, 671–760 (2020).
  * [13] J. Cen, F. Correa, and A. Fring, Degenerate multi-solitons in the sine-Gordon equation, J. Phys. A Math. Theor. 50, 435201 (2017).
  * [14] J. Cen and A. Fring, Asymptotic and scattering behaviour for degenerate multi-solitons in the Hirota equation, Physica D 397, 17–24, (2019).
  * [15] S. Cuccagna, R. Jenkins, On the Asymptotic Stability of Soliton Solutions of the Defocusing Nonlinear Schrödinger Equation. Commun. Math. Phys. 343, 921–969 (2016).
  * [16] G. Dean, T. Klotz, B. Prinari and F. Vitale, Dark–dark and dark–bright soliton interactions in the two-component defocusing nonlinear Schrödinger equation, Appl. Anal. 92, 379–397 (2013).
  * [17] P. A. Deift, A. R. Its, and X. Zhou, Long-time asymptotics for integrable nonlinear wave equations, Important developments in soliton theory, Springer Series in Nonlinear Dynamics, Springer, Berlin, 1993.
  * [18] P. A. Deift and X. Zhou, A steepest descent method for oscillatory Riemann-Hilbert problems. Asymptotics for the MKdV equation, Ann. Math. 37, 295–368 (1993).
  * [19] G. A. El, A. M. Kamchatnov, M. V. Pavlov, and S. A. Zykov, Kinetic Equation for a Soliton Gas and Its Hydrodynamic Reductions, J. Nonlinear Sci. 21, 151–191, (2011).
  * [20] G. El and A. Tovbis, Spectral theory of soliton and breather gases for the focusing nonlinear Schrödinger equation, Phys. Rev. E 101, 052207 (2020).
  * [21] B.-F. Feng, General $N$-soliton solution to a vector nonlinear Schrödinger equation, J. Phys. A: Math. Theor. 47, 355203 (2014).
  * [22] I. P. Gabitov, V. E. Zakharov, and A. V. Mikhaǐlov, Superfluorescence pulse shape, JETP Lett. 37, 279–282 (1983).
  * [23] I. P. Gabitov, V. E. Zakharov, and A. V. Mikhaǐlov, Non-linear theory of superfluorescence, JETP 59, 703–709 (1984).
  * [24] I. R. Gabitov, V. E. Zakharov, and A. V. Mikhaǐlov, Maxwell-Bloch equation and the inverse scattering method, Theor. Math. Phys. 63, 328–343 (1985).
  * [25] A. Gelash, D. Agafontsev, V. Zakharov, G. El, S. Randoux, and P. Suret, Bound State Soliton Gas Dynamics Underlying the Spontaneous Modulational Instability, Phys. Rev. Lett. 123, 234102 (2019).
  * [26] A. A. Gelash and V. E. Zakharov, Superregular solitonic solutions: a novel scenario for the nonlinear stage of modulation instability. Nonlinearity 27, R1–38 (2014).
  * [27] M. Girotti, T. Grava, R. Jenkins, and K. D. T.-R. McLaughlin, Rigorous Asymptotics of a KdV Soliton Gas, Commun. Math. Phys. 384, 733–784 (2021).
  * [28] M. Girotti, T. Grava, R. Jenkins, K.T.-R. McLaughlin, and A. Minakov, Soliton versus the gas: Fredholm determinants, analysis, and the rapid oscillations behind the kinetic equation, Comm. Pure Appl. Math. 76, 3233–3299 (2023).
  * [29] A. Gkogkou and B. Prinari, Soliton interactions in certain square matrix nonlinear Schrödinger systems, Eur. Phys. J. Plus 135, 609 (2020).
  * [30] R. Guo, B. Tian, X. Lu, H. Zhang, and W. Liu, Darboux Transformation and Soliton Solutions for the Generalized Coupled Variable Coefficient Nonlinear Schrödinger-Maxwell-Bloch System with Symbolic Computation, Comput. Math. Math. Phys. 52, 565–577 (2012).
  * [31] R. Guo, B. Tian and L. Wang, Soliton solutions for the reduced Maxwell–Bloch system in nonlinear optics via the N-fold Darboux transformation, Nonlinear Dyn. 69, 2009–2020, (2012).
  * [32] L. Huang and J. Lenells, Construction of solutions and asymptotics for the sine-Gordon equation in the quarter plane, J. Integrable Syst. 3, xyy014, (2018).
  * [33] S. Kamvissis, K. D. T-R McLaughlin, and P. D. Miller, Semiclassical Soliton Ensembles for the Focusing Nonlinear Schrodinger Equation, Princeton University Press, (2003).
  * [34] V. I. Karpman and V. V. Solov’ev, A perturbational approach to the 2-soliton systems, Phys. D 3, 487–502 (1981).
  * [35] V. Kotlyarov and O. Minakov, Maxwell–Bloch equations without spectral broadening: the long-time asymptotics of an input pulse in a long two-level laser amplifier, Nonlinearity 36, 5007 (2023).
  * [36] G. L. Lamb, Analytical descriptions of ultrashort optical pulse propagation in a resonant medium, Rev. Mod. Phys. 43, 99–124 (1971).
  * [37] S. Li and P. D. Miller, On the Maxwell-Bloch System in the Sharp-Line Limit Without Solitons, Commun. Pure Appl. Math., 77, 457–542 (2024).
  * [38] S. Li and G. Biondini, Soliton interactions and degenerate soliton complexes for the focusing nonlinear Schrödinger equation with nonzero background, Eur. Phys. J. Plus 133, 400 (2018).
  * [39] S. Li, G. Biondini, G. Kovačič and I. Gabitov, Resonant optical pulses on a continuous wave back-ground in two-level active media, EPL 121, 20001 (2018).
  * [40] S. Li, G. Biondini and C. Schiebold, On the degenerate soliton solutions of the focusing nonlinear Schrödinger equation, J. Math. Phys. 58, 033507 (2017).
  * [41] L. Ling and X. Zhang, Large and infinite-order solitons of the coupled nonlinear Schrödinger equation, Physica D 457, 133981 (2024).
  * [42] C. Liu, S.-C. Chen, and N. Akhmediev, Fundamental and Second-Order Superregular Breathers in Vector Fields, Phys. Rev. Lett. 132, 027201 (2024).
  * [43] S. V. Manakov, On the theory of two-dimensional stationary self-focusing of electromagnetic waves, Sov. Phys. JETP 38, 248 (1974).
  * [44] S. L. McCall and E. L. Hahn, Self-induced transparency by pulsed coherent light, Phys. Rev. Lett. 18, 908–911 (1967).
  * [45] S. L. McCall and E. L. Hahn, Self-induced transparency, Phys. Rev. 183, 457–485 (1969).
  * [46] S. L. McCall and E. L. Hahn, Pulse-area-pulse-energy description of a traveling-wave laser amplifier, Phys. Rev. A 2, 861–870 (1970).
  * [47] J. Rao, D. Mihalache, J. He, and F. Zhou, Degenerate and non-degenerate vector solitons and their interactions in the two-component long-wave–short-wave model of Newell type, Chaos Solit. Fractals 166, 112963 (2023).
  * [48] C. Schiebold, Asymptotics for the multiple pole solutions of the nonlinear Schrödinger equation, Nonlinearity 30, 2930 (2017).
  * [49] Y. S. Xue, B. Tian, W. B. Ai, F. H. Qi, R. Guo, and B. Qin, Soliton interactions in a generalized inhomogeneous coupled Hirota–Maxwell–Bloch system, Nonlinear Dyn. 67, 2799–2806 (2012).
  * [50] N. J. Zabusky and M. D. Kruskal, Interaction of “Solitons” in a Collisionless Plasma and the Recurrence of Initial States, Phys. Rev. Lett. 15, 240 (1965).
  * [51] V. E. Zakharov, Propagation of an amplifying pulse in a two-level medium, JETP Lett. 32, 589–593 (1980).
  * [52] V. E. Zakharov, Turbulence in Integrable Systems, Stud. Appl. Math. 122, 219–234 (2009).
  * [53] V.E. Zakharov, A.B. Shabat, Exact Theory of Two-Dimensional Self-Focusing and One-Dimensional Self-Modulation of Waves in Nonlinear Media, Sov. Phys. JETP 34, 62 (1972).
  * [54] J. H. Zhang, L. Wang, and C. Liu, Superregular breathers, characteristics of nonlinear stage of modulation instability induced by higher-order effects, Proc. R. Soc. A 473, 20160681 (2017).
  * [55] G. Zhang, Z. Yan, Focusing and defocusing mKdV equations with nonzero boundary conditions: Inverse scattering transforms and soliton interactions, Physica D 410, 132521 (2020)
  * [56] X. Zhou, The Riemann-Hilbert problem and inverse scattering, SIAM J. Math. Anal. 20, 966–986 (1989).
  * [57] X. Zhou, The L2-Sobolev space bijectivity of the scattering and inverse scattering transforms. Comm. Pure Appl. Math. 51, 697–731 (1998).


# NT @ UW-21-13
Confinement in Two-Dimensional QCD and the Infinitely Long Pion

Colin M. Weller and Gerald A. Miller  University of Washington, Seattle
Seattle, WA 98195-1560

###### Abstract

Three current models of QCD in (1+1) dimensions are examined and extended in
light-front coordinates. A pion of high momentum is found to have an infinite
extent along its direction of motion.

There has been a recent revival of interest in theories involving two-
dimensional(one-space-one time) treatments of QCD De Téramond and Brodsky
(2021); Li and Vary (2021); Ahmady  _et al._ (2021a, b). This stems from the
need to include the effects of non-vanishing quark masses and also to enlarge
the number of space-time variables of light-front holographic QCD to four. In
the original treatments (see the review Brodsky  _et al._ (2015)) the chiral
limit is used and the longitudinal light-front momentum fraction, $x$, is
frozen Sheckler and Miller (2021), so that effectively the only degrees of
freedom are light-front time and transverse position. The first effort aimed
at including the effects of mass was contained in Ref. Chabysheva and Hiller
(2013). The present epistle aims to unify the earlier treatments and to
exhibit the confining aspects of the approaches in three-spatial dimensions
using the spatial coordinate, ${\tilde{z}}$, that is canonically conjugate to
the variable $x$ Miller and Brodsky (2020).

We begin by briefly summarizing light-front holographic QCD (LFHQCD) following
Brodsky  _et al._ (2015). Light-front (LF) quantization is a relativistic,
frame-independent approach to describing the constituent structure of hadrons.
The simple structure of the light-front vacuum allows an unambiguous
definition of the partonic content of a hadron in QCD and of hadronic light-
front wave functions. The QCD light-front Hamiltonian is constructed from the
energy-momentum tensor of QCD. The spectrum and light-front wave functions of
relativistic bound states are obtained from the resulting eigenvalue equation,
an infinite set of coupled integral equations for the LF components in a
complete basis of noninteracting $n$-particle states. This provides a quantum-
mechanical probabilistic interpretation of the structure of hadronic states in
terms of their constituents at the same light-front time $x^{+}=x^{0}+x^{3}$.

To a first semiclassical approximation, where quantum loops and quark masses
are not included, the relativistic bound-state equation for light hadrons can
be reduced to an effective LF Schrödinger equation. In conjugate coordinate
space, the relevant dynamical variable is an invariant impact kinematical
variable $\mbox{\boldmath$\zeta$}={\bf b}\sqrt{x(1-x)}$, where ${\bf b}$ is
the transverse separation of the constituents. (Boldface notation specifies
vectors of the two-dimensional transverse space.) The complexities of the
strong interaction dynamics are then hidden in an effective potential
$U(\zeta)$. Remarkably, the resulting light-front Hamiltonian has a structure
that matches exactly the eigenvalue equations in anti-deSitter space Brodsky
_et al._ (2015). This offers the possibility to explicitly connect the AdS
wave function $\Phi(z)$ to the internal constituent structure of hadrons.
Moreover, one can obtain the AdS wave equations by starting from the
semiclassical approximation to light-front QCD in physical spacetime. This
connection yields a relation between the coordinate $z$ of AdS space with the
impact LF variable $\zeta$, thus giving the holographic variable $z$ a precise
definition and intuitive meaning in light-front QCD Brodsky and de Teramond
(2009).

An effective light-front Schrödinger equation for the quark
($m_{1})$-antiquark ($m_{2})$, wave function $\psi(x,{\bf k}),$ of a meson is
given by

$\left[\frac{m_{1}^{2}}{x}+\frac{m_{2}^{2}}{1-x}+\frac{{\bf
k}^{2}}{x(1-x)}+V_{\rm eff}\right]\psi=M_{h}^{2}\psi,$ (1)

where ${\bf k}$ is the transverse relative momentum, the first three terms are
the light-front energy of two non-interacting partons, $M_{h}$ is the
invariant mass of the hadron, and $V_{\rm eff}$ is an effective potential that
acts in three-dimensional space.

Note that the kinetic energy term depends only on the two-dimensional vector
${\bf k}/\sqrt{x(1-x)}$, if the chiral limit ($m_{1,2}=0)$ is taken. Then,
using coordinate space, and taking $V_{\rm eff}$ to depend on
$\mbox{\boldmath$\zeta$},\,V_{\rm eff}=U_{\perp}(\mbox{\boldmath$\zeta$})$ one
obtains

$\displaystyle\left(-\frac{d^{2}}{d\zeta^{2}}+\frac{L^{2}}{\zeta^{2}}+U_{\perp}(\mbox{\boldmath$\zeta$})\right)\varphi(\mbox{\boldmath$\zeta$})=M^{2}\varphi(\mbox{\boldmath$\zeta$}),$
(2)

which is the same as the equation of motion in the soft-wall AdS model Karch
_et al._ (2006) if $L$ is taken to be angular momentum, $\zeta$ is identified
with the fifth dimension $z$ in AdS space, and the wave function
$\varphi(\mbox{\boldmath$\zeta$})$ is identified with the string modes
$\Phi(z)$.

Assuming duality with AdS5 suggests models for
$U_{\perp}(\mbox{\boldmath$\zeta$})$, through a correspondence between the
transverse Schrödinger equation (1) and the equation of motion for a spin-$J$
field in AdS5 Brodsky  _et al._ (2015). For the soft-wall model Karch  _et
al._ (2006), the effective potential reduces to an oscillator potential
$U_{\perp}(\zeta)=\kappa^{4}\zeta^{2}+2\kappa^{2}(J-1),$ where $\kappa$ is the
strength of the holographic confinement and $J$ is the total angular momentum.
For this potential, the spectrum of masses is
$M_{\perp}^{2}=4\kappa^{2}\left(n+(J+|L_{z}|)/2\right)$, with $n$ the radial
quantum number, and the transverse wave functions are the two-dimensional
oscillator eigenfunctions. The spectrum of the model provides for a linear
Regge trajectory and a good fit to light meson masses. Using this potential
leads to a massless pion ($n=0,L_{z}=0,J=0$) in the chiral limit.

In Eq. (2) the variable $x$ is held fixed, leading to a $(2+1)$ space-time
description that provides an excellent representation of the hadronic
spectrum, but is manifestly incomplete because of the missing dynamics of the
longitudinal direction. One incorporates De Téramond and Brodsky (2021); Li
and Vary (2021); Ahmady  _et al._ (2021a, b); Chabysheva and Hiller (2013)
these dynamics by asserting that the total potential $V_{\rm eff}$ of Eq. (1)
is the sum of two terms: $V_{\rm
eff}=U_{\perp}(\mbox{\boldmath$\zeta$})+V_{\parallel}(x)$, with

$\displaystyle\left[\frac{m_{1}^{2}}{x}+\frac{m_{2}^{2}}{1-x}+V_{\parallel}\right]X_{n}(x)=M_{\parallel}^{2}X_{n}(x).$
(3)

Under these assumptions the full wave function $\psi$ is given by the product

$\displaystyle\psi(x,{\bf b})=\varphi(\mbox{\boldmath$\zeta$})X_{n}(x)$ (4)

and $M_{h}^{2}=M_{\perp}^{2}+M_{\parallel}^{2}$. Here the normalization
convention Callan  _et al._ (1976)

$\displaystyle\int_{0}^{1}\frac{|X_{n}(x)|^{2}}{x(1-x)}\,dx=1$ (5)

is used. It is helpful to define $\chi_{n}(x)=X_{n}(x)/\sqrt{x(1-x)}$.

One expects that the correct QCD-potential is not a sum of two independent
terms. In that case the space of product wave functions forms a useful,
complete, relativistic basis, such as advocated in Chabysheva and Hiller
(2013) and implemented in Li  _et al._ (2016). Nevertheless, Refs. De Téramond
and Brodsky (2021); Li and Vary (2021); Ahmady  _et al._ (2021a, b) compare
the resulting values of $M_{h}^{2}$ with measured mesonic spectra.

Next, we consider the equation of the form
$H_{\parallel}\chi_{n}=M_{n}^{2}\chi_{n}$, with $H_{\parallel}$ Hermitian.
This is the form of the wave equation used by both ’t Hooft and Callen, Coote
& Gross Callan  _et al._ (1976) (and many others). Note that the $V_{||}$ of
Refs. Li and Vary (2021) and ’t Hooft (1974) are not diagonal in $x$. This is
because confining potentials must have an explicit dependence on the
coordinate-space variable, ${\tilde{z}}$, that is canonically conjugate to
$x$. We shall see that none of the models of $V_{\parallel}$ in current use is
of the form $V_{||}(x)$.

With the fundamental longitudinal wave equation the usual orthonormal equation
applies

$\displaystyle\int\,dx\chi^{*}_{n}(x)\chi_{m}(x)=\delta_{nm},$ (6)

as explicitly stated in Ref. Callan  _et al._ (1976). Next, consider the
matrix element.

$\displaystyle\langle n|H|m\rangle$ $\displaystyle=\int
dx\,dy\,\chi^{*}_{n}(x)H(x,y)\chi_{m}(y).$ (7)

The Hamiltonian can be expressed as
$H(x,y)=\frac{m^{2}}{x(1-x)}\delta(x-y)+V_{||}(x,y),\,$so that:

$\displaystyle\langle n|H|m\rangle$ $\displaystyle=\int
dx\frac{m^{2}}{x(1-x)}\chi_{n}^{*}(x)\chi_{m}(x)$
$\displaystyle+\int\,dx\,dy\chi^{*}_{n}(x)V_{||}(x,y)\chi_{m}(y).$ (8)

The longitudinal potential $V_{\parallel}$ must be chosen to determine the
function $\chi_{n}(x)$. There are two choices in the literature. The first two
(LV) Li and Vary (2021) and the ’t Hooft model (tH) ’t Hooft (1974), are given
by

$\displaystyle V_{\rm
LV}(x)\chi_{n}(x)=-\sigma^{2}\partial_{x}\,x(1-x)\partial_{x}\chi_{n}(x)$
$\displaystyle(V_{\rm
tH}\,\chi_{n})(x)={g^{2}\over\pi}P\int_{0}^{1}dy{\chi_{n}(x)-\chi_{n}(y)\over(x-y)^{2}},$
(9)

with the principal value is defined Brower  _et al._ (1979) as
$P{f(x,y)\over(x-y)^{2}}\equiv{1\over
2}[{f(x,y)\over(x-y+i\epsilon)^{2}}+{f(x,y)\over(x-y-i\epsilon)^{2}}]$ in the
limit $\epsilon\to 0$. The function $V_{\rm LV}(x)$ is used in Li and Vary
(2021), and has the advantage that exact solutions are available in terms of
Jacobi polynomials. The ’t Hooft model ’t Hooft (1974), obtained in the
large-$N$ limit of two-dimensional QCD, is used in Ahmady  _et al._ (2021a,
b); Chabysheva and Hiller (2013). This is the natural choice for a confining
potential in one spatial dimension. When QCD is quantized in light-cone gauge,
such a potential appears automatically as an instantaneous Coulomb-like
interaction, $V_{\rm
tH}({\tilde{z}})=g^{2}|{\tilde{z}}|e^{-\epsilon|{\tilde{z}}|}$, between quark
currents Brower  _et al._ (1979), with ${\tilde{z}}$ as the longitudinal
position operator Miller and Brodsky (2020) that is canonically conjugate to
the longitudinal momentum variable $x$. Taking the Fourier transform of
$V_{\rm tH}({\tilde{z}})$ using the transformation $\langle
x|{\tilde{z}}\rangle=e^{-ix{\tilde{z}}}/\sqrt{2\pi}$ and including effects of
the quark self-energy via the term $\chi(x)$ term of the principle-value
integral leads to the expression appearing in Eq. (9). Note that in the ’t
Hooft model, the masses $m_{1,2}$ are explicitly current quark masses. The ’t
Hooft model was extensively studied during the 1970’s; see the review Ellis
(1977).

The third approach De Téramond and Brodsky (2021) is to assert that $\chi(x)$
is a Gaussian in the invariant mass-squared: $\chi(x)={\cal
N}\exp{[-1/(2\kappa^{2})(-m_{1}^{2}/x+m_{2}^{2}/(1-x))]}$, where $\cal N$ is a
normalization constant Brodsky  _et al._ (2015). This model is termed the
invariant mass wave function (IMWF).

To obtain an expression for the Hamiltonian using the $X$ normalization:
define $X_{n}(x)\equiv\sqrt{x(1-x)}\chi_{n}(x)$ so that

$\displaystyle\langle n|H|m\rangle$
$\displaystyle=\int\frac{dx}{x(1-x)}{\tilde{X}}_{n}(x)\frac{m^{2}}{x(1-x)}X_{m}(x)+\langle
n|V_{\parallel}|m\rangle,$ (10)

with

$\displaystyle\langle n|V_{||}|m\rangle$
$\displaystyle=\int\frac{dx\,dy}{x(1-x)}{\tilde{X}}_{n}(x){\sqrt{x(1-x)}\over\sqrt{y(1-y)}}V(x,y)X_{m}(y).$
(11)

The effective potential, ${\sqrt{x(1-x)}\over\sqrt{y(1-y)}}V(x,y)$ appearing
between the two wave functions is not Hermitian. The function $\tilde{X}_{n}$
are obtained using the effective potential
${\sqrt{y(1-y)}\over\sqrt{x(1-x)}}V(x,y)$.

The consequence of this lack of hermiticity is that, for example, the
resulting wave equation for the t’ Hooft model would be of the form:

$\displaystyle M_{n}^{2}X_{n}(x)=$
$\displaystyle\frac{m^{2}}{x(1-x)}X_{n}(x)-\frac{g^{2}}{\pi}\mathcal{P}\int\,dy\frac{(X_{n}(x)-\frac{\sqrt{x(1-x)}}{\sqrt{y(1-y)}}X_{n}(y))}{(x-y)^{2}}.$
(12)

This equation was solved in Ref. Chabysheva and Hiller (2013). The non-
Hermitian nature of the effective potential seems to represent a difficulty in
extending light-front holographic QCD in the conformal limit to include non-
zero quark masses and excitation of longitudinal modes. Instead, the
normalization given in Eq. (6) is used for the rest of the paper.

We now begin to discuss the differences between the two potentials of Eq. (9).
In the chiral limit, the ground state wave function of both Li and Vary
(2021); Callan  _et al._ (1976) can be seen from the differential equation for
$M_{\parallel}^{2}$ to be simply ( dropping the subscript) $\chi(x)=1$, with
an eigenvalue $M_{\parallel}^{2}=0$ . Ref. De Téramond and Brodsky (2021) also
have $\chi(x)=1$ in the chiral limit.

The LFHQCD formalism Brodsky  _et al._ (2015) achieved excellent reproduction
of the hadronic spectrum. Using $\chi(x)=1$ along with $M_{\parallel}^{2}=0$
preserves that ground-state spectrum.

The useful identity Callan  _et al._ (1976):

$\displaystyle
M_{\parallel}^{2}\int_{0}^{1}\,dx\,\chi(x)=\int_{0}^{1}\,dx\,\chi(x)\left[\frac{m_{1}^{2}}{x}+\frac{m_{2}^{2}}{1-x}\right]$
(13)

holds in both the Li-Vary and ’t Hooft models. The need to have finite results
for both sides of this equation signifies that $\chi(x)$ vanishes at the end-
points, $x=0,1$.

Given the similarities between the Li-Vary and ’t Hooft models, it is natural
to initially expect that any differences between the models might be thought
to be small. This is far from the case. Using coordinate space helps in
understanding that there indeed are significant differences between the
potentials. This is done by evaluating $V_{\rm LV}(x)$ in coordinate space.

$\displaystyle\langle{\tilde{z}}|V_{\rm
LV}(x)|{\tilde{z}}^{\prime}\rangle=\sigma^{2}{\tilde{z}}{\tilde{z}}^{\prime}\int_{0}^{1}{dx\over
2\pi}e^{ix({\tilde{z}}^{\prime}-{\tilde{z}})}x(1-x)$
$\displaystyle=\frac{\sigma^{2}}{2\pi}{\tilde{z}}{\tilde{z}}^{\prime}e^{i\,{({\tilde{z}}^{\prime}-{\tilde{z}})\over
2}}{j_{1}({{\tilde{z}}^{\prime}-{\tilde{z}}\over
2})\over{\tilde{z}}^{\prime}-{\tilde{z}}},$ (14)

where $j_{1}$ is a spherical Bessel function. The diagonal elements of this
operator are given by $\sigma^{2}{\tilde{z}}^{2}/6$ which is in stark contrast
to the linear behavior of the potential in the ’t Hooft model. Moreover, there
is a significant coordinate-space non-locality in $V_{\rm LV}(x)$ that does
not occur for $V_{\rm tH}$. Another difference is that the potential $V_{\rm
LV}(x)$ is complex although Hermitian.

The most salient difference between the two potentials is in the high-energy
spectrum: $M^{2}_{\parallel}({\rm tH)}\approx g^{2}k,$ Callan  _et al._ (1976)
where $k$ is an integer. The high energy spectrum of the ’t Hooft model is
much more compact than that of the Li-Vary model. For large values of $k$,
$\chi_{\rm tH}(x)=\sqrt{2}\sin(\pi kx+\delta(k,x))$. The phase shift
$\delta(k,x)$ is given in Ref. Brower  _et al._ (1979).

The spectrum corresponding to the Brodsky-de Teramond IMWF model De Téramond
and Brodsky (2021) has not been provided before. The lore is that this model
does not come from a wave equation. Here we develop a wave equation of the
Sturm-Liouville form that does yield the IMWF and discuss the properties of
the solutions. The first step is asserting that the IMWF is the solution of
some wave equation. The presentation is simplified by using $m_{1,2}=m$ and
also by employing the variable $y>0$, with $y^{2}=1/(x(1-x)).$ Then using
$\phi(y)=e^{-m^{2}y^{2}\over 2\kappa^{2}}$, consider the differential
equation:

$\displaystyle-{d\over dy}f(y){d\phi(y)\over
dy}+m^{2}y^{2}\phi(y)=M^{2}\phi(y),$ (15)

with $f(y)$ to be determined. The term $m^{2}y^{2}$ is the familiar kinetic
energy term. Using the expression for $\phi(y)$ and equating the left- and
right-sides of the equation leads to the results $f=\kappa^{4}/m^{2}$ and
$M^{2}=\kappa^{2}$. A peculiar feature is that the ground-state mass is
independent of $m$. The wave function that leads to the ground state IMWF is
then given by

$\displaystyle-{\kappa^{4}\over
m^{2}}\phi^{\prime\prime}+m^{2}y^{2}\phi=M^{2}\phi,$ (16)

which is of a familiar harmonic oscillator form. This can be converted to an
equation in the variable $x$ by using ${dx\over dy}={1\over
x-1/2}x^{3/2}(1-x)^{3/2}\equiv g(x),\,$ so that Eq. (16) becomes

$\displaystyle-g(x){\kappa^{4}\over m^{2}}{d\over dx}g(x){d\phi(x)\over
dx}+{m^{2}\phi(x)\over x(1-x)}={M^{2}}\phi(x),$ (17)

which is of the Sturm-Liouville form for $0<x<1/2$ and $1/2<x<1$. The ground
state energy is $\kappa^{2}$. Each excitation in energy increases $M^{2}$ by
$2\kappa^{2}$. The feature that the spectrum does not depend on quark masses
does not seem reasonable and we therefore discard Eq. (16).

The next step is to compare two remaining approaches Li and Vary (2021); ’t
Hooft (1974) in the limit that the quark masses are non-zero, but small
compared to the strength parameter of the model. For simplicity we remain with
the case that $m_{1,2}=m$. The Li-Vary result is that

$\displaystyle M_{\parallel}^{2}({\rm LV})=2\sigma m+4m^{2}.$ (18)

Their model uses the parameters $m=15$ MeV, and $\sigma=620$ MeV. ’t Hooft
found that the ground-state wave function $\chi(x)$ is well approximated by
the form $\chi_{0}(x)=x^{\beta}(1-x)^{\beta}$. Analysis of the behavior for
the solution of the wave equation for $x\to 0\,\rm or\,1$ shows that $\beta$
should satisfy the transcendental equation
$\frac{m^{2}\pi}{g^{2}}-1+\pi\beta\cot\pi\beta=0.$ For small values of $m$,
the quantity $\beta\pi$ must also be small so that
$\beta=\sqrt{3\over\pi}{m\over g}$. Using $\chi_{0}$ on both sides of Eq. (13)
yields that $M_{\parallel}^{2}=m^{2}(\frac{2}{\beta}+4),$ or equivalently

$\displaystyle M_{\parallel}^{2}({\rm tH})=2\sqrt{\pi\over 3}gm+4m^{2}.$ (19)

Using the average value of the $u$ and $d$ quark masses, $m=3.5$ MeV Zyla  _et
al._ (2020) with $M_{L}^{2}=140$ MeV, tells us that $g=2700$ MeV, and
$\beta=0.00126.$ Thus, these two models contain a (1+1)-dimensional version of
the Gell-Mann-Oakes-Renner Gell-Mann  _et al._ (1968) relation in which the
squared mass of the ground state is proportional to the current quark mass.

(a)

(b)

Figure 1: Plots of the resulting density, $\rho(\tilde{z})$, and
$\tilde{z}^{2}\rho(\tilde{z})$ are given for all three models. The discrepancy
between the models in $\rho(\tilde{z})$ is small but finite, and the
$1/\tilde{z}^{2}$ asymptotic behavior is confirmed.

It is worthwhile to note that the use of current-quark masses in the ’t Hooft
model causes the two different models to obtain very different masses of the
first excited state. Ref Li and Vary (2021) obtains excitation energies of
approximately 1 GeV and associates these values with excited states of the
pion. In the present work, using modern values of the current quark masse Zyla
_et al._ (2020) gives the lowest excited state a mass on the order of $g$, or
about 3 GeV. This is high enough into the continuum of states with large
widths to be unobservable. Thus, the version of the ’t Hooft model used here
preserves the spectra produced by LFHQCD.

The confining aspects of the ’t Hooft model have been well-studied long ago ’t
Hooft (1974); Callan  _et al._ (1976), using a momentum-space ($x$) dependence
approach based on studying the cancellation of infrared singularities.
Another, possibly more intuitive approach, may be obtaining by examining
coordinate-space wave functions that depend on the canonically conjugate
spatial variable, ${\tilde{z}}$. An intuitive way to think about this variable
is that it is the separation between the quark and anti-quark in the direction
of motion of a pion moving with high momentum.

Coordinate-space wave functions are obtained using the transformation

$\displaystyle\chi({\tilde{z}})\equiv\int_{0}^{1}{dx\over\sqrt{2\pi}}e^{ix{\tilde{z}}}\chi(x).$
(20)

It is useful to examine the density
$\rho({\tilde{z}})\equiv|\chi({\tilde{z}})|^{2}$, a real-valued quantity, for
the three models are shown in Fig. 1(a). These results seem very similar
because of the relatively small quark masses of the three models. The
densities seem to go to 0 for large absolute values of ${\tilde{z}}.$ They do,
but there in an interesting feature seen by plotting
${\tilde{z}}^{2}\rho({\tilde{z}})$ in Fig. 1(b): the densities fall as
$1/{\tilde{z}}^{2}$.

This behavior may be understood by making a asymptotic expansion, obtained by
using $e^{i{\tilde{z}}x}=1/(i{\tilde{z}})\partial_{x}e^{i{\tilde{z}}x}$ and
the feature that $\chi(x)$ vanishes at the end points. Then

$\displaystyle\lim_{|{\tilde{z}}|\to\infty}\chi({\tilde{z}})={i\over{\tilde{z}}}\int_{0}^{1}{dx\over\sqrt{2\pi}}e^{ix{\tilde{z}}}\partial_{x}\chi(x).$
(21)

Squaring this quantity leads to the stated $1/{\tilde{z}}^{2}$ dependence.

The similarity of the behaviors for small masses suggests that the chiral
limit should be examined. In this case, $\chi(x)=1$ for all three models. Then
a simple closed-form expression for $\chi({\tilde{z}})$ can be obtained. The
resulting density, $\rho_{\chi}({\tilde{z}})$ is given by

$\displaystyle\rho_{\chi}({\tilde{z}})={2\over\pi}{\sin^{2}{{\tilde{z}}/2}\over{\tilde{z}}^{2}},$
(22)

an expression that accounts explicitly for the oscillatory behavior as well as
the $1/{\tilde{z}}^{2}$ asymptotic behavior. One may also examine the spatial
extent of the pion wave function by considering the mean-square value of
${\tilde{z}}^{2}$, that is given by the ground-state expectation value:

$\displaystyle\langle\chi|{\tilde{z}}^{2}|\chi\rangle={2\over\pi}\int_{-\infty}^{\infty}\,dz\,\sin^{2}({\tilde{z}}/2)=\infty.$
(23)

In the chiral limit the pion has an infinite spatial extent, true for all
three models.

It is necessary to see how or if this infinite size conflicts with current
understanding. First note that the elastic pion form factor of LFHQCD has
already been computed using $\chi(x)=1$; see e.g. Ref. de Teramond  _et al._
(2018). The infinite extent is buried within the integrals needed to compute
the elastic form factor, a consequence the ability to probe with only
transverse momentum transfers; see e.g. the review Miller (2010).

The relevance of the infinite longitudinal extent can be understood in analogy
with the familiar Ioffe-time argument Gribov  _et al._ (1965); Ioffe (1969);
Kovchegov and Levin (2012); Braun and Müller (2008) for deep inelastic
scattering at small values of Bjorken $x$. The idea is that an incident
virtual photon fluctuates into a $q\bar{q}$ pair. The energy difference,
$\Delta E$, between the two states is very small, so that the fluctuation
lives for a long time $\Delta t>1/\Delta E$. For virtual photons of high
energy, the extent of the fluctuation is $c\Delta t$ is very large. Now
consider a high-momentum pion incident on a target. The value of $\Delta E$
goes to 0 if both the pion and its constituents are massless. Thus the
infinite extent of a mass-less pion is part of standard lore.

We summarize. The similarities and differences between the three models of
Refs. De Téramond and Brodsky (2021); Li and Vary (2021); ’t Hooft (1974) are
exhibited. Very significant differences in the excitation spectra, at both low
and high energies, are obtained, even though all of the ground-state wave
functions are the same in the chiral limit. If the ’t Hooft model is used
along with current values of light-quark masses, the original spectrum
calculations reviewed in Ref. Brodsky  _et al._ (2015) are preserved. Finally,
the pion is shown to have an infinite extent in the longitudinal direction.

###### Acknowledgements.

This work was supported by the U.S. Department of Energy Office of Science,
Office of Nuclear Physics under Award No. DE-FG02-97ER-41014. We thank J.
Hiller, Y. Li, J. P. Vary, R. Sandapen, G.F. de Téramond and S.J. Brodsky for
helpful comments.

## References

  * De Téramond and Brodsky (2021) Guy F. De Téramond and Stanley J. Brodsky, “Longitudinal dynamics and chiral symmetry breaking in holographic light-front QCD,” (2021), arXiv:2103.10950 [hep-ph] .
  * Li and Vary (2021) Yang Li and James P. Vary, “Light-front holography with chiral symmetry breaking,” (2021), arXiv:2103.09993 [hep-ph] .
  * Ahmady  _et al._ (2021a) Mohammad Ahmady, Harleen Dahiya, Satvir Kaur, Chandan Mondal, Ruben Sandapen, and Neetika Sharma, “Extending light-front holographic QCD using the ’t Hooft Equation,” Phys. Lett. B 823, 136754 (2021a), arXiv:2105.01018 [hep-ph] .
  * Ahmady  _et al._ (2021b) Mohammad Ahmady, Satvir Kaur, Sugee Lee MacKay, Chandan Mondal, and Ruben Sandapen, “Hadron spectroscopy using the light-front holographic Schrödinger equation and the ’t Hooft equation,” Phys. Rev. D 104, 074013 (2021b), arXiv:2108.03482 [hep-ph] .
  * Brodsky  _et al._ (2015) Stanley J. Brodsky, Guy F. de Teramond, Hans Gunter Dosch, and Joshua Erlich, “Light-Front Holographic QCD and Emerging Confinement,” Phys. Rept. 584, 1–105 (2015), arXiv:1407.8131 [hep-ph] .
  * Sheckler and Miller (2021) Aiden B. Sheckler and Gerald A. Miller, “Mystery of Bloom-Gilman duality: A light-front holographic QCD perspective,” Phys. Rev. D 103, 096018 (2021), arXiv:2101.00100 [hep-ph] .
  * Chabysheva and Hiller (2013) S. S. Chabysheva and John R. Hiller, “Dynamical model for longitudinal wave functions in light-front holographic QCD,” Annals Phys. 337, 143–152 (2013), arXiv:1207.7128 [hep-ph] .
  * Miller and Brodsky (2020) Gerald A. Miller and Stanley J. Brodsky, “Frame-independent spatial coordinate $\tilde{z}$: Implications for light-front wave functions, deep inelastic scattering, light-front holography, and lattice QCD calculations,” Phys. Rev. C 102, 022201 (2020), arXiv:1912.08911 [hep-ph] .
  * Brodsky and de Teramond (2009) Stanley J. Brodsky and Guy F. de Teramond, “Light-Front Holography and AdS/QCD: A First Approximation to QCD,” eCONF C0906083, 17 (2009), arXiv:0910.0625 [hep-ph] .
  * Karch  _et al._ (2006) Andreas Karch, Emanuel Katz, Dam T. Son, and Mikhail A. Stephanov, “Linear confinement and AdS/QCD,” Phys. Rev. D 74, 015005 (2006), arXiv:hep-ph/0602229 .
  * Callan  _et al._ (1976) Curtis G. Callan, Jr., Nigel Coote, and David J. Gross, “Two-Dimensional Yang-Mills Theory: A Model of Quark Confinement,” Phys. Rev. D 13, 1649 (1976).
  * Li  _et al._ (2016) Yang Li, Pieter Maris, Xingbo Zhao, and James P. Vary, “Heavy Quarkonium in a Holographic Basis,” Phys. Lett. B 758, 118–124 (2016), arXiv:1509.07212 [hep-ph] .
  * ’t Hooft (1974) Gerard ’t Hooft, “A Two-Dimensional Model for Mesons,” Nucl. Phys. B 75, 461–470 (1974).
  * Brower  _et al._ (1979) R. C. Brower, W. L. Spence, and J. H. Weis, “Bound states and asymptotic limits for quantum chromodynamics in two dimensions,” Phys. Rev. D 19, 3024–3049 (1979).
  * Ellis (1977) John R. Ellis, “Applications of Two-Dimensional QCD,” Acta Phys. Polon. B 8, 1019–1059 (1977).
  * Zyla  _et al._ (2020) P.A. Zyla _et al._ (Particle Data Group), “Review of Particle Physics,” PTEP 2020, 083C01 (2020).
  * Gell-Mann  _et al._ (1968) Murray Gell-Mann, R. J. Oakes, and B. Renner, “Behavior of current divergences under SU(3) x SU(3),” Phys. Rev. 175, 2195–2199 (1968).
  * de Teramond  _et al._ (2018) Guy F. de Teramond, Tianbo Liu, Raza Sabbir Sufian, Hans Günter Dosch, Stanley J. Brodsky, and Alexandre Deur (HLFHS), “Universality of Generalized Parton Distributions in Light-Front Holographic QCD,” Phys. Rev. Lett. 120, 182001 (2018), arXiv:1801.09154 [hep-ph] .
  * Miller (2010) Gerald A. Miller, “Transverse Charge Densities,” Ann. Rev. Nucl. Part. Sci. 60, 1–25 (2010), arXiv:1002.0355 [nucl-th] .
  * Gribov  _et al._ (1965) V. N. Gribov, B. L. Ioffe, and I. Ya. Pomeranchuk, “What is the range of interactions at high-energies,” Yad. Fiz. 2, 768–776 (1965).
  * Ioffe (1969) B. L. Ioffe, “Space-time picture of photon and neutrino scattering and electroproduction cross-section asymptotics,” Phys. Lett. B 30, 123–125 (1969).
  * Kovchegov and Levin (2012) Yuri V. Kovchegov and Eugene Levin,  _Quantum chromodynamics at high energy_, Vol. 33 (Cambridge University Press, 2012).
  * Braun and Müller (2008) V. Braun and Dieter Müller, “Exclusive processes in position space and the pion distribution amplitude,” Eur. Phys. J. C 55, 349–361 (2008), arXiv:0709.1348 [hep-ph] .