# Spin-Dependent Interactions and Heavy-Quark Transport in the QGP

Zhanduo Tang  Ralf Rapp  Cyclotron Institute and Department of Physics and
Astronomy, Texas A&M University, College Station, TX 77843-3366, USA

###### Abstract

We extend a previously constructed $T$-matrix approach to the quark-gluon
plasma (QGP) to include the effects of spin-dependent interactions between
partons. Following earlier work within the relativistic quark model, the spin-
dependent interactions figure as relativistic corrections to the Cornell
potential. When applied to the vacuum spectroscopy of quarkonia, in particular
their mass splittings in $S$\- and $P$-wave states, the issue of the Lorentz
structure of the confining potential arises. We confirm that a significant
admixture of a vector interaction (to the previously assumed scalar
interaction) improves the description of the experimental mass splittings. The
temperature corrections to the in-medium potential are constrained by results
from thermal lattice-QCD for the equation of state (EoS) and heavy-quark (HQ)
free energy in a selfconsistent set-up for heavy- and light-parton spectral
functions in the QGP. We then deploy the refined in-medium heavy-light
$T$-matrix to compute the charm-quark transport coefficients in the QGP. The
vector component of the confining potential, through its relativistic
corrections, enhances the friction coefficient for charm quarks in the QGP
over previous calculations by tens of percent at low momenta and temperatures,
and more at higher momenta. Our results are promising for improving the
current phenomenology of open heavy-flavor observables at Relativistic Heavy
Ion Collider (RHIC) and the Large Hadron Collider (LHC).

$T$-matrix, heavy-quarkonium spectroscopy, spin-dependent interaction, mixed
confining potential, heavy-quark transport

##  I Introduction

The exploration of hadron properties in vacuum and the properties of the
quark-gluon plasma (QGP) are usually regarded as rather independent areas in
the study of Quantum Chromodynamics (QCD). However, in both areas the basic
building block are soft parton interactions rooted in the non-perturbative
sector of the theory, albeit in different environments. Of particular interest
are heavy quarks: heavy-quarkonium spectroscopy in vacuum has provided deep
insights into potential between a heavy (charm or bottom) quark ($Q=c,b$) and
its antiquark ($\bar{Q}$). The Cornell potential and its refinements remain a
phenomenologically successful tool in the description of the pertinent bound
states, taking advantage of expansion in the inverse heavy-quark (HQ) mass,
$1/M_{Q}$ Brambilla et al. (2004). The long-range (linear) part of the
potential, which by now is also well established in lattice QCD (lQCD), is
arguably one of the most direct manifestations of the confining force in QCD.
In the context of high-temperature QCD and its study in ultrarelativistic
heavy-ion collisions (URHICs), this led to the idea of utilizing quarkonium
production as a probe of deconfinement, although the originally proposed
suppression signature has evolved considerably over the last three decades
Rapp et al. (2010); Kluberg and Satz (2010); Braun-Munzinger and Stachel
(2010). Specifically, transport approaches have been developed toward the more
general objective of deducing the in-medium QCD force from quarkonium
observables, by implementing it into the transport coefficients that govern
both suppression and regeneration reactions in the evolving fireball of a
heavy-ion collision, see, e.g., Ref. Du et al. (2019). This effort is
critically aided by ample information from lQCD on the in-medium properties of
quarkonia through HQ free energies and Euclidean correlation functions
Petreczky and Teaney (2006); Petreczky et al. (2006); Petreczky (2009); Ding
et al. (2012), which constrain calculations of spectral functions that can
serve as an interface to phenomenological applications Rapp and van Hees
(2010); Mocsy et al. (2013); He et al. (2022).

Open heavy-flavor (HF) particles have emerged as an excellent probe of the
transport properties of the QCD medium in URHICs Prino and Rapp (2016);
Beraudo et al. (2018); Dong and Greco (2019). Produced in initial hard
processes, low-momentum heavy quarks exert a Brownian motion through the QGP
characterized by a spatial diffusion coefficient, hadronize in different HF
hadrons and subsequently are further transported through the hadronic medium.
The large HQ mass implies the dominance of elastic interactions with small
energy transfer amenable to potential approximations, and the final HF baryon
spectra carry a memory of their interaction history due to a thermalization
time being comparable or larger than the fireball lifetime.

The present work builds on previous efforts to develop a quantum many-body
theory to describe the spectral and transport properties of open and hidden HF
particles in a strongly coupled QGP Cabrera and Rapp (2007); Riek and Rapp
(2010), including the 1- and 2-body Green’s functions of thermal partons for
obtaining the equation of state (EoS) in a selfconsistent Brueckner scheme
Mannarelli and Rapp (2005); Liu et al. (2019). The basic ingredient to this
framework is the 2-body interaction kernel for the in-medium $T$-matrix for
which we employ an ansatz using a Cornell potential, whose temperature
corrections are constrained by lQCD data for the HQ free energy. A salient
feature of this approach is that it recovers basic features of vacuum
spectroscopy (such as masses of quarkonia, $D$ and $B$ mesons, and non-
Goldstone light hadrons), providing a baseline for the calculation of medium
effects. In the spirit of a $1/M_{Q}$ expansion, spin-orbit and spin-spin
interactions were not included thus far. In the present paper, we take the
next step by including the latter by benchmarking them against the hyper-/fine
mass splittings of quarkonia in vacuum.

Our study raises the question of the Lorentz structure of the confining
potential. Historically, a default assumption of a purely scalar interaction
has been employed Mur et al. (1994); Lucha et al. (1991), implying a vanishing
long-range magnetic contribution that nevertheless could reproduce the
empirical fine structure for heavy quarkonium Buchmuller (1982). However,
studies of the Wilson loop suggest that the confining potential cannot be a
purely scalar kernel Brambilla and Vairo (1997a, b), and the latter also
causes problems in constructing a stable vacuum of QCD Szczepaniak and Swanson
(1997). In the relativistic quark model Ebert et al. (1998, 2003) a mixing of
scalar and vector structures in the confining potential has been found to
yield a quarkonium spectroscopy in good overall agreement with experimental
data. The approach employed in our work is close in spirit to these works, and
therefore we will incorporate the possibility of a mixed confining Lorentz
structure in the $T$-matrix kernel; the pertinent relativistic corrections
will turn out to have significant ramifications for the HQ diffusion
coefficient.

This article is organized as follows. In Sec. II, we briefly recollect the
main elements of the thermodynamic $T$-matrix approach. In Sec. III we
implement spin-dependent interactions as well as a vector component of the
confining force into the potential. In Sec. IV we compute heavy-quarkonium
spectral functions from the $T$-matrix and discuss the charmonium and
bottomonium spectroscopy in vacuum. In Sec. V we lay out our constraints on
the in-medium corrections to the potential using lQCD data for static HQ free
energies (Sec. V.1) and the QGP equation of state (Sec. V.2), and discuss the
pertinent numerical results (Sec. V.3). In Sec. VI we outline the calculation
of the HQ transport coefficients and highlight the implications of the vector
component in the confining interaction on the numerical results for charm
quarks. We summarize and conclude in Sec. VII.

##  II T-matrix Approach

The thermodynamic $T$-matrix is a 2-particle irreducible (PI) quantum many-
body scheme that selfconsistently solves the 1- and 2-body Green’s functions
and is thus suitable for strongly interacting systems. In Refs. Mannarelli and
Rapp (2005); Riek and Rapp (2010) it has been initially developed to study the
properties of HF flavor particles in the QGP allowing for a reduction of the
4-dimensional (4D) Bethe-Salpeter 2-body scattering equation to a 3D one which
allows for tractable numerical solutions. Subsequently, it has also been
extended to the light-parton sector Liu and Rapp (2018), based on the notion
that the effective masses of the QGP’s constituents are typically large
compared to temperatures not too far above the pseudocritcal one of $T_{\rm
pc}\simeq 160$ MeV. The starting point can be formulated in terms of an
effective Hamiltonian with a relativistic potential,

$\displaystyle H=$
$\displaystyle\sum\varepsilon_{i}(\mathbf{p})\psi_{i}^{\dagger}(\mathbf{p})\psi_{i}(\mathbf{p})+\frac{1}{2}\psi_{i}^{\dagger}\left(\frac{\mathbf{P}}{2}-\mathbf{p}\right)$
(1)
$\displaystyle\times\psi_{j}^{\dagger}\left(\frac{\mathbf{P}}{2}+\mathbf{p}\right)V_{ij}^{a}\psi_{j}\left(\frac{\mathbf{P}}{2}+\mathbf{p}^{\prime}\right)\psi_{i}\left(\frac{\mathbf{P}}{2}-\mathbf{p}^{\prime}\right)$

which emphasizes the implementation of unitarity through resummations of the
propagators (also referred to as a Dyson-Schwinger set-up). Here, $\mathbf{p}$
and $\mathbf{p^{\prime}}$ denote the relative momentum of the incoming and
outgoing states, and $\mathbf{P}$ the total momentum of the two-body system.
Furthermore, $\varepsilon_{i}(\mathbf{p})=\sqrt{M_{i}^{2}+\mathbf{p}^{2}}$ is
the dispersion relation of a parton with mass $M_{i}$, and the $V_{ij}^{a}$
are the potentials between particles $i$ and $j$ in a color channel $a$. The
summation includes momentum, spin, color and flavor for quarks and gluons. The
infinite series of ladder diagrams generated by the Hamiltonian in Eq. (1)
straightforwardly results in the $T$-matrix equation, depicted in Fig. 1.

Figure 1: $T$-matrices resummation for ladder diagrams.

In the center-of-mass (CM) frame, one has

$\displaystyle\ T_{ij}^{a}\left(z,\mathbf{p},\mathbf{{p}^{\prime}}\right)$
$\displaystyle=$ $\displaystyle
V_{ij}^{a}\left(\mathbf{p},\mathbf{{p}^{\prime}}\right)+\int_{-\infty}^{\infty}\frac{d^{3}\mathbf{k}}{\left(2\pi\right)^{3}}V_{ij}^{a}\left(\mathbf{p},\mathbf{k}\right)$
(2) $\displaystyle\times
G_{ij}^{0}\left(z,\mathbf{k}\right)T_{ij}^{a}\left(z,\mathbf{k},\mathbf{{p}^{\prime}}\right)\
,$

where $G_{ij}^{0}$ is the two-body propagator, $z=E\pm i\epsilon$ the
analytical energy variable, and $\mathbf{p}$ and $\mathbf{p^{\prime}}$ are the
incoming and outgoing 3-momenta in the CM frame, respectively. The reduction
scheme from 4D to 3D is not unique Brockmann and Machleidt (1999) but its
specific choice has minor impact on the results; we choose the Thompson scheme
following our previous studies Riek and Rapp (2010); Liu and Rapp (2018). In
this scheme, the two-body propagator in spectral representation can be written
as

$\displaystyle G_{ij}^{0}(z,\mathbf{k})$ $\displaystyle=$
$\displaystyle\int_{-\infty}^{\infty}d\omega_{1}d\omega_{2}\frac{\left[1\pm
n_{i}\left(\omega_{1}\right)\pm
n_{j}\left(\omega_{2}\right)\right]}{z-\omega_{1}-\omega_{2}}$ (3)
$\displaystyle\times\rho_{i}\left(\omega_{1},\mathbf{k}\right)\rho_{j}\left(\omega_{2},\mathbf{k}\right)\
,$

with the single-particle propagator

$G_{i}(z)=\frac{1}{\left[G_{i}^{0}(z,k)\right]^{-1}-\Sigma_{i}(z,k)}=\frac{1}{z-\varepsilon_{i}(k)-\Sigma_{i}(z,k)}$
(4)

and the single-particle spectral function

$\rho_{i}\left(\omega_{,}\mathbf{k}\right)=-\frac{1}{\pi}\operatorname{Im}G_{i}(\omega+i\epsilon)\
.$ (5)

The ± signs in Eq. (3) correspond to bosons (upper) or fermions (lower)111The
convention that upper/lower signs denote bosons/fermions is applied throughout
this work., and $n_{i}$ is the Bose or Fermi distribution function for parton
$i$. In quasi-particle approximation Eq. (3) reduces to222This differs from
Refs. Brockmann and Machleidt (1999); Riek and Rapp (2010) by a factor of
$m_{ij}(\mathbf{k})=\frac{M_{i}M_{j}}{\varepsilon_{i}(\mathbf{k})\varepsilon_{j}(\mathbf{k})}$;
here we keep the convention of Ref. Liu and Rapp (2018) where
$m_{ij}(\mathbf{k})$ is absorbed into the relativistic corrections to the
potential which will be elaborated in Sec. III.2.

$\displaystyle G_{ij}^{0}(z,\mathbf{k})$ $\displaystyle=$
$\displaystyle\frac{1}{z-\varepsilon_{i}(\mathbf{k})-\varepsilon_{j}(\mathbf{k})-\Sigma_{i}(\mathbf{k})-\Sigma_{j}(\mathbf{k})}\
.$ (6)

The single-particle selfenergies in the QGP, $\Sigma_{i}(\mathbf{k})$, are
obtained by closing the $T$-matrix with an in-medium single-parton propagator
from the heat bath; its spectral representation is

$\displaystyle\Sigma_{i}\left(z,\mathbf{p}_{1}\right)$ $\displaystyle=$
$\displaystyle\frac{1}{d_{i}}\int\frac{d^{3}\mathbf{p}_{2}}{(2\pi)^{3}}\int_{-\infty}^{\infty}d\omega_{2}\frac{dE}{\pi}\frac{-1}{z+\omega_{2}-E}$
(7)
$\displaystyle\times\sum_{a,j}d_{s}^{ij}d_{a}^{ij}\operatorname{Im}T_{ij}^{a}\left(E,\mathbf{p}_{1},\mathbf{p}_{2}\mid\mathbf{p}_{1},\mathbf{p}_{2}\right)$
$\displaystyle\times\rho_{j}\left(\omega_{2},\mathbf{p}_{2}\right)\left[n_{j}\left(\omega_{2}\right)\mp
n_{ij}(E)\right],$

with
$T\left(E,\mathbf{p}_{1},\mathbf{p}_{2}\mid\mathbf{p}_{1},\mathbf{p}_{2}\right)$
the forward-scattering $T$-matrix, i.e.,
$\mathbf{p}_{1}^{\prime}=\mathbf{p}_{1}$ and
$\mathbf{p}_{2}^{\prime}=\mathbf{p}_{2}$, where $\mathbf{p}_{1,2}$ and
$\mathbf{p}_{1,2}^{\prime}$ are the incoming and outgoing momenta for particle
1 and 2 respectively, defined in the thermal frame. The $n_{ij}$ refers to the
thermal distribution for the two-body state $ij$, while $\mp$ refers to the
bosonic/fermionic single-parton state $i$. The $d_{a,s}^{ij}$ are color and
spin degeneracies of the two-body system, and $d_{i}$ is the spin-color
degeneracy of the single parton $i$. We also need to add the purely real
thermal Fock term Fetter and Walecka (2012),

$\displaystyle\Sigma_{i}\left(\mathbf{p}_{1}\right)$ $\displaystyle=$
$\displaystyle\mp\int\frac{d^{3}\mathbf{p}_{2}}{(2\pi)^{3}}\int_{-\infty}^{\infty}d\omega_{2}V_{i\bar{i}}^{a=1}\left(\mathbf{p}_{1}-\mathbf{p}_{2}\right)\rho_{i}\left(\omega_{2},\mathbf{p}_{2}\right)$
(8) $\displaystyle\times n_{i}\left(\omega_{2}\right)\ ,$

which is not part of the selfenergy in Eq. (7). The $V_{i\bar{i}}^{a=1}$
refers to the color-singlet potential between particle and antiparticle. The
selfenergy can be solved selfconsistently by iterating Eqs. (2), (7) and (8)
numerically. In doing so, the $T$-matrix in the thermal frame,
$T_{ij}^{a}\left(\omega_{1}+\omega_{2},\mathbf{p}_{1},\mathbf{p}_{2}\mid\mathbf{p}_{1}^{\prime},\mathbf{p}_{2}^{\prime}\right)$
needs to be transformed into the the CM frame,
$T_{ij}^{a}\left(E_{\mathrm{cm}},p_{\mathrm{cm}},p_{\mathrm{cm}}^{\prime},\cos(\theta_{\mathrm{cm}})\right)$.
This is accomplished by

$\displaystyle E_{\mathrm{cm}}$ $\displaystyle=$
$\displaystyle\sqrt{\left(\omega_{1}+\omega_{2}\right)^{2}-\left(\mathbf{p}_{1}+\mathbf{p}_{2}\right)^{2}}$
$\displaystyle s_{\mathrm{on}}$ $\displaystyle=$
$\displaystyle\left(\varepsilon_{1}\left(\mathbf{p}_{1}\right)+\varepsilon_{2}\left(\mathbf{p}_{2}\right)\right)^{2}-\left(\mathbf{p}_{1}+\mathbf{p}_{2}\right)^{2}$
$\displaystyle p_{\mathrm{cm}}$ $\displaystyle=$
$\displaystyle\sqrt{\frac{\left(s_{\mathrm{on}}-M_{i}^{2}-M_{j}^{2}\right)^{2}-4M_{i}^{2}M_{j}^{2}}{4s_{\mathrm{on}}}}$
$\displaystyle\cos\left(\theta_{\mathrm{cm}}\right)$ $\displaystyle=$
$\displaystyle\frac{\mathbf{p}_{\mathrm{cm}}\cdot\mathbf{p}_{\mathrm{cm}}^{\prime}}{p_{\mathrm{cm}}p_{\mathrm{cm}}^{\prime}},$
(9)

where $\cos(\theta_{\mathrm{cm}})$ is the angle between the incoming and
outgoing momenta in the CM frame, and $p_{\mathrm{cm}}^{\prime}$ can be
obtained by substituting $s_{\mathrm{on}}(\mathbf{p}_{1},\mathbf{p}_{2})$ with
$s_{\mathrm{on}}(\mathbf{p}_{1}^{\prime},\mathbf{p}_{2}^{\prime})$. As
discussed in Ref. Liu and Rapp (2018), the reason for using the on-shell
value, $s_{on}$, for $p_{\mathrm{cm}}$ is to preserve the analytical
properties of the $T$-matrix after the transformation into the CM frame.

The 3D $T$-matrix integral equation can be further reduced to a 1D one by
applying the partial-wave expansion in the CM frame (from hereon the subscript
“cm” is suppressed for simplicity),

$X(\mathbf{p},\mathbf{p^{\prime}})=4\pi\sum_{L}(2L+1)X^{L}(p,p^{\prime})P_{L}[\cos(\theta)]\
,$ (10)

where $X$ denotes $V$ or $T$, $L$ the angular-momentum quantum number, and $p$
and $p^{\prime}$ are the moduli of $\mathbf{p}$ and $\mathbf{p^{\prime}}$. The
1D $T$-matrix equation then takes the form

$\displaystyle\ T_{ij}^{L,a}(z,p,p^{\prime})$ $\displaystyle=$ $\displaystyle
V_{ij}^{{L,a}}(p,p^{\prime})+\frac{2}{\pi}\int_{-\infty}^{\infty}k^{2}dkV_{ij}^{{L,a}}(p,k)$
(11) $\displaystyle\times G_{ij}^{0}(z,k)T_{ij}^{L,a}(z,k,p^{\prime}).$

Equation (11) can be solved by discretizing the momenta to convert it into a
matrix equation and solve it by matrix inversion.

##  III Two-Body Potentials in Vacuum

In this section we first discuss the static potentials in Sec. III.1, then
introduce the relativistic corrections to the potential between particles $i$
and $j$, and construct the confining potential with mixed Lorentz structures
in Sec. III.2. The potential is generalized to different color channels at the
end of this section. For simplicity, we suppress the color factors indices
until the end of this section.

###  III.1 Static Potential

The kernel of the $T$-matrix equation (2) is based on the Cornell potential,
with a color-Coulomb potential, $V_{C}$, plus a confining potential (“string”
term), $V_{S}$. In coordinate space the common ansatz is

$\widetilde{V}(r)={\widetilde{V}}_{\mathcal{C}}(r)+{\widetilde{V}_{\mathcal{S}}}(r)=-\frac{4}{3}\frac{\alpha_{s}}{r}+\sigma
r\ ,$ (12)

where $\alpha_{s}$ and $\sigma$ are the perturbative coupling constant and
nonperturbative string tension, respectively. To obtain the momentum-space
potentials, ${V}_{\mathcal{C}/\mathcal{S}}(\mathbf{k})$, depending on the
momentum transfer $\mathbf{k}=\mathbf{p-p^{\prime}}$, we use the subtracted
quantities
${V}_{\mathcal{C}/\mathcal{S}}(r)={\widetilde{V}}_{\mathcal{C}/\mathcal{S}}(r)-{\widetilde{V}}_{\mathcal{C}/\mathcal{S}}(\infty)$
to ensure the convergence of the Fourier transforms. A running coupling is
implemented in the Coulomb potential for off-shell scattering in momentum
space as Riek and Rapp (2010)
$F_{run}(p,p^{\prime})=\ln\left[\frac{\Delta^{2}}{\Lambda^{2}}\right]/\ln\left[\frac{(p-p^{\prime})^{2}+\Delta^{2}}{\Lambda^{2}}\right]$.
For the confining potential, we enforce a flat potential above a string
breaking scale of about $r_{SB}=1$ fm to account for string breaking. We
employ the same potential parameters as in previous studies Liu and Rapp
(2018), i.e., $\alpha_{s}=0.27$, $\sigma=0.225$ $\textup{GeV}^{2}$, $\Delta=1$
GeV and $\Lambda=0.2$ GeV, which are fitted to the lQCD data of the vacuum
free energy Cheng et al. (2008); Petreczky (2010); Mocsy et al. (2013);
Bazavov et al. (2014a); Bazavov et al. (2014b, 2018) (equivalent to the vacuum
potential ) as shown in Fig. 2.

Figure 2: The fitted vacuum potential versus the lQCD data. The blue line
denotes the fitted vacuum potential, the colored dots the lQCD data from Refs.
Cheng et al. (2008); Petreczky (2010); Mocsy et al. (2013); Bazavov et al.
(2014a); Bazavov et al. (2014b, 2018).

###  III.2 Relativistic Corrections and Spin-Dependent Interactions

Relativistic effects in the one-gluon exchange amplitude are well known,
containing spin-independent and spin-dependent ones. For the pertinent vector
potential (denoted as $V^{vec}$), the spin-independent correction amounts to
multiplying the momentum-space potential by a factor

$\displaystyle\mathcal{R}_{ij}$ $\displaystyle=$
$\displaystyle\sqrt{\frac{1}{m_{ij}(p)}}\sqrt{1+\frac{p^{2}}{\varepsilon_{i}(p)\varepsilon_{j}(p)}}$
(13)
$\displaystyle\times\sqrt{\frac{1}{m_{ij}(p^{\prime})}}\sqrt{1+\frac{p^{\prime
2}}{\varepsilon_{i}(p^{\prime})\varepsilon_{j}(p^{\prime})}}\ ,$

which is known as the Breit correction, representing magnetic effects. For
scalar potentials, denoted as $V^{sca}$, no relativistic correction arises to
leading order in $1/M_{Q}$, see Ref. Riek and Rapp (2010). We write the total
spin-independent potential in momentum space as

$\displaystyle V_{ij}\left(\mathbf{p},\mathbf{p}^{\prime}\right)$
$\displaystyle=$
$\displaystyle\mathcal{R}_{ij}V^{vec}\left(\mathbf{p}-\mathbf{p}^{\prime}\right)+V^{sca}\left(\mathbf{p}-\mathbf{p}^{\prime}\right).$
(14)

To implement spin-dependent interactions, including spin-orbit ($V^{LS}$),
spin-spin ($V^{SS}$) and tensor ($V^{T}$) channels, we follow Ref. Lucha and
Schoberl (1995) where the detailed procedure to derive the Fermi-Breit
Hamiltonian is laid out. The pertinent corrections for vector and scalar
potentials between two partons with equal masses ($M_{i}=M_{j}\equiv M$) in
coordinate space are given by

$\displaystyle{V^{LS}}(r)$ $\displaystyle=$
$\displaystyle\frac{1}{2M^{2}r}\left\langle\mathbf{L\cdot
S}\right\rangle\left(3\frac{d}{dr}V^{vec}(r)-\frac{d}{dr}V^{sca}(r)\right),$
$\displaystyle{V^{SS}}(r)$ $\displaystyle=$
$\displaystyle\frac{2}{3M^{2}}\left\langle\mathbf{S_{1}\cdot
S_{2}}\right\rangle\Delta V^{vec}(r),$ $\displaystyle{V}^{T}(r)$
$\displaystyle=$ $\displaystyle\frac{1}{12M^{2}}\left\langle
S_{12}\right\rangle\left(\frac{1}{r}\frac{d}{dr}V^{vec}(r)-\frac{d^{2}}{dr^{2}}V^{vec}(r)\right)\
,$

where $\Delta\equiv\nabla^{2}$ in the $SS$ interaction is the Laplace
operator. Note that the scalar interactions do not contribute to the $SS$ and
$T$ corrections. We note that the vector potential, $V^{vec}$, in the spin-
dependent potentials above do not receive the spin-independent Breit
correction, $\mathcal{R}$, introduced in Eq. (13). Following Ref. Hong et al.
(2022), we smear the Dirac delta function $\delta$ in the $SS$ part by a
Gaussian,
$\tilde{\delta}(r)=\left(\frac{b}{\sqrt{\pi}}\right)^{3}e^{-b^{2}r^{2}}$, to
avoid the singularity. We take $b=10$ in this work and have checked that for
$b>$10 the $SS$ interaction saturates in the quarkonium spectroscopy.

The expectation values take the standard form (with $L$, $S$ and $J$ denoting
the orbital, spin and total angular-momentum quantum numbers, respectively):
$\left\langle\mathbf{L\cdot
S}\right\rangle=\frac{1}{2}[J(J+1)-L(L+1)-S(S+1)]$,
$\left\langle\mathbf{S_{1}\cdot
S_{2}}\right\rangle=\frac{1}{2}[S(S+1)-\frac{3}{2}]$, and $\left\langle
S_{12}\right\rangle=\frac{4}{(2L+3)(2L-1)}\left[S(S+1)J(J+1)-\frac{3}{2}\left\langle\mathbf{L\cdot
S}\right\rangle-3(\left\langle\mathbf{L\cdot S}\right\rangle)^{2}\right]$ for
$L\neq 0$ and $S=1$, but $\left\langle S_{12}\right\rangle$ vanishes for
either $L=0$ or $S=0$. The total potential (with relativistic corrections)
between, e.g., a heavy quark and anti-quark in momentum-space reads

$\displaystyle
V_{Q\bar{Q}}\left(\mathbf{p},\mathbf{p}^{\prime}\right)=\mathcal{R}_{Q\bar{Q}}V^{vec}\left(\mathbf{p}-\mathbf{p}^{\prime}\right)+V^{sca}(\mathbf{p}-\mathbf{p^{\prime}})$
$\displaystyle\quad+V^{LS}(\mathbf{p}-\mathbf{p^{\prime}})+V^{SS}(\mathbf{p}-\mathbf{p^{\prime}})+V^{T}(\mathbf{p}-\mathbf{p^{\prime}})\
,$ (16)

where the spin-dependent terms in momentum-space are obtained through Fourier
transform, $V^{a}(\mathbf{k}=\mathbf{p}-\mathbf{p^{\prime}})=\int
d^{3}\mathbf{r}e^{-i\mathbf{k}\cdot\mathbf{r}}V^{a}(\mathbf{r})$ with
$a=LS,SS,T$. We absorb $m_{ij}(\mathbf{k})$ in the two-body propagator into
the relativistic corrections for the potentials to keep the same convention as
in Ref. Liu and Rapp (2018), thus Eqs. (14) become

$\displaystyle
V_{ij}\left(\mathbf{p},\mathbf{p}^{\prime}\right)\rightarrow\sqrt{{m_{ij}(p)}}\sqrt{{m_{ij}(p^{\prime})}}V_{ij}\left(\mathbf{p},\mathbf{p}^{\prime}\right)$
$\displaystyle\qquad=\mathcal{R}_{ij}^{vec}V^{vec}\left(\mathbf{p}-\mathbf{p}^{\prime}\right)+\mathcal{R}_{ij}^{sca}V^{sca}\left(\mathbf{p}-\mathbf{p}^{\prime}\right),$
(17)

with

$\displaystyle\mathcal{R}_{ij}^{vec}$ $\displaystyle\equiv$
$\displaystyle\sqrt{{m_{ij}(p)}}\sqrt{{m_{ij}(p^{\prime})}}\mathcal{R}_{ij}$
$\displaystyle=$
$\displaystyle\sqrt{1+\frac{p^{2}}{\varepsilon_{i}(p)\varepsilon{j}(p)}}\sqrt{1+\frac{p^{\prime
2}}{\varepsilon_{i}(p^{\prime})\varepsilon_{j}(p^{\prime})}},$
$\displaystyle\mathcal{R}_{ij}^{sca}$ $\displaystyle\equiv$
$\displaystyle\sqrt{{m_{ij}(p)}}\sqrt{{m_{ij}(p^{\prime})}}$ (18)
$\displaystyle=$
$\displaystyle\sqrt{\frac{M_{i}M_{j}}{\varepsilon_{i}(p)\varepsilon_{j}(p)}}\sqrt{\frac{M_{i}M_{j}}{\varepsilon_{i}(p^{\prime})\varepsilon_{j}(p^{\prime})}}\
.$

and Eq. (16) becomes

$\displaystyle
V_{Q\bar{Q}}\left(\mathbf{p},\mathbf{p}^{\prime}\right)\rightarrow\sqrt{{m_{ij}(p)}}\sqrt{{m_{ij}(p^{\prime})}}V_{Q\bar{Q}}\left(\mathbf{p},\mathbf{p}^{\prime}\right)\qquad\qquad$
$\displaystyle=\mathcal{R}_{Q\bar{Q}}^{vec}V^{vec}\left(\mathbf{p}-\mathbf{p}^{\prime}\right)+\mathcal{R}_{Q\bar{Q}}^{sca}V^{sca}\left(\mathbf{p}-\mathbf{p}^{\prime}\right)\qquad$
$\displaystyle\qquad+\mathcal{R}_{Q\bar{Q}}^{spin}[V^{LS}(\mathbf{p}-\mathbf{p^{\prime}})+V^{SS}(\mathbf{p}-\mathbf{p^{\prime}})+V^{T}(\mathbf{p}-\mathbf{p^{\prime}})]$
(19)

with

$\displaystyle\mathcal{R}_{ij}^{spin}$ $\displaystyle\equiv$
$\displaystyle\sqrt{{m_{ij}(p)}}\sqrt{{m_{ij}(p^{\prime})}}$ (20)
$\displaystyle=$
$\displaystyle\sqrt{\frac{M_{i}M_{j}}{\varepsilon_{i}(p)\varepsilon_{j}(p)}}\sqrt{\frac{M_{i}M_{j}}{\varepsilon_{i}(p^{\prime})\varepsilon_{j}(p^{\prime})}}\
.$

The Lorentz structure for Coulomb potential is entirely vector, and a common
assumption for the confining one is to be entirely scalar, i.e.,
$V^{vec}=V_{\mathcal{C}}$ and $V^{sca}=V_{\mathcal{S}}$. As was mentioned in
the introduction, there are reasons to believe that the confining potential is
not a purely scalar one but a mixture of vector and scalar Lorentz structures,
i.e., $V^{vec}=V_{\mathcal{C}}+(1-\chi)V_{\mathcal{S}}$ and $V^{sca}=\chi
V_{\mathcal{S}}$. The key parameter is the mixing coefficient, $\chi$, defined
such that for $\chi=1$ the interaction reduces to the case with a purely
scalar confining potential, while values below one characterize a vector
admixture.

The potentials in the various color channels,
$V_{ij}^{a}\left(\mathbf{p},\mathbf{{p}^{\prime}}\right)$, are obtained by the
substitutions
$V_{\mathcal{C}}\left(\mathbf{p}-\mathbf{p}^{\prime}\right)\rightarrow\mathcal{F}_{a}^{\mathcal{C}}V_{\mathcal{C}}\left(\mathbf{p}-\mathbf{p}^{\prime}\right)$
and
$V_{\mathcal{S}}\left(\mathbf{p}-\mathbf{p}^{\prime}\right)\rightarrow\mathcal{F}_{a}^{\mathcal{S}}V_{\mathcal{S}}\left(\mathbf{p}-\mathbf{p}^{\prime}\right)$.
For the Coulomb interaction, $\mathcal{F}_{a}^{\mathcal{C}}$ are the standard
Casimir coefficients listed in Tab. 1 (together with the pertinent degeneracy
factors), and we take the absolute values of the Casimir coefficients for the
string interaction, $\mathcal{F}_{a}^{\mathcal{S}}$, to ensure a positive
definite string tension Liu and Rapp (2018).

Table 1: Casimir and degeneracy factors for different color channels (Casimir factor, degeneracy). $qq$ | $q\bar{q}$ | $(q/\bar{q})g$ | $gg$
---|---|---|---
$(1/2,3)$ | $(1,1)$ | $(9/8,3)$ | $(9/4,1)$
$(-1/4,6)$ | $(-1/8,8)$ | $(3/8,6)$ | $(9/8,16)$
|  | $(-3/8,15)$ | $(-3/4,27)$

The parton masses are also related to the potential introduced above. The
constituent masses of the heavy quarks, $M_{Q}$, receives two contributions,
the first one is calculated by the selfenergy from the color-singlet ($a=1$)
potential (including the relativistic factors) and the second one is a “bare
mass”, $M^{0}_{Q}$, which is associated with condensate contributions that we
do not calculate explicitly in the present framework,

$\displaystyle
M_{Q}=-\frac{1}{2}\int\frac{d^{3}\mathbf{p}}{(2\pi)^{3}}V_{Q\bar{Q}}^{a=1}(\mathbf{p})+M^{0}_{Q},$
(21)

##  IV Heavy-Quarkonium Spectroscopy in Vacuum

In this section we introduce the correlation and spectral functions including
their non-relativistic classifications in angular momentum (Sec. IV.1) and
discuss our fits to the vacuum spectra including the spin-related interactions
(Sec. IV.2).

###  IV.1 Correlation and Spectral Functions

To evaluate the quarkonia spectra in both charm and bottom sectors, we compute
the quark-antiquark spectral functions for different mesonic quantum-numbers
channels using the pertinent $T$-matrices as described in Sec. II. The bound-
state masses are then determined from the peak values of corresponding mesonic
spectral functions. In the vacuum, we introduce a small width in the single-
quark propagators which allows us to numerically resolve the bound-state mass
while not distorting their masses. Since we only account for the $Q\bar{Q}$
channels (off-shell) couplings to intermediate two-meson states (e.g.,
$DD^{*}$ channels) are not accounted for, which could affect the masses near
the $D\bar{D}$ threshold somewhat.

Table 2 lists the $L$, $S$ and $J$ assignments in the scalar (S), pseudoscalar
(PS), vector (V), axial-vector (AV) and tensor (T) mesonic channels. In
practice, a cutoff $r_{c}=0.01$ fm is introduced in the Fourier transform for
the spin-dependent potentials to avoid ultraviolet divergences; we have
checked that the results are not sensitive to variations in $r_{c}$ by $\pm
50\%$.

Table 2: Non-relativistic classification of angular-momentum quantum numbers in different mesonic channels. Channels |  $L$ |  $S$ |  $J$
---|---|---|---
S |  $1$ |  $1$ |  $0$
PS |  $0$ |  $0$ |  $0$
V |  $0$ |  $1$ |  $1$
AV1 |  $1$ |  $0$ |  $1$
AV2 |  $1$ |  $1$ |  $1$
T |  $1$ |  $1$ |  $2$

Figure 3: Diagrammatic representation of the $Q\bar{Q}$ correlation function.
The dots denote meson current operators $\Gamma_{M}$.

The spectral functions are obtained form the correlation functions, $G$, of
the meson currents. The latter are obtained by closing the two incoming and
outgoing legs of $T$-matrix (plus a non-interacting contribution) with the
corresponding projection operator for the different quarkonium channels, see
Fig. 3. Writing

$G=G_{0}+\Delta G\ ,$ (22)

the non-interacting part of correlation function in the CM frame is given by

$\displaystyle\ G_{0}(E,T)$ $\displaystyle=$ $\displaystyle
N_{f}N_{c}\int\frac{d^{3}p}{(2\pi)^{3}}\mathcal{R}_{Q\bar{Q}}^{sca}$
$\displaystyle\times\mathrm{Tr}\left\\{\Gamma_{\alpha}\Lambda_{+}(\mathbf{p})\Gamma_{\alpha}\Lambda_{-}(-\mathbf{p})\right\\}G_{Q\bar{Q}}^{0}(E,p),$

where
$\Lambda_{\pm}(\mathbf{p})=[\varepsilon_{Q}(p)\gamma^{0}-(\mathbf{p}\cdot\bm{\gamma})\pm
M_{Q}]/2M_{Q}$ are the positive/negative energy projectors for quark and
antiquark, respectively, and
$\Gamma_{\alpha}\in(1,i\gamma_{5},\gamma^{\mu},\gamma^{\mu}\gamma_{5},\frac{i}{2}[\gamma^{\mu},\gamma^{\nu}])$
the vertex operators for the mesonic S, PS, V, AV and T channels,
respectively; the pertinent traces are listed in Tab. 3. Furthermore,
$N_{f}=1$ and $N_{c}=3$ are the numbers of flavor and color for the heavy
quark.

Table 3: Values for the trace coefficients in different mesonic channels for the non-interacting part, $G_{0}$, of the correlation functions. $\Gamma_{\alpha}$ |  $G_{0}$ trace
---|---
S |  $2\frac{p^{2}}{M_{Q}^{2}}$
PS |  $2(1+\frac{p^{2}}{M_{Q}^{2}})$
V |  $6+4\frac{p^{2}}{M_{Q}^{2}}$
AV |  $4\frac{p^{2}}{M_{Q}^{2}}$
T |  $4\frac{p^{2}}{M_{Q}^{2}}$

The interacting part of correlation function in the CM frame is given by

$\displaystyle\Delta G(E,T)$ $\displaystyle=$
$\displaystyle\frac{N_{f}N_{c}}{8\pi^{4}}\int dpp^{2}G_{Q\bar{Q}}^{0}(E,p)\int
dp^{\prime}p^{\prime 2}G_{Q\bar{Q}}^{0}(E,p^{\prime})$ (24)
$\displaystyle\times\mathcal{R}_{Q\bar{Q}}^{sca}\mathcal{T}(\Gamma_{\alpha};E,p,p^{\prime}),$

with the scattering amplitude

$\displaystyle\begin{aligned}
\mathcal{T}(\Gamma_{\alpha};E,p,p^{\prime})&=\int
d(\cos\theta)\textup{Tr}(\Gamma_{\alpha};p,p^{\prime},\theta)T_{Q\bar{Q}}(E,\mathbf{p},\mathbf{{p}^{\prime}})\\\
&=8\pi[a_{0}(p,p^{\prime})T_{Q\bar{Q}}^{0}+a_{1}(p,p^{\prime})T_{Q\bar{Q}}^{1}]\
.\end{aligned}$ (25)

The $T$-matrix, $T_{Q\bar{Q}}^{L}$, is calculated from Eq. (11) with the
interaction kernel in Eq. (19). The $a_{L}$ denote the coefficients of the
orbital-angular momentum, $L$, in the partial-wave expansion of the trace,

$\displaystyle\textup{Tr}(\Gamma_{\alpha};p,p^{\prime},\theta)$
$\displaystyle=$
$\displaystyle\textup{Tr}\left\\{\Lambda_{+}(\mathbf{p})\Gamma_{\alpha}\Lambda_{-}(\mathbf{-p})\Lambda_{-}(\mathbf{-p^{\prime}})\right.$
$\displaystyle\left.\times\Gamma_{\alpha}\Lambda_{+}(\mathbf{p^{\prime}})\right\\}$
$\displaystyle=$ $\displaystyle
a_{0}(p,p^{\prime})P_{0}(\cos\theta)+a_{1}(p,p^{\prime})P_{1}(\cos\theta)\ ;$

they are listed in Tab. 4. For the evaluation of the traces in the V, AV and T
channels we focus on the spatial components. The correlation functions defined
in Refs. Cabrera and Rapp (2007) and Riek and Rapp (2010) do not have the
$\mathcal{R}_{Q\bar{Q}}^{sca}$ factor due to the different definitions of
potentials and two-body propagators, but they are equivalent to those in Refs.
Liu and Rapp (2018) and in this work.

Table 4: Coefficients of orbital angular momentum (up to $L=1$) in different mesonic channels for the interacting part, $G$, of the correlation functions. $\Gamma_{\alpha}$ | $a_{0}(p,p^{\prime})$ | $a_{1}(p,p^{\prime})$
---|---|---
S | $-\frac{p^{2}p^{\prime 2}}{M_{Q}^{4}}$ | $(1+\frac{\varepsilon(p)\varepsilon_{Q}(p^{\prime})}{M_{Q}^{2}})\frac{pp^{\prime}}{M_{Q}^{2}}$
PS |  | $1+\frac{p^{2}+p^{\prime 2}}{M_{Q}^{2}}$
---
$+\frac{\varepsilon_{Q}(p)\varepsilon_{Q}(p^{\prime})}{M_{Q}^{2}}+\frac{p^{2}p^{\prime
2}}{M_{Q}^{4}}$
$-\frac{\varepsilon_{Q}(p)\varepsilon_{Q}(p^{\prime})pp^{\prime}}{M_{Q}^{4}}$
V |  | $3(1+\frac{\varepsilon_{Q}(p)\varepsilon_{Q}(p^{\prime})}{M_{Q}^{2}})$
---
$+2\frac{p^{2}+p^{\prime 2}}{M_{Q}^{2}}+\frac{4}{3}\frac{p^{2}p^{\prime
2}}{M_{Q}^{4}}$
$-(1+2\frac{\varepsilon_{Q}(p)\varepsilon_{Q}(p^{\prime})}{M_{Q}^{2}})\frac{pp^{\prime}}{M_{Q}^{2}}$
AV | $-\frac{4}{3}\frac{p^{2}p^{\prime 2}}{M_{Q}^{4}}$ | $2(1+\frac{\varepsilon_{Q}(p)\varepsilon_{Q}(p^{\prime})}{M_{Q}^{2}})\frac{pp^{\prime}}{M_{Q}^{2}}$
T | $-\frac{2}{3}\frac{p^{2}p^{\prime 2}}{M_{Q}^{4}}$ | $2(1+\frac{\varepsilon_{Q}(p)\varepsilon_{Q}(p^{\prime})}{M_{Q}^{2}})\frac{pp^{\prime}}{M_{Q}^{2}}$

At higher orders in the $1/M_{Q}$ expansion, the partial-wave expansion leads
to a mixing between $S$\- and $P$\- wave components in a given meson channel.
However, to keep with the (non-relativistic) classification of the meson
channels with definite quantum numbers of $L$, $S$ and $J$, we terminate the
expansion when the “unnatural” partial waves admix Cabrera and Rapp (2007).
The pertinent leading orders are collected in Tab. 5.

Table 5: Leading orders of the angular-momentum coefficients, $a_{L}$, in evaluating the traces for the non-interacting (second column) and interacting (third and fourth columns) part of the different mesonic channels specified in the first column. $\Gamma_{\alpha}$ |  $G_{0}$ trace | $a_{0}(p,p^{\prime})$ | $a_{1}(p,p^{\prime})$
---|---|---|---
S | $2\frac{p^{2}}{M_{Q}^{2}}$ | $\mathcal{O}(p^{4}/M_{Q}^{4})$ | $-2\frac{pp^{\prime}}{M_{Q}^{2}}+\mathcal{O}(p^{4}/M_{Q}^{4})$
PS | $2+\mathcal{O}(p^{2}/M_{Q}^{2})$ | $-2+\mathcal{O}(p^{2}/M_{Q}^{2})$ | $\mathcal{O}(p^{2}/M_{Q}^{2})$
V | $6+\mathcal{O}(p^{2}/M_{Q}^{2})$ | $-6+\mathcal{O}(p^{2}/M_{Q}^{2})$ | $\mathcal{O}(p^{2}/M_{Q}^{2})$
AV | $4\frac{p^{2}}{M_{Q}^{2}}$ | $\mathcal{O}(p^{4}/M_{Q}^{4})$ | $-4\frac{pp^{\prime}}{M_{Q}^{2}}+\mathcal{O}(p^{4}/M_{Q}^{4})$
T | $4\frac{p^{2}}{M_{Q}^{2}}$ | $\mathcal{O}(p^{4}/M_{Q}^{4})$ | $-4\frac{pp^{\prime}}{M_{Q}^{2}}+\mathcal{O}(p^{4}/M_{Q}^{4})$

From the correlation functions, the mesonic spectral functions follow from the
imaginary part,

$\sigma_{\alpha}(E,T)=-\frac{1}{\pi}\operatorname{Im}G_{\alpha}(E+i\epsilon,T),$
(27)

where the subscript $\alpha$ denotes the different meson channels.

###  IV.2 Heavy-Quarkonium Spectra in Vacuum

We are now in position to investigate the quantitative consequences of the
spin-dependent interactions and the scalar-vector mixing effect in the
confining potential on the charmonium and bottomonium spectroscopy in vacuum.
In practice, we adopt a value for the HQ width of $\Gamma_{Q}=20$ MeV which is
small enough to not affect the vacuum masses, but large enough to allow for
straightforward numerical computations and plotting. Our fit procedure is as
follows: For a given mixing coefficient, $\chi$, we adjust the bare-quark
masses $M_{c}^{0}$ ($M_{b}^{0}$) to find the best fit for all the masses of
charmonium (bottomonium) states given by the Particle Data Group (PDG)
Tanabashi et al. (2018) using a $\chi^{2}$ statistical test. In principle we
could optimize the value for the mixing coefficient, $\chi$, by strictly
minimizing the variance of the fit. However, in practice we found that
$\chi=$0.6 already provides most of the improvement in the mass splittings
compared to $\chi=1$, while for still smaller values the constituent quark
masses become rather large producing uncomfortably large binding energies; in
addition, a more precise evaluation of the quarkonium masses near the open
heavy-flavor threshold would also require the inclusion of hadronic loop
corrections. In particular, we do not pursue here more extreme scenarios, as
proposed, e.g., in Refs. Ebert et al. (1998, 2003) where optimized fits with
$\chi=-1$ were found (implying a negative string tension for the scalar term,
counteracted by a twice larger vector component).

We start by displaying the comparison of charmonium spectral functions between
purely scalar ($\chi=1$) and mixed ($\chi=0.6$) confining potential for all
states below the $D\bar{D}$ threshold in Fig. 4 (we only plot the interacting
parts of spectral functions since the free part does not affect the bound-
state locations). The various peaks in each quantum-number channel are readily
assigned as, S: $\chi_{c0}(1P)$; PS: $\eta_{c}(1S)$ and $\eta_{c}(2S)$; V:
$J/\Psi(1S)$ and $\Psi(2S)$; AV1: $h_{c}(1P)$; AV2: $\chi_{c1}(1P)$; T:
$\chi_{c2}(1P)$. The masses extracted from the pole positions are listed in
Tab. 6 and compared to the experimental values. The potential with mixing
effect generates more attraction from the additional relativistic corrections
(recall Eq. (13)), thus requiring a larger HQ massfor $\chi=0.6$ as quoted in
Tab. 6. For $\chi=1$, neither the $S$\- nor $P$-wave mass splittings are well
reproduced; the results are much improved by introducing the mixing effect
with $\chi=0.6$.

Figure 4: The vacuum charmonium spectral functions (interacting parts) in S, PS, V, AV and T channels with mixing coefficient $\chi=1$ (upper panel) and $\chi=0.6$ (lower panel). Table 6:  Experimental and theoretical values for charmonium masses (in GeV) with the mixing coefficient $\chi=1$ and $\chi=0.6$; $M_{c}^{0}$ and $M_{c}$ are the bare and constituent charm-quark mass, respectively. Channel | Particle | Exp. |  | Th.
---
$\chi=1$
$M_{c}^{0}=1.352$
$M_{c}=1.872$
| Th.
---
$\chi=0.6$
$M_{c}^{0}=1.359$
$M_{c}=1.916$
S | $\chi_{c0}(1P)$ | 3.415 | 3.498 | 3.448
PS | $\eta_{c}(1S)$ | 2.984 | 3.079 | 3.022
$\eta_{c}(2S)$ | 3.638 | 3.624 | 3.600
V | $J/\Psi(1S)$ | 3.097 | 3.120 | 3.104
$\Psi(2S)$ | 3.686 | 3.650 | 3.650
AV1 | $h_{c}(1P)$ | 3.525 | 3.518 | 3.500
AV2 | $\chi_{c1}(1P)$ | 3.511 | 3.519 | 3.499
T | $\chi_{c2}(1P)$ | 3.556 | 3.519 | 3.544

To explain better understand the impact of the mixing effect, it is useful to
summarize the expectation values for $\left\langle\mathbf{L\cdot
S}\right\rangle$, $\left\langle\mathbf{S_{1}\cdot S_{2}}\right\rangle$ and
$\left\langle S_{12}\right\rangle$ in Tab. 7. We take the mass splitting
between PS and V channels, where only spin-spin interaction are operative, as
an example. For simplicity, we will use the Cornell potential in Eq. (12) to
make the argument. According to Eq. (LABEL:eq_VLSSST), the spin-spin
interaction is ${V^{SS}}(r)\sim\left\langle\mathbf{S_{1}\cdot
S_{2}}\right\rangle\Delta V^{vec}(r)=\left\langle\mathbf{S_{1}\cdot
S_{2}}\right\rangle[\frac{16\pi\alpha_{s}}{3}\delta(r)+(1-\chi)\frac{2\sigma}{r}]$.
Since the quantity in the bracket is positive, the negative (positive) value
for $\left\langle\mathbf{S_{1}\cdot S_{2}}\right\rangle$ gives a more
attractive (repulsive) interactions, which is the origin of the
$J/\Psi$-$\eta_{c}$ (i.e., V-PS) splitting. The mixing with $\chi<1$ obviously
enhances this effect. Similar arguments can be made for the spin-orbit,
$V^{LS}(r)\sim\left\langle\mathbf{L\cdot
S}\right\rangle\left[\frac{4\alpha_{s}}{r^{3}}+(3-4\chi)\frac{\sigma}{r}\right]$,
and tensor, $V^{LS}(r)\sim\left\langle
S_{12}\right\rangle\left[\frac{4\alpha_{s}}{r^{3}}+(1-\chi)\frac{\sigma}{r}\right]$,
interactions; i.e., the strength of $V^{LS}$, $V^{SS}$ and $V^{T}$ are all
enhanced by introducing a vector component in the confining potential, thereby
improving the splittings in comparison to experiment.

Table 7: Couplings of $LS$, $SS$ and $T$ interactions for different quarkonium channels. Channel | $\left\langle\mathbf{L\cdot S}\right\rangle$ | $\left\langle\mathbf{S_{1}\cdot S_{2}}\right\rangle$ | $\left\langle S_{12}\right\rangle$
---|---|---|---
S | $-2$ | $1/4$ | $-4$
PS | $0$ | $-3/4$ | $0$
V | $0$ | $1/4$ | $0$
AV1 | $0$ | $-3/4$ | $0$
AV2 | $-1$ | $1/4$ | $2$
T | $1$ | $1/4$ | $-2/5$

Figure 5: The vacuum bottomonium spectral functions (interacting parts) in S,
PS, V, AV and T channels with mixing coefficient $\chi=1$ (upper panel) and
$\chi=0.6$ (lower panel). The spectral functions in S, AV and T channels are
multiplied by a factor of 4 for better visibility.

We have carried out a similar analysis for bottomonium spectral functions. In
Fig. 5 bound-state spectral functions between $\chi=1$ and 0.6 are compared.
We identify the peaks in each channel as follows (not all of which have an
experimental counterpart (yet)): S: $\chi_{b0}(1P)$ and $\chi_{b0}(2P)$; PS:
$\eta_{b}(1S)$, $\eta_{b}(2S)$ and $\eta_{b}(3S)$; V: $\Upsilon(1S)$,
$\Upsilon(2S)$ and $\Upsilon(3S)$; AV1: $h_{b}(1P)$, $h_{b}(2P)$ and
$h_{b}(3P)$; AV2: $\chi_{b1}(1P)$, $\chi_{b1}(2P)$ and $\chi_{b1}(3P)$; T:
$\chi_{b2}(1P)$, $\chi_{b2}(2P)$ and $\chi_{b2}(3P)$. The comparison between
the experimental values and the masses extracted from spectral functions is
compiled in Tab. 8. Also here an improvement in the mass splittings is found
by introducing the vector confining potential, but it is not as significant as
the charmonium sector, primarily due to the larger $1/M_{b}$ suppression for
the spin-induced forces.

Table 8: Experimental and theoretical bottomonium spectroscopy with the mixing coefficient $\chi=1$ and $\chi=0.6$ (in GeV). $M_{b}^{0}$ and $M_{b}$ are the bare and constituent masses for bottom quark, respectively. Channel | Particle | Exp. |  | Th.
---
$\chi=1$
$M_{b}^{0}=4.681$
$M_{b}=5.247$
| Th.
---
$\chi=0.6$
$M_{b}^{0}=4.681$
$M_{b}=5.266$
S | $\chi_{b0}(1P)$ | 9.859 | 9.871 | 9.864
$\chi_{b0}(2P)$ | 10.233 | 10.227 | 10.220
PS | $\eta_{b}(1S)$ | 9.399 | 9.496 | 9.470
V | $\Upsilon(1S)$ | 9.460 | 9.520 | 9.500
$\Upsilon(2S)$ | 10.023 | 9.999 | 9.994
$\Upsilon(3S)$ | 10.355 | 10.345 | 10.324
AV1 | $h_{b}(1P)$ | 9.899 | 9.896 | 9.893
AV2 | $\chi_{b1}(1P)$ | 9.893 | 9.894 | 9.877
$\chi_{b1}(2P)$ | 10.255 | 10.248 | 10.243
$\chi_{b1}(3P)$ | 10.513 | 10.520 | 10.500
T | $\chi_{b2}(1P)$ | 9.912 | 9.897 | 9.899
$\chi_{b2}(2P)$ | 10.269 | 10.249 | 10.249

Finally, we have evaluated the spin-induced interactions in the heavy-light
sector, which is the key ingredient to calculating the heavy-quark transport
coefficients discussed in Sec. VI. Specifically, in the $S$-wave color-singlet
$D$-meson channel, the mass splitting between the pseudoscalar $D$-meson and
the vector $D^{*}$-meson improves from 30 MeV for $\chi$=1 to 120 MeV for
$\chi$=0.6.

##  V In-Medium Potential and Selfconsistent QGP

In this section, we briefly introduce (and carry out) the framework for
determining the medium modifications to the potential and its application to
the EoS and spectral functions of the QGP within a selfconsistent quantum
many-body approach Liu and Rapp (2018). In a nutshell the procedure consists
of 2 selfconsistency loops as follows. First, the in-medium potential will be
constrained through calculating the HQ free energies from the $T$-matrix and
fitting it to pertinent lQCD data. The key fit parameters in this step are the
screening masses, $m_{d}$ and $m_{s}$ of the color-Coulomb and string
interactions. The in-medium potentials are then applied in a selfconsistent
2-PI scheme to compute the EoS of the QGP and fit those results to pertinent
lQCD data as well. The main parameters in this step are the in-medium light-
quark and gluon masses, but the EoS is computed including the full off-shell
properties of the one-body spectral functions and two-body scattering
amplitudes. Since the parton selfenergies are computed from their
$T$-matrices, this forms a selfconsistency problem which is solved by
numerical iteration. However, the calculation of the HQ free energy also
requires the spectral functions (selfenergies) of the heavy quarks, calculated
from heavy-light $T$-matrices closed off with thermal parton spectral
functions. Thus, after constraining the light sector with the EoS, the in-
medium HQ spectral functions are re-calculated and inserted into the
computation of the HQ free energies. Refitting the screening masses to the
lQCD data, a refinement of the in-medium two-body potential is obtained which
is then reprocessed in a new fit to the EoS. This constitutes the second
(“outer”) iteration loop which is also iterated numerically.

In the remainder of this section, we first introduce the the main equations to
compute the HQ free energy (Sec. V.1) and the EoS (Sec. V.2), and then discuss
the numerical results with the updated in-medium potential (Sec. V.3).

###  V.1 Static HQ Free Energy

Our starting point is an ansatz for the medium modifications of the potential;
following previous studies Liu and Rapp (2018) we employ

$\displaystyle{\widetilde{V}}_{\mathcal{C}}(r)$ $\displaystyle=$
$\displaystyle-\frac{4}{3}\alpha_{s}\frac{e^{-m_{d}r}}{r}-\frac{4}{3}\alpha_{s}m_{d}$
$\displaystyle{\widetilde{V}}_{\mathcal{S}}(r)$ $\displaystyle=$
$\displaystyle-\frac{\sigma
e^{-m_{s}r-\left(c_{b}m_{s}r\right)^{2}}}{m_{s}}+\frac{\sigma}{m_{s}}\ ,$ (28)

where $m_{d}$ and $m_{s}$ are the respective Debye screening masses for
Coulomb and confining potentials, related by
$m_{s}=\left(c_{s}m_{d}^{2}\sigma/\alpha_{s}\right)^{1/4}$ Liu and Rapp
(2018). The quadratic term in the exponential,
$-\left(c_{b}m_{s}r\right)^{2}$, accelerates the suppression of the long-range
part of the confining potential to simulate string breaking. In the limit of
vanishing screening masses the vacuum potential of Eq. (12) is recovered.

The HQ free energy, $F_{Q\bar{Q}}(r,T)$, is defined as the difference between
the free energies of the QGP without and with a static quark and antiquark
(not counting their infinite masses) separated by a distance $r$ (see, e.g.,
Ref. Beraudo et al. (2008)),

$F_{Q\bar{Q}}(r,T)=-\frac{1}{\beta}\ln\left[G_{Q\bar{Q}}^{>}(-i\beta,r)\right],$
(29)

where $G_{Q\bar{Q}}^{>}(-i\tau,r)$ is the Euclidean time Green function and
$\beta=1/T$ the inverse temperature. In the vacuum, this simply corresponds to
the potential between $Q$ and $\bar{Q}$, cf. Sec. III. In medium, the free
energy and the potential are no longer identical to each other due to the
presence of entropy contributions resulting from medium effects encoded in the
HQ selfenergies (which we calculate from the in-medium heavy-light $T$-matrix)
and the potential. In Ref. Liu and Rapp (2018) a compact form of the free
energy has been derived as

$\displaystyle F_{Q\bar{Q}}(r,T)=$
$\displaystyle-\frac{1}{\beta}\ln\left[-\int_{-\infty}^{\infty}\frac{dE}{\pi}e^{-\beta
E}\right.$ (30)
$\displaystyle\left.\times\operatorname{Im}\left[\frac{1}{E+i\epsilon-\widetilde{V}(r)-\Sigma_{Q\bar{Q}}(E+i\epsilon)}\right]\right],$

with the color-singlet potential $\widetilde{V}(r)$ (color-flavor indices are
suppressed for simplicity) from Eq. (28). The relationship between the two-
body selfenergy, $\Sigma_{Q\bar{Q}}(z)$, and the two-body propagator,
$G_{Q\bar{Q}}^{0}(z)$, is Liu and Rapp (2018)

$\left[G_{Q\bar{Q}}^{0}(z)\right]^{-1}=z-2\Delta M_{Q}-\Sigma_{Q\bar{Q}}(z)\
,$ (31)

with a Fock mass term for each quark, $\Delta M_{Q}=\widetilde{V}(\infty)/2$.
In the static limit, $G_{Q\bar{Q}}^{0}(z)$ reduces to

$G_{Q\bar{Q}}^{0}(z)=\int_{-\infty}^{\infty}d\omega_{1}d\omega_{2}\frac{\rho_{Q}\left(\omega_{1}\right)\rho_{\bar{Q}}\left(\omega_{2}\right)}{z-\omega_{1}-\omega_{2}}\
,$ (32)

where
$\rho_{Q/\bar{Q}}(\omega)=-\frac{1}{\pi}\operatorname{Im}G_{Q/\bar{Q}}(\omega+i\epsilon)$
are the single-particle spectral functions with propagators
$G_{Q/\bar{Q}}(z)=\frac{1}{z-M_{Q/\bar{Q}}-\Sigma_{Q/\bar{Q}}(z)}$ in the
static limit. Then the single-particle selfenergy, $\Sigma_{Q}(z)$, can be
solved selfconsistently by combining the $T$-matrix and the selfenergy
equations. By taking the heavy-light $T$-matrix from Eq. (2) in the “half-
static” limit, where the $\mathbf{p_{1}}$ dependence is suppressed due to the
infinite static-quark mass, Eq. (7) takes the form

$\displaystyle\Sigma_{Q}(z)=$
$\displaystyle\int\frac{d^{3}\mathbf{p}_{2}}{(2\pi)^{3}}\int_{-\infty}^{\infty}d\omega_{2}\frac{dE}{\pi}\frac{-1}{z+\omega_{2}-E}\frac{1}{d_{Q}}\sum_{a,j}d_{s}^{Qj}d_{a}^{Qj}$
(33) $\displaystyle\times
T_{Qj}^{a}\left(E,\mathbf{p}_{2}\mid\mathbf{p}_{2}\right)\rho_{j}\left(\omega_{2},\mathbf{p}_{2}\right)n_{j}(\omega_{2})\
.$

The CM transformation in Eq. (9) reduces to

$E_{\mathrm{cm}}=\omega_{1}+\omega_{2},\quad
p_{\mathrm{cm}}=p_{2},\quad\cos\left(\theta_{\mathrm{cm}}\right)=\cos(\theta),$
(34)

with $\omega_{1}+\omega_{2}\gg|\mathbf{p}_{1}+\mathbf{p}_{2}|$. The resulting
single-particle selfenergy, $\Sigma_{Q}(z)$, is inserted into Eq. (32) to
obtain the $Q\bar{Q}$ propagator, and Eq. (31) yields the two-body selfenergy,
$\Sigma_{Q\bar{Q}}(z)$.

Interference effects lead to a suppression of the imaginary part of the two-
body selfenergy (relative to the sum of the single-particle absorptive parts),
sometimes referred to as “imaginary part of the potential” (it is, in fact, a
suppression of the latter) Laine et al. (2007). In the $T$-matrix formalism
this amounts to a 3-body contribution which is rather challenging to compute
explicitly Liu and Rapp (2018). Instead, an $r$-dependent suppression factor
of a functional form guided by perturbative results has been implemented Liu
and Rapp (2018) with a factorizable ansatz,
$\Sigma_{Q\bar{Q}}(z,r)=\Sigma_{Q\bar{Q}}(z)\phi(r)$, where the function
$\phi(r)$ will be part of the constraints from the lQCD data for static HQ
free energies at each temperature. The interference effect is mostly relevant
for deeply bound heavy quarkonia, where, in the color singlet channel, the
imaginary part should vanish in the limit of $r\to 0$ (corresponding to a
color-neutral object).

###  V.2 Equation of State

The equation of state of a many-body system is encoded in the pressure,
$P(T,\mu)$, as a function of temperature and chemical potential. The EoS is
driven by the dominant degrees of freedom in the medium, and is therefore
sensitive to their spectral properties, including their masses. For a
homogeneous grand canonical ensemble, the relationship between the EoS and the
grand potential per unit volume is given by $\Omega=-P$. We adopt the
Luttinger-Ward-Baym (LWB) formalism which provides a diagrammatic and
thermodynamically consistent quantum approach that allows to incorporate the
off-shell dynamics of the one- and two-body correlation functions. Quantum
effects are expected to become particularly important for a strongly coupled
system with large scattering rates (widths) Luttinger and Ward (1960); Baym
and Kadanoff (1961); Baym (1962). One has

$\Omega=\mp\frac{-1}{\beta}\sum_{n}\operatorname{Tr}\left\\{\ln\left(-G^{-1}\right)+\left[\left(G^{0}\right)^{-1}-G^{-1}\right]G\right\\}\pm\Phi,$
(35)

where “$\operatorname{Tr}$” denotes the trace over spin, color, flavor and
3-momentum, $\sum_{n}$ the Matsubara frequency sum, and the $G^{0}$ and $G$
are the free and fully dressed single-particle Green’s function. The two-body
interaction contribution is encoded in the Luttinger-Ward functional (LWF),
$\Phi=\sum_{v=1}^{\infty}\Phi_{v}$, where

$\Phi_{v}=\frac{-1}{\beta}\sum_{n}\operatorname{Tr}\left\\{\frac{1}{2v}\left(\frac{-1}{\beta}\right)^{v}\left[(-\beta)^{v}\Sigma_{v}(G)\right]G\right\\}$
(36)

with

$\Sigma_{v}(G)=\int d\tilde{p}\left[VG_{(2)}^{0}VG_{(2)}^{0}\cdots V\right]G\
,$ (37)

using the notation $\int d\tilde{p}\equiv-\beta^{-1}\sum_{n}\int
d^{3}\mathbf{p}/(2\pi)^{3}$ with
$\tilde{p}\equiv\left(i\omega_{n},\mathbf{p}\right)$. The $\phi_{\nu}$
correspond to the “skeleton diagrams” of the $\nu^{\rm th}$ order in the
potential expansion. To account for possible bound-states formation and their
contribution to the pressure, one has to resum the skeleton series. For non-
separable interactions this has been achieved through a matrix-logarithm
resummation technique in Refs. Liu and Rapp (2018, 2017, 2020a), resulting in
a structure similar to the $T$-matrix resummation in Eq. (7):

$\displaystyle\Omega=$ $\displaystyle\sum_{j}\mp d_{j}\int
d\tilde{p}\left\\{\ln\left(-G_{j}(\tilde{p})^{-1}\right)\right.$ (38)
$\displaystyle\left.+\left[\Sigma_{j}(\tilde{p})-\frac{1}{2}\log\Sigma_{j}(\tilde{p})\right]G_{j}(\tilde{p})\right\\}$

with

$\displaystyle\log\Sigma_{i}\left(z,\mathbf{p}_{1}\right)$
$\displaystyle=\int\frac{d^{3}\mathbf{p}_{2}}{(2\pi)^{3}}\int_{-\infty}^{\infty}d\omega_{2}\frac{dE}{\pi}\frac{-1}{z+\omega_{2}-E}$
(39)
$\displaystyle\times\frac{1}{d_{i}}\sum_{a,j}d_{s}^{ij}d_{a}^{ij}\operatorname{Im}\left[\log
T_{ij}^{a}\left(E,\mathbf{p}_{1},\mathbf{p}_{2}\mid\mathbf{p}_{1},\mathbf{p}_{2}\right)\right]$
$\displaystyle\times\rho_{j}\left(\omega_{2},\mathbf{p}_{2}\right)\left[n_{j}\left(\omega_{2}\right)\mp
n_{ij}(E)\right].$

The transformation of the $T$-matrices in Eq. (39) between the thermal and CM
frame is given by Eq. (9). The grand potential can then be obtained after
carrying out the sum over Matsubara frequencies in Eq. (38).

###  V.3 Selfconsistent in-Medium Results

We now turn to the selfconsistent in-medium results at four temperatures,
$T=0.194$ GeV, 0.258 GeV, 0.320 GeV and 0.400 GeV, constrained by the lQCD
data for static HQ free energies (Sec. V.3.1) and QGP EoS (Sec. V.3.2). All
in-medium calculations are carried out with the mixing coefficient for
$\chi=0.6$ and 1 in this study, but we do not yet incorporate the spin-
dependent corrections. In particular in the light sector, i.e., for the QGP
EoS, their effect can be rather significant and deserves a separate study
(some compensatory effects are expected due to both attractive and repulsive
contributions). Nevertheless, we want to ensure that the medium within which
the heavy quarks are embedded satisfies basic constraints from lQCD.

####  V.3.1 Static HQ Free Energies

Recalling Eq. (30), the HQ free energy, $F_{Q\bar{Q}}(r,T)$, is a functional
of the potential, $\widetilde{V}(r)$, and the two-body selfenergy,
$\Sigma_{Q\bar{Q}}(E+i\epsilon)$. Note that $F_{Q\bar{Q}}(r,T)$ increases with
increasing $\widetilde{V}(r)$ but with decreasing $|\Sigma_{Q}(E+i\epsilon)|$.
A larger Debye screening mass, $m_{d}$ and/or $m_{s}$, suppresses
$\widetilde{V}(r)$ so that the partons become more weakly coupled, which in
turn lowers $F_{Q\bar{Q}}(r,T)$; at the same time, a larger $m_{d}$ reduces
$|\Sigma_{Q}(E+i\epsilon)|$ in medium, which in turn enhances
$F_{Q\bar{Q}}(r,T)$. It is this competition between $\widetilde{V}(r)$ and
$|\Sigma_{Q}(E+i\epsilon)|$ that leads to a non-monotonic behavior of
$F_{Q\bar{Q}}(r,T,m_{d})$ with $m_{d}$. Since $m_{d}$ is directly related to
the free energy at infinite distance (cf. Sec. V.1), we define a function
$F_{trial}(r\rightarrow\infty,T,m_{d})$ calculated from our many-body approach
which we require to be equal to the lQCD value,
$F_{lQCD}(r\rightarrow\infty,T)$. We typically find two solutions for $m_{d}$
for any fixed parameter set provided the maximum of the trial free energy lies
above the lQCD value. We denote the solutions with the smaller and the larger
$m_{d}$ as strongly coupled solution (SCS) and weakly coupled solution (WCS),
respectively (in analogy to the two solutions which were found in Ref. Liu and
Rapp (2018)). Here, we only focus on the SCS which results in transport
parameters in much better agreement with heavy-ion phenomenology (i.e., a
liquid-like behavior with interaction energies comparable to the parton
masses, as well as HQ transport parameters) than the WCS Liu and Rapp (2018).

The resulting potentials and fits to lQCD results Bazavov and Petreczky (2013)
for the HQ free energies with $c_{b}=1.55$ (1.3) and $c_{s}=0.06$ (0.01) for
$\chi=0.6$ (1) at different temperatures are shown in Fig. 6. As found in
earlier studies, large HQ widths lead to a substantial enhancement of the
potential over the free energies; in particular, at the lowest temperature of
0.194 GeV, the potential is close to the vacuum one, but becomes notably
suppressed at higher $T$.

Figure 6: Left and middle panels: The in-medium HQ free energies (blue) and
potentials (orange) for $\chi=0.6$ (solid) and 1 (dashed) at different
temperatures in comparison to lQCD data for the HQ free energies from Ref.
Bazavov and Petreczky (2013) for $N_{f}=2+1$ light flavors (black dots). The
right panel shows the temperature dependence of the screening masses, $m_{d}$
(blue) and $m_{s}$ (orange), for $\chi=0.6$ (solid) and 1 (dashed) as
resulting from our fit. The $\chi=1$ results are taken from Ref. Liu and Rapp
(2018).

Consequently, the screening masses for Coulomb ($m_{d}$) and confining
($m_{s}$) potentials, shown in the right panel of Fig. 6, have rather small
values at low $T$, with the string interaction exhibiting a weaker screening
with increasing $T$. This implies that remnants of the confining force survive
in the QGP well above the critical region. The potential with mixed confining
interaction ($\chi=0.6$) is enhanced by the extra relativistic corrections
(cf. Eq. (13)), requiring a stronger screening to fit the lQCD free-energy
data. Therefore, the screening masses for confining potential for $\chi=0.6$
are larger than those for $\chi=1$, see the right panel of Fig. 6. However,
note that the $\chi=0.6$ solution generates a stronger force at relatively
small distances, a feature that will figure importantly in the QGP structure
and HQ transport properties.

####  V.3.2 Equation of State

In Fig. 7 we display the pressure together with the fitted light-parton masses
for $\chi=0.6$ and 1, which allow for a good reproduction of the lQCD data.
However, while the parton masses are effective in achieving a quantitative
agreement with the lQCD results, the underlying quark and gluon spectral
functions for $\chi=0.6$ and 1 both feature large selfenergies, especially
imaginary parts which, at low momentum and temperatures, are comparable to, or
even larger, than the parton masses, cf. the spectral function widths in Fig.
8 (first and second rows). The large scattering rates are mostly driven by
dynamical resonance formation in the underlying $T$-matrices (which in turn
are generated by resumming the strong potential). These resonances contribute
through the resummed LWF functional $\Phi\sim 1/2\textup{log}\Sigma G$
introduced in Sec. V.2, whose contribution for $\chi=0.6$ and 1 is displayed
in Fig. 7. The increasing proportion of LWF contribution with decreasing
temperature indicates the onset of a change in the degrees of freedom.
Specifically, the LWF parts make up more than 70(50)% of the pressure at
$T=0.194$ GeV for $\chi=0.6(1)$. While the spectral functions for $\chi=0.6$
generally share the main features with those for $\chi=1$ at low momenta, a
notable quantitative difference is that the widths for $\chi=0.6$ do not fall
off with momentum as much as those for $\chi=1$. In the former case, this is a
consequence of the 3-momentum dependence of the confining interaction whose
vector component, through relativistic effects, generates additional
interaction strength and thus larger scattering rates at larger momenta
relative to the $\chi=1$ case.

Figure 7: The pressure (normalized by $T^{4}$) in comparison to the lQCD data
(black dots) from Ref. Bazavov et al. (2014b) (upper panel), and the in-medium
light-quark and gluon masses as a function of temperature (lower panel). The
$\chi=1$ results are taken from Ref. Liu and Rapp (2018).

Figure 8: Single-parton spectral functions for light quarks (first row),
gluons (second row) and charm quarks (third row) with $\chi=0.6$ (solid) and 1
(dashed) as a function of energy for various 3-momenta in each panel. From
left to right, the four columns correspond to temperatures of $T=194$, 258,
320 and 400 MeV, respectively. The $\chi=1$ results are taken from Ref. Liu
and Rapp (2018).

##  VI Charm-Quark Transport Coefficients

With the parameters of the interaction potential and parton masses fixed with
the aid of lQCD data, we can now investigate the effect of the mixed potential
on charm-quark transport properties in the QGP. As elaborated in Ref. Liu et
al. (2019) it is important to account for the off-shell properties of both
charm quarks and thermal partons in the evaluation of the transport
coefficient, especially due to the formation of near-threshold bound states
which only provide limited phase for quasiparticle (on-shell) scattering. This
point is further corroborated upon inspecting the equilibrium spectral
functions, $\rho_{q,g,c}$, of the partons displayed in Fig. 8, exhibiting
large widths of $\sim$0.6 GeV or so at low momentum. As already mentioned in
Sec. V.3.2, the main difference between $\chi=0.6$ and 1 is that the widths
for $\chi=0.6$ do not fall off with momentum as much as those for $\chi=1$.
This feature persists in the heavy-light scattering amplitudes, which are the
main ingredient to the charm-quark transport coefficients discussed below, see
Fig. 9. At $T=0.194$ GeV, the peak value of the imaginary part of the $S$-wave
color-singlet heavy-light scattering amplitude for $\chi=0.6$ still shows a
rather marked decrease with increasing center-of-mass mass momentum of the
colliding partons, but it is significantly weaker than for $\chi=1$ with a
purely scalar confining potential; e.g., the peak reduction from the
$p_{cm}=0$ to $p_{cm}=0.5$ GeV case is almost a factor 3 for the latter but
only $\sim$1.5 for $\chi=0.6$. This trend continues to higher temperatures; at
$T=0.400$ GeV, the peak reduction from $p_{cm}=0$ to $p_{cm}=0.6$ GeV is
essentially absent for $\chi=0.6$, while it is still a factor of 1.6 for the
purely scalar confining potential. The $T$-matrix amplitudes for $\chi=0.6$
are smaller than that for $\chi=1$ at low momenta due to its stronger
screening in confining potential (recall its larger Debye screening masses in
the right panel of Fig. 6); however, they exceed the ones for $\chi=1$ for
$p_{cm}$$\,\stackrel{{\scriptstyle\textstyle>}}{{\sim}}\,$0.5 GeV due to the
harder 3-momentum dependence of confining potential through relativistic
effects. The $T$-matrices for other partial waves and color channels share
similar features, and thus we do not reproduce therm here.

Figure 9: The imaginary part of the $S$-wave charm-light $T$-matrices in the
color-singlet channel at different temperatures. The $T$-matrix is displayed
for four different values of the CM momentum ($p_{cm}$) in each panel. The
$\chi=1$ results are taken from Ref. Liu and Rapp (2018).

Turning now to the HQ transport coefficients in the QGP, we adopt their
standard definition through a Fokker-Planck equation where they amount to the
first and second momentum of the momentum transfer of the heavy-light
scattering amplitude squared, integrated over the thermal-parton
distributions. At this level, off-shell effects can be readily implemented by
an additional energy convolution over the light-parton spectral functions.
However, since also charm quarks acquire widths which are not small, the
inclusion of their spectral width is also warranted. This has been worked out
in Ref. Liu et al. (2019) employing the Kadanoff-Baym equations, resulting in
the following expression for the HQ friction coefficient (or relaxation rate):

$\displaystyle A(p)=$
$\displaystyle\sum_{i}\frac{1}{2\varepsilon_{c}(p)}\int\frac{d\omega^{\prime}d^{3}\mathbf{p}^{\prime}}{(2\pi)^{3}2\varepsilon_{c}\left(p^{\prime}\right)}\frac{d\nu
d^{3}\mathbf{q}}{(2\pi)^{3}2\varepsilon_{i}(q)}\frac{d\nu^{\prime}d^{3}\mathbf{q}^{\prime}}{(2\pi)^{3}2\varepsilon_{i}\left(q^{\prime}\right)}$
(40)
$\displaystyle\times\delta^{(4)}\frac{(2\pi)^{4}}{d_{c}}\sum_{a,l,s}|M|^{2}\rho_{c}\left(\omega^{\prime},p^{\prime}\right)\rho_{i}(\nu,q)\rho_{i}\left(\nu^{\prime},q^{\prime}\right)$
$\displaystyle\times\left[1-n_{c}\left(\omega^{\prime}\right)\right]n_{i}(\nu)\left[1\pm
n_{i}\left(\nu^{\prime}\right)\right](1-\frac{\mathbf{p}\cdot\mathbf{p^{\prime}}}{\mathbf{p}^{2}})\
.$

As before (recall Sec. II) , $\varepsilon_{i/c}$, $\rho_{i/c}$ and $n_{i/c}$
are the dispersion relations, spectral and thermal-distribution functions for
partons $i/c$, respectively, $\delta^{(4)}$ is a short-hand notation for the
energy-momentum conserving $\delta$-function in the 2$\to$2 scattering
process, and $d_{c}=6$ the spin-color degeneracy of charm quarks. The
summation $\sum_{i}$ is over all light-flavor quarks and gluons,
$i=u,\bar{u},d,\bar{d},s,\bar{s},g$, where the masses for light and strange
quarks are assumed to be the same. In the above expression, the quasiparticle
approximation is only applied to the incoming charm quark by assigning it a
sharp energy $\varepsilon_{c}(\mathbf{p})$ at momentum $\mathbf{p}$, while all
other partons are treated via off-shell integrations. The heavy-light
scattering matrix elements, $|M|^{2}$, in Eq. (40) are related to the
$T$-matrix in the CM frame by

$\displaystyle\sum_{a,L,s}\left|M^{2}\right|=16\varepsilon_{c}\left(p_{\mathrm{cm}}\right)\varepsilon_{i}\left(p_{\mathrm{cm}}\right)\varepsilon_{c}\left(p_{\mathrm{cm}}^{\prime}\right)\varepsilon_{i}\left(p_{\mathrm{cm}}^{\prime}\right)d_{s}^{ci}$
(41)
$\displaystyle\times\sum_{a}d_{a}^{ci}\left|4\pi\sum_{L}(2L+1)T_{ci}^{a,L}\left(E_{\mathrm{cm}},p_{\mathrm{cm}},p_{\mathrm{cm}}^{\prime}\right)P_{L}(x)\right|^{2}$

with the color and spin degeneracies of the two-body system, $d_{a,s}^{ci}$.
The heavy-light $T$-matrix,
$T_{ci}^{a,L}(E_{\mathrm{cm}},p_{\mathrm{cm}},p_{\mathrm{cm}}^{\prime})$, is
calculated in the CM frame in all possible two-body color channels, $a$, and
partial-wave channels, $L$ (expanded up to $L=8$ to ensure convergence at high
momenta). The CM energy $E_{\mathrm{cm}}$, incoming CM momentum
$p_{\mathrm{cm}}$, outgoing CM momentum $p^{\prime}_{\mathrm{cm}}$, and
scattering angle, $x=\cos\theta_{\mathrm{cm}}$, are expressed as functions of
$E,\mathbf{p},\mathbf{q},\mathbf{p}^{\prime},\mathbf{q}^{\prime}$, through the
transformation in Eq. (9). Instead of only the moduli of $p_{\mathrm{cm}}$ and
$p^{\prime}_{\mathrm{cm}}$, their explicit vector form is required here Liu et
al. (2019):

$\displaystyle
p_{\mathrm{cm}\|}=\frac{\varepsilon_{\mathrm{p}_{2}}p_{1\|}-\varepsilon_{\mathrm{p}_{1}}p_{2\|}}{\sqrt{s_{\mathrm{on}}}},$
$\displaystyle\mathbf{p}_{\mathrm{cm}\perp}=\frac{\mathbf{p}_{1}p_{2\|}-\mathbf{p}_{2}p_{1\|}}{\left|\mathbf{p}_{1}+\mathbf{p}_{2}\right|},$
(42)

with $\|$ and $\perp$ indicating parallel and perpendicular to the relative
velocity, respectively, and likewise for the outgoing (primed) momenta.

Figure 10: The charm-quark friction coefficient at different temperatures for
$\chi=0.6$ (solid) and 1 (dashed). The $\chi=1$ results are taken from Ref.
Liu et al. (2019).

In Fig. 10 we plot our results for the friction coefficient $A(p)$ with the
mixed confining potential ($\chi=0.6$) in comparison to the results with a
purely scalar confining potential ($\chi=1$). We stipulate that both
calculations are carried out for a thermal QGP medium which satisfies the
constraints from the EoS and HQ free energy. With the vector component in the
confining potential, the low-momentum values of the relaxation rate are
enhanced by several tens of percent, but the more significant effect is the
increase at higher momenta, for the same reasons as discussed above in the
context of the single-parton spectral functions and their scattering
amplitudes. For example, for a charm-quark momentum of 4 GeV, the enhancement
is about a factor 2.6, while at momenta of 10 GeV it reaches an even larger
factor of $\sim$3.5 at the lowest temperature. However, at the latter
momentum, radiative contributions are expected to be large. At first sight it
might be surprising that the enhancement due to the vector component in the
confining potential also transpires at low momenta although the pertinent
$T$-matrix amplitudes are smaller than those with purely scalar confining
potential at low CM momenta (cf. Fig. 9). To some extent this can be
understood due to the fact that even at vanishing HQ momentum the thermal
motion of the surrounding medium partons creates a finite momentum in the CMS,
but there is also a non-trivial interference effect in the expression (41)
that plays a role.

To scrutinize different contributions, we take the charm-light contribution
($c\bar{q}$) for $A(p=0)$ at $T=194$ MeV as an example and collect in Tab. 9
partial-wave components of the collision rate (obtained by replacing
$(1-\frac{\mathbf{p}\cdot\mathbf{p^{\prime}}}{\mathbf{p}^{2}})$ by 1 in Eq.
(40)) and relaxation rate up to angular momenta of 2 (note that for the
collision rate the interference contributions should vanish, whereby the small
numerical values quoted in the table, which are of the order of 1-2 permille
of the total, are an indication of our numerical accuracy) . We denote by
“$LL^{\prime}$” the terms $\sim
T^{L}(T^{L^{\prime}})^{*}=\textup{Re}T^{L}\textup{Re}T^{L^{\prime}}+\textup{Im}T^{L}\textup{Im}T^{L^{\prime}}+i(-\textup{Re}T^{L}\textup{Im}T^{L^{\prime}}+\textup{Im}T^{L}\textup{Re}T^{L^{\prime}})$
in Eq. (41), where the imaginary part vanishes by definition since
$T^{L}(T^{L^{\prime}})^{*}+(T^{L})^{*}T^{L^{\prime}}=2(\textup{Re}T^{L}\textup{Re}T^{L^{\prime}}+\textup{Im}T^{L}\textup{Im}T^{L^{\prime}})$.
In accord with the $T$-matrix behavior at low momenta in Fig. 9, the collision
rate for $\chi=1$ is larger than that for $\chi=0.6$ for each partial-wave
component, leading to a larger total collision rate at low momentum. The
situation is more involved for the relaxation rate: the diagonal partial-wave
components ($L=L^{\prime}$) for $\chi=1$ are smaller than those for
$\chi=0.6$, and one also notices the relatively more important role of the
higher partial waves compared to the collision rate (which is dominated by the
$S$-wave contribution). In addition, the presence of the
${\mathbf{p}\cdot\mathbf{p^{\prime}}}$ term, together with the Legendre
polynomials, causes negative interference components ($L\neq L^{\prime}$), and
their absolute values are larger for $\chi=1$. Upon adding the diagonal and
interference components the relaxation rate for $\chi=0.6$ becomes larger.

The widely discussed spatial diffusion coefficient, $D_{s}=T/(M_{c}A(p=0))$,
is related to the relaxation time, $\tau_{c}=1/A(p=0)$, at vanishing
3-momentum of the heavy quark. It is commonly scaled by the inverse thermal
wavelength, $2\pi T$, to render a dimensionless quantity for which we display
our results in Fig. 11 as a function of temperature. The $\chi=0.6$ result
shows a mild reduction relative to the $\chi=1$ one, which is again caused by
the larger average momenta of the thermal partons probed by the charm quark.

Table 9: The contributions of various partial-wave components (specified in the first row) of the $c\bar{q}$ scattering amplitude to the $c$-quark collision rate (lines 2 and 3) and relaxation rate (lines 4 and 5) for $p=0$ at $T=194$ MeV (in units of $\operatorname{fm^{-1}}$). $LL^{\prime}$ | $00$ | $11$ | $22$ |  $01$=$10$ |  $12$=$21$
---|---|---|---|---|---
$\chi=1$ | $0.8058$ | $0.3321$ | $0.0678$ | $0.0038$ | $0.0016$
$\chi=0.6$ | $0.7308$ | $0.2943$ | $0.0533$ | $0.0022$ | $0.0014$
$\chi=1$ | $0.1501$ | $0.1138$ | $0.0354$ | $-0.0521$ | $-0.0306$
$\chi=0.6$ | $0.1522$ | $0.1313$ | $0.0359$ | $-0.0322$ | $-0.0305$

Figure 11: The charm-quark spatial diffusion coefficient for $\chi=0.6$
(solid) and 1 (dashed). The $\chi=1$ results are taken from Ref. Liu et al.
(2019).

The increase in the elastic charm-quark friction coefficient, and in
particular its harder 3-momentum dependence, found here could have significant
ramifications for the phenomenology of open HF probes in URHICs. In a recent
work He and Rapp (2020) a good description of $D$, $D_{s}$ and $\Lambda_{c}$
observables in heavy-ion collisions has been achieved using the $T$-matrix
based transport coefficients from Refs. Riek and Rapp (2010); Huggins and Rapp
(2012), which are based on the internal-energy ($U$) as a potential proxy but
with an extra $K$ factor of about 1.6. The pertinent results for $A(p,T)$
(with $K=1.6$) are slightly larger than the ones from the SCS with $\chi=1$ at
low momentum, but much larger at higher momenta. However, with our new
$\chi=0.6$ results, the low-momentum deficit can be overcome, while they still
fall below the high-momentum results of the $U$-potential with $K=1.6$. Yet,
the inclusion of radiative processes, as computed within the $T$-matrix
approach in Ref. Liu and Rapp (2020b) could result in a total transport
coefficients that are quite comparable to the one employed in Ref. He and Rapp
(2020), without the need of any phenomenological adjustments.

##  VII Conclusions

We have augmented the thermodynamic $T$-matrix approach to include the effects
of spin-dependent interactions between heavy quarks, including spin-orbital,
spin-spin and tensor contributions, as part of the more general objective to
assess 1/$M_{Q}$ corrections. Toward this end we have utilized the Breit-Fermi
Hamiltonian to derive these interactions in the context of the Cornell
potential as the two-body interaction kernel for the $T$-matrix equation. When
benchmarking these interactions using the experimentally observed splittings
in vacuum quarkonium spectroscopy, we have found that, in accordance with
previous studies, a moderate admixture of a Lorentz-vector component in the
confining potential allows for a much improved description especially in the
charmonium sector. We have then implemented the amended interaction kernel
into our quantum many-body approach for the QGP. While the spin-dependent
interactions themselves are expected to be of minor importance (and therefore
have been neglected), the vector component of the confining potential turns
out to be rather significant. After selfconsistently refitting the in-medium
HQ free energies and the QGP EoS under the inclusion of the vector component,
quantitative modifications of the QGP properties toward shorter distances
(larger momenta) were found. A strong broadening of the thermal-parton
spectral functions persists to higher 3-momenta as a consequence of an
increased interaction strength in the thermodynamic 2-body scattering
amplitudes at larger momenta. The harder amplitudes and spectral functions are
a consequence of the relativistic corrections induced by the vector part of
the confining interaction, as opposed to a purely scalar interaction. This
suggests that the nature of the confining force in the QCD vacuum has an
impact on the properties of the strongly-coupled QGP, with liquid properties
that extend to higher resolution scales compared to a purely scalar confining
force. Finally, we have applied the modified set-up to calculate the friction
coefficient of charm quarks. As compared to the results with a purely scalar
string potential, a slightly larger relaxation rate is found at small momentum
(and a pertinent decrease in the diffusion coefficient), but a much larger
increase of a factor of $\sim$2-3 (or more) at momenta of around 5 GeV (and
above). These are promising features to make a significant step forward in
achieving a quantitative description of HF diffusion in heavy-ion collisions
at RHIC and the LHC based on microscopically and non-perturbatively calculated
transport coefficients. Work in this direction is in progress.

###### Acknowledgements.

This work has been supported by the U.S. National Science Foundation under
grant nos. PHY-1913286 and PHY-2209335, and by the U.S. Department of Energy,
Office of Science, Office of Nuclear Physics through the Topical Collaboration
in Nuclear Theory on Heavy-Flavor Theory (HEFTY) for QCD Matter under award
no. DE-SC0023547.

## References

  * Brambilla et al. (2004) N. Brambilla et al. (Quarkonium Working Group) (2004), eprint hep-ph/0412158.
  * Rapp et al. (2010) R. Rapp, D. Blaschke, and P. Crochet, Prog. Part. Nucl. Phys. 65, 209 (2010), eprint 0807.2470.
  * Kluberg and Satz (2010) L. Kluberg and H. Satz, _Color Deconfinement and Charmonium Production in Nuclear Collisions_ (2010), eprint 0901.3831.
  * Braun-Munzinger and Stachel (2010) P. Braun-Munzinger and J. Stachel, Landolt-Bornstein 23, 424 (2010), eprint 0901.2500.
  * Du et al. (2019) X. Du, S. Y. F. Liu, and R. Rapp, Phys. Lett. B 796, 20 (2019), eprint 1904.00113.
  * Petreczky and Teaney (2006) P. Petreczky and D. Teaney, Phys. Rev. D 73, 014508 (2006), eprint hep-ph/0507318.
  * Petreczky et al. (2006) P. Petreczky, K. Petrov, D. Teaney, and A. Velytsky, PoS LAT2005, 185 (2006), eprint hep-lat/0510021.
  * Petreczky (2009) P. Petreczky, Eur. Phys. J. C 62, 85 (2009), eprint 0810.0258.
  * Ding et al. (2012) H. T. Ding, A. Francis, O. Kaczmarek, F. Karsch, H. Satz, and W. Soeldner, Phys. Rev. D 86, 014509 (2012), eprint 1204.4945.
  * Rapp and van Hees (2010) R. Rapp and H. van Hees (2010), pp. 111–206, eprint 0903.1096.
  * Mocsy et al. (2013) A. Mocsy, P. Petreczky, and M. Strickland, Int. J. Mod. Phys. A 28, 1340012 (2013), eprint 1302.2180.
  * He et al. (2022) M. He, H. van Hees, and R. Rapp (2022), eprint 2204.09299.
  * Prino and Rapp (2016) F. Prino and R. Rapp, J. Phys. G 43, 093002 (2016), eprint 1603.00529.
  * Beraudo et al. (2018) A. Beraudo et al., Nucl. Phys. A 979, 21 (2018), eprint 1803.03824.
  * Dong and Greco (2019) X. Dong and V. Greco, Prog. Part. Nucl. Phys. 104, 97 (2019).
  * Cabrera and Rapp (2007) D. Cabrera and R. Rapp, Phys. Rev. D 76, 114506 (2007), eprint hep-ph/0611134.
  * Riek and Rapp (2010) F. Riek and R. Rapp, Phys. Rev. C 82, 035201 (2010), eprint 1005.0769.
  * Mannarelli and Rapp (2005) M. Mannarelli and R. Rapp, Phys. Rev. C 72, 064905 (2005), eprint hep-ph/0505080.
  * Liu et al. (2019) S. Y. F. Liu, M. He, and R. Rapp, Phys. Rev. C 99, 055201 (2019), eprint 1806.05669.
  * Mur et al. (1994) V. D. Mur, V. S. Popov, Y. A. Simonov, and V. P. Yurov, J. Exp. Theor. Phys. 78, 1 (1994), eprint hep-ph/9401203.
  * Lucha et al. (1991) W. Lucha, F. F. Schoberl, and D. Gromes, Phys. Rept. 200, 127 (1991).
  * Buchmuller (1982) W. Buchmuller, Phys. Lett. B 112, 479 (1982).
  * Brambilla and Vairo (1997a) N. Brambilla and A. Vairo, Phys. Rev. D 55, 3974 (1997a), eprint hep-ph/9606344.
  * Brambilla and Vairo (1997b) N. Brambilla and A. Vairo, Phys. Lett. B 407, 167 (1997b), eprint hep-ph/9703377.
  * Szczepaniak and Swanson (1997) A. P. Szczepaniak and E. S. Swanson, Phys. Rev. D 55, 3987 (1997), eprint hep-ph/9611310.
  * Ebert et al. (1998) D. Ebert, V. O. Galkin, and R. N. Faustov, Phys. Rev. D 57, 5663 (1998), [Erratum: Phys.Rev.D 59, 019902 (1999)], eprint hep-ph/9712318.
  * Ebert et al. (2003) D. Ebert, R. N. Faustov, and V. O. Galkin, Phys. Rev. D 67, 014027 (2003), eprint hep-ph/0210381.
  * Liu and Rapp (2018) S. Y. F. Liu and R. Rapp, Phys. Rev. C 97, 034918 (2018), eprint 1711.03282.
  * Brockmann and Machleidt (1999) R. Brockmann and R. Machleidt, Int. Rev. Nucl. Phys. 8, 121 (1999), eprint nucl-th/9612004.
  * Fetter and Walecka (2012) A. L. Fetter and J. D. Walecka, _Quantum theory of many-particle systems_ (Courier Corporation, 2012).
  * Cheng et al. (2008) M. Cheng et al., Phys. Rev. D 77, 014511 (2008), eprint 0710.0354.
  * Petreczky (2010) P. Petreczky, J. Phys. G 37, 094009 (2010), eprint 1001.5284.
  * Bazavov et al. (2014a) A. Bazavov, Y. Burnier, and P. Petreczky, Nucl. Phys. A 932, 117 (2014a), eprint 1404.4267.
  * Bazavov et al. (2014b) A. Bazavov et al. (HotQCD), Phys. Rev. D 90, 094503 (2014b), eprint 1407.6387.
  * Bazavov et al. (2018) A. Bazavov, N. Brambilla, P. Petreczky, A. Vairo, and J. H. Weber (TUMQCD), Phys. Rev. D 98, 054511 (2018), eprint 1804.10600.
  * Lucha and Schoberl (1995) W. Lucha and F. F. Schoberl, in _International Summer School for Students on Development in Nuclear Theory and Particle Physics_ (1995), eprint hep-ph/9601263.
  * Hong et al. (2022) K.-H. Hong, H.-C. Kim, and U. Yakhshiev, PTEP 2022, 103D02 (2022), eprint 2208.01851.
  * Tanabashi et al. (2018) M. Tanabashi et al. (Particle Data Group), Phys. Rev. D 98, 030001 (2018).
  * Beraudo et al. (2008) A. Beraudo, J. P. Blaizot, and C. Ratti, Nucl. Phys. A 806, 312 (2008), eprint 0712.4394.
  * Laine et al. (2007) M. Laine, O. Philipsen, P. Romatschke, and M. Tassler, JHEP 03, 054 (2007), eprint hep-ph/0611300.
  * Luttinger and Ward (1960) J. M. Luttinger and J. C. Ward, Phys. Rev. 118, 1417 (1960).
  * Baym and Kadanoff (1961) G. Baym and L. P. Kadanoff, Phys. Rev. 124, 287 (1961).
  * Baym (1962) G. Baym, Phys. Rev. 127, 1391 (1962).
  * Liu and Rapp (2017) S. Y. F. Liu and R. Rapp, J. Phys. Conf. Ser. 779, 012034 (2017), eprint 1609.04877.
  * Liu and Rapp (2020a) S. Y. F. Liu and R. Rapp, Eur. Phys. J. A 56, 44 (2020a), eprint 1612.09138.
  * Bazavov and Petreczky (2013) A. Bazavov and P. Petreczky, J. Phys. Conf. Ser. 432, 012003 (2013), eprint 1211.5638.
  * He and Rapp (2020) M. He and R. Rapp, Phys. Rev. Lett. 124, 042301 (2020), eprint 1905.09216.
  * Huggins and Rapp (2012) K. Huggins and R. Rapp, Nucl. Phys. A 896, 24 (2012), eprint 1206.6537.
  * Liu and Rapp (2020b) S. Y. F. Liu and R. Rapp, JHEP 08, 168 (2020b), eprint 2003.12536.