# Longitudinal dynamics for mesons on the light cone

Yang Li  Department of Modern Physics, University of Science and Technology of
China, Hefei, Anhui 230026, China  James P. Vary  Department of Physics and
Astronomy, Iowa State University, Ames, IA 50011, USA

###### Abstract

We survey a set of proposed light-front models for confinement in the valence
quark sector of the mesons and portray similarities as well as differences. We
present the spectroscopies for the light mesons that result from a selection
of longitudinal confinement forms. We note that the Sturm-Liouville theory
provides a unifying framework for many elements of this comparison.

##  I Introduction

The light-front Schrödinger wave equation (LFSWE) provides a relativistic
semiclassical first approximation to QCD for mesons represented as quark and
antiquark bound states deTeramond:2008ht ,

$\Big{[}\frac{\vec{k}^{2}_{\perp}+m_{q}^{2}}{x}+\frac{\vec{k}^{2}_{\perp}+m^{2}_{\bar{q}}}{1-x}+V\Big{]}\psi(x,\vec{k}_{\perp})=M^{2}\psi(x,\vec{k}_{\perp}).$
(1)

Here, the effective potential $V$ encodes the non-perturbative dynamics of QCD
at low energy resolution. Remarkably, Brodsky and de Téramond discovered that
the LFSWE in the chiral limit ($m_{q}=m_{\bar{q}}=0$) can be identified with
the equation of motion of strings in the fifth dimension of Anti-de Sitter
(AdS5) space thus establishing an extraordinary connection between light-front
QCD and AdS/QCD, a bottom-up approach known as light-front holography (LFH,
Brodsky:2006uqa ). LFH has been successfully applied to describing the hadron
spectroscopy, including the spectra of mesons Brodsky:2006uqa ; Dosch:2015bca
; Dosch:2016zdv ; Nielsen:2018ytt ; Zou:2019tpo , the infrared behavior of the
strong coupling Deur:2016tte , hadron form factors Brodsky:2006uqa ;
Brodsky:2007hb ; Brodsky:2008pf ; Sufian:2016hwn  and parton distributions
deTeramond:2018ecg ; Liu:2019vsn . An excellent review of this approach is
Ref. Brodsky:2014yha .

Light-front holography addresses the chiral limit where only the dynamics of
the transverse degrees of freedom (d.o.f.) are considered Li:2022ytx .
Recently, attention has focused on the role of longitudinal dynamics
Chabysheva:2012fe ; Glazek:2013jba ; Trawinski:2014msa ; Li:2015zda ;
Li:2021jqb ; DeTeramond:2021jnn ; Ahmady:2021lsh ; Ahmady:2021yzh ;
Shuryak:2021hng ; Shuryak:2021mlh ; Shuryak:2022thi . The longitudinal
dynamics is important for incorporating finite quark masses, chiral dynamics
and longitudinal excitations as well as for identifying the physical states in
the excitation spectrum Gutsche:2012ez . Unlike the transverse d.o.f., the
longitudinal dynamics in the presence of transverse dynamics of the LFSWE has
not been posed from underlying principles.

In light of various proposals in the literature, it is tempting to compare
these longitudinal confining interactions. Ahmady et al. compared the hadron
spectra employing the ’t Hooft model and our model Ahmady:2021lsh . Shuryak &
Zahed compared the effective light-front Hamiltonian obtained from the Nambu-
Goto string having massive ends with a phenomenological model we (with others)
introduced Shuryak:2021hng . Weller & Miller further compared the longitudinal
confining potentials of these three models and the corresponding wave
functions Weller:2021wog . In this work, we discuss the role of longitudinal
dynamics in LFSWE. In particular, we focus on the endpoint behavior of the
wave functions and the scaling of the eigenvalues. We show that this behavior
is closely related to the singularities of the corresponding LFSWE. This
observation enables us to consider a broader class of longitudinal confining
potentials based on the Sturm-Liouville theory. We present several concrete
examples.

This perspective is consistent with the traditional quantum many-body approach
to self-bound systems, e.g., atoms, molecules and nuclei. There, first
approximations are employed to chart out a set of “orbits” or “shells” of the
system. Then, the wave functions of these orbits are adopted as the basis
within which the full Hamiltonian operator is diagonalized to obtain the
“exact” results. In light of this link to the quantum many-body treatment,
basis light-front quantization (BLFQ Vary:2009gt ), we referred to our recent
LFSWE approach for light mesons as BLFQ0. In the heavy quark systems, the non-
relativistic potential model is a valid first approximation Godfrey:1985xj .
This model is augmented to the Hamiltonian QCD in Coulomb gauge to address the
full QCD dynamics Szczepaniak:1995cw  (see Ref. Reinhardt:2017pyr  for a
recent review).In the light hadron sector, light-front holography plays an
analogous role with its remarkable phenomenological successes. Recently, it
was argued that the enpoint asymptotics of the longitudinal dynamics is
essential for implementing the Gell-Mann-Oakes-Renner (GMOR) relation – a
direct consequence of chiral symmetry breaking in the light sector
Gutsche:2012ez ; Li:2021jqb .

The reminder of this work is organized as follows. Sect II introduces the
LFSWE within the separation of variables ansatz and its application in QCD,
the light-front holographic QCD. We then discuss the need for a longitudinal
dynamics and compare various proposals from the literature in Sect. III. In
the next section, Sect. IV, we construct a general class of longitudinal
potentials based on the Sturm-Liouville theory. Several specific cases are
investigated. Finally, we conclude in Sect. V.

##  II Formalism

Formally, the LFSWE is a low energy effective theory of QCD. Therefore, it may
be obtained from a Hamiltonian renormalization group evolution, either the
Okubo-Suzuki-Lee type Okubo:1954 ; Suzuki:1980 ; Suzuki:1982a ; Suzuki:1982b ;
Suzuki:1983 ; Suzuki:1994  or the flow equation type Glazek:1993rc ;
Glazek:1994qc ; Wegner:1994 . Alternatively, it can be viewed as the
truncation up to the valence Fock sector hence the Tamm-Dancoff type of
approximation is applied Lepage:1980fj ; Perry:1990mz . These formal
connections to QCD make it possible to systematically improve the LFSWE, as
shown in Refs. Wilson:1994fk ; Glazek:2012qj ; Glazek:2017rwe ;
Serafin:2018aih . In either approach, one can write down the wave equation,

$\Big{[}\frac{\vec{k}_{\perp}^{2}+m^{2}_{q}+\Sigma_{s}(x,\vec{k}_{\perp})}{x}+\frac{\vec{k}_{\perp}^{2}+m^{2}_{\bar{q}}+\Sigma_{\bar{s}}(1-x,-\vec{k}_{\perp})}{1-x}\Big{]}\psi_{s\bar{s}/h}^{(J,m_{J})}(x,\vec{k}_{\perp})\\\
+\sum_{s^{\prime},\bar{s}^{\prime}}\int\frac{{\mathrm{d}}x^{\prime}}{2x^{\prime}(1-x^{\prime})}\int\frac{{\mathrm{d}}^{2}k_{\perp}^{\prime}}{(2\pi)^{3}}V_{s\bar{s}s^{\prime}\bar{s}^{\prime}}(x,\vec{k}_{\perp},x^{\prime},\vec{k}^{\prime}_{\perp})\psi_{s^{\prime}\bar{s}^{\prime}/h}^{(J,m_{J})}(x^{\prime},\vec{k}^{\prime}_{\perp})=M^{2}_{h}\psi_{s\bar{s}/h}^{(J,m_{J})}(x,\vec{k}_{\perp}).$
(2)

$\Sigma_{\sigma}$ is the quark self-energy and needs to be solved from the gap
equation. To the first approximation, we can absorb the self-energy in an
effective quark mass. $J$ is the total angular momentum and $m_{J}$ is its
magnetic projection. $\psi$ is the wave function which is normalized according
to,

$\sum_{s\bar{s}}\int_{0}^{1}\frac{{\mathrm{d}}x}{2x(1-x)}\int\frac{{\mathrm{d}}^{2}k_{\perp}}{(2\pi)^{3}}\bar{\psi}_{s\bar{s}/h}^{(Jm_{J})}(x,\vec{k}_{\perp})\psi_{s\bar{s}/h^{\prime}}^{(J^{\prime}m_{J}^{\prime})}(x,\vec{k}_{\perp})=\delta_{hh^{\prime}}.$
(3)

Subscript $h$ indicates the hadron species. This wave equation corresponds to
the effective Hamiltonian
$H_{\mathrm{eff}}=({\vec{k}_{\perp}^{2}+m^{2}_{q}})/{x}+({\vec{k}_{\perp}^{2}+m^{2}_{\bar{q}}})/({1-x})+V$.
The normalization adopted here stems from the Lorentz invariant phase space
element,

$\int\frac{{\mathrm{d}}^{4}p}{(2\pi)^{4}}2\pi\delta(p^{2}-m^{2})=\int\frac{{\mathrm{d}}^{3}p}{(2\pi)^{3}2p^{0}}\theta(p^{0})=\int\frac{{\mathrm{d}}p^{+}{\mathrm{d}}^{2}p_{\perp}}{(2\pi)^{3}2p^{+}}\theta(p^{+})$
(4)

To associate the eigenvalues and eigenfunctions to hadrons, we need to compute
the discrete symmetries $P$ and $C$. They are represented as Brodsky:2006ez ,

$\displaystyle C=\;$
$\displaystyle\int_{0}^{1}\frac{{\mathrm{d}}x}{2x(1-x)}\int\frac{{\mathrm{d}}^{2}k_{\perp}}{(2\pi)^{3}}\psi^{*}_{s\bar{s}}(x,\vec{k}_{\perp})\psi_{\bar{s}s}(1-x,-\vec{k}_{\perp}),$
(5) $\displaystyle(-1)^{J}P=\;$
$\displaystyle\int_{0}^{1}\frac{{\mathrm{d}}x}{2x(1-x)}\int\frac{{\mathrm{d}}^{2}k_{\perp}}{(2\pi)^{3}}\psi^{*}_{s\bar{s}}(x,\vec{k}_{\perp})\psi_{-s-\bar{s}}(x,\tilde{k}_{\perp}),$
(6)

where $\tilde{k}_{\perp}=(-k_{x},k_{y})$. The mirror parity $(-1)^{J}P$ is
employed here since it is a more convenient observable on the light front
Brodsky:2006ez . Since the mirror parity flips the sign of the spin, we can
use the $m_{J}=0$ state to compute the mirror parity.

###  II.1 Separation of variables

(a)  (b)

Figure 1: Representative light-front QCD interactions contributed to the (a)
transverse and (b) longitudinal effective interactions, respectively.

In the massless limit, the kinetic energy becomes
${\vec{k}_{\perp}^{2}}/{x(1-x)}$, which depends only on a 2D variable
$\vec{\kappa}_{\perp}\equiv{\vec{k}_{\perp}}/{\sqrt{x(1-x)}}$. It is natural
to assume the system in the chiral limit also depends only on
$\vec{\kappa}_{\perp}$ and its conjugate coordinate
$\vec{\zeta}_{\perp}=i\nabla_{\kappa_{\perp}}=\sqrt{x(1-x)}\vec{r}_{\perp}$
Brodsky:2006uqa . In the more general case, especially when the quark masses
are not zero, we can assume that the interaction is separable:
$V=V_{\perp}+V_{\|}$, an ansatz first explicitly introduced by Chabysheva and
Hiller Chabysheva:2012fe  (cf. Ref. Brodsky:2008pf ). We note that the
separation of the transverse and longitudinal dynamics is natural in light-
front QCD Brodsky:1997de . Fig. 1 shows two representative light-front QCD
interactions. Fig. 1a involves the exchange of the purely transverse gluons
while Fig. 1b involves the instantaneous interaction in the longitudinal
direction. The latter survives in 1+1D as the Schwinger-’t Hooft type of
interaction Schwinger:1962tp ; tHooft:1974pnl ; Mo:1992sv .

Of course, the true effective interaction between a quark and an antiquark may
not be separable, for instance, the one-gluon-exchange interaction. A proper
metric of the separability is the entanglement entropy $S_{\perp}$ between the
transverse and the longitudinal d.o.f.’s, which is defined as the von Neumann
entropy for the reduced density matrix,
$S_{\perp}=-\mathrm{tr}\big{[}\rho_{\perp}\log\rho_{\perp}\big{]}$ Eisert:2010
. Fig. 2 shows the entanglement entropy between the transverse and the
longitudinal d.o.f.’s for charmonium states obtained from a non-separable
interaction Li:2017mlw . The effective Hamiltonian adopts a separable
confining interaction (holographic confinement plus the Li-Maris-Zhao-Vary
longitudinal confinement) along with a non-separable one-gluon-exchange
interaction. As one can see, $S_{\perp}$ is generally small for the ground
states – indicating that these states are physically separable within this
model. In this work, we adopt the separation of variable ansatz, which should
provide a good description for the low-lying states.

Figure 2: Entanglement entropy between the transverse and the longitudinal
d.o.f.’s for charmonium states obtained from BLFQ. The effective Hamiltonian
adopts a separable confining (holographic confinement plus the LMZV
longitudinal confinement) along with a non-separable one-gluon-exchange
interaction Li:2017mlw .

With this ansatz, the mass eigenvalues $M^{2}=M_{\perp}^{2}+M^{2}_{\|}$ and
the wave function are also separable,

$\psi(x,\vec{k}_{\perp})=\sqrt{4\pi}\phi(\vec{k}_{\perp}/\sqrt{x(1-x})\chi(x).$
(7)

We obtain two LFSWEs Chabysheva:2012fe ,

$\displaystyle\Big{[}{\vec{\kappa}^{2}_{\perp}}+V(\vec{\zeta}_{\perp})\Big{]}\phi(\vec{\kappa}_{\perp})=M^{2}_{\perp}\phi(\vec{\kappa}_{\perp}),$
(8)
$\displaystyle\Big{[}\frac{m_{q}^{2}}{x}+\frac{m^{2}_{\bar{q}}}{1-x}+V_{\|}(\tilde{z})\Big{]}\chi(x)=M^{2}_{\|}\chi(x).$
(9)

Here, the longitudinal effective potential $V_{\|}$ depends on the boost
invariant longitudinal coordinate $\tilde{z}={\frac{1}{2}}p^{+}x^{-}$ (or
other related longitudinal coordinate variables, e.g. (24)), introduced by
Miller and Brodsky Miller:2019ysh . We note that even in the massless limit,
Eq. (9) is still required to have non-trivial solutions. Hereafter we consider
specific spin configurations (the leading-twist spin components) and the
explicit spin indices will be suppressed, except for the purpose of state
identification. N.B. the separation of variables ansatz requires that the
variables $\vec{\zeta}_{\perp}$ and $\tilde{z}$ (or other related longitudinal
coordinate variables) are independent. Since
$\vec{\zeta}_{\perp}=\sqrt{x(1-x)}\vec{r}_{\perp}$ depends on $x$, $\tilde{z}$
should be defined as a differential operator against $\vec{\zeta}_{\perp}$,
viz. $\tilde{z}=i\frac{\partial}{\partial x}\big{|}_{\vec{\zeta}_{\perp}}$.

It is convenient to normalize the transverse and longitudinal wave functions,

$\int\frac{{\mathrm{d}}^{2}\kappa_{\perp}}{(2\pi)^{2}}\big{|}\phi(\vec{\kappa}_{\perp})\big{|}^{2}=1,\quad\int_{0}^{1}{\mathrm{d}}x\big{|}\chi(x)\big{|}^{2}=1.$
(10)

It is also useful to introduce a longitudinal wave function
$X(x)=\sqrt{x(1-x)}\chi(x)$ with the normalization,

$\int_{0}^{1}\frac{{\mathrm{d}}x}{x(1-x)}\big{|}X(x)\big{|}^{2}=1.$ (11)

In our convention, $X(x)$ is proportional to the distribution amplitude. This
is related to the alternative normalization convention of the wave function
$\widetilde{\Psi}(x,\vec{\zeta}_{\perp})$ or its transverse Fourier conjugate,

$\Psi(x,\vec{\kappa}_{\perp})=\int\frac{{\mathrm{d}}^{2}\kappa_{\perp}}{(2\pi)^{2}}\widetilde{\Psi}(x,\vec{\zeta}_{\perp})e^{-i\vec{\kappa}_{\perp}\cdot\vec{\zeta}_{\perp}}=\psi(x,\vec{k}_{\perp})/\sqrt{x(1-x)}$
(12)

Under this convention, the normalizations of the wave functions are (recall
$\vec{\kappa}_{\perp}=\vec{k}_{\perp}/\sqrt{x(1-x)}$,
$\vec{\zeta}_{\perp}=\sqrt{x(1-x)}\vec{r}_{\perp}$),

$\int_{0}^{1}{\mathrm{d}}x\int\frac{{\mathrm{d}}^{2}k_{\perp}}{16\pi^{3}}\big{|}\Psi(x,\vec{k}_{\perp}/\sqrt{x(1-x)})\big{|}^{2}=1,\quad\frac{1}{4\pi}\int_{0}^{1}\frac{{\mathrm{d}}x}{x(1-x)}\int{\mathrm{d}}^{2}\zeta_{\perp}\big{|}\widetilde{\Psi}(x,\vec{r}_{\perp})\big{|}^{2}=1.$
(13)

With this convention, the separable wave function takes the form,
$\widetilde{\Psi}(x,\vec{\zeta}_{\perp})=\sqrt{4\pi}\widetilde{\phi}(\zeta_{\perp})X(x)$,
where $\tilde{\phi}$ is the Fourier transform of $\phi$,

$\tilde{\phi}(\zeta_{\perp})=\int\frac{{\mathrm{d}}^{2}\kappa_{\perp}}{(2\pi)^{2}}e^{i\vec{\kappa}_{\perp}\cdot\vec{\zeta}_{\perp}}\phi(\vec{\kappa}_{\perp}).$
(14)

As an example, the pion wave function in LFH in the massless limit is
$\chi(x)=1$ and $X(x)=\sqrt{x(1-x)}$.

###  II.2 Light-front holography

In LFH, the transverse effective potential is determined by the dilaton
profile $\Phi(z)$ in AdS/QCD,
$V_{\perp}=(1/2)\Phi^{\prime\prime}+(1/4)\Phi^{\prime
2}+\big{[}(2J-3)/2z\big{]}\Phi^{\prime}$ Brodsky:2008pg . For the soft-wall
profile, $V_{\perp}=\kappa^{4}\zeta^{2}_{\perp}+2\kappa^{2}(J-1)$. Here,
$\kappa$ is the strength of the confinement in mass dimension. $J$ is the
total angular momentum. The parton transverse separation
$\zeta_{\perp}=\sqrt{x(1-x)}r_{\perp}$ is mapped to the fifth coordinate $z$
in Anti-de Sitter space. The mass eigenvalues are,

$M^{2}_{nmJ}=4\kappa^{2}\Big{(}n+\frac{J+|m|}{2}\Big{)}.$ (15)

where $n,m$ are the radial and angular quantum numbers in the transverse plane
($\vec{\zeta}_{\perp}$). The corresponding wave functions are 2D harmonic
oscillator functions,

$\phi_{nm}(\vec{\zeta}_{\perp})=\kappa\sqrt{\frac{n!}{\pi(n+|m|)!}}e^{-\frac{\kappa^{2}\zeta_{\perp}^{2}}{2}}(\kappa\zeta)^{|m|}_{\perp}L_{n}^{|m|}(\kappa^{2}\zeta_{\perp}^{2}),$
(16)

where $n$ and $m$ are the transverse radial and angular quantum numbers. $N$
is a normalization constant. The ground state $n=0,m=0,J=0$ describes the
pion, whose mass is predicted to vanish. Its wave function is Gaussian,

$\psi_{\pi}(x,\vec{k}_{\perp})=\frac{4\pi}{\kappa}\exp\Big{[}-\frac{\vec{k}^{2}_{\perp}}{2\kappa^{2}x(1-x)}\Big{]}.$
(17)

To incorporate the finite quark masses, Brodsky and de Téramond employed a
longitudinal wave function without explicitly introducing the longitudinal
dynamics Brodsky:2008pf . This ansatz, known as the invariant mass ansatz
(IMA), is based on the observation that with the presence of finite quark
mass, the light-front kinetic energy becomes,

$\frac{\vec{k}^{2}_{\perp}}{x(1-x)}\;\to\;\frac{\vec{k}^{2}_{\perp}+m_{q}^{2}}{x}+\frac{\vec{k}^{2}_{\perp}+m_{\bar{q}}^{2}}{1-x}.$
(18)

Hence it is instructive to make the substitution in the wave functions. Then
the pion wave function becomes,

$\psi_{\pi}(x,\vec{k}_{\perp})=N\frac{4\pi}{\kappa}\exp\Big{[}-\frac{\vec{k}^{2}_{\perp}+m_{q}^{2}}{2\kappa^{2}x(1-x)}\Big{]}.$
(19)

In other words, the longitudinal wave function becomes,

$\chi_{\textsc{ima}}(x)=N\exp\Big{[}-\frac{m_{q}^{2}}{2\kappa^{2}x(1-x)}\Big{]}.$
(20)

The pion mass is also shifted to a non-zero value by,

$M_{\pi}^{2}=\int_{0}^{1}\frac{{\mathrm{d}}x}{2x(1-x)}\int\frac{{\mathrm{d}}^{2}k_{\perp}}{(2\pi)^{3}}\big{|}\psi_{\pi}(x,\vec{k}_{\perp})\big{|}^{2}\frac{m_{q}^{2}}{x(1-x)}.$
(21)

In the vicinity of the chiral limit, the theory predicts a near quadratic
quark mass dependence of the pion mass, $M_{\pi}^{2}\approx
2m_{q}^{2}(\ln\kappa^{2}/m_{q}^{2}-\gamma_{\textsc{e}})$, where
$\gamma_{\textsc{e}}\approx 0.577216$ is Euler’s constant Li:2021jqb . The
scaling of the pion mass $M_{\pi}$ as a function of the quark mass $m_{q}$ is
inconsistent with the results from chiral symmetry breaking, which predicts a
linear dependence instead, a result also known as the Gell-Mann-Oakes-Renner
(GMOR) relation Gell-Mann:1968hlm ,

$f^{2}_{\pi}M_{\pi}^{2}=2m_{q}\langle 0|\bar{q}q|0\rangle+O(m_{q}^{2}).$ (22)

Gutsche et al. suggested replacing the longitudinal wave function (20) from
IMA by a power-law like function to generate the GMOR relation Gutsche:2012ez
.

Another issue with IMA is the lack of longitudinal excitations, which leads to
inconsistencies in state identification Li:2015zda ; Li:2021jqb . For example,
in LFH the mass of $\rho(1D)$ is degenerate with the mass of $\rho(2S)$, and
is 400 MeV lower than the mass of $\rho(1700)$, which is usually identified as
$\rho(1D)$. Indeed, in LFH, $\rho(1700)$ is identified as $\rho(3S)$ instead.
Similarly, $K_{1}(1400)$ in LFH is identified as $2{}^{3}P_{1}$ as opposed to
the conventional identification $1{}^{1}P_{1}$. The resulting mass is higher
than the experimental measurements.

For states identified with quantum numbers $(n,m,S,J)$ in LFHQCD with IMA,
$(-1)^{J}P=(-1)^{m+S+1}$, $C=(-1)^{m+S}$. From these relations, quantum number
assignments of the scalars $a_{0}$ ($0^{++}$), the axial vectors $b_{1}$
($1^{+-}$), and the axial tensors $a_{2}$ ($2^{++}$) etc. are not compatible
with the discrete symmetries. The source of the discrepancy is the lack of
longitudinal degrees of freedom (d.o.f.) in LFH with IMA. Indeed, in the basis
function approach to the same quantum numbers, these states are excited in the
longitudinal direction Li:2017mlw .

A possible alternative is to invoke the total orbital angular momentum $L$ and
to assume that the LFH wave functions with IMA are only valid for states
$L_{z}=|m|=L$. Note that $L$ is not a good quantum number in relativistic
quantum mechanics. The obtained wave function may not be in leading twist
($L_{z}=0$), either. By contrast, incorporating the longitudinal dynamics will
provide the complete set of wave functions and the assignment of $P$ and $C$
does not rely upon these interpretations.

##  III Longitudinal confinement

The issues encountered in LFH point to the need for longitudinal dynamics.
Indeed, as we mentioned, the decoupling of the transverse and longitudinal
d.o.f. in Eqs. (8–9) does not necessarily imply the absence of the
longitudinal dynamics. The presence of longitudinal confinement naturally
leads to longitudinal excitations which ensure the correct state
identification, similar to the non-relativistic cases but with exact discrete
quantum numbers on the light front. Now, the new set of quantum numbers
becomes $(n,m,l,S,J)$, where $l$ counts the number of nodes in the
longitudinal direction. Under charge conjugation, $x\leftrightarrow(1-x)$
generates an extra sign $(-1)^{l}$. Hence, the charge conjugation quantum
number should become $C=(-1)^{m+l+S}$. One can see that this resolves
discrepancies associated with state identification in LFHQCD for states with
longitudinal excitations. For non-separable states involving more than one
eigenstate of the effective Hamiltonian, one can always go back Eqs. (5–6).
Hence, a longitudinal confinement is required to maintain the 3D structure of
hadrons.

Of course, the longitudinal dynamics is also needed to incorporate the finite
quark masses dynamically. This was first done by Chabysheva and Hiller
utilizing the ’t Hooft model Chabysheva:2012fe ,

$\big{(}V_{\text{tH}}\circ\chi\big{)}(x)=\frac{g^{2}}{\pi}\mathrm{P}\int{\mathrm{d}}y\frac{\chi(x)-\chi(y)}{(x-y)^{2}},$
(23)

where the principle value prescription
$\mathrm{P}({1}/{x^{2}})\equiv({1}/{2})\big{[}{1}/{(x^{2}+i\epsilon)}+{1}/{(x^{2}-i\epsilon)}\big{]}$
is applied to the pole in the integrand. Chabysheva and Hiller focused on a
specific solution with $m_{q}=g/\sqrt{\pi}$. The ’t Hooft model as the
longitudinal confinement was recently revisited by Ahmady et al.
Ahmady:2021lsh ; Ahmady:2021yzh .

Another longitudinal confinement based on a construction of the harmonic
oscillator potential was proposed by Głazek and Trawiński (GT) Glazek:2013jba
; Trawinski:2014msa . They introduce a third momentum
$\kappa_{\|}\equiv(m_{q}+m_{\bar{q}})(x-\hat{m}_{q})/\sqrt{x(1-x)}$, where
$\hat{m}_{q}=m_{q}/(m_{q}+m_{\bar{q}})$. The longitudinal confining potential
is defined as the harmonic oscillator potential of the conjugate coordinate
$\zeta_{\|}$ of $\kappa_{\|}$, viz.

$V_{\text{GT}}=\kappa^{4}\zeta_{\|}^{2}.$ (24)

The resulting mass squared eigenvalues obey the linear Regge trajectory,
$M^{2}_{l}=\kappa^{2}(l+\frac{1}{2})$. The ground-state wave function is
Gaussian,

$\chi_{\text{GT}}(x)=N\exp\Big{[}-(m_{q}+m_{\bar{q}})^{2}\frac{(x-\hat{m}_{q})^{2}}{2\kappa^{2}x(1-x)}\Big{]}=N^{\prime}\exp\Big{[}-\frac{1}{2\kappa^{2}}\Big{(}\frac{m_{q}^{2}}{x}+\frac{m_{\bar{q}}^{2}}{1-x}\Big{)}\Big{]}$
(25)

Therefore, this model can be viewed as a generalization of the IMA.

The quark masses adopted in Refs. Chabysheva:2012fe ; Glazek:2013jba ;
Trawinski:2014msa  are of the order of the constituent quark mass ($\sim$ 300
MeV). In Ref. Li:2015zda , Li, Maris, Zhao and Vary proposed an analytically
solvable longitudinal confining potential and applied it to the quarkonium in
the framework of BLFQ. In Ref. Li:2021jqb , we further show that in the
vicinity of the chiral limit, the LMZV model reproduces the GMOR relation as
well as the signature power-law wave function. This is sometimes termed BLFQ0
to distinguish from the applications to the heavy flavors where the
longitudinal coupling $\sigma$ is chosen to match to the transverse confining
strength in the NR limit. The longitudinal confinement takes the form,

$\big{(}V_{\text{LMZV}}\circ\chi\big{)}(x)=-\sigma^{2}\frac{{\mathrm{d}}}{{\mathrm{d}}x}\big{(}x(1-x)\frac{{\mathrm{d}}}{{\mathrm{d}}x}\chi(x)\big{)},$
(26)

where $\sigma$ is the strength of the longitudinal confinement. The same
longitudinal confinement is adopted by de Téramond and Brodsky in Ref.
DeTeramond:2021jnn  to incorporate chiral symmetry breaking. Both works note
the close relation between the LMZV/BLFQ0 model and the ’t Hooft model. The
latter is known to exhibit the chiral symmetry breaking via the
Berezinskii–Kosterlitz–Thouless (BKT) mechanism Ji:2020bby . From the ground
state masses in the vicinity of the chiral limit
$M^{2}_{\textsc{lmzv}}=\sigma(m_{q}+m_{\bar{q}})+O(m_{q,\bar{q}}^{2})$, and
$M^{2}_{\text{tH}}=g\sqrt{\frac{\pi}{3}}(m_{q}+m_{\bar{q}})+O(m_{q,\bar{q}}^{2})$,
we can identify
$\sigma=g\sqrt{\frac{\pi}{3}}=|\langle\bar{q}q\rangle|/f^{2}_{\pi}$. The
longitudinal wave function of the LMZV/BLFQ0 model is,

$\chi_{\textsc{lmzv}}(x)=Nx^{\alpha_{1}}(1-x)^{\alpha_{2}}$ (27)

where $\alpha_{i}=m_{i}/\sigma$. The ’t Hooft wave function is,

$\chi_{\text{tH}}(x)\dot{=}Nx^{\beta_{1}}(1-x)^{\beta_{2}}.$ (28)

The exponents $\beta_{1,2}$ satisfy $\pi
m_{i}^{2}/g^{2}=1-\pi\beta_{i}\cot\pi\beta_{i}$. The solution in the vicinity
of the chiral limit is, $\beta_{i}=(m_{i}/g)(\sqrt{3/\pi})=\alpha_{i}$,
identical to those of the LMZV/BLFQ0 model. In other words, the ground-state
mass eigenvalue and eigenfunctions of the both models are identical. These two
models differ in the excited state masses. For large excitation, the mass
eigenvalues of the ’t Hooft model obey the linear Regge trajectory,
$M^{2}_{n}\propto n$, whereas those of the LMZV/BLFQ0 model are quadratic,
$M^{2}_{n}\propto n(n+1)$. A detailed comparison of the hadron spectra using
the ’t Hooft model and the LMZV/BLFQ0 model with the original predictions from
LFH with IMA is reported by Ahmady et al. Ahmady:2021yzh .

Note that Ahmady et al. imposed a constraint $l\geq n+m$ on the hadron spectra
and abandoned low-lying states without the required longitudinal excitations.
For example, the scalar ($0^{++}$) $a_{0}(980)$ is identified as
$(n,m,S)=(0,1,1)$ in LFH with IMA, with a predicted mass
$M_{\textsc{ima}}=0.76\,\text{GeV}$. In the LMZV/BLFQ0 model, this particle is
identified as $(n,m,l,S)=(0,0,1,1)$ with a predicted mass
$M_{\textsc{lmzv}}=0.91\,\text{GeV}$. Both models are consistent with the
quark model identification $1{}^{3}\\!P_{0}$ while the latter improves the
theoretical prediction. By contrast, in Ahmady et al., this state is
identified as a tetra-quark with $(n_{T},m_{T},l_{T},S_{T})=(0,0,7,0)$. The
resulting mass is $M_{\textsc{amkms}}=1.39\,\text{GeV}$. In fact, in Ref.
Ahmady:2021yzh , all scalars, including $\chi_{c0}(1P)$ and $\chi_{b0}(1P)$,
are identified as tetra-quarks.

The pion wave functions of the LMZV/BLFQ0 model are similar to those adopted
by Gutsche et al. Gutsche:2012ez . The pion mass is exactly two times of
Gutsche et al.’s result. The difference is caused by the presence of the
confining interaction in LMZV/BLFQ0 model, which contributes to the other half
of the pion mass.

Figure 3 compares the mass spectra of unflavored light mesons and kaons as
predicted by LFH with IMA Brodsky:2014yha , LFH with with the power-law like
wave function by Gutsche et al.Gutsche:2012ez , LFH with LMZV/BLFQ0 Li:2021jqb
and LFH with ’t Hooft potential by Ahmady et al. Ahmady:2021yzh . Overall,
LMZV/BLFQ0 provides significant improvement of LFH with IMA towards the
experimental data.

Figure 3: Comparison of unflavored light meson spectra and kaons as predicted
by LFH with IMA Brodsky:2014yha , with LMZV/BLFQ0 Li:2021jqb  and with ’t
Hooft potential as implemented by Ahmady et al. Ahmady:2021yzh . Results from
Gutsche et al. which adopt LFH with a longitudinal wave function similar to
LMZV are also attached for comparison Gutsche:2012ez .

In the non-relativistic (NR) limit, rotational symmetry requires that the
combination of the transverse and longitudinal confining interactions is a
rotational invariant. The transverse confinement reduces to a harmonic
oscillator potential in this limit
$\kappa^{4}\zeta_{\perp}^{2}\to({\kappa^{4}}/{4})\vec{r}_{\perp}^{2}$. It is
expected that the longitudinal confinement also reduces to the quadratic form.
The GT model satisfies rotational symmetry in the NR limit by construction.
The LMZV model also satisfies this symmetry with the identification
$\sigma=\kappa^{2}/(m_{q}+m_{\bar{q}})$ in the NR limit. It is curious to note
that heavy quark effective theory requires
$\kappa\propto\sqrt{M_{H}}\approx\sqrt{m_{q}+m_{\bar{q}}}$. Hence, the
confining strength $\sigma$ approaches to a constant in the heavy quark limit.
From the fits to the physical quark masses, $\sigma_{\textsc{hq}}\approx
0.24\,\text{GeV}$, where ”HQ” stands for ”Heavy Quark” DeTeramond:2021jnn .
Physically, this value is determined by the gluon condensate $\langle
G_{\mu\nu}G^{\mu\nu}\rangle$, which is in general different from the values in
the light sector $\sigma=\langle\bar{q}q\rangle/f_{\pi}^{2}\approx
0.6\,\text{GeV}$ determined by the quark condensate DeTeramond:2021jnn .

In coordinate space, the ’t Hooft interaction is (recall
$\tilde{z}=\frac{1}{2}P^{+}x^{-}$),

$V_{\text{tH}}=g^{2}\big{|}\tilde{z}\big{|}$ (29)

In the NR limit, $\tilde{z}\to-2Mr_{z}$. The rotational symmetry is not
restored with the quadratic transverse confinement from LFH. Note that in Ref.
Ahmady:2021yzh , Ahmady et al. argued that rotational symmetry can be restored
with the matching $g=\kappa$, as the transverse LFH confinement and the ’t
Hooft interaction are equivalent to the c.m. instant form potential
$V_{\perp}=(\kappa^{2}/2)r_{\perp}$ and $V_{\|}=(g^{2}/2)r_{z}$, respectively,
based on an relation
$U_{\text{LF}}=V_{\text{IF}}+2(m_{q}+m_{\bar{q}})V_{\text{IF}}$ proposed by
Trawiński et al. Trawinski:2014msa . Even if this relation is held in the NR
limit, the corresponding potential
$V=V_{\perp}+V_{\|}=(\kappa^{2}/2)(\sqrt{r_{x}^{2}+r_{y}^{2}}+r_{z})$ is still
not rotationally invariant. Shuryak and Zahed proposed an alternative
confining interaction that generalizes the ’t Hooft model to 3+1D while
combining the light-front holographic interaction Shuryak:2021hng ;
Shuryak:2021mlh .

$V_{\text{SZ}}=2\sigma_{T}\sqrt{\tilde{z}^{2}+M^{2}r^{2}_{\perp}}.$ (30)

The rotational symmetry is restored in the NR limit in this model.

Figure 4: Comparison of the longitudinal wave functions of ground states from
the IMA and from the LMZV model.

Figure 4 compares the ground-state longitudinal wave functions
$X(x)=\sqrt{x(1-x)}\chi(x)$ from the IMA and from the LMZV model. The
parameters of the IMA are adopted from LFHQCD Brodsky:2014yha ;
Nielsen:2018ytt . For $q\bar{q}$ and $s\bar{s}$, the parameters of the LMZV
model are adopted from Ref. Li:2021jqb . For $c\bar{c},b\bar{b}$, the
parameters are matched to the IMA parameters with the NR reduction,
$\sigma=\kappa^{2}/(m_{q}+m_{\bar{q}})$. As we can see, the main difference
between these two confining potentials is the endpoint asymptotics. Despite
the striking similarities of the longitudinal wave functions from IMA and the
LMZV model shown in Fig. 4, we would like to emphasize that the endpoint
behaviors are dramatically different and this fact has important physical
consequences. In the ’t Hooft and LMZV models, the endpoint asymptotics is a
direct consequence of chiral symmetry breaking and the exponent is related to
the chiral condensate. The endpoint behavior also has a dramatic impact on
hadronic observables in high energy collisions as hard kernels $T_{H}$ are
sensitive to the endpoint singularities Lepage:1980fj .

As an example, one may consider the pion form factor Li:2021jqb . In light-
front holography, the pion form factor is Brodsky:2007hb ; Brodsky:2014yha ,

$F_{\pi}(Q^{2})=\int_{0}^{1}\frac{{\mathrm{d}}x}{x(1-x)}\big{|}X(x)\big{|}^{2}\exp\Big{(}-\frac{1-x}{x}\frac{Q^{2}}{4\kappa^{2}}\Big{)}.$
(31)

The large $Q^{2}$ behavior of $F_{\pi}(Q^{2})\sim 1/Q^{1+2\delta}$ is directly
related to the endpoint asymptotics of the (normalized) distribution amplitude
$X\sim(1-x)^{\delta}$ Lepage:1980fj . In ’t Hooft and LMZV models,
$\delta=1/2+m_{q}/\sigma$. Hence, $Q^{2}F_{\pi}(Q^{2})\sim 1/\log Q^{2}$ in
the vicinity of the chiral limit. This is in contrast to the LFH prediction
$Q^{2}F_{\pi}(Q^{2})\to\text{const}$ Brodsky:2007hb  and to the prediction
from LFH with IMA, $Q^{2}F_{\pi}(Q^{2})\to e^{-cQ^{2}}$. These results are
shown in Fig. (5), where we adopt parameters from the original LFH
Brodsky:2007hb . Of course, a comprehensive investigation of the pion form
factor in holography requires the dressed holographic current Brodsky:2014yha
, which is beyond the scope of the present work. We only point out that the
confined current has similar dependence on the endpoint behavior with the bare
current at large $Q^{2}$ Brodsky:2011xx .

A similar example is the pion radiative transition form factor111Here we have
neglected the evolution of the pion distribution amplitude since the evolution
is extremely slow for large but finite $Q^{2}$, say $1\;\text{GeV}^{2}\lesssim
Q^{2}\lesssim 1000\;\text{GeV}^{2}$. Brodsky:2011xx ; Swarnkar:2015osa ,

$Q^{2}F_{\pi\gamma}(Q^{2})=\frac{4}{\sqrt{3}}\int_{0}^{1}{\mathrm{d}}x\frac{\phi_{\pi}(x,\tilde{Q})}{1-x}.\quad(\tilde{Q}=(1-x)Q)$
(32)

Its value at large $Q^{2}$ is predicted to approach different asymptotic
values with different asymptotics at the endpoints.

Figure 5: Large-$Q^{2}$ asymptotics of the pion form factor
$Q^{2}F_{\pi}(Q^{2})$ from LFH (solid), LFH with IMA (dotted), and LMZV/BLFQ0
(dashed), for the original parameters of LFH with
$\kappa=0.375\,\mathrm{GeV}$. The quark masses are
$m_{\\{u,d\\}}=50\,\mathrm{MeV}$ for LFH with IMA and
$m_{\\{u,d\\}}=12\,\mathrm{MeV}$ for LMZV/BLFQ0.

##  IV Sturm-Liouville Theory

After a comparison of various models, a natural question to ask is whether
there is freedom for further proposals. This task is not as easy as it seems.
Because the construction of a valid Hermitian operator in the functional space
may be jeopardized by the singularities in the light-front kinetic energy
term. Weller & Miller unified the previous longitudinal potentials with Miller
and Brodsky’s spatial coordinate $\tilde{z}$ Weller:2021wog . An alternative
approach is the momentum representation,

$(V_{\text{F}}\circ\chi)(x)=\int_{0}^{1}{\mathrm{d}}yK(x,y)\chi(y)$ (33)

The kernel $K(x,y)$ may be related to the coordinate space ($\tilde{z}$)
potential.

$K(x-x^{\prime})=\int{\mathrm{d}}\tilde{z}\,e^{ix\tilde{z}}V_{\|}(\tilde{z}).$
(34)

The ’t Hooft model $K(x-x^{\prime})\propto 1/(x-x^{\prime})^{2}$ corresponds
to the linear potential $V_{\text{tH}}=|\tilde{z}|$, as mentioned. In general,
the power-law potential $|\tilde{z}|^{p}$ leads to a kernel
$K(x-x^{\prime})=\frac{1}{|x-x^{\prime}|^{p+1}}$ ($p\notin 2\mathbb{N}$). The
kernel for the quadratic potential $\tilde{z}^{2}$ is the derivative of a
Dirac-$\delta$, which is more convenient to be expressed as differential
operators. Due to the low dimensionality, these constructions may be
associated with severe singularities.

Another useful tool is the Sturm-Liouville (SL) theory. The corresponding
interaction takes the form,

$V_{\text{SL}}=-\partial_{x}p(x)\partial_{x}+s(x).$ (35)

Clearly, the LMZV interaction belongs to this category. The longitudinal
Schrödinger equation becomes (9),

$\Big{[}\frac{m_{q}^{2}}{x}+\frac{m^{2}_{\bar{q}}}{1-x}+V_{\text{SL}}\Big{]}\chi(x)=M^{2}_{\|}\chi(x).$
(36)

It may be written as a Sturm-Liouville problem

$-\big{(}p(x)\chi^{\prime}(x)\big{)}^{\prime}+q(x)\chi(x)=\lambda\chi(x),\quad(\chi(0)=\chi(1)=0)$
(37)

where $q(x)=m_{q}^{2}/x+m^{2}_{\bar{q}}/(1-x)+s(x)$. Note that the derivative
$\partial_{x}$ is taken with respect to $\vec{\zeta}_{\perp}$ not
$\vec{r}_{\perp}$, viz $\partial_{x}=(\partial/\partial
x)_{\vec{\zeta}_{\perp}}$. In the NR limit, rotational symmetry can be
restored if $p(x=\hat{m}_{q})=\kappa^{4}/[4(m_{q}+m_{\bar{q}})^{2}]$.

We shall call this class of potentials as the Sturm-Liouville (SL) type, since
the longitudinal LFSWE forms a Sturm-Liouville eigenvalue problem (SLEVP). The
SL longitudinal confinement is naturally compatible with the transverse
holographic confinement in the NR limit. The advantage of these extensions is
that the known properties of the SLEVP, e.g. the distribution of eigenvalues,
can be applied to the longitudinal LFSWE.

To the first approximation, function $p(x)$, stemming from the inter-quark
potential, should be symmetric wrt the interchange $x\leftrightarrow 1-x$.
Therefore, it can be written as a function of $x(1-x)$ and
$(x-\frac{1}{2})^{2}$. Unless the quark masses are equal,
$(x-\frac{1}{2})^{2}$ is not consistent with the NR reduction
$x\approx\hat{m}_{q}$ and $x\approx\hat{m}_{\bar{q}}$. Therefore, the natural
choice of $p(x)$ is $p(x)=\sigma^{2}[x(1-x)]^{\gamma}$. A typical source of
$s(x)$ is the self-energy correction. A dynamically generated mass may have
the asymptotics $\sim B/x(1-x)$, which can be absorbed into the quark mass. In
general, we can assume that $s(x)$ is no more singular than the light-front
kinetic energy. A convenient choice is $s(x)=0$.

A comment on the normalization is in order. If, instead, the following
normalization is adopted,

$\int_{0}^{1}\frac{{\mathrm{d}}x}{x(1-x)}\big{|}X(x)\big{|}^{2}=1,$ (38)

where a non-trivial weight function $w(x)=1/[x(1-x)]$ is introduced, the
Sturm-Liouville problem for $X(x)$ should be,

$-\big{(}P(x)X^{\prime}(x)\big{)}^{\prime}+Q(x)X(x)=\lambda
w(x)X(x),\quad(X(0)=X(1)=0).$ (39)

where

$\displaystyle P(x)=p(x)w(x)=\frac{p(x)}{x(1-x)},$ (40) $\displaystyle
Q(x)=q(x)w(x)-\frac{pw^{\prime
2}}{4w}+\frac{1}{2}(pw^{\prime})^{\prime}=\frac{q(x)}{x(1-x)}-\frac{8x^{2}-8x+3}{4x^{3}(1-x)^{3}}p(x)+\frac{2x-1}{2x^{2}(1-x)^{2}}p^{\prime}(x).$
(41)

For example, the corresponding interaction for the LMZV model reads,

$U_{\|}=-\sigma^{2}\partial^{2}_{x}+\frac{\sigma^{2}}{4x^{2}(1-x)^{2}}$ (42)

Now let us turn to the general properties of the SLEVP, which is covered by
the Sturm-Liouville theorem Krall:1986 ; Zettl:2005 . The SL problem (37) is
called singular if either $p(x)=0$ at the endpoints or $q(x)$ has
singularities. Otherwise, the SL problem is called regular. Clearly, the
light-front kinetic energy term in $q(x)$ leads to a singular SLEVP. We remind
that the Sturm-Liouville theorem is only applicable to the regular SLEVP or
the singular SLEVP in limit circle non-oscillation (LCNO) Krall:1986 ;
Zettl:2005 . For other singular SLEVPs, including the limit circle oscillation
(LCO) case and the limit point (LP) cases, the properties of the eigenvalues
and eigenfunctions are not known in general. For example, the eigenvalues may
not be bounded from below; the eigen-solutions may have infinite oscillations
near the endpoints. Such pathological solutions are known in finite
truncations of relativistic bound-state equations, such as Bethe-Salpeter
equations. What is worse, the classification of singular SLEVP generally
relies on the solution. Mathematically, only a finite number of cases are
classified. A catalog of the singular SLEVP is given by Ref. Everitt:2000 .

In fact, the most famous SLEVP of the LCNO case on the finite interval
$[0,\,1]$, is the familiar Legendre function, corresponding to $p(x)=x(1-x)$,
and $q(x)=0$. Here, we have applied a shift $x=(t+1)/2$. The eigenvalues
$\lambda_{\ell}=\ell(\ell+1)$ and the eigen-functions are the Legendre
polynomials $P_{\ell}(2x-1)$. The longitudinal confinement
$V_{\|}=-\sigma^{2}\partial_{x}\big{(}x(1-x)\partial_{x}\big{)}$ covers a
number of well-known special functions, the Legendre function, the Chebyshev
function and the Gegenbauer function. The most general case is the Jacobi
function which is applicable to the general unequal quark mass case,

$\Big{[}\frac{m_{q}^{2}}{x}+\frac{m^{2}_{\bar{q}}}{1-x}-\sigma^{2}\frac{{\mathrm{d}}}{{\mathrm{d}}x}x(1-x)\frac{{\mathrm{d}}}{{\mathrm{d}}x}\Big{]}\chi(x)=M^{2}_{\|}\chi(x).$
(43)

The associated SL problem is,

$-\big{(}x(1-x)\chi^{\prime}\big{)}^{\prime}+\Big{(}\frac{\beta^{2}}{4x}+\frac{\alpha^{2}}{4(1-x)}\Big{)}\chi=\lambda\chi,\quad(\chi(0)=\chi(1)=0)$
(44)

where $\alpha=2m_{\bar{q}}/\sigma,\beta=2m_{q}/\sigma$. The solution is the
Jacobi polynomials:

$\chi_{n}(x)=N_{n}x^{\frac{\beta}{2}}(1-x)^{\frac{\alpha}{2}}P_{n}^{(\alpha,\beta)}(2x-1).$
(45)

The eigenvalues are $\lambda_{n}=(n+\alpha+\beta)(n+1+\alpha+\beta)\propto
n^{2}$. The mass eigenvalue
$M^{2}_{\|n}=(m_{q}+m_{\bar{q}})^{2}+\sigma(m_{q}+m_{\bar{q}})(2n+1)+\sigma^{2}n(n+1)$.
The ground state mass obeys the GMOR relation,
$M^{2}=(m_{q}+m_{\bar{q}})^{2}+\sigma(m_{q}+m_{\bar{q}})$. The ground state
wave function is power-law like.

The Jacobi function belongs to an even broader class of special functions, the
hypergeometric function AS ; DLMF . It is expected that a broader class of
longitudinal confining potentials may be constructed there with hypergeometric
functions as the solutions. Unfortunately, the classification of the general
hypergeometric differential equation on the interval $[-1,1]$ is not known.

A similar situation occurs for another large class of second order
differential equations on $[0,1]$, the Heun equation Everitt:2000 ;
Slavyanov:2000 . There are at least two special cases of the Heun equation
that are related to the LMZV interaction. The first one stems from the
azimuthal part of the three-body Schrödinger equation that describes the
hydrogen-molecule ion $H^{+}_{2}$ (two protons plus one electron)
Slavyanov:2000 . It corresponds to a longitudinal potential
$V_{\|}=-\sigma^{2}\big{[}\partial_{x}\big{(}x(1-x)\partial_{x}\big{)}+\eta(2x-1)^{2}\big{]}$.
The eigensolutions are shown to acquire additional exponential factors:
$\exp(\pm c\big{|}2x-1\big{|})$ Mitin:2015 . The second case is the Teukolsky
equation that describes the stability of the Kerr black hole. It corresponds
to a longitudinal potential
$V_{\|}=-\sigma^{2}\big{[}\partial_{x}\big{(}x(1-x)\partial_{x}\big{)}+a^{2}(2x-1)^{2}+4a(2x-1)\big{]}$
and a mass function $M=m+b(2x-1)$ Slavyanov:2000 . The mass correction can
also be absorbed into the potential.

The scaling of the eigenvalues is given by the Atkinson–Mingarelli theorem
Atkinson:1987 . For large $n$, $\lambda_{n}$ scales quadratically,

$\lambda_{n}\sim\frac{n^{2}\pi^{2}}{\Big{(}\int_{a}^{b}\frac{1}{\sqrt{p(x)}}\Big{)}^{2}}$
(46)

From this theorem, the eigenvalues of the SL confining potentials are
generally quadratic, unless the integral in the denominator diverges.

A natural generalization of the LMZV model is to consider
$p(x)=\sigma^{2}x^{\gamma}(1-x)^{\gamma}$. Except for $\gamma=1$, the obtained
SLEVPs are not well studied in mathematics. However, as we will see, at least
some of these problems (e.g., $\gamma=0$) are physically well defined. We will
investigate the longitudinal potentials with $\gamma=0,1,2,3$. The results are
summarized in Table 1. Note that for $\gamma\geq 2,\gamma\in\mathbb{Z}$, the
estimate Eq. (46) is not applicable and other distributions of the eigenvalues
are possible.

###  IV.1 Case I: $p(x)=\sigma^{2}$

In the first case, we take $p(x)=\sigma^{2},s(x)=0$. This is in fact one of
the simplest cases. The longitudinal confining potential
$V_{\|}=\sigma^{2}\tilde{z}^{2}$ is quadratic in the longitudinal coordinate
space. Then the longitudinal LFSWE becomes,

$\Big{[}\frac{m_{q}^{2}}{x}+\frac{m^{2}_{\bar{q}}}{1-x}-\sigma^{2}\frac{{\mathrm{d}}^{2}}{{\mathrm{d}}x^{2}}\Big{]}\chi(x)=M^{2}_{\|}\chi(x).$
(47)

The associated Sturm-Liouville problem, for the equal mass case
$m_{q}=m_{\bar{q}}$, is,

$-\chi^{\prime\prime}(x)+\frac{\mu^{2}}{x(1-x)}\chi(x)=\lambda\chi(x),\quad(\chi(0)=\chi(1)=0)$
(48)

where $\mu=m_{q}/\sigma$, and $\lambda=M^{2}/\sigma^{2}$.

This is a 1D Schrödinger equation with a potential $V\propto 1/(x(1-x))$. In
the chiral limit ($\mu=0$), the potential well becomes an infinite square
potential well. The solution in this case is the trigonometry function
$\chi_{n}(x)=\sqrt{2}\sin n\pi x$, and the $n$-th eigenvalue
$\lambda_{n}=n^{2}\pi^{2}$, ($n=1,2,\cdots$).222The massless solution
$\lambda=0$, $\chi(x)=a+bx$, is eliminated by the boundary condition. The
ground state energy is the zero-point energy, $M^{2}=\sigma^{2}\pi^{2}$.

For $m_{q}>0$, no analytic solution is known so far, although in special
cases, the solution is related to the (regularized) hypergeometric function
${}_{2}F_{1}$ AS ; DLMF . Numerical solutions are shown in Fig. 6. The $n$-th
eigenvalue, again, scales as $\lambda_{n}\propto n^{2}$. At the endpoints
$x\to 0$, the wave function $\chi(x)\sim x$. Since there is no singularity in
the integral of the kinetic energy, the quark mass contribution can be
estimated as, $\Delta M^{2}=\big{[}\gamma_{\textsc{e}}+\ln
2\pi-\mathrm{Ci}(2\pi)\big{]}m_{q}^{2}+O(m_{q}^{4})$, where,
$\gamma_{\textsc{e}}$ is Euler’s constant,
$\mathrm{Ci}(2\pi)\approx-0.0225607$ is the cosine integral. AS ; DLMF

Figure 6: Solutions of the SL problem (48).

###  IV.2 Case II: $p(x)=\sigma^{2}x(1-x)$

This case is the LMZV model discussed above.

###  IV.3 Case III: $p(x)=\sigma^{2}x^{2}(1-x)^{2}$

In this case, we choose $p(x)=\sigma^{2}x^{2}(1-x)^{2}$. Then the longitudinal
LFSWE becomes,

$\Big{[}\frac{m_{q}^{2}}{x}+\frac{m^{2}_{\bar{q}}}{1-x}-\sigma^{2}\frac{{\mathrm{d}}}{{\mathrm{d}}x}x^{2}(1-x)^{2}\frac{{\mathrm{d}}}{{\mathrm{d}}x}\Big{]}\chi(x)=M^{2}_{\|}\chi(x).$
(49)

The associated SL problem is,

$-\big{(}x^{2}(1-x)^{2}\chi^{\prime}\big{)}^{\prime}+\Big{(}\frac{\beta^{2}}{4x}+\frac{\alpha^{2}}{4(1-x)}\Big{)}\chi=\lambda\chi,\quad(\chi(0)=\chi(1)=0)$
(50)

The analytical solutions are not known. By matching to the singularities, it
can be shown that at the endpoints, the wave functions scales as
$\chi(x)\sim\exp\big{[}-\frac{2}{\sigma}\big{(}\frac{m_{q}}{\sqrt{x}}+\frac{m_{\bar{q}}}{\sqrt{1-x}}\big{)}\big{]}$.
Numerical solutions show the eigenvalues $\lambda_{n}\propto n$ for low $n$
and $\lambda_{n}\propto n^{2}$ for large $n$.

###  IV.4 Case IV: $p(x)=\sigma^{2}x^{3}(1-x)^{3}$

In this case, we choose $p(x)=\sigma^{2}x^{3}(1-x)^{3}$. Then the longitudinal
LFSWE becomes,

$\Big{[}\frac{m_{q}^{2}}{x}+\frac{m^{2}_{\bar{q}}}{1-x}-\sigma^{2}\frac{{\mathrm{d}}}{{\mathrm{d}}x}x^{3}(1-x)^{3}\frac{{\mathrm{d}}}{{\mathrm{d}}x}\Big{]}\chi(x)=M^{2}_{\|}\chi(x).$
(51)

The associated SL problem is,

$-\big{(}x^{3}(1-x)^{3}\chi^{\prime}\big{)}^{\prime}+\Big{(}\frac{\beta^{2}}{4x}+\frac{\alpha^{2}}{4(1-x)}\Big{)}\chi=\lambda\chi,\quad(\chi(0)=\chi(1)=0)$
(52)

The analytical solutions are not known. By matching to the singularities, it
can be shown that at the endpoints, the wave functions scale as
$\chi(x)\sim\exp\big{[}-\frac{1}{\sigma}\big{(}\frac{m_{q}}{x}+\frac{m_{\bar{q}}}{1-x}\big{)}\big{]}$.
Numerical calculations show that the eigenvalues $\lambda_{n}\propto n$.

This case is closely related to the Głazek-Trawiński potential introduced
above. To see this, we can express the longitudinal coordinate
$\zeta_{\|}=i\nabla_{\kappa_{\|}}$ of Głazek and Trawiński in terms of
$\tilde{z}=i\partial_{x}$. Recall, in the equal mass case,
$\kappa_{\|}=m_{q}\frac{2x-1}{\sqrt{x(1-x)}}$, or equivalently,
$x=(1/2)(1+\kappa_{\|}/\sqrt{4m_{q}^{2}+\kappa^{2}_{\|}})$. Then, we can
obtain $\zeta_{\|}=(2/m_{q})\sqrt{x^{3}(1-x)^{3}}i\partial_{x}$. Therefore,
the longitudinal confinement in this model becomes,

$V_{\|}^{(\text{GT})}=-\frac{4\lambda^{2}}{m_{q}^{2}}\partial_{x}x^{3}(1-x)^{3}\partial_{x}.$
(53)

The confining strength $\sigma=\frac{2\lambda}{m_{q}}$. With this
identification, the wave function at the endpoints becomes,
$\chi(x)\sim\exp\big{(}-\frac{m_{q}^{2}}{2\lambda x(1-x)}\big{)}$, in
agreement with Głazek and Trawiński’s wave function.

There is a caveat here since $\kappa_{\perp}=\vec{k}_{\perp}/\sqrt{x(1-x)}$
contains the longitudinal coordinate $x$. Therefore,
$\vec{k}_{\perp}=\vec{\kappa}_{\perp}/\sqrt{4+\kappa^{2}_{\|}/m_{q}^{2}}$ also
depending on $\kappa_{\|}$. To obtain the 3rd coordinate $\zeta_{\|}$, we can
apply the chain rule:

$\displaystyle\zeta_{\|}\equiv\,$ $\displaystyle
i\frac{\partial}{\partial\kappa_{\|}}$ (54) $\displaystyle=\,$
$\displaystyle\frac{\partial\vec{k}_{\perp}}{\partial\kappa_{\|}}\cdot\vec{r}_{\perp}+\frac{\partial
x}{\partial\kappa_{\|}}P^{+}r_{+}$ (55) $\displaystyle=\,$
$\displaystyle-{\frac{1}{2}}(x-{\frac{1}{2}})x^{{\frac{1}{2}}}(1-x)^{{\frac{1}{2}}}m_{q}^{-1}\vec{k}_{\perp}\cdot\vec{r}_{\perp}+2x^{{\frac{3}{2}}}(1-x)^{{\frac{3}{2}}}m_{q}^{-1}\tilde{z}.$
(56)

Here $\vec{r}_{\perp}=i\partial/\partial{\vec{k}_{\perp}}$,
$r_{+}=i\partial/\partial k^{+}$, and
$P^{+}r_{+}\equiv\tilde{z}=i\partial_{x}$ is Brodsky and Miller’s frame-
independent longitudinal coordinate. Since
$\partial\vec{k}_{\perp}/\partial\kappa_{\|}\neq 0$, in general $\zeta_{\|}$
depends on $\vec{r}_{\perp}$. In other words, $\vec{\kappa}$ and $\vec{\zeta}$
are not independent conjugate pairs, and the longitudinal interaction
constructed here may not be exactly the same as the Głazek-Trawiński model.

Table 1: Summary of several longitudinal interactions discussed in this work. interaction |  $\tilde{z}$-space | $M_{\text{gs}}^{2}$ | $M^{2}_{n}$ | $\chi(x)$ |
---|---|---|---|---|---
$-\sigma^{2}\partial_{x}^{2}$ | $\sigma^{2}\tilde{z}^{2}$ |  | $\sim n^{2}$ | $\sim x(1-x)$ |
$-\sigma^{2}\partial_{x}x(1-x)\partial_{x}$ |  | $(m_{q}+m_{\bar{q}})^{2}+\sigma(m_{q}+m_{\bar{q}})$ | $\sim n^{2}$ | $\sim x^{\frac{m_{q}}{\sigma}}(1-x)^{\frac{m_{\bar{q}}}{\sigma}}$ |
$-\sigma^{2}\partial_{x}x^{2}(1-x)^{2}\partial_{x}$ |  |  | $\sim n^{2}$ | $\sim\exp\big{(}-\frac{2}{\sigma}\sqrt{\frac{m_{q}^{2}}{x}+\frac{m_{\bar{q}}^{2}}{1-x}}\big{)}$ |
$-\sigma^{2}\partial_{x}x^{3}(1-x)^{3}\partial_{x}$ |  |  | $\sim n$ | $\sim\exp\big{(}-\frac{1}{\sigma}\big{(}{\frac{m_{q}}{x}+\frac{m_{\bar{q}}}{1-x}}\big{)}\big{)}$ |
$-g^{2}\frac{1}{(x-x^{\prime})^{2}}$ | $g^{2}|\tilde{z}|$ | $\sqrt{\frac{\pi}{3}}g(m_{q}+m_{\bar{q}})$ | $\sim n$ | $\sim x^{\beta}(1-x)^{\beta^{\prime}}$ |

##  V Summary and outlooks

In this work, we discussed the role of longitudinal dynamics in semi-classical
light-front Schrödinger wave equations. Based on a separation of variables
ansatz, we relate the longitudinal confining interaction to the endpoint
asymptotics of the resulting wave function as well as to the scaling of the
eigenvalues. This enables us to compare longitudinal interactions proposed in
the recent literature as well as a broader class of longitudinal interactions
based on the Sturm-Liouville theory. We analyzed several concrete cases
$V_{\|}\sim-\partial_{x}x^{\gamma}(1-x)^{\gamma}\partial_{x}$
($\gamma=0,1,2,3$). We also pointed out the challenges for constructing valid
longitudinal interactions due to the presence of the endpoint singularities in
the light front kinetic energy.

Our motivation is to search for a plausible first approximation to mesons in
conjunction with light-front holography. Both the ’t Hooft interaction
$\propto\big{|}\tilde{z}\big{|}$ and the LMZV interaction
$\propto\partial_{x}\big{(}x(1-x)\big{)}\partial_{x}$ implement the chiral
symmetry breaking for the pion. They both produce power-law-like wave
functions at the endpoints. On the other hand, interactions of quadratic form
naturally match the transverse holographic confining potential in the
nonrelativistic limit. We further expect that the resulting wave functions can
be used as a basis for the full quantum many-body calculations. In this
regard, the LMZV interaction is useful for its computational convenience.
Indeed, the longitudinal wave functions from this interaction, the power-law
weighted Jacobi polynomials, form a standard basis for wave functions with
power-law behavior, such as the solution of the ’t Hooft model. Other
longitudinal interactions may be useful in different contexts.

We focused on the light mesons. The reason is that the hyperfine structure of
the heavy systems, including the heavy-light mesons are dominated by the one-
gluon-exchange interaction Li:2015zda ; Godfrey:1985xj . A single confining
interaction is not sufficient. Indeed, the works on heavy systems have to
introduce additional parameters. Moreover, the restoration of the 3D
rotational symmetry in the non-relativistic limit puts strong constraints on
the confining interactions.333Here we only retain those longitudinal confining
forms that reduce to the harmonic oscillator in the non-relativistic limit due
to another symmetry constraint: we require factorization of the CM motion from
the internal motion. As such, different forms of the longitudinal confinement
become nearly identical.

One of the interesting questions is how to generalize the present
investigation to the multi-parton sectors, including notably the baryons.
There have been several proposals in the literature, including the pairwise
longitudinal interactions Mondal:2019jdg ; Xu:2021wwj , the quark-diquark
model Yu:2019 ; Ahmady:2021yzh , and the orthogonal basis over a $k$-simplex
(triangle, tetrahedron, …) Chabysheva:2013oka ; Chabysheva:2014rra ;
Burkardt:2016ffk ; Shuryak:2022thi .

## Acknowledgements

We wish to thank S.J. Brodsky and G.F. de Téramond for enlightening
discussions and for their critical reading of an earlier version of the
manuscript. We also gratefully acknowledge fruitful discussions with S.J.
Brodsky, G.F. de Téramond, G.A. Miller, X. Zhao, C. Mondal, J. Hiller, R.
Sandapen and Q. Wang. This work is supported in part by the Department of
Energy under Grants No. DE-FG02-87ER40371, and No. DE-SC0018223
(SciDAC4/NUCLEI). Y.L. is supported by the New faculty start-up fund of the
University of Science and Technology of China.

## References

  * (1) G. F. de Teramond and S. J. Brodsky, Phys. Rev. Lett. 102, 081601 (2009) [arXiv:0809.4899 [hep-ph]].
  * (2) S. J. Brodsky and G. F. de Teramond, Phys. Rev. Lett. 96, 201601 (2006) [arXiv:hep-ph/0602252 [hep-ph]].
  * (3) H. G. Dosch, G. F. de Teramond and S. J. Brodsky, Phys. Rev. D 92, no.7, 074010 (2015); [arXiv:1504.05112 [hep-ph]].
  * (4) H. G. Dosch, G. F. de Teramond and S. J. Brodsky, Phys. Rev. D 95, no.3, 034016 (2017); [arXiv:1612.02370 [hep-ph]].
  * (5) M. Nielsen, S. J. Brodsky, G. F. de Téramond, H. G. Dosch, F. S. Navarra and L. Zou, Phys. Rev. D 98, no.3, 034002 (2018); [arXiv:1805.11567 [hep-ph]].
  * (6) L. Zou, H. G. Dosch, G. F. De Téramond and S. J. Brodsky, Phys. Rev. D 99, no.11, 114024 (2019); [arXiv:1901.11205 [hep-ph]].
  * (7) A. Deur, S. J. Brodsky and G. F. de Teramond, Nucl. Phys. 90, 1 (2016); [arXiv:1604.08082 [hep-ph]].
  * (8) S. J. Brodsky and G. F. de Teramond, Phys. Rev. D 77, 056007 (2008) [arXiv:0707.3859 [hep-ph]].
  * (9) S. J. Brodsky and G. F. de Teramond, Phys. Rev. D 78, 025032 (2008) [arXiv:0804.0452 [hep-ph]].
  * (10) R. S. Sufian, G. F. de Téramond, S. J. Brodsky, A. Deur and H. G. Dosch, Phys. Rev. D 95, no.1, 014011 (2017); [arXiv:1609.06688 [hep-ph]].
  * (11) G. F. de Teramond et al. [HLFHS], Phys. Rev. Lett. 120, no.18, 182001 (2018); [arXiv:1801.09154 [hep-ph]].
  * (12) T. Liu, R. S. Sufian, G. F. de Téramond, H. G. Dosch, S. J. Brodsky and A. Deur, Phys. Rev. Lett. 124, no.8, 082003 (2020); [arXiv:1909.13818 [hep-ph]].
  * (13) S. J. Brodsky, G. F. de Teramond, H. G. Dosch and J. Erlich, Phys. Rept. 584, 1 (2015); [arXiv:1407.8131 [hep-ph]].
  * (14) Y. Li, P. Maris and J. P. Vary, “Chiral sum rule on the light front,” [arXiv:2203.14447 [hep-th]].
  * (15) S. S. Chabysheva and J. R. Hiller, Annals Phys. 337, 143-152 (2013); [arXiv:1207.7128 [hep-ph]].
  * (16) S. D. Glazek and A. P. Trawiński, Phys. Rev. D 88, no.10, 105025 (2013) [arXiv:1307.2059 [hep-ph]].
  * (17) A. P. Trawiński, S. D. Głazek, S. J. Brodsky, G. F. de Téramond and H. G. Dosch, Phys. Rev. D 90, no.7, 074017 (2014) [arXiv:1403.5651 [hep-ph]].
  * (18) Y. Li, P. Maris, X. Zhao and J. P. Vary, Phys. Lett. B 758, 118 (2016); [arXiv:1509.07212 [hep-ph]].
  * (19) Y. Li and J. P. Vary, Phys. Lett. B 825, 136860 (2022) [arXiv:2103.09993 [hep-ph]].
  * (20) G. F. De Téramond and S. J. Brodsky, [arXiv:2103.10950 [hep-ph]].
  * (21) M. Ahmady, H. Dahiya, S. Kaur, C. Mondal, R. Sandapen and N. Sharma, [arXiv:2105.01018 [hep-ph]].
  * (22) M. Ahmady, S. Kaur, S. L. MacKay, C. Mondal and R. Sandapen, Phys. Rev. D 104, no.7, 074013 (2021) [arXiv:2108.03482 [hep-ph]].
  * (23) E. Shuryak and I. Zahed, [arXiv:2111.01775 [hep-ph]].
  * (24) E. Shuryak and I. Zahed, [arXiv:2112.15586 [hep-ph]].
  * (25) E. Shuryak and I. Zahed, [arXiv:2202.00167 [hep-ph]].
  * (26) T. Gutsche, V. E. Lyubovitskij, I. Schmidt and A. Vega, Phys. Rev. D 87, no. 5, 056001 (2013); [arXiv:1212.5196 [hep-ph]].
  * (27) C. M. Weller and G. A. Miller, [arXiv:2111.03194 [hep-ph]].
  * (28) J. P. Vary, H. Honkanen, J. Li, P. Maris, S. J. Brodsky, A. Harindranath, G. F. de Teramond, P. Sternberg, E. G. Ng and C. Yang, Phys. Rev. C 81, 035205 (2010); [arXiv:0905.1411 [nucl-th]].
  * (29) S. Godfrey and N. Isgur, Phys. Rev. D 32, 189-231 (1985)
  * (30) A. Szczepaniak, E. S. Swanson, C. R. Ji and S. R. Cotanch, Phys. Rev. Lett. 76, 2011-2014 (1996) [arXiv:hep-ph/9511422 [hep-ph]].
  * (31) H. Reinhardt, G. Burgio, D. Campagnari, E. Ebadati, J. Heffner, M. Quandt, P. Vastag and H. Vogt, Adv. High Energy Phys. 2018, 2312498 (2018) [arXiv:1706.02702 [hep-th]].
  * (32) S. Okubo, Progr. Theor. Phys. 12 (1954) 603.
  * (33) K. Suzuki, S.Y. Lee, Progr. Theor. Phys. 64 (1980) 2091.
  * (34) K. Suzuki, Progr. Theor. Phys. 68 (1982) 246.
  * (35) K. Suzuki, Progr. Theor. Phys. 68 (1982) 1999.
  * (36) K. Suzuki, R. Okamoto, Progr. Theor. Phys. 70 (1983) 439.
  * (37) K. Suzuki, R. Okamoto, Progr. Theor. Phys. 92 (1994) 1045.
  * (38) S. D. Glazek and K. G. Wilson, Phys. Rev. D 48, 5863-5872 (1993)
  * (39) S. D. Glazek and K. G. Wilson, Phys. Rev. D 49, 4214-4218 (1994)
  * (40) F. Wegner, Ann. Phys. 506, 77 (1994).
  * (41) G. P. Lepage and S. J. Brodsky, Phys. Rev. D 22, 2157 (1980)
  * (42) R. J. Perry, A. Harindranath and K. G. Wilson, Phys. Rev. Lett. 65, 2959-2962 (1990)
  * (43) K. G. Wilson, T. S. Walhout, A. Harindranath, W. M. Zhang, R. J. Perry and S. D. Glazek, Phys. Rev. D 49, 6720-6766 (1994) [arXiv:hep-th/9401153 [hep-th]].
  * (44) S. D. Glazek, Acta Phys. Polon. B 43, 1843-1862 (2012); [arXiv:1204.4760 [hep-th]].
  * (45) S. D. Głazek, M. Gómez-Rocha, J. More and K. Serafin, Phys. Lett. B 773, 172-178 (2017); [arXiv:1705.07629 [hep-ph]].
  * (46) K. Serafin, M. Gómez-Rocha, J. More and S. D. Głazek, Eur. Phys. J. C 78, no.11, 964 (2018); [arXiv:1805.03436 [hep-ph]].
  * (47) S. J. Brodsky, S. Gardner and D. S. Hwang, Phys. Rev. D 73, 036007 (2006) [arXiv:hep-ph/0601037 [hep-ph]].
  * (48) S. J. Brodsky, H. C. Pauli and S. S. Pinsky, Phys. Rept. 301, 299-486 (1998) [arXiv:hep-ph/9705477 [hep-ph]].
  * (49) J. S. Schwinger, Phys. Rev. 128, 2425-2429 (1962)
  * (50) G. ’t Hooft, Nucl. Phys. B 75, 461-470 (1974)
  * (51) Y. z. Mo and R. J. Perry, J. Comput. Phys. 108, 159-174 (1993).
  * (52) J. Eisert, M. Cramer, M.B. Plenio, Rev. Mod. Phys. 82, 277 (2010).
  * (53) Y. Li, P. Maris and J. P. Vary, Phys. Rev. D 96, no.1, 016022 (2017); [arXiv:1704.06968 [hep-ph]].
  * (54) G. A. Miller and S. J. Brodsky, Phys. Rev. C 102, no.2, 022201 (2020); [arXiv:1912.08911 [hep-ph]].
  * (55) S. J. Brodsky and G. F. de Teramond, Subnucl. Ser. 45, 139 (2009); [arXiv:0802.0514 [hep-ph]].
  * (56) M. Gell-Mann, R. J. Oakes and B. Renner, Phys. Rev. 175, 2195-2199 (1968)
  * (57) S. J. Brodsky, F. G. Cao and G. F. de Teramond, Phys. Rev. D 84, 075012 (2011) [arXiv:1105.3999 [hep-ph]].
  * (58) R. Swarnkar and D. Chakrabarti, Phys. Rev. D 92, no.7, 074023 (2015) [arXiv:1507.01568 [hep-ph]].
  * (59) X. Ji, Y. Liu and I. Zahed, Phys. Rev. D 103, no.7, 074002 (2021) [arXiv:2010.06665 [hep-ph]].
  * (60) A.M. Krall, “Applied analysis”, D. Reidel Publishing Co. 1986.
  * (61) A. Zettl, “Sturm-Liouville theory”, AMS, 2005.
  * (62) W. N. Everitt, “A catalogue of Sturm-Liouville differential equations. Sturm-Liouville theory: past and present.” Eds. W. O. Amrein, A. M. Hinz, and D. B. Pearson. Birkhäuser Verlag. ISBN: 3-7643-7066-1.
  * (63) Abramowitz, M.; Stegun, I. A., eds. (1983), “Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables”. 9th reprint with additional corrections of tenth original printing with corrections. New York. ISBN 978-0-486-61272-0.
  * (64) NIST Digital Library of Mathematical Functions. http://dlmf.nist.gov/, Release 1.1.4 of 2021-01-15. F. W. J. Olver, A. B. Olde Daalhuis, D. W. Lozier, B. I. Schneider, R. F. Boisvert, C. W. Clark, B. R. Miller, B. V. Saunders, H. S. Cohl, and M. A. McClain, eds.
  * (65) S. Yu. Slavyanov, W. Lay, “Special functions: A unified theory based on singularities.” Oxford University Press (2000).
  * (66) A.V. Mitin, arXiv:1508.01359v2 [quant-ph]
  * (67) F. V. Atkinson, A. B. Mingarelli, J. für die Reine und Ang. Math. (Crelle), 375/376 (1987), 380–393.
  * (68) C. Mondal, S. Xu, J. Lan, X. Zhao, Y. Li, D. Chakrabarti and J. P. Vary, Phys. Rev. D 102, no.1, 016008 (2020) [arXiv:1911.10913 [hep-ph]].
  * (69) S. Xu et al. [BLFQ], Phys. Rev. D 104, no.9, 094036 (2021) [arXiv:2108.03909 [hep-ph]].
  * (70) A. Yu, https://dr.lib.iastate.edu/handle/20.500.12876/31997
  * (71) S. S. Chabysheva, B. Elliott and J. R. Hiller, Phys. Rev. E 88, no.6, 063307 (2013) [arXiv:1307.4964 [physics.comp-ph]].
  * (72) S. S. Chabysheva and J. R. Hiller, Phys. Rev. E 90, no.6, 063310 (2014) [arXiv:1409.6333 [hep-ph]].
  * (73) M. Burkardt, S. S. Chabysheva and J. R. Hiller, Phys. Rev. D 94, no.6, 065006 (2016) [arXiv:1607.00026 [hep-th]].


C\sup_{b\in\mathcal{B}}\int_{\mathbb{R}^{d+p}\setminus\mathcal{B}_{0,\kappa}^{d+p}}(1+\|(s,e)\|)^{2}\,d(P^{S}\otimes
P^{\varepsilon})(s,e)$ $\displaystyle\leq\frac{C^{\prime}}{n},$

for some constant $C^{\prime}>0$ depending on
$d,\alpha,\beta,\sigma_{1},\sigma_{2},\gamma_{1},\gamma_{2}$, where we used
the fact that $P^{(S,\varepsilon)}\sim(\rho,\alpha\wedge\beta)\textup{-sub-
Weibull}$ and Lemma 26 in the final inequality. It therefore remains to
control

$G:=\sup_{b\in\mathcal{B}}\Bigl{|}\int_{E_{b}}\psi_{b}^{*}(v)\,d(P^{A(b)}-P_{n}^{A(b)})(v)\Bigr{|}=\sup_{b\in\mathcal{B}}\Bigl{|}\int_{\mathcal{B}_{0,\kappa}^{d+p}}\psi_{b}^{*}\circ
T_{b}(s,e)\,d(P^{S}\otimes P^{\varepsilon}-P^{S}_{n}\otimes
P^{\varepsilon}_{n})(s,e)\Bigr{|}.$

To simplify the notation, define the centered function

$\bar{\psi}_{b}^{*}(s,e):={\psi}_{b}^{*}\circ
T_{b}(s,e)\operatorname{\mathds{1}}\\{\|(s,e)\|\leq\kappa\\}-\operatorname{\mathbb{E}}[{\psi}_{b}^{*}\circ
T_{b}(S,\varepsilon)\operatorname{\mathds{1}}\\{\|(S,\varepsilon)\|\leq\kappa\\}],$

then it follows that $\|\bar{\psi}_{b}^{*}\|_{\infty}\leq 2C(\kappa+1)^{2}\leq
C(\log m)^{\frac{1}{2\wedge\alpha\wedge\beta}}$. In this notation, we have
$G=\sup_{b\in\mathcal{B}}|n^{-1}\sum_{i\in[n]}\bar{\psi}_{b}^{*}(S_{i},\varepsilon_{i})|$.
By Markov’s inequality, we then have

$\displaystyle\operatorname{\mathbb{E}}(G)$
$\displaystyle=\int_{0}^{+\infty}\operatorname{\mathbb{P}}(G\geq t)\,dt\leq
n^{-1/2}+C\int_{n^{-1/2}}^{+\infty}\frac{(\log
m)^{\frac{2}{2\wedge\alpha\wedge\beta}}}{nt^{2}}dt\lesssim\frac{(\log
m)^{\frac{2}{2\wedge\alpha\wedge\beta}}}{\sqrt{n}}.$ (70)

We now claim that $G$, when viewed as a function of
$(s_{1},e_{1}),\ldots,(s_{n},e_{n})$, satisfies the bounded difference
property [73, see e.g.]. By symmetry, it suffices to consider a perturbation
on $(s_{1},e_{1})$. Define
$v=(v_{i})_{i=1}^{n},v^{\prime}=(v_{i}^{\prime})_{i=1}^{n}$ where each
$v_{i}=(s_{i},e_{i}),v_{i}^{\prime}=(s_{i}^{\prime},e_{i}^{\prime})\in\mathbb{R}^{d+p}$,
such that $v_{i}=v_{i}^{\prime}$ for any $i\not=1$. We have

$\displaystyle\sup_{b\in\mathcal{B}}\biggl{|}\frac{1}{n}\sum_{i=1}^{n}\bar{\psi}_{b}^{*}(v_{i})\biggr{|}-\sup_{b\in\mathcal{B}}\biggl{|}\frac{1}{n}\sum_{i=1}^{n}\bar{\psi}_{b}^{*}(v_{i}^{\prime})\biggr{|}$
$\displaystyle\leq\sup_{b\in\mathcal{B}}\biggl{\\{}\biggl{|}\frac{1}{n}\sum_{i=1}^{n}\bar{\psi}_{b}^{*}(v_{i})\biggr{|}-\biggl{|}\frac{1}{n}\sum_{i=1}^{n}\bar{\psi}_{b}^{*}(v_{i}^{\prime})\biggr{|}\biggr{\\}}$
$\displaystyle\leq\frac{1}{n}\sup_{b\in\mathcal{B}}\Bigl{|}\bar{\psi}_{b}^{*}(v_{1})-\bar{\psi}_{b}^{*}(v_{1}^{\prime})\Bigr{|}\leq\frac{2C(\log
m)^{\frac{1}{2\wedge\alpha\wedge\beta}}}{n},$

estbalishing, the bounded difference property for $G$. Thus by McDiarmid’s
inequality [73, see e.g.], we obtain that the event

$\Lambda_{1}:=\Bigl{\\{}G\leq\mathbb{E}G+\frac{\sqrt{2}C(\log
m)^{\frac{2}{2\wedge\alpha\wedge\beta}}}{\sqrt{n}}\Bigr{\\}},$

occurs with probability at least $1-1/m$.

Thus, working on the event $\Omega_{1}\cap\Lambda_{1}$, we deduce from (70)
that

$\displaystyle\sup_{b\in\mathcal{B}}\Bigl{|}\int\psi_{b}^{*}(v)d(P^{A(b)}-P_{n}^{A(b)})(v)\Bigr{|}\leq\frac{C(\log
m)^{\frac{2}{2\wedge\alpha\wedge\beta}}}{\sqrt{n}}+\frac{C^{\prime}}{n}\leq\frac{C(\log
m)^{\frac{2}{2\wedge\alpha\wedge\beta}}}{\sqrt{n}},$

for some constant $C>1$ depends on
$d,\alpha,\beta,\gamma_{1},\gamma_{2},\sigma_{1},\sigma_{2}$, which completes
the first claim of the lemma.

For the second claim, in order to bound
$\int\psi_{b}(u)\,d(P^{U}-P_{m}^{U})(u)$, we notice that the anti-
concentration property of $P^{U}$ holds due to the Gaussian assumption. Thus
[68, Theorem 11 ] implies that
$\|\nabla\psi_{b}(u)\|\leq\tilde{C}(1+\|u\|)^{\frac{2}{\alpha\wedge\beta}}$
for some $\tilde{C}>0$ depending on $d,\sigma_{1},\sigma_{2},\alpha,\beta$,
and it follows that
$|\psi_{b}(u)|\leq\tilde{C}(1+\|u\|)^{\frac{2}{\alpha\wedge\beta}+1}$.

Define $\Omega_{2}:=\\{\max_{1\leq i\leq m}\|U_{i}\|\leq\gamma\\}$. From the
proof of Lemma 18 again, we have $\mathbb{P}(\Omega_{2})\geq 1-2/n$. Working
on $\Omega_{2}$, we have

$\displaystyle\sup_{b\in\mathcal{B}}\biggl{|}\int_{\mathbb{R}^{d}\setminus\mathcal{B}_{0,\gamma}}\psi_{b}(u)d(P^{U}-P_{m}^{U})(u)\biggr{|}$
$\displaystyle=\sup_{b\in\mathcal{B}}\biggl{|}\int_{\mathbb{R}^{d}\setminus\mathcal{B}_{0,\gamma}}\psi_{b}(u)dP^{U}(u)\biggr{|}$
$\displaystyle\leq\tilde{C}\int_{\mathbb{R}^{d}\setminus\mathcal{B}_{0,\gamma}}(1+\|u\|)^{\frac{2}{\alpha\wedge\beta}+1}dP^{U}(u)\leq\frac{\tilde{C}}{m},$

for some constant $\tilde{C}>0$ depending on
$d,\alpha,\beta,\sigma_{1},\sigma_{2}$. Now, defining
$\tilde{G}:=\sup_{b\in\mathcal{B}}\bigl{|}\int_{\mathcal{B}_{0,\gamma}}\psi_{b}d(P^{U}-P_{m}^{U})\bigr{|}$,
by the same argument as in the proof of the first part of this lemma, there is
an event $\Lambda_{2}$ with probability at least $1-m^{-1}$, such that on
$\Omega_{2}\cap\Lambda_{2}$, we have

$\displaystyle\sup_{b\in\mathcal{B}}\Bigl{|}\int\psi_{b}(u)d(P^{U}-P_{m}^{U})(u)\Bigr{|}\leq\tilde{G}+\frac{\tilde{C}^{\prime}}{m}\leq\mathbb{E}\tilde{G}+\frac{\bar{C}(\log
m)^{\frac{2}{2\wedge\alpha\wedge\beta}}}{\sqrt{m}}+\frac{\tilde{C}}{m}\leq\frac{\bar{C}(\log
m)^{\frac{2}{2\wedge\alpha\wedge\beta}}}{\sqrt{m}},$ (71)

for $\bar{C}>0$ depending only on
$d,\alpha,\beta,\sigma_{1},\sigma_{2},\gamma_{2}$. ∎

######  Lemma 21.

There exists $C>0$ depending only on
$d,\alpha,\beta,\sigma_{1},\sigma_{2},\gamma_{1},\gamma_{2}$, such that with
probability at least $1-5/n$, we have

$\displaystyle\sup_{b\in\mathcal{B}}\biggl{|}\int\|v\|^{2}\,d\bigl{(}P_{n}^{A(b)}-P^{A(b)}\bigr{)}(v)\biggr{|}\leq
C(\log m)^{\frac{2}{2\wedge\alpha\wedge\beta}}n^{-1/2},$
$\displaystyle\biggl{|}\int\|u\|^{2}\,d\bigl{(}P^{U}-P_{m}^{U}\bigr{)}(u)\biggr{|}\leq
C\sqrt{\frac{\log m}{m}}.$

###### Proof.

Observe that the only property of $\psi_{b}^{*}$ that we used in the first
part of the proof of Lemma 20 is that $\|\nabla\psi_{b}^{*}(v)\|\leq
C(\|v\|+1)$ for all $b\in\mathcal{B}$ and $v\in\mathbb{R}^{d}$. The same
property is satisfied by the function $v\mapsto\|v\|^{2}$. Hence, a very
similar proof to that of Lemma 20 will establish the first claim here.

As for the second inequality, since $U_{i}\stackrel{{\scriptstyle
i.i.d.}}{{\sim}}\mathcal{N}(0,I_{d})$, we have
$\sum_{i=1}^{m}\|U_{i}\|^{2}\sim\chi^{2}_{md}$. By [67, Lemma 1 ] we deduce
that

$\displaystyle\operatorname{\mathbb{P}}\biggl{(}\Bigl{|}\frac{1}{m}\sum_{i=1}^{m}\|U_{i}\|^{2}-\operatorname{\mathbb{E}}\|U\|^{2}\Bigr{|}\geq\sqrt{\frac{2d\log
m}{m}}+\frac{2\log m}{m}\biggr{)}\leq\frac{2}{m},$

which implies the second claim. ∎

###### Proof.

of Theorem 8 Recalling that in the regime of (6) event $\Theta$ holds with
probability at least $1-4(\log n)^{-1}$, and working on $\Theta$ we have
$\hat{b}\in\mathcal{B}$. Thus there exists $M>0$ depending only on $d$,
$\alpha$, $\beta$, $\sigma_{1}$, $\sigma_{2}$, $\gamma_{1}$, $\gamma_{2}$ such
that with probability at least $1-33(\log n)^{-1}$, we have

$\displaystyle\mathcal{L}(\hat{b})-\mathcal{L}(b^{*})$ $\displaystyle\leq
2\sup_{b\in\mathcal{B}}|\mathcal{L}(b)-\mathcal{L}_{n,m}(b)|$
$\displaystyle\leq\Bigl{|}\frac{1}{m}\sum_{i=1}^{m}\|U_{i}\|^{2}-\operatorname{\mathbb{E}}\|U\|^{2}\Bigr{|}+\sup_{b\in\mathcal{B}}\Bigl{|}\frac{1}{n}\sum_{i=1}^{n}\|T_{b}(S_{i},\varepsilon_{i})\|^{2}-\operatorname{\mathbb{E}}\|T_{b}(S_{i},\varepsilon_{i})\|^{2}\Bigr{|}$
$\displaystyle\hskip
128.0374pt+\sup_{b\in\mathcal{B}}\Bigl{|}\mathcal{W}_{2}^{2}(P^{A(b)},P^{U})-\mathcal{W}_{2}^{2}(P_{n}^{A(b)},P_{m}^{U})\Bigr{|}$
$\displaystyle\leq M(\log
m)^{\frac{8}{2\wedge\alpha\wedge\beta}}\biggl{(}\sqrt{\frac{p}{n}}+\frac{1}{n^{2/d}}\biggr{)},$
(72)

where the second inequality uses the definition of
$\langle\\!\langle\cdot,\cdot\rangle\\!\rangle_{\mathcal{W}_{2}}$ and in the
final inequality, we used Lemma 21 to control the first two terms and
Proposition 17 for the last term.

On the other hand, by the lower bound developed in (31) and Lemma 22 we have
for $r:=\langle\\!\langle
P^{\varepsilon},P^{U}\rangle\\!\rangle_{\mathcal{W}_{2}}$ that

$\displaystyle\mathcal{L}(\hat{b})-\mathcal{L}(b^{*})\geq\sqrt{r^{2}+\|b^{*}-\hat{b}\|_{\Sigma}^{2}}-r\geq\frac{1}{2}(1+r^{2})^{-1/2}\|b^{*}-\hat{b}\|_{\Sigma}^{2}.$
(73)

Combining (72) with (73), we obtain that

$\displaystyle\|b^{*}-\hat{b}\|_{\Sigma}\leq M(\log
m)^{\frac{4}{2\wedge\alpha\wedge\beta}}\biggl{\\{}\biggl{(}\frac{p}{n}\biggr{)}^{1/4}+\frac{1}{n^{1/d}}\biggr{\\}},$
(74)

with probability at least $1-33(\log n)^{-1}$. Here we close the proof. ∎

##  Appendix B Auxillary results

######  Lemma 22.

For any $a\geq 0$, we have inequality

$\sqrt{a+x^{2}}\leq\begin{cases}\frac{x^{2}}{2\sqrt{a}}+\sqrt{a}&,\text{if
$0\leq x\leq 1$},\\\ (x-1)+\frac{1}{2\sqrt{a}}+\sqrt{a}&,\text{if
$x>1$}.\end{cases},$

and

$\sqrt{a+x^{2}}\geq\begin{cases}\frac{x^{2}}{2\sqrt{a+1}}+\sqrt{a}&,\text{if
$0\leq x\leq 1$},\\\
\frac{x-1}{\sqrt{a+1}}+\frac{1}{2\sqrt{a+1}}+\sqrt{a}&,\text{if
$x>1$}.\end{cases}$

###### Proof.

Write

$\displaystyle\sqrt{a+x^{2}}=\int_{0}^{x}\frac{t}{\sqrt{a+t^{2}}}dt+\sqrt{a}.$

Thus the first inequality can be obtained by utilizing $t/\sqrt{a+t^{2}}\leq
t/\sqrt{a}$ and $t/\sqrt{a+t^{2}}\leq 1$ in the case of $0\leq t\leq 1$ and
$t\geq 1$ respectively. The second inequality follows by noting that
$t/\sqrt{a+t^{2}}\geq t/\sqrt{a+1}$ when $0\leq t\leq 1$ and
$t/\sqrt{a+t^{2}}\geq 1/\sqrt{a+1}$ when $t\geq 1$. ∎

######  Lemma 23.

There exist independent random vectors $Z$ and $\varepsilon$ such that
$P^{Z},P^{\varepsilon}\in\mathcal{P}_{2}(\mathbb{R}^{d})\cap\mathcal{P}_{\mathrm{ac}}(\mathbb{R}^{d})$
such that $\langle\\!\langle
Z+\varepsilon,U\rangle\\!\rangle_{\mathcal{W}_{2}}^{2}=\langle\\!\langle
Z,U\rangle\\!\rangle_{\mathcal{W}_{2}}^{2}+\langle\\!\langle\varepsilon,U\rangle\\!\rangle_{\mathcal{W}_{2}}^{2}$.

###### Proof.

Consider independent random vectors $Z\sim\mathcal{N}(0,\Sigma)$ and
$\varepsilon\sim\mathcal{N}(0,\Gamma)$. By the same argument as in (30), we
have

$\displaystyle\langle\\!\langle
Z+\varepsilon,U\rangle\\!\rangle_{\mathcal{W}_{2}}$
$\displaystyle=\operatorname{\mathrm{Tr}}((\Sigma+\Gamma)^{1/2})$
$\displaystyle\langle\\!\langle Z,U\rangle\\!\rangle_{\mathcal{W}_{2}}$
$\displaystyle=\operatorname{\mathrm{Tr}}(\Sigma^{1/2})$
$\displaystyle\langle\\!\langle\varepsilon,U\rangle\\!\rangle_{\mathcal{W}_{2}}$
$\displaystyle=\operatorname{\mathrm{Tr}}(\Gamma^{1/2}).$

Hence, the desired result hold if we take $\Sigma=\sigma^{2}I_{d}$ and
$\Gamma=\gamma^{2}I_{d}$. ∎

######  Proposition 24.

Let $X$ be a random vecotor. Then the following properties are equivalent:

  1. (i)

There exists $\sigma>0$ such that $\operatorname{\mathbb{P}}(\|X\|\geq x)\leq
2e^{-\frac{1}{2}(x/\sigma)^{\beta}}$ for all $x\geq 0$.

  2. (ii)

There exists $K_{\sigma}>0$ such that $\\{\mathbb{E}\|X\|^{k}\\}^{1/k}\leq
K_{\sigma}k^{1/\beta}$.

  3. (iii)

There exists $K_{\sigma}^{\prime}>0$ such that
$\operatorname{\mathbb{E}}\exp{\bigl{(}(\lambda\|X\|)^{\beta}\bigr{)}}\leq\exp{\bigl{(}(\lambda
K_{\sigma}^{\prime})^{\beta}\bigr{)}}$ for all $|\lambda|\leq
1/K_{\sigma}^{\prime}.$

  4. (iv)

$X$ follows the $(\sigma,\beta)$-sub-weibull distribution.

The proof follows by [72, Theorem 2.1 ].

######  Proposition 25.

For $p_{1},p_{2}\in\mathbb{N}$, let
$X\in\mathbb{R}^{p_{1}},Y\in\mathbb{R}^{p_{2}}$ be two independent sub-Weibull
random vectors with parameter $(\sigma_{1},\alpha)$ and $(\sigma_{2},\beta)$
respectively. Then the following statements holds:

  1. (i)

For matrices $A\in\mathbb{R}^{d\times{p_{1}}}$ and
$B\in\mathbb{R}^{d\times{p_{2}}}$, there exists $\sigma>0$ depending only on
$\sigma_{1}$, $\sigma_{2}$, $\|A\|_{\mathrm{op}}$, $\|B\|_{\mathrm{op}}$ such
that $AX+BY\sim(\sigma,\alpha\wedge\beta)\textup{-sub-Weibull}$.

  2. (ii)

There exists $\sigma>0$ depending only on $\sigma_{1},\sigma_{2}$ such that
the concatenation of two random vectors $Z:=(X,Y)\in\mathbb{R}^{p_{1}+p_{2}}$
is a sub-Weibull random vector with parameter $(\sigma,\alpha\wedge\beta)$.

###### Proof.

(i) Suppose $K_{\sigma_{1}}$ and $K_{\sigma_{2}}$ are the induced constants of
$X$ and $Y$ by the part (ii) of Proposition 24. Then it follows that

$\displaystyle\bigl{(}\operatorname{\mathbb{E}}\|AX+BY\|^{k}\bigr{)}^{1/k}$
$\displaystyle\leq(\operatorname{\mathbb{E}}\|AX\|^{k})^{1/k}+(\operatorname{\mathbb{E}}\|BY\|^{k})^{1/k}$
$\displaystyle\leq\|A\|_{\mathrm{op}}(\operatorname{\mathbb{E}}\|X\|^{k})^{1/k}+\|B\|_{\mathrm{op}}(\operatorname{\mathbb{E}}\|Y\|^{k})^{1/k}$
$\displaystyle\leq\|A\|_{\mathrm{op}}\vee\|B\|_{\mathrm{op}}\cdot\bigl{(}K_{\sigma_{1}}k^{1/\alpha}+K_{\sigma_{2}}k^{1/\beta}\bigr{)}$
$\displaystyle\leq
2(\|A\|_{\mathrm{op}}\vee\|B\|_{\mathrm{op}})\cdot(K_{\sigma_{1}}\vee
K_{\sigma_{2}})k^{1/(\alpha\wedge\beta)}.$

This proves that $AX+BY$ satisfies part (ii) in the Proposition 24 thus the
conclusion follows by the equivalence of part (ii) and (iv).

(ii) For any integer $k\geq 1$, we have

$\displaystyle\bigl{(}\operatorname{\mathbb{E}}\|(X,Y)\|^{k}\bigr{)}^{1/k}$
$\displaystyle\leq\bigl{(}\operatorname{\mathbb{E}}(\|X\|+\|Y\|)^{k}\bigr{)}^{1/k}$
$\displaystyle\leq\bigl{(}\operatorname{\mathbb{E}}\|X\|^{k}\bigr{)}^{1/k}+\bigl{(}\operatorname{\mathbb{E}}\|Y\|^{k}\bigr{)}^{1/k}\leq(K_{\sigma_{1}}\vee
K_{\sigma_{2}})k^{1/(\alpha\wedge\beta)},$

where the sub-Weibull assumption on $X$ and $Y$ have been exploited. The
conclusion follows by employing Proposition 24. ∎

######  Lemma 26.

If $X$ is a $(\sigma,\beta)$-sub-Weibull random vector as defined in (11),
then for any $s>0$, there exists $C>0$, depending on $s,\sigma,\beta$, such
that $\mathbb{E}\bigl{(}\|X\|^{s}\operatorname{\mathds{1}}\\{\|X\|\geq
t\\}\bigr{)}\leq Ce^{-\frac{1}{4}(t/\sigma)^{\beta}}$.

###### Proof.

We have

$\displaystyle\mathbb{E}\bigl{(}\|X\|^{s}\operatorname{\mathds{1}}\\{\|X\|\geq
t\\}\bigr{)}$
$\displaystyle=\mathbb{E}\Bigl{[}\|X\|^{s}\operatorname{\mathds{1}}\Bigl{\\{}e^{\frac{1}{4}(\|X\|/\sigma)^{\beta}}\geq
e^{\frac{1}{4}(t/\sigma)^{\beta}}\Bigr{\\}}\Bigr{]}$
$\displaystyle\leq\operatorname{\mathbb{E}}\bigl{\\{}\|X\|^{s}e^{\frac{1}{4}(\|X\|/\sigma)^{\beta}}e^{-\frac{1}{4}(t/\sigma)^{\beta}}\bigr{\\}}$
$\displaystyle\leq
e^{-\frac{1}{4}(t/\sigma)^{\beta}}\bigl{\\{}\operatorname{\mathbb{E}}\|X\|^{2s}\bigr{\\}}^{1/2}\Bigl{\\{}\operatorname{\mathbb{E}}e^{\frac{1}{2}(\|X\|/\sigma)^{\beta}}\Bigr{\\}}^{1/2}$
$\displaystyle\leq
2^{1/2}e^{-\frac{1}{4}(t/\sigma)^{\beta}}\bigl{\\{}\operatorname{\mathbb{E}}\|X\|^{2s}\bigr{\\}}^{1/2},$

where we used the definition of $X$ being $(\sigma,\beta)$-sub-Weibull in
final step. The desired bound follows since by Proposition 24, we have
$\mathbb{E}\|X\|^{2s}\leq C$ for some constant $C$ that depends on
$s,\sigma,\beta$. ∎

######  Lemma 27.

Suppose $X,Y$ are independent $d$-dimensional random vectors with finite
second moment. If $X$ follows an absolutely continuous distribution with
density function $f_{X}$ which admits the following anti-concentration
inequality for some constant $\gamma_{1},\gamma_{2}>0$:

$f_{X}(x)\geq\gamma_{1}\exp\bigl{(}-\gamma_{2}\|x\|^{2}\bigr{)},\quad\forall\
\|x\|\geq\operatorname{\mathbb{E}}\|Y\|^{2}.$

Then there exists a constant $K_{1}$ depends on $\gamma_{1}$ and $\gamma_{2}$
such that the density function of $V:=X+Y$, write as $f_{V}$, satisfying

$\displaystyle f_{V}(v)\geq K_{1}\exp{(-2\gamma_{2}\|v\|^{2})},\quad\forall\
\|v\|\geq 2\operatorname{\mathbb{E}}\|Y\|^{2}.$

###### Proof.

Write $M_{2}:=\operatorname{\mathbb{E}}\|Y\|^{2}<+\infty$. For all $\|v\|\geq
2M_{2}$, we have

$\displaystyle f_{V}(v)$ $\displaystyle=\int
f_{X}(v-y)f_{Y}(y)dy\geq\int_{\|y\|\leq
M_{2}}\gamma_{1}\exp{\bigl{(}-\gamma_{2}\|v-y\|^{2}\bigr{)}}f_{Y}(y)dy$
$\displaystyle\geq\int_{\|y\|\leq
M_{2}}\gamma_{1}^{\prime}\exp{\bigl{(}-2\gamma_{2}\|v\|^{2}\bigr{)}}f_{Y}(y)dy\geq\gamma_{1}^{\prime}\bigl{(}1-\frac{1}{M_{2}}\bigr{)}\exp{\bigl{(}-2\gamma_{2}\|v\|^{2}\bigr{)}},$

where $\gamma_{1}^{\prime}=\gamma_{1}\exp(-2\gamma_{2}M_{2}^{2})$ and the last
inequality is followed by the Markov inequality. Thus the result holds by
letting $K_{1}=\gamma_{1}^{\prime}\bigl{(}1-\frac{1}{M_{2}}\bigr{)}$. ∎

######  Lemma 28.

Let $\mathcal{X},\mathcal{Y}\subseteq\mathbb{R}^{d}$ are Borel sets such that
$L_{2}$ is bounded on $\mathcal{X}\times\mathcal{Y}$, i.e.
$\|L_{2}\|_{\infty}:=\sup_{(x,y)\in\mathcal{X}\times\mathcal{Y}}L_{2}(x,y)<+\infty$.
Then for any $\mu\in\mathcal{P}_{2}(\mathcal{X})$ and
$\nu\in\mathcal{P}_{2}(\mathcal{Y})$ we have

$\displaystyle\inf$
$\displaystyle\big{\\{}J_{\mu,\nu}(\varphi,\psi):(\varphi,\psi)\in\tilde{\Phi}\\}$
$\displaystyle=\inf\big{\\{}J_{\mu,\nu}(\varphi,\psi):(\varphi,\psi)\in\tilde{\Phi},-\|L_{2}\|_{\infty}\leq\varphi-\|\cdot\|^{2}/2\leq
0,\ 0\leq\psi-\|\cdot\|^{2}/2\leq\|L_{2}\|_{\infty}\\},$

where $\tilde{\Phi}:=\\{(\varphi,\psi)\in L^{1}(\mathcal{X})\times
L^{1}(\mathcal{Y}):\varphi(x)+\psi(y)\geq x^{T}y,\
\forall(x,y)\in\mathcal{X}\times\mathcal{Y}\\}$.

###### Proof.

Note that by the argument same as (17) we have

$\displaystyle\inf$
$\displaystyle\big{\\{}J_{\mu,\nu}(\varphi,\psi):(\varphi,\psi)\in\tilde{\Phi}\\}$
$\displaystyle=\int_{\mathcal{X}}\frac{\|x\|^{2}}{2}d\mu(x)+\int_{\mathcal{Y}}\frac{\|y\|^{2}}{2}d\nu(y)-\sup\bigl{\\{}J_{\mu,\nu}(\varphi,\psi):(\varphi,\psi)\in\Phi_{2}\bigr{\\}},$
(75)

where $\Phi_{2}:=\\{(\varphi,\psi)\in L^{1}(\mathcal{X})\times
L^{1}(\mathcal{Y}):\varphi(x)+\psi(y)\leq L_{2}(x,y),\
\forall(x,y)\in\mathbb{R}^{d}\times\mathbb{R}^{d}\\}$. Note by [71, Remark
1.13 ], we may restrict the supremum in the right-hand side of (75) over some
bounded functions:

$\displaystyle\sup$
$\displaystyle\bigl{\\{}J_{\mu,\nu}(\varphi,\psi):(\varphi,\psi)\in\Phi_{2}\bigr{\\}}$
$\displaystyle=\sup\bigl{\\{}J_{\mu,\nu}(\varphi,\psi):(\varphi,\psi)\in\Phi_{2},0\leq\varphi\leq\|L_{2}\|_{\infty},\,-\|L_{2}\|_{\infty}\leq\psi\leq
0\bigr{\\}}.$ (76)

By [70, Theorem 5.10 ] we may further impose that $\varphi$ be c-concave and
$\psi=\varphi^{c}$. Suppose $(\varphi_{0},\varphi_{0}^{c})$ be a solution to
the right-hand side of (76). Define
$\tilde{\varphi}:=\|\cdot\|^{2}/2-\varphi_{0}$,
$\tilde{\psi}:=\|\cdot\|^{2}/2-\varphi_{0}^{c}$. Then by (75) we have

$\displaystyle\inf\big{\\{}J_{\mu,\nu}(\varphi,\psi):(\varphi,\psi)\in\tilde{\Phi}\\}=\int_{\mathcal{X}}\tilde{\varphi}(x)d\mu(x)+\int_{\mathcal{Y}}\tilde{\psi}(y)d\nu(y).$
(77)

Moreover, note

$\displaystyle\tilde{\varphi}(x)=\|x\|^{2}/2-\varphi_{0}(x)=\|x\|^{2}/2-\inf_{y\in\mathcal{Y}}\\{c(x,y)-\varphi_{0}^{c}(y)\\}=\sup_{y\in\mathcal{Y}}\\{x^{T}y-(\|y\|^{2}/2-\varphi_{0}^{c}(y))\\},$

which implies that $(\tilde{\varphi},\tilde{\psi})\in\tilde{\Phi}$. Combine
this with (77), we proved that $(\tilde{\varphi},\tilde{\psi})$ is an optimal
solution to the left-hand side of (75). Finally, by the boundedness of
$\varphi_{0}$ and $\varphi_{0}^{c}$, we have

$\displaystyle 0\leq\|x\|^{2}/2-\tilde{\varphi}(x)\leq\|L_{2}\|_{\infty}\
\text{and}\ -\|L_{2}\|_{\infty}\leq\|y\|^{2}/2-\tilde{\psi}(y)\leq 0,$

as desired. ∎

######  Theorem 29.

[62, Theorem 1] Let $X\sim P^{X}$ be a probability measure on $\mathbb{R}^{d}$
such that $M_{\ell}:=\operatorname{\mathbb{E}}\|X\|^{\ell}<+\infty$ with
$\ell\in(2,+\infty)$. If $P_{n}^{X}$ is the corresponding empirical
distribution, then there exists a constant $C>0$ depending only on $d$ and
$\ell$ such that for all $n\geq 1$,

$\displaystyle\operatorname{\mathbb{E}}\left[\mathcal{W}_{2}^{2}(P^{X},P_{n}^{X})\right]\leq
CM_{\ell}^{2/\ell}\tau_{n}(d,\ell),$ (78)

where

$\displaystyle\tau_{n}(d,\ell):=\begin{cases}n^{-\frac{1}{2}}&\text{ if
$d<4$}\\\ n^{-\frac{1}{2}}\log(1+n)&\text{if $d=4$}\\\
n^{-\frac{2}{d}}&\text{if
$d>4$}\end{cases}+\begin{cases}n^{-\frac{1}{d}}&\text{if $\ell>4$}\\\
n^{-\frac{1}{2}}\log(1+n)&\text{if $\ell=4$}\\\
n^{\frac{2-\ell}{\ell}}&\text{if $2<\ell<4$}.\end{cases}$

##  Appendix C Spatial reference distribution

In this section, we derive the MCQR loss function under reference distribution
$U[-1,1]$, which may provide an intuitive example for the verification of
Proposition 3. In one dimension, the traditional rank and quantile can be
understood as a pair of optimal transport maps between the distribution of
interest $X\sim P$ and the uniform distribution $U\sim U[0,1]$. When $P$ does
not assign mass to sets with Hausdorff dimension $0$, the corresponding
distribution function $F$ and its inverse map $Q:=F^{-1}$ serve as the
corresponding optimal transport map. This concept can be generalized to other
reference distributions, for instance, $U[-1,1]$. In this case, the spatial
distribution function $F_{\mathrm{sp}}(\cdot):=2F(\cdot)-1$ takes on the role
of $F$ in the previous case. Moreover, the corresponding check function needs
to be modified as

$\displaystyle\rho_{\tau}^{\mathrm{sp}}(X-\theta):=(1+\tau)(X-\theta)-2(X-\theta)\operatorname{\mathds{1}}\\{X-\theta<0\\},\quad\forall\tau\in[-1,1].$

Suppose $V\sim U[-1,1]$, then the composite quantile regression optimization
becomes

$\displaystyle\operatorname{\mathbb{E}}\int_{-1}^{1}\rho_{\tau}^{\mathrm{sp}}\bigl{(}Y-\beta^{\top}X$
$\displaystyle-q(\tau)\bigr{)}\cdot\frac{1}{2}\,d\tau=\operatorname{\mathbb{E}}\int_{-1}^{1}\bigl{(}Y-\beta^{\top}X-q(\tau)\bigr{)}^{-}\,d\tau+\int_{-1}^{1}\int_{\tau}^{-1}\frac{1}{2}q(\tau)\,dtd\tau$
$\displaystyle=\operatorname{\mathbb{E}}\max_{t\in[-1,1]}\int_{t}^{1}-\bigl{(}Y-bX-q(\tau)\bigr{)}\,d\tau+\int_{-1}^{1}\int_{t}^{1}\frac{1}{2}q(\tau)\,d\tau
dt$
$\displaystyle=\operatorname{\mathbb{E}}\max_{t\in[-1,1]}\bigl{(}-(1-t)(Y-bX)+\phi(t)\bigr{)}+\operatorname{\mathbb{E}}\phi(V)$
$\displaystyle=\operatorname{\mathbb{E}}\max_{t\in[-1,1]}\bigl{(}t(Y-bX)+\phi(t)\bigr{)}+\operatorname{\mathbb{E}}\phi(V),$

where $\phi(t)=\int_{t}^{1}q(\tau)d\tau.$ Thus, applying the same argument as
Lemma 2 we can see that the composition quantile regression estimator of
$b^{*}$ is once again

$\displaystyle b^{*}=\operatorname*{arg\,min}\langle\\!\langle
P^{Y-bX},P^{V}\rangle\\!\rangle_{\mathcal{W}_{2}}.$

This gives some intuition on Proposition 3. However, if choose the standard
normal distribution as the reference distribution, we may not be able to find
a straightforward optimal transport map as $F$ or $F_{\mathrm{sp}}$, but
Proposition 3 demonstrates the validity of this extension.

##  Appendix D Spatial quantile

The concept of the spatial (or geometric) quantile was initially introduced by
[59]. Uniquely characterizing the underlying probability distribution as a
special case of M-quantile, as demonstrated in [65, Theorem 2.5], this
quantile permits a seamless extension to the regression framework [58] and
functional quantile regression [57, 60]. A more recent development involves an
extension to the hypersphere, as explored by [66].

The definition of Spatial quantile starts from rewriting the check function
$\rho_{\tau}(\cdot)$ as

$\displaystyle\rho_{\tau}(z)$
$\displaystyle=\frac{1}{2}\bigl{(}|z|+(2\tau-1)z\bigr{)}=\frac{1}{2}(|z|+vz),\
\text{for any $z\in\mathbb{R}$},$

with $v=2\tau-1$. Thus a natural extension of the check function to the multi-
dimensional case is by substituting the absolute value function by the
$L_{1}$-loss function:

$\displaystyle\Phi_{v}(z):=\frac{1}{2}(\|z\|+v^{\top}z),$

where $v=\tau u$, and $u\in\mathcal{S}^{d-1}$. This extension of the check
function immediately leads to the following definition of spatial quantile:

######  Definition 30.

Suppose $Y\sim\operatorname{\mathbb{P}}^{Y}$ is a random variable on
$\mathbb{R}^{d}$ $(d\geq 1)$. Then for any $\tau\in[0,1]$ and
$u\in\mathcal{S}^{d-1}$, the $\tau u$-spatial quantile of $P^{Y}$ is defined
as

$\displaystyle\mathrm{Q}_{\tau
u}=\operatorname*{arg\,min}_{y\in\mathbb{R}^{d}}\operatorname{\mathbb{E}}\Phi_{\tau
u}(Y-y).$ (79)

Note the solution of (79) are such that

$\displaystyle\operatorname{\mathbb{E}}\bigl{(}\frac{Y-Q_{\tau u}}{\|Y-Q_{\tau
u}\|}\bigr{)}=-\tau u.$

Intuitively speaking, this indicates that $Q_{\tau u}$ defines a point in
$\mathbb{R}^{d}$ such that the average unit vector from it to other random
samples should be $\tau u$.

The generalization to quantile regression setting is simply by applying the
spatial quantile definition to $Y-b^{*}X-a$, where $X\in\mathbb{R}^{p}$ is the
covariate vector, $b^{*}\in\mathbb{R}^{d\times p}$ is the regression
coefficient and $a$ is the intercept term. Specifically, for fixed
$\tau\in[0,1]$ and $u\in\mathcal{S}^{d-1}$,

$\displaystyle(a_{\tau
u},b_{\tau,u})=\operatorname*{arg\,min}_{b\in\mathbb{R}^{d\times p},\
a\in\mathbb{R}^{d}}\operatorname{\mathbb{E}}\Phi_{\tau u}(Y-bX-a).$

Therefore, given observations $(Y_{1},X_{1}),\ldots,(Y_{n},X_{n})$ satisfying
equations

$Y_{i}=b^{*}X_{i}+\varepsilon_{i},\ i=1,\ldots,n,$

for some random residue terms $\varepsilon_{i}$’s that are independent with
with $X_{i}$’s, the spatial quantile estimator of $b^{*}$ can be obtained by

$\displaystyle(\hat{b}^{(\mathrm{sp})},\hat{a}_{\tau
u}^{(\mathrm{sp})})\operatorname*{arg\,min}\limits_{b\in\mathbb{R}^{d\times
p},\ a\in\mathbb{R}^{p}}\frac{1}{n}\sum_{i=1}^{n}\Phi_{\tau
u}(Y_{i}-bX_{i}-a).$

Therefore, the optimizer can be obtained by applying classical convex
optimization algorithms.