††thanks: Corresponding author

# Bayesian Inference of the Specific Shear and Bulk Viscosities of the Quark-
Gluon Plasma at Crossover from $\phi$ and $\Omega$ Observables

Zhidong Yang  School of Physics and Astronomy, Shanghai Key Laboratory for
Particle Physics and Cosmology, and Key Laboratory for Particle Astrophysics
and Cosmology (MOE), Shanghai Jiao Tong University, Shanghai 200240, China
Lie-Wen Chen<EMAIL_ADDRESS>School of Physics and Astronomy, Shanghai Key
Laboratory for Particle Physics and Cosmology, and Key Laboratory for Particle
Astrophysics and Cosmology (MOE), Shanghai Jiao Tong University, Shanghai
200240, China

###### Abstract

Due to their weak final state interactions, the $\phi$ meson and $\Omega$
baryon provide unique probes of the properties of the quark-gluon plasma (QGP)
formed in relativistic heavy-ion collisions. Using the quark recombination
model with the quark phase-space information parameterized in a viscous
blastwave, we study the transverse-momentum spectra and elliptic flows of
$\phi$ and $\Omega$ in Au+Au collisions at $\sqrt{s_{\rm NN}}=200$ GeV and
Pb+Pb collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV. The viscous blastwave
includes non-equilibrium deformations of thermal distributions due to shear
and bulk stresses and thus carries information on the specific shear viscosity
$\eta/s$ and the specific bulk viscosity $\zeta/s$ of the QGP. We perform a
model-to-data comparison with Bayesian inference and simultaneously obtain
$\eta/s=(2.08^{+1.10}_{-1.09})/4\pi$ and $\zeta/s=0.06^{+0.04}_{-0.04}$ at
$90\%$ C.L. for the baryon-free QGP at crossover temperature of about $160$
MeV. Our work provides a novel approach to simultaneously determine the
$\eta/s$ and $\zeta/s$ of the QGP at hadronization.

Heavy Ion Collisions, Quark Recombination, Strange quark, Shear viscosity,
bulk viscosity

###### pacs:

24.85.+p,25.75.Ag,25.75.Ld

##  I Introduction

A new state of matter that consists of deconfined quarks and gluons, called
the quark-gluon plasma (QGP), is believed to have been created in relativistic
heavy-ion collisions (HICs) at the Relativistic Heavy Ion Collider (RHIC) and
the Large Hadron Collider (LHC). The QGP is expected to fill the early
universe at about $10^{-6}$ s after the Big Bang. Studies show that the QGP
behaves like a near-perfect fluid and has very small specific shear viscosity
$\eta/s$ Romatschke and Romatschke (2007); Song and Heinz (2008), i.e., the
ratio of shear viscosity $\eta$ to entropy density $s$, close to the universal
lower bound $\eta/s=1/4\pi$ based on the AdS/CFT correspondence Kovtun  _et
al._ (2005). It turns out that the $\eta/s$ directly characterizes the
behavior of QGP fluid and has significant impacts on the final observables of
HICs. Studies also indicate the bulk viscosity ($\zeta$) has non-negligible
effects on HIC observables Ryu  _et al._ (2015, 2018).

Determining the temperature and baryon density dependence of QGP’s $\eta/s$
and $\zeta/s$ is of fundamental importance and remains a big challenge in the
field Shen (2021). For zero baryon density or bayron chemical potential
($\mu_{B}$=0), the first principle lattice QCD (LQCD) predicts that the
transition from hadrons to QGP is not a real phase transition but a smooth
crossover at a pseudo-critical temperature $T_{\rm pc}\approx 160$ MeV
Bhattacharya  _et al._ (2014); Borsanyi  _et al._ (2020). Theoretical
calculations on the temperature dependence of the viscosities of the baryon-
free QGP have been widely explored in various approaches in the past years,
including LQCD Meyer (2007, 2009); Mages  _et al._ (2015); Astrakhantsev  _et
al._ (2017); Borsányi  _et al._ (2018), Functional Renormalization Group (FRG)
method Christiansen  _et al._ (2015), T-matrix Liu and Rapp (2020) and
perturbative QCD (pQCD) Ghiglieri  _et al._ (2018) for the $\eta/s$, as well
as LQCD Karsch  _et al._ (2008); Meyer (2008); Astrakhantsev  _et al._ (2018)
and holographic model Rougemont  _et al._ (2017) for the $\zeta/s$. Overall,
the theoretical predictions remain largely uncertain even at $T_{\rm pc}$ and
the first principle calculation is still challenging Moore (2020).
Alternately, many efforts have focused on constraining the viscosities of QGP
with experimental observables.

Early works using viscous hydrodynamics have given semiquantitative estimates
on the $\eta/s$ of QGP Romatschke and Romatschke (2007); Schenke  _et al._
(2011); Gale  _et al._ (2013); Niemi  _et al._ (2016). Usually they assumed a
constant $\eta/s$ over the entire evolution and found $\eta/s=1\sim
2.5/(4\pi)$. Recently, results using multi-stage approaches which combine
initial conditions, viscous hydrodynamics and hadronic transport are reported
in Refs. Bernhard  _et al._ (2015, 2016, 2019); Everett  _et al._ (2021a, b);
Nijs  _et al._ (2021a, b); Nijs and van der Schee (2022); Parkkila  _et al._
(2021a, b); Heffernan  _et al._ (2023). A variety of experimental observables,
often focused on protons, kaons and pions, and a number of parameters are used
in these works. Although they used very similar methods, the obtained results
tend to be different, depending on the details of the model and observables.
For example, for the baryon-free QGP at $T_{\rm pc}$, JETSCAPE Collaboration
reported $\eta/s\approx 1.6/(4\pi)$ and $\zeta/s\approx 0.09$ Everett  _et
al._ (2021a), while Trajectum group reported $\eta/s\approx 0.9/(4\pi)$ and
$\zeta/s\approx 0.004$ Nijs  _et al._ (2021a), displaying a significant
discrepancy. Very recently Trajectum group updated their results to
$\eta/s\approx 2.3/(4\pi)$ and $\zeta/s\approx 0.02$ at $T_{\rm pc}$ after
considering the recent measured nucleon nucleon cross section Nijs and van der
Schee (2022). Meanwhile, they reported $\eta/s(T)$ decreases as temperature
increases Nijs and van der Schee (2022), which is opposite to the results of
Duke group Bernhard  _et al._ (2019).

In this work, we present a novel approach to constrain the $\eta/s$ and
$\zeta/s$ of the baryon-free QGP at $T_{\rm pc}$ by using the measured
transverse-momentum spectra and elliptic flows of $\phi$ and $\Omega$ in Au+Au
collisions at $\sqrt{s_{\rm NN}}=200$ GeV and Pb+Pb collisions at
$\sqrt{s_{\rm NN}}=2.76$ TeV where the QGP has a negligible baryon density and
its transition to hadrons is a smooth crossover. In our approach, hadrons are
produced through quark recombination Greco  _et al._ (2003a, b); Fries  _et
al._ (2003a, b, 2008); He  _et al._ (2010) with the phase-space distribution
of quarks at hadronization parameterized in a viscous blastwave Teaney (2003);
Jaiswal and Koch (2015); Yang and Fries (2017, 2020, 2022) which includes non-
equilibrium deformations of thermal distributions due to shear and bulk
stresses. The viscous corrections on the QGP at hadronization are then
imported into $\phi$ and $\Omega$ through the recombination process. Since the
$\phi$ and $\Omega$ have weak hadronic interactions Shor (1985), they thus
carry direct information of QGP at hadronization with negligible hadronic
effects Shor (1985); van Hecke  _et al._ (1998); Chen and Ko (2006); Chen  _et
al._ (2008); Auvinen  _et al._ (2017); Hwa and Zhu (2017); Ye  _et al._
(2017); Pu  _et al._ (2018); Song  _et al._ (2020), making it possible to
determine the $\eta/s$ and $\zeta/s$ of the QGP at crossover from the $\phi$
and $\Omega$ observables.

##  II Model and methods

###  II.1 Theoretical model

Quark recombination or coalescence models were first proposed to explain the
baryon-over-meson enhancement and valence quark number scaling in RHIC Au+Au
collisions Greco  _et al._ (2003a, b); Fries  _et al._ (2003a, b). In these
models, valence quarks are assumed to be abundant in the phase space and
recombine into hadrons through quark recombination. The hadron formation
process is usually assumed to be instantaneous and takes a very thin
hypersurface ($\Delta\tau\approx 0$). Following Refs. Fries  _et al._ (2003b,
a, 2008), the momentum distribution of baryons is given by

$\displaystyle E\frac{dN_{B}}{d^{3}\bf p}$ $\displaystyle=$ $\displaystyle
C_{B}\int_{\Sigma}\frac{p^{\mu}\cdot
d\mathbf{\sigma_{\mu}}}{(2\pi)^{3}}\int_{0}^{1}dx_{1}dx_{2}dx_{3}\Phi_{B}(x_{1},x_{2},x_{3})$
(1) $\displaystyle\times f_{a}({\bf r},x_{1}\mathbf{p})f_{b}({\bf
r},x_{2}\mathbf{p})f_{c}({\bf r},x_{3}\mathbf{p}).\ $

where $C_{B}$ is the spin degeneracy factor of a given baryon species,
$\Sigma$ is the hypersurface of hadronization, $\varPhi_{B}$ is the effective
wave function squared of baryons, $x_{1,2,3}$ are light cone coordinates
defined as ${\bf p}_{1,2,3}=x_{1,2,3}\bf p$, and $f_{a,b,c}$ are the parton
phase-space distributions. Here we make the assumption that the partons are
essentially ultra-relativistic and traveling along the light cone, which
should be valid for the LHC and the top RHIC energies. The $\varPhi_{B}$ is
parameterized as polynomial Fries  _et al._ (2003b), and here we use a
Gaussian ansatz
$\Phi_{B}\sim\exp(-\frac{(x_{1}-x_{a})^{2}+(x_{2}-x_{b})^{2}+(x_{3}-x_{c})^{2}}{\sigma_{B}^{2}})\delta(x_{1}+x_{2}+x_{3}-1)$,
where $x_{a,b,c}=m_{1,2,3}/(m_{1}+m_{2}+m_{3})$ are the peak values, and
$m_{1,2,3}$ are the masses of the constitute partons. Similar expression can
be derived for mesons.

The quark phase-space distribution is parameterized in a viscous blastwave
Yang and Fries (2017, 2020, 2022), based on Retiere and Lisa (RL) blastwave
Retiere and Lisa (2004). The quark distribution is given by

$f(r,p)=f_{0}(r,p)+\delta f_{\rm{shear}}(r,p)+\delta f_{\rm{bulk}}(r,p)$ (2)

where $f_{0}$ is the equilibrium Bose/Fermi distribution, $\delta
f_{\rm{shear}}$ and $\delta f_{\rm{bulk}}$ are corrections from the shear and
bulk viscosities, respectively. For the shear viscous corrections, we use the
Grad’s method Song and Heinz (2008); Damodaran  _et al._ (2017)

$\delta
f_{\rm{shear}}=\frac{1}{2s}\frac{p^{\mu}p^{\nu}}{T^{3}}\pi^{\mu\nu}f_{0}(1\pm
f_{0})$ (3)

where $\pi^{\mu\nu}$ is the shear stress tensor and $+(-)$ for bosons
(fermions). In the Navier-Stokes approximation
$\pi^{\mu\nu}=2\eta\sigma^{\mu\nu}$ where $\sigma^{\mu\nu}$ is the shear
gradient tensor defined as
$\sigma^{\mu\nu}=\frac{1}{2}\left(\nabla^{\mu}u^{\nu}+\nabla^{\nu}u^{\mu}\right)-\frac{1}{3}\Delta^{\mu\nu}\nabla_{\lambda}u^{\lambda}$
with flow field $u^{\mu}$, $\nabla^{\mu}=\Delta^{\mu\nu}\partial_{\nu}$ and
$\Delta^{\mu\nu}=g^{\mu\nu}-u^{\mu}u^{\nu}$. For the bulk viscous corrections,
we use 14-moment approximation Denicol  _et al._ (2014); Paquet  _et al._
(2016)

$\delta f_{\rm bulk}=-f_{0}(1\pm
f_{0})\Pi\frac{\tau_{\Pi}}{\zeta}\left[\frac{1}{3}\frac{m^{2}}{T}\frac{1}{p^{\mu}u_{\mu}}+\frac{p^{\mu}u_{\mu}}{T}\left(c_{s}^{2}-\frac{1}{3}\right)\right]$
(4)

where $\Pi$ is the bulk viscous pressure and $\tau_{\Pi}$ is the bulk
relaxation time. At the first order approximation, one has
$\Pi=-\zeta\partial_{\mu}u^{\mu}$. The expression for
$\frac{\tau_{\Pi}}{\zeta}$ is given in Denicol  _et al._ (2014).

Now we present the blastwave parameterization for the flow field $u^{\mu}$. In
the following, $R_{x,y}$ are the semi axes of the fireball at freeze-out,
$\rho=\sqrt{x^{2}/R_{x}^{2}+y^{2}/R_{y}^{2}}$ is the reduced radius,
$\eta_{s}=\frac{1}{2}\log\frac{t+z}{t-z}$ is the space-time rapidity. The
hypersurface is assumed to be constant $\tau=\sqrt{t^{2}-z^{2}}$. The flow
field is parameterized as

$\displaystyle u^{\mu}$ $\displaystyle=$
$\displaystyle(\cosh\eta_{s}\cosh\eta_{T},\sinh\eta_{T}\cos\phi_{b},$ (5)
$\displaystyle\qquad\sinh\eta_{T}\sin\phi_{b},\sinh\eta_{s}\cosh\eta_{T})$

where $\eta_{T}$ is the transverse flow rapidity and $\phi_{b}$ is the
azimuthal angle of $u^{\mu}$ in the transverse plane. Here we assume boost
invariance and set longitudinal flow rapidity to be equal to the space-time
rapidity $\eta_{s}$. $\eta_{T}$ is given by the transverse velocity
$v_{T}=\tanh\eta_{T}$ with

$v_{T}=\rho^{n}(\alpha_{0}+\alpha_{2}\cos 2\phi_{b})$ (6)

where $\alpha_{0}$ is the average surface velocity, $\alpha_{2}$ is an
elliptic deformation of the flow field and $n$ is a power term. In this work,
we use a linear expression and set $n=1$. We choose the transverse flow vector
perpendicular to the elliptic surface at $\rho=1$, i.e.,
$\tan\phi_{b}=R_{x}^{2}/R_{y}^{2}\tan\phi$, where $\phi=\arctan y/x$ is the
azimuthal angle of the position. In terms of the geometric parameters, it
turns out the ratio $R_{y}/R_{x}$ has a large influence on elliptical flow, so
we choose $R_{y}/R_{x}$ as a fit parameter and constrain $R_{x}$, $R_{y}$ and
$\tau$ by fitting the ratio $R_{y}/R_{x}$ and by adding the simple geometric
estimate $R_{x}\approx(R_{0}-b/2)+\tau c_{\tau}(\alpha_{0}+\alpha_{2})$ where
$R_{0}$ is the radius of the colliding nucleus, $b$ is the impact parameter
and $c_{\tau}=\bar{\alpha}_{0}/\alpha_{0}$ relates the time-averaged surface
velocity. In this work we use $R_{0}=7.1$ fm (7.0 fm) for Pb (Au) based on
formula in Nerlo-Pomorska and Pomorski (1994) and $c_{\tau}=0.65$. The values
of $b$ for each centrality bin are taken from Glauber Monte Carlo simulations
used by related experiments Abelev  _et al._ (2013); Adamczewski-Musch  _et
al._ (2018).

It should be noted that the viscous blastwave parameterizes the flow field and
freeze-out hypersurface with a simple ansatza, which can be regarded as an
approximate snapshot of a viscous hydrodynamic system at a fixed time Teaney
(2003); Jaiswal and Koch (2015); Yang and Fries (2017, 2020, 2022). Therefore,
the viscous blastwave carries information on the viscosities of the fluid at a
fixed time, e.g., the QGP at hadronization in our present work.

###  II.2 Experimental data and fit parameters

To compare our model to experiment, we utilize the transverse-momentum
($p_{T}$) spectra and elliptic flows $v_{2}$ of $\phi$ and $\Omega$ as our
observables. The data are taken from the STAR Collaboration for Au+Au
collisions at $\sqrt{s_{\rm NN}}=200$ GeV Adams  _et al._ (2007); Abelev  _et
al._ (2007); Adamczyk  _et al._ (2016) and the ALICE Collaboration for Pb+Pb
collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV Abelev  _et al._ (2014, 2015); Adam
_et al._ (2017) at mid-rapidity. For $v_{2}$ of $\phi$ and $\Omega$, we use
the data in centralities 0-30%, 30-80% from STAR and 10-20%, 20-30%, 30-40%,
40-50% from ALICE. Due to the limit of available data, the centrality bins for
$p_{T}$ spectra are slightly different from that of $v_{2}$. For the $\phi$
spectra, we use the data in 10-20%, 40-50% from STAR and 10-20%, 20-30%,
30-40%, 40-50% from ALICE. For the $\Omega$ spectra, we use the data in 0-5%
and 40-60% from STAR and 10-20%, 20-40%, 40-60% from ALICE. When the
centrality of spectra are not matched to that of $v_{2}$, we choose the
spectra that belongs to the neighbor centrality.

The fit ranges are $p_{T}<$ 1.8 GeV/$c$ for $\phi$ and $p_{T}<$ 3.2 GeV/$c$
for $\Omega$. We do a combined analysis for both STAR and ALICE data and have
116 data points in total. We note that the correction $\delta f(r,p)$ is small
for both $\phi$ and $\Omega$, i..e., less than $20\%$ of $f_{0}$ for the
majority of points (very few points up to $40\%$ of $f_{0}$), which is much
smaller than the usually adopted upper bound $\sim 1$ Everett  _et al._
(2021b); Teaney (2003), guaranteeing the applicability of the viscous
corrections.

The parameters entering the model are ($\tau$, $T$, $\alpha_{0}$,
$\alpha_{2}$, $R_{y}/R_{x}$, $\eta/s$, $\zeta/s$) from blastwave and
($\sigma_{B}$, $\sigma_{M}$) from the quark recombination. Here
$\sigma_{B}$($\sigma_{M}$) is the wave function width for baryons (mesons). We
find the values of $\sigma_{M,B}$ mainly affect the particle yields and have
little impact on $v_{2}$. To describe the yields of $\phi$ and $\Omega$, we
set $\sigma_{M}=0.3$, $\sigma_{B}=0.09$. We note that our conclusion changes
little by varying the values of $\sigma_{M}$ and $\sigma_{B}$. The freeze-out
or recombination hypersurface is related to the fireball volume and can be
determined from the particle yield.

For the other constants, we use sound speed squared $c^{2}_{s}=0.15$ (see Eq.
(4)) for the baryon-free QGP at $T_{\rm pc}$ Huovinen and Petreczky (2010),
quark mass $m_{s}=0.5$ GeV, spin degeneracy factor $C_{M}=3$ for
$\mathrm{\phi}$ and $C_{B}=4$ for $\Omega^{-}$ (or $\bar{\Omega}^{+}$). To fit
the yield of $\phi$ and $\Omega$ as well as determine the value of freeze-out
time $\tau$, we introduce a fugacity factor $\gamma_{s,\bar{s}}$ for strange
quarks and set $\gamma_{s,\bar{s}}=0.8$ for all centrality bins which is close
to the value found in Pu  _et al._ (2018). Now the parameters left in our
model are ($T,\alpha_{0},\alpha_{2},R_{y}/R_{x},\eta/s,\zeta/s$). For each
centrality bin of STAR (0-30%, 30-80%) and ALICE (10-20%, 20-30%, 30-40%,
40-50%), the fluid has unique values for $\alpha_{0}$, $\alpha_{2}$ and
$R_{y}/R_{x}$ with the same shared values for $T$, $\eta/s$ and $\zeta/s$, and
thus we have 21 parameters.

Figure 1:  Posterior distributions and correlations of selected parameters:
$T$, $\eta/s$ and $\zeta/s$. The numbers indicate the median values with the
90%-credibility range.

###  II.3 Bayesian method

To determine the above parameters, we use the Bayesian analysis package from
the Models and Data Analysis Initiative (MADAI) project Bass  _et al._ (2016);
Bernhard  _et al._ (2016). The MADAI package includes a Gaussian process
emulator and a Bayesian analysis tool. According to Bayes’ theorem, for model
parameters $\rm x=(x_{1},x_{2},x_{3},...)$ and experimental observables $\rm
y=(y_{1},y_{2},y_{3},...)$, the probability for the true parameters $\rm
x_{\star}$ is

$P({\rm x}_{\star}|X,Y,{\rm y}_{\rm{exp}})\propto P(X,Y,{\rm
y}_{\rm{exp}}|{\rm x}_{\star})P({\rm x}_{\star}).$ (7)

The left-hand side is the _posterior_ probability of ${\rm x}_{\star}$ given
the design ($X$, $Y$) and the experimental data ${\rm y}_{\rm{exp}}$. On the
right-hand side, $P({\rm x}_{\star})$ is the _prior_ probability and
$P(X,Y,{\rm y}_{\rm{exp}}|{\rm x}_{\star})$ is the likelihood, i. e. , the
probability of the model describing the data ${\rm y}_{\rm{exp}}$ at ${\rm
x}_{\star}$, given by

$P(X,Y,{\rm y}_{\rm{exp}}|{\rm x}_{\star})\propto\exp(-\frac{1}{2}\Delta{\rm
y}^{\top}\Lambda^{-1}\Delta{\rm y})$ (8)

where $\Delta{\rm y}={\rm y}_{\star}-{\rm y}_{\rm{exp}}$ is the difference
between the measurement and the prediction, and $\Lambda$ is the covariance
matrix including the experimental and model uncertainties. More details can be
found in Bernhard  _et al._ (2016); Everett  _et al._ (2021b).

All model parameters have _uniform_ prior distribution except for the
recombination temperature $T$, which has a _Gaussian_ prior. For example, we
use prior ranges 0-4.4/($4\pi$) for $\eta/s$, 0-0.16 for $\zeta/s$ and
0.72-0.8$c$ for $\alpha_{0}$, 0.02-0.036$c$ for $\alpha_{2}$, 1.1-1.18 for
$R_{y}/R_{x}$ in 0-30% centrality bin of STAR. For $T$, we use Gaussian
distribution with a mean value of 160 MeV and a standard deviation of 5 MeV to
reflect our prior knowledge of quark recombination temperature of a continuous
crossover at $T_{\rm pc}=158\pm 15$ MeV for zero net baryon density Borsanyi
_et al._ (2020).

The Bayesian analysis package works as follows. First, we set prior ranges for
each parameter and generate a set of training points in the parameter space.
Second, we calculate all fitted observables at each training point. The
package then builds a Gaussian process emulator, which can estimate the
observables for random parameter values. Finally a Markov Chain Monte Carlo
(MCMC) provides a likelihood analysis and gives the maximum likelihood or
best-fit parameters. Here we have 21 parameters and use N=500 training points.
To verify the Bayesian analysis works properly, we perform a validation for
the Bayesian inference framework and find the Bayesian framework correctly
reproduces model parameters within reasonable uncertainties. The likelihood
analysis has used ${\rm N}_{\star}=10^{7}$ predicted points to search for the
best-fit parameters, which is sufficient for MCMC to converge.

##  III Results and discussions

Using the data and prior parameter ranges, we perform a model-to-data
comparison with MADAI package. Fig. 1 shows the univariate posterior
distributions for the recombination temperature $T$, $\eta/s$ and $\zeta/s$.
We obtain $T=159\pm 8$ MeV, $\eta/s=(2.08^{+1.10}_{-1.09})/4\pi$ and
$\zeta/s=0.06^{+0.04}_{-0.04}$ at $90\%$ confidence level (C.L.). Fig. 1 also
shows the correlations between $T$, $\eta/s$ and $\zeta/s$. We find there is a
strong correlation between $\eta/s$ and $\zeta/s$. This can be qualitatively
understood by the different effects of $\eta/s$ and $\zeta/s$ on the particle
momentum distribution, i.e., the $\eta/s$ decreases the azimuthal momentum
asymmetry of particles, while the $\zeta/s$ reduces their momenta
isotropically Everett  _et al._ (2021a). As a result, the $\zeta/s$ tends to
increase $v_{2}$ while the $\eta/s$ prefers to reduces $v_{2}$ at higher
$p_{T}$. We also find there is a slight anti-correlation between $T$ and
$\eta/s$, and a moderate correlation between $T$ and $\zeta/s$.

Table 1:  The best-fit values of the viscous blastwave parameters for different centrality bins with $T=159$ MeV, $\eta/s=2.08/(4\pi)$ and $\zeta/s=0.06$. Centrality |  $\tau$ (fm/$c$) | $\alpha_{0}(c)$ | $\alpha_{2}(c)$ | $R_{y}/R_{x}$
---|---|---|---|---
Au+Au 200 GeV
0-30% | 7.8 | 0.76 | 0.027 | 1.13
30-80% | 5.1 | 0.70 | 0.039 | 1.31
Pb+Pb 2.76 TeV
10-20% | 9.4 | 0.84 | 0.028 | 1.18
20-30% | 9 | 0.81 | 0.037 | 1.18
30-40% | 7.8 | 0.81 | 0.036 | 1.25
40-50% | 6.8 | 0.79 | 0.031 | 1.36

For the other parameters, we show the best-fit values, which are defined as
the mean value given by maximum likelihood analysis, as given in Tab. 1. Using
the best-fit parameters, we calculate the $p_{T}$ spectra and $v_{2}$ of
$\phi$ and $\Omega$ and compare them with data. Fig. 2 shows the results of
$p_{T}$ spectra and $v_{2}$ in different centrality bins for Au+Au collisions
at $\sqrt{s_{\rm NN}}=200$ GeV and Pb+Pb at $\sqrt{s_{\rm NN}}=2.76$ TeV,
respectively. As expected, our calculations describe the data rather well.

Figure 2:  Transverse-momentum spectra and elliptic flows $v_{2}$ of $\phi$
and $\Omega^{-}$+$\bar{\Omega}^{+}$ at mid-rapidity in Au+Au collisions at
$\sqrt{s_{\rm NN}}=200$ GeV and Pb+Pb collisions at $\sqrt{s_{\rm NN}}=2.76$
TeV from the recombination model using the best-fit parameters given by
Baysian analysis. The data from STAR Adams  _et al._ (2007); Abelev  _et al._
(2007); Adamczyk  _et al._ (2016) and ALICE Abelev  _et al._ (2014, 2015);
Adam  _et al._ (2017) are included for comparison.

Figure 3 shows our Bayesian inference of the $\eta/s$ and $\zeta/s$ for the
baryon-free QGP at $T_{\rm pc}\approx 160$ MeV (with 68.3% and 90% C.L.,
indicated by SJTU). For comparison, we include in Fig. 3(a) the temperature
dependence of $\eta/s$ for the baryon-free QGP/hadronic matter from other
approaches, i.e., the multi-stage methods (with 90% C.L.) from Duke Bernhard
_et al._ (2019), JETSCAPE Everett  _et al._ (2021a) and Trajectum Nijs and van
der Schee (2022), the viscous blastwave (BW) Yang and Fries (2022), Chapman-
Enskog (Chap-Ensk) method Dash  _et al._ (2019), hadron resonance gas (HRG)
model with Hagedorn states (HS) Noronha-Hostler  _et al._ (2009), FRG
Christiansen  _et al._ (2015), LQCD (Lattice1 Meyer (2007, 2009), Lattice2
Mages  _et al._ (2015), Lattice3 Astrakhantsev  _et al._ (2017) and Lattice4
Borsányi  _et al._ (2018), with 68.3% C.L.), T-matrix Liu and Rapp (2020), and
next-to-leading order pQCD Ghiglieri  _et al._ (2018). Similarly, Fig. 3(b)
includes the corresponding results for $\zeta/s$ from the multi-stage methods
from Duke Bernhard  _et al._ (2019), JETSCAPE Everett  _et al._ (2021a) and
Trajectum Nijs and van der Schee (2022), LQCD (Lattice1 Meyer (2008), Lattice3
Astrakhantsev  _et al._ (2018) and Lattice5 Karsch  _et al._ (2008), with
68.3% C.L.), hybrid model (McGill) Ryu  _et al._ (2015), holographic model
(Holo) Rougemont  _et al._ (2017), HRG model with HS Noronha-Hostler  _et al._
(2009) and SMASH transport model Rose  _et al._ (2020).

Figure 3:  Temperature dependence of $\eta/s$ (a) and $\zeta/s$ (b) for
QGP/hadronic matter at $\mu_{B}$=0 (see text for details).

It is interesting to see from Fig. 3(a) that different methods give roughly
consistent results for $\eta/s$, indicating a trend that the $\eta/s$
approaches to a minimum at the crossover temperature $T_{\rm pc}$, with a
rather steep rise towards lower temperatures and a slow rise towards higher
temperatures. The behavior with minimum $\eta/s$ at $T_{\rm pc}$ is also
observed in recent work based on 2 + 1D hydrodynamical model with the Eskola-
Kajantie-Ruuskanen-Tuominen (EKRT) initial state Auvinen  _et al._ (2020) and
the quasiparticle model prediction Mykhaylova  _et al._ (2019). We note our
present constraint on the $\eta/s$ at $T_{\rm pc}$ is compatible with those of
Duke, JETSCAPE and Trajectum from the multi-stage methods, as well as the
Chap-Ensk, BW, FRG and Lattice3.

On the other side, one sees from Fig. 3(b) that the results of $\zeta/s$
exhibit significant discrepancy around $T_{\rm pc}$. For example, at $T_{\rm
pc}$, a rather large value of $\zeta/s\approx 0.35$ is obtained from the LQCD
(Lattice3 and Lattice5) and McGill group, $\zeta/s\approx 0.09$ from JETSCAPE,
$\zeta/s\approx 0.06$ from our present work, $\zeta/s\approx 0.03$ from Duke
group, and $\zeta/s\approx 0.02$ from Trajectum. In general, our result on
$\zeta/s$ is in good agreement with those from holographic model, SMASH
transport model and JETSCAPE, and has minor overlap with that from Duke and
Trajectum. The discrepancies between Duke, JETSCAPE and Trajectum mainly come
from different treatments of observations and parameters Nijs  _et al._
(2021a) although they use the similar multi-stage methods. The $\zeta/s$ is
believed to have a peak around $T_{\rm pc}$ and go to zero at sufficiently
high temperatures Karsch  _et al._ (2008); Ryu  _et al._ (2015); Everett  _et
al._ (2021b), while the magnitude of the peak value remains unknown and
depends on the model Shen (2021).

We note that using higher $p_{T}$ ranges of experimental data has a minor
influence on our results. For example, using $p_{T}<2.3$ GeV/$c$ for $\phi$
mesons and $p_{T}<3.6$ GeV/$c$ for $\Omega$ baryons leads to $T=160\pm 8$ MeV,
$4\pi\eta/s=2.14\pm 0.44$ and $\zeta/s=0.04\pm 0.03$ at $90\%$ C.L.. As
expected, the higher $p_{T}$ ranges have more data points and put stronger
constraints on $\eta/s$ and $\zeta/s$. On the other hand, the higher $p_{T}$
ranges will cause larger viscous corrections for $f_{0}$, and the threshold of
applicable $p_{T}$ is not known precisely. Since using higher $p_{T}$ ranges
is at risk to violate the applicability of the model, we make a conservative
choice and use lower ranges in this work. We look forward to using more
precise data from future experimental measurements to reduce the
uncertainties.

We now discuss the uncertainties in our analysis. In general, we can group
them into three categories: (a) Uncertainties from assumptions made in
blastwave parameterization, e.g. the simple ansatz for the flow field and the
recombinte hypersurface, and the Navier-Stokes approximation for shear stress.
(b) Uncertainties from assumptions made in quark recombination formulism, e.g.
the instantaneous hadronization process and simple form for wave function
squared of hadrons. (c) Uncertainties from the errors in experimental data and
the quality of the Gaussian emulator. To quantify uncertainty (a), one can
compare hydro simulations with blastwave, which is what we did in Yang and
Fries (2020). We find the uncertainties ($\approx 0.1/4\pi$) from
approximations used in blastwave is rather small compared with uncertainties
from Bayesian analysis ($\approx 1/4\pi$). To give a rough estimate for
uncertainty (b), we vary the wave function width of mesons and baryons by
increasing their values by 33%, i.e., from the default $\sigma_{M}=0.3$ and
$\sigma_{B}=0.09$ to $\sigma_{M}=0.4$ and $\sigma_{B}=0.12$, and we obtain
$\eta/s=2.0/4\pi$ and $\zeta/s=0.07$, which correspond to minor changes, i.e.,
a decrease of $\eta/s$ by 4% and a increase of $\zeta/s$ by 8%. We note that
for uncertainty (c), one could use more refined recombination formula, such as
formula in Pu  _et al._ (2018) to decrease the uncertainty. Uncertainty (c) is
provided by the MADAI code and is shown in our final results. In addition,
there could be uncertainties from the assumption that $\phi$ and $\Omega$ have
weak hadronic interactions. In principle, a quantitative estimate of the
uncertainties could be obtained by including hadronic transport simulations
after hadronization, but this is beyond the scope of the present work and it
would be interesting to pursue in future.

Finally, we would like to mention that in the present work, the quark phase-
space distribution is parameterized by a viscous blastwave and the $\eta/s$
and $\zeta/s$ of the QGP is only constrained at crossover. An alternative and
perhaps more realistic way is to take the phase-space distribution of quarks
from viscous hydrodynamic simulations, which can include parameterized
$\eta/s(T)$ and $\zeta/s(T)$ of the QGP. By doing this, one may extract
information on the temperature dependence of the $\eta/s$ and $\zeta/s$ of the
QGP from the $\phi$ and $\Omega$ observables. Such hybrid approach combining
hydrodynamic model and quark recombination has been recently proposed in Ref.
Zhao  _et al._ (2020), and in the future, we may follow the similar approach
to investigate the $\eta/s$ and $\zeta/s$ of the QGP from the $\phi$ and
$\Omega$ observables.

##  IV Conclusions

Within the quark recombination model with the quark phase-space distribution
parameterized in a viscous blastwave, we have demonstrated that the $\phi$ and
$\Omega$ produced in relativistic heavy-ion collisions can be used to
constrain the specific shear viscosity $\eta/s$ and the specific bulk
viscosity $\zeta/s$ of the QGP at hadronization. By preforming Bayesian
analyses on the measured transverse-momentum spectra and elliptic flows of
$\phi$ and $\Omega$ in Au+Au collisions at $\sqrt{s_{\rm NN}}=200$ GeV and
Pb+Pb collisions at $\sqrt{s_{\rm NN}}=2.76$ TeV, we obtain
$\eta/s=(2.08^{+1.10}_{-1.09})/4\pi$ and $\zeta/s=0.06^{+0.04}_{-0.04}$ at
$90\%$ C.L. for the baryon-free QGP at crossover temperature $T_{\rm
pc}\approx 160$ MeV. Our work suggests that high quality data of $\phi$ and
$\Omega$ in heavy-ion collisions at various energies may provide important
information on the temperature and baryon density dependence of QGP’s $\eta/s$
and $\zeta/s$.

## acknowledgments

The authors would like to thank Rainer J. Fries, Guang-You Qin and Yifeng Sun
for useful discussions. This work was supported by the National Natural
Science Foundation of China under Grant Nos. 12205182, 12235010 and 11625521,
and the National SKA Program of China No. 2020SKA0120300.

## References

  * Romatschke and Romatschke (2007) P. Romatschke and U. Romatschke, Phys. Rev. Lett. 99, 172301 (2007), arXiv:0706.1522 [nucl-th] .
  * Song and Heinz (2008) H. Song and U. W. Heinz, Phys. Rev. C 77, 064901 (2008), arXiv:0712.3715 [nucl-th] .
  * Kovtun  _et al._ (2005) P. Kovtun, D. T. Son, and A. O. Starinets, Phys. Rev. Lett. 94, 111601 (2005), arXiv:hep-th/0405231 [hep-th] .
  * Ryu  _et al._ (2015) S. Ryu, J. F. Paquet, C. Shen, G. S. Denicol, B. Schenke, S. Jeon, and C. Gale, Phys. Rev. Lett. 115, 132301 (2015), arXiv:1502.01675 [nucl-th] .
  * Ryu  _et al._ (2018) S. Ryu, J.-F. Paquet, C. Shen, G. Denicol, B. Schenke, S. Jeon, and C. Gale, Phys. Rev. C 97, 034910 (2018), arXiv:1704.04216 [nucl-th] .
  * Shen (2021) C. Shen, Nucl. Phys. A 1005, 121788 (2021), arXiv:2001.11858 [nucl-th] .
  * Bhattacharya  _et al._ (2014) T. Bhattacharya _et al._ , Phys. Rev. Lett. 113, 082001 (2014), arXiv:1402.5175 [hep-lat] .
  * Borsanyi  _et al._ (2020) S. Borsanyi, Z. Fodor, J. N. Guenther, R. Kara, S. D. Katz, P. Parotto, A. Pasztor, C. Ratti, and K. K. Szabo, Phys. Rev. Lett. 125, 052001 (2020), arXiv:2002.02821 [hep-lat] .
  * Meyer (2007) H. B. Meyer, Phys. Rev. D 76, 101701 (2007), arXiv:0704.1801 [hep-lat] .
  * Meyer (2009) H. B. Meyer, Nucl. Phys. A 830, 641C (2009), arXiv:0907.4095 [hep-lat] .
  * Mages  _et al._ (2015) S. W. Mages, S. Borsányi, Z. Fodor, A. Schäfer, and K. Szabó, PoS LATTICE2014, 232 (2015).
  * Astrakhantsev  _et al._ (2017) N. Astrakhantsev, V. Braguta, and A. Kotov, JHEP 04, 101 (2017), arXiv:1701.02266 [hep-lat] .
  * Borsányi  _et al._ (2018) S. Borsányi, Z. Fodor, M. Giordano, S. D. Katz, A. Pasztor, C. Ratti, A. Schafer, K. K. Szabo, and B. C. Tóth, Phys. Rev. D 98, 014512 (2018), arXiv:1802.07718 [hep-lat] .
  * Christiansen  _et al._ (2015) N. Christiansen, M. Haas, J. M. Pawlowski, and N. Strodthoff, Phys. Rev. Lett. 115, 112002 (2015), arXiv:1411.7986 [hep-ph] .
  * Liu and Rapp (2020) S. Y. F. Liu and R. Rapp, Eur. Phys. J. A 56, 44 (2020), arXiv:1612.09138 [nucl-th] .
  * Ghiglieri  _et al._ (2018) J. Ghiglieri, G. D. Moore, and D. Teaney, JHEP 03, 179 (2018), arXiv:1802.09535 [hep-ph] .
  * Karsch  _et al._ (2008) F. Karsch, D. Kharzeev, and K. Tuchin, Phys. Lett. B 663, 217 (2008), arXiv:0711.0914 [hep-ph] .
  * Meyer (2008) H. B. Meyer, Phys. Rev. Lett. 100, 162001 (2008), arXiv:0710.3717 [hep-lat] .
  * Astrakhantsev  _et al._ (2018) N. Y. Astrakhantsev, V. V. Braguta, and A. Y. Kotov, Phys. Rev. D 98, 054515 (2018), arXiv:1804.02382 [hep-lat] .
  * Rougemont  _et al._ (2017) R. Rougemont, R. Critelli, J. Noronha-Hostler, J. Noronha, and C. Ratti, Phys. Rev. D 96, 014032 (2017), arXiv:1704.05558 [hep-ph] .
  * Moore (2020) G. D. Moore, in  _Criticality in QCD and the Hadron Resonance Gas_ (2020) arXiv:2010.15704 [hep-ph] .
  * Schenke  _et al._ (2011) B. Schenke, S. Jeon, and C. Gale, Phys. Rev. Lett. 106, 042301 (2011), arXiv:1009.3244 [hep-ph] .
  * Gale  _et al._ (2013) C. Gale, S. Jeon, B. Schenke, P. Tribedy, and R. Venugopalan, Phys. Rev. Lett. 110, 012302 (2013), arXiv:1209.6330 [nucl-th] .
  * Niemi  _et al._ (2016) H. Niemi, K. J. Eskola, and R. Paatelainen, Phys. Rev. C 93, 024907 (2016), arXiv:1505.02677 [hep-ph] .
  * Bernhard  _et al._ (2015) J. E. Bernhard, P. W. Marcy, C. E. Coleman-Smith, S. Huzurbazar, R. L. Wolpert, and S. A. Bass, Phys. Rev. C 91, 054910 (2015), arXiv:1502.00339 [nucl-th] .
  * Bernhard  _et al._ (2016) J. E. Bernhard, J. S. Moreland, S. A. Bass, J. Liu, and U. Heinz, Phys. Rev. C 94, 024907 (2016), arXiv:1605.03954 [nucl-th] .
  * Bernhard  _et al._ (2019) J. E. Bernhard, J. S. Moreland, and S. A. Bass, Nature Phys. 15, 1113 (2019).
  * Everett  _et al._ (2021a) D. Everett _et al._ (JETSCAPE), Phys. Rev. Lett. 126, 242301 (2021a), arXiv:2010.03928 [hep-ph] .
  * Everett  _et al._ (2021b) D. Everett _et al._ (JETSCAPE), Phys. Rev. C 103, 054904 (2021b), arXiv:2011.01430 [hep-ph] .
  * Nijs  _et al._ (2021a) G. Nijs, W. van der Schee, U. Gürsoy, and R. Snellings, Phys. Rev. Lett. 126, 202301 (2021a), arXiv:2010.15130 [nucl-th] .
  * Nijs  _et al._ (2021b) G. Nijs, W. van der Schee, U. Gürsoy, and R. Snellings, Phys. Rev. C 103, 054909 (2021b), arXiv:2010.15134 [nucl-th] .
  * Nijs and van der Schee (2022) G. Nijs and W. van der Schee, Phys. Rev. Lett. 129, 232301 (2022), arXiv:2206.13522 [nucl-th] .
  * Parkkila  _et al._ (2021a) J. E. Parkkila, A. Onnerstad, and D. J. Kim, Phys. Rev. C 104, 054904 (2021a), arXiv:2106.05019 [hep-ph] .
  * Parkkila  _et al._ (2021b) J. E. Parkkila, A. Onnerstad, F. Taghavi, C. Mordasini, A. Bilandzic, and D. J. Kim, (2021b), arXiv:2111.08145 [hep-ph] .
  * Heffernan  _et al._ (2023) M. R. Heffernan, C. Gale, S. Jeon, and J.-F. Paquet, (2023), arXiv:2302.09478 [nucl-th] .
  * Greco  _et al._ (2003a) V. Greco, C. M. Ko, and P. Levai, Phys. Rev. Lett. 90, 202302 (2003a), arXiv:0301093 .
  * Greco  _et al._ (2003b) V. Greco, C. M. Ko, and P. Levai, Phys. Rev. C 68, 034904 (2003b), arXiv:0305024 [nucl-th] .
  * Fries  _et al._ (2003a) R. J. Fries, B. Muller, C. Nonaka, and S. A. Bass, Phys. Rev. Lett. 90, 202303 (2003a), arXiv:0301087 [nucl-th] .
  * Fries  _et al._ (2003b) R. J. Fries, B. Muller, C. Nonaka, and S. A. Bass, Phys. Rev. C 68, 044902 (2003b), arXiv:0306027 [nucl-th] .
  * Fries  _et al._ (2008) R. J. Fries, V. Greco, and P. Sorensen, Ann. Rev. Nucl. Part. Sci. 58, 177 (2008), arXiv:0807.4939 [nucl-th] .
  * He  _et al._ (2010) M. He, R. J. Fries, and R. Rapp, Phys. Rev. C 82, 034907 (2010), arXiv:1006.1111 [nucl-th] .
  * Teaney (2003) D. Teaney, Phys. Rev. C 68, 034913 (2003), arXiv:nucl-th/0301099 .
  * Jaiswal and Koch (2015) A. Jaiswal and V. Koch, (2015), arXiv:1508.05878 [nucl-th] .
  * Yang and Fries (2017) Z. Yang and R. J. Fries, J. Phys. Conf. Ser. 832, 012056 (2017), arXiv:1612.05629 [nucl-th] .
  * Yang and Fries (2020) Z. Yang and R. J. Fries, (2020), arXiv:2007.11777 [nucl-th] .
  * Yang and Fries (2022) Z. Yang and R. J. Fries, Phys. Rev. C 105, 014910 (2022), arXiv:1807.03410 [nucl-th] .
  * Shor (1985) A. Shor, Phys. Rev. Lett. 54, 1122 (1985).
  * van Hecke  _et al._ (1998) H. van Hecke, H. Sorge, and N. Xu, Phys. Rev. Lett. 81, 5764 (1998), arXiv:9804035 .
  * Chen and Ko (2006) L.-W. Chen and C. M. Ko, Phys. Rev. C 73, 044903 (2006), arXiv:nucl-th/0602025 .
  * Chen  _et al._ (2008) J. H. Chen, F. Jin, D. Gangadharan, X. Z. Cai, H. Z. Huang, and Y. G. Ma, Phys. Rev. C 78, 034907 (2008), arXiv:0801.2265 [nucl-th] .
  * Auvinen  _et al._ (2017) J. Auvinen, K. Redlich, and S. A. Bass, J. Phys. Conf. Ser. 779, 012045 (2017), arXiv:1610.00589 [nucl-th] .
  * Hwa and Zhu (2017) R. C. Hwa and L. Zhu, J. Phys. Conf. Ser. 779, 012049 (2017), arXiv:1609.06365 [nucl-ex] .
  * Ye  _et al._ (2017) Y. J. Ye, J. H. Chen, Y. G. Ma, S. Zhang, and C. Zhong, Chin. Phys. C 41, 084101 (2017), arXiv:1706.00894 [nucl-th] .
  * Pu  _et al._ (2018) J. Pu, K.-J. Sun, and L.-W. Chen, Phys. Rev. C 98, 064905 (2018), arXiv:1808.04053 [nucl-th] .
  * Song  _et al._ (2020) J. Song, F.-l. Shao, and Z.-t. Liang, Phys. Rev. C 102, 014911 (2020), arXiv:1911.01152 [nucl-th] .
  * Retiere and Lisa (2004) F. Retiere and M. A. Lisa, Phys. Rev. C 70, 044907 (2004), arXiv:0312024 [nucl-th] .
  * Damodaran  _et al._ (2017) M. Damodaran, D. Molnar, G. G. Barnaföldi, D. Berényi, and M. Ferenc Nagy-Egri, (2017), arXiv:1707.00793 [nucl-th] .
  * Denicol  _et al._ (2014) G. S. Denicol, S. Jeon, and C. Gale, Phys. Rev. C 90, 024912 (2014), arXiv:1403.0962 [nucl-th] .
  * Paquet  _et al._ (2016) J.-F. Paquet, C. Shen, G. S. Denicol, M. Luzum, B. Schenke, S. Jeon, and C. Gale, Phys. Rev. C 93, 044906 (2016), arXiv:1509.06738 [hep-ph] .
  * Nerlo-Pomorska and Pomorski (1994) B. Nerlo-Pomorska and K. Pomorski, Z. Phys. A 348, 169 (1994), arXiv:nucl-th/9401015 .
  * Abelev  _et al._ (2013) B. Abelev _et al._ (ALICE), Phys. Rev. C 88, 044909 (2013), arXiv:1301.4361 [nucl-ex] .
  * Adamczewski-Musch  _et al._ (2018) J. Adamczewski-Musch _et al._ (HADES), Eur. Phys. J. A 54, 85 (2018), arXiv:1712.07993 [nucl-ex] .
  * Adams  _et al._ (2007) J. Adams _et al._ (STAR), Phys. Rev. Lett. 98, 062301 (2007), arXiv:nucl-ex/0606014 .
  * Abelev  _et al._ (2007) B. I. Abelev _et al._ (STAR), Phys. Rev. Lett. 99, 112301 (2007), arXiv:nucl-ex/0703033 .
  * Adamczyk  _et al._ (2016) L. Adamczyk _et al._ (STAR), Phys. Rev. Lett. 116, 062301 (2016), arXiv:1507.05247 [nucl-ex] .
  * Abelev  _et al._ (2014) B. B. Abelev _et al._ (ALICE), Phys. Lett. B 728, 216 (2014), [Erratum: Phys.Lett.B 734, 409–410 (2014)], arXiv:1307.5543 [nucl-ex] .
  * Abelev  _et al._ (2015) B. B. Abelev _et al._ (ALICE), JHEP 06, 190 (2015), arXiv:1405.4632 [nucl-ex] .
  * Adam  _et al._ (2017) J. Adam _et al._ (ALICE), Phys. Rev. C 95, 064606 (2017), arXiv:1702.00555 [nucl-ex] .
  * Huovinen and Petreczky (2010) P. Huovinen and P. Petreczky, Nucl. Phys. A 837, 26 (2010), arXiv:0912.2541 [hep-ph] .
  * Bass  _et al._ (2016) S. A. Bass _et al._ , “MADAI collaboration: http://madai.phy.duke.edu/,” (accessed June, 2016).
  * Dash  _et al._ (2019) A. Dash, S. Samanta, and B. Mohanty, Phys. Rev. D 100, 014025 (2019), arXiv:1905.07130 [nucl-th] .
  * Noronha-Hostler  _et al._ (2009) J. Noronha-Hostler, J. Noronha, and C. Greiner, Phys. Rev. Lett. 103, 172302 (2009), arXiv:0811.1571 [nucl-th] .
  * Rose  _et al._ (2020) J. B. Rose, J. M. Torres-Rincon, and H. Elfner, J. Phys. G 48, 015005 (2020), arXiv:2005.03647 [hep-ph] .
  * Auvinen  _et al._ (2020) J. Auvinen, K. J. Eskola, P. Huovinen, H. Niemi, R. Paatelainen, and P. Petreczky, Phys. Rev. C 102, 044911 (2020), arXiv:2006.12499 [nucl-th] .
  * Mykhaylova  _et al._ (2019) V. Mykhaylova, M. Bluhm, K. Redlich, and C. Sasaki, Phys. Rev. D 100, 034002 (2019), arXiv:1906.01697 [hep-ph] .
  * Zhao  _et al._ (2020) W. Zhao, C. M. Ko, Y.-X. Liu, G.-Y. Qin, and H. Song, Phys. Rev. Lett. 125, 072301 (2020), arXiv:1911.00826 [nucl-th] .


Method

[width=0.6]images/3d_tokens.png

Figure : The illustration of 3D tokens and shifted windows mechanism. At
first, each frame is splitted into patches. The patches are then grouped into
windows. In the two consecutive attention layers, the windows configuration
are then shifted. In this way, the attention can happen across windows while
keeping the computation cheap, because the attention are only calculated
within window (not calculated globally)

.

In this paper, we propose a novel method to estimate EF from a cardiac
ultrasound video. Our method uses a deep learning model, based on hierarchical
vision Transformers. We named the model as UltraSwin. UltraSwin adopts
Transformers [Vaswani2017], a popular deep learning model in Natural Language
Processing (NLP) and its derivative work on Computer Vision (CV)
[Dosovitskiy2020, Liu2022].

##  Model Architecture

The architecture of UltraSwin is described in figure LABEL:fig:ultraswin-
overall. The model is received a video as an input, specifically ultrasound
video containing short time cardiac recording. The output of the model is the
estimation of EF for the cardiac in the ultrasound video. UltraSwin has 2 main
modules, Transformers Encoder (TE) and EF Regressor. The TE module acts as
feature extractor while EF Regressor as regressor head. The TE module is used
to learn representation from input video and then output the feature maps.
They are then processed by EF Regressor and transformed into scalar value.
This value is then used as an EF estimation for the input ultrasound video.
Instead of treat a frame of the video as input token like in UVT [101007],
UltraSwin uses 3D video patches as input token following work of Liu et.
al.[Liu2021].

##  Pre-processing

In this research, we use ultrasound video from EchoNet-Dynamic dataset
[Ouyang2019]. This dataset contains echocardiography videos with variety of
frame length and contains at least one heart-beat. Although, it can have more
than one heart-beat per video, ES (End Systolic) and ED (End Dyastolic)
ground-truth are only for one heart-beat. Each frame in the video have spatial
dimension of $112\times 112$ pixels. The frame’s width and height must satisfy
$2^{n}$, so it can be processed into patches. Each input for Transformers must
have same length, so we cut the video into fixed length of $128$ frames. We
choose $128$, because it is the closest $2^{n}$ value from $112$. We select ES
and ED frames and any frames between them then cut them. The Echonet Dynamic
dataset that we used in this research, contains varied length (total frames),
frame rate and image quality [Ouyang2019]. Therefore, if the total frames in
the cut out video are more than $128$, we subsample the frames between ES and
ED. Otherwise, we repeat or mirroring the frames between ES and ED and place
them after ED to get the total $128$ frame length. Suppose the sequence of
frames $F=[m_{ES},m_{b_{1}},\cdots{},m_{b_{n}},m_{ED}]$, we repeat the frames
between ED and ES to create new sequences
$\hat{F}=[m_{ES},m_{b_{1}},\cdots{},m_{b_{n}},m_{ED},m_{b_{1}},\cdots{}]$. We
choose this technique based on the research by Reynaud et al., where the
mirroring technique gives better results than the random sampling technique in
terms of fitting the number of total frames to 128 [101007]. After that, we
pad the frames with blank pixels, so its dimension becomes $128\times 128$
pixels. We also tried to augment the video with standard augmentations like
horizontal flip, vertical flip, random rotation and others. Suprisingly, we
found that augmentations lead to worse performance. This result indicates that
ultrasound video dataset are sensitive to augmentation operations.

##  Transformers Encoder

This module contains 4 stages. Unlike Vision Transformers (ViT)
[Dosovitskiy2020] that has fixed patch size along the stages, UltraSwin use
hierarchical architecture following Swin Transformers [Liu2022] in the spatial
spaces. At every stage, the patch size is downsampled into half of the patch
size in the previous stages. To make the model learn temporal information,
UltraSwin follows Video Swin Transformers [Liu2021] to process the video input
in the shape the 3D patches. TE module contains two main components.

###  3D Patch Partition

Suppose an input video has dimension of $T\times H\times W\times 3$, where
$T,H,W$ and $3$ represent number of frames, frame’s height, frame’s width and
number of channels respectively. The video input is then partitioned into 3D
patch with dimension $2\times 4\times 4\times 3$. In the Transformers world,
this 3D patch is called token. Each token contains embedding features with
length $96$. Actually, we can use any number other than $96$, but greater
number can significantly affects the computation cost. We use $96$ following
Video Swin Transformers[Li2021] as it gives a good performance. This process
yields $\frac{T}{2}\times\frac{H}{4}\times\frac{W}{4}$ tokens in total. The
tokens are then flattened into sequences before processed by the Transformers.
Figure  illustrates the 3D tokens. The features of each token are then
transformed by a linear layer into an arbritrary $C$ dimensions. So, the
dimension of the tokens is now
$\frac{T}{2}\times\frac{H}{4}\times\frac{W}{4}\times C$. This number is hyper-
parameter and we can use arbitrary number for $C$.

###  Block Swin Transformers

Transformers [Vaswani2017] and Vision Transformers (ViT) [Dosovitskiy2020] use
global self-attention (SA) and compute softmax score between each tokens, thus
making the computation and memory resources grow quadratically with token
length. This approach is efficient enough for single input image. On the other
side, the video have multiple image frames, so the approach are not suitable
for video related tasks like video classification, video segmentation and
others. UltraSwin use local window self-attention following Swin Transformers
[Liu2022] that is proven more efficient in video related task than global
self-attention [Li2021]. While efficient, local window self-attention lacks of
connection accross window. This can cause performance degradation on the
model. To solve this issue, UltraSwin shifts the window partition in two
consecutive Swin Transformers Block as illustrated in figure . As Transformers
can have multiple layers of blocks, UltraSwin shifted the window configuration
in every two consecutive blocks. This design is proven to be effective in
image recognition task [Liu2022]. The main reason why it is effective is
because it enables the connections between non-overlapping windows with their
neighbours. Suppose a sequence of 3D tokens with size $T^{{}^{\prime}}\times
H^{{}^{\prime}}\times W^{{}^{\prime}}\times 3$. In the first layer, these
tokens are then arranged into regular non-overlapping window of size $P\times
M\times M$, thus resulting
$\lceil{}\frac{T^{{}^{\prime}}}{P}\rceil{}\times\lceil{}\frac{H^{{}^{\prime}}}{M}\rceil{}\times\lceil{}\frac{W^{{}^{\prime}}}{M}\rceil{}$
non-overlapping 3D windows in total. In the second layer, configuration of
every window is shifted within width, height and temporal axes by
$(\frac{P}{2}\times\frac{M}{2}\times\frac{M}{2})$. The self-attention
mechanism is applied multiple times in parallel. This is called heads. In
multi-head scenario, the output from each self-attention are concatenated. In
first layer, we called multi-head self-attention (MSA) and shifted window
multi-head self-attention (SW-MSA) for second layer. Formally, we stated MSA
as $[\text{SA}_{1},\text{SA}_{2},\cdots{},\text{SA}_{n}]$ and SW-MSA as
$[\text{SW-SA}_{1},\text{SW-SA}_{2},\cdots{},\text{SW-SA}_{n}]$, where
$\text{SA}_{i}$ and $\text{SW-SA}_{i}$ refer to self-attention in layer-$i$
and shifted-window self-attention in layer-$i$ respectivelly. The SA itself
can be formulated as follow: SA(Q,K,V) &= SoftMax(QKTd + B)V where
$K,V,Q\in\mathbb{R}^{PM^{2}\times d}$ are matrices for _key_ , _value_ and
_query_ respectively, while $d$ is _query_ and _key_ dimension, $PM^{2}$ is
the number of tokens in 3D window, and $B\in\mathbb{R}^{P^{2}\times
M^{2}\times M^{2}}$ is matrix of relative position bias.

Table : UltraSwin variants

Parameter | UltraSwin-base | UltraSwin-small
---|---|---
embedding dimension | 128 | 96
number of head | 4, 8, 16, 321 | 3, 6, 12, 241
layer depth | 2, 2, 18, 21 | 2, 2, 18, 21
Total parameter | 88.2M | 49.7M

These values are for stage 1, 2, 3 and 4 respectivelly

The self-attention blocks are then followed by feed forward networks, which is
2 layers MLP with GELU [Hendrycks2016] non-linearity in between. Layer
Normalization (LN) [Ba2016] is applied before self-attention module and before
the MLP. Residual connection [7780459] is then applied after self-attention
and after MLP. Two consecutives of Swin Transformers blocks in layer-$l$ and
layer-$l+1$ can be formulated as follow: ^z^l &= MSA(LN(z^l-1)) + z^l-1
z^l = MLP(GELU(LN(^z^l))) + ^z^l
^z^l+1 = SW-MSA(LN(z^l)) + z^l
z^l+1 = MLP(GELU(LN(^z^l+1))) + ^z^l+1

##  EF Regressor

The EF regressor take the output of the TE module as input. The input is a
features map with dimension
$\frac{T}{2}\times\frac{H}{32}\times\frac{W}{32}\times 8C$. The temporal axes
are then reduced from the map, resulting new dimension
$\frac{H}{32}\times\frac{W}{32}\times 8C$. A linier layer is then applied to
reduce the last axes of feature maps from $8C$ into $4C$. A Layer
Normalization (LN) is then applied followed by linier layer to reduce the
feature axes into $1$ dimension. Spatial reduction is then applied to the map
resulting $1\times 1$ dimesion scalar. This scalar value is then used as EF
estimation.

##  Model Variants

We propose two variants of UltraSwin: UltraSwin-base and UltraSwin-small.
Table  summarizes the two variants. Number of head and layer depth values in
table  refer to configurations on stage $1,2,3$ and $4$ respectivelly. The
total parameter for UltraSwin-small is almost half from total parameter for
UltraSwin-base.

##  Loss Function

Both variants are trained to minimize MSE (Mean Squared Error). We use MSE
because it is commonly used in regression task and gives the best
performances. MSE is defined as follow: L(y,^y) &= 1N ∑_i=1^N(y_i - ^y_i)^2
where $y_{i}$ and $\hat{y}_{i}$ refer to EF ground-truth and EF prediction
from the model respectivelly.


# The strongly intensive observable in pp collisions at LHC energies in the
string fusion model

V V Vechernin1 S N Belokurova1 1 Physics Faculty, Saint-Petersburg State
University, 7-9 Universitetskaya emb., 199034 St.Petersburg, Russia
<EMAIL_ADDRESS>
###### Abstract

The properties of the strongly intensive variable characterizing correlations
between the number of particles in two separated rapidity interval in pp
interactions at LHC energies are studied in the framework of the string fusion
model. We perform the MC simulations of string distributions in the impact
parameter plane to take into account the experimental conditions of pp
collisions. We account the string fusion processes, leading to the formation
of string clusters, embedding a finite lattice (a grid) in the impact
parameter plane. As a result, we found the dependence of this variable both on
the distance between the centers of the observation windows and their width
for the minbias pp collisions at several initial energies. Analyzing these
dependencies we can extract the important information on the properties of
string clusters. We show that in pp collisions at LHC energies the string
fusion effects have a significant impact on the behavior of this strongly
intensive variable. The role of these effects is increasing with the initial
energy and centrality of collisions. In particular, we found that the increase
of this variable with initial energy takes place due to the growth of the
portion of the fused string clusters in string configurations arising in pp
interactions.

##  1 Introduction

By now, the string (color flux tubes) model [1, 2] has become a standard
approach for the description of the soft part of hadronic interactions at high
energies. Different versions of the string model are applied in the various MC
event generators: PYTHIA, VENUS, HIJNG, AMPT, EPOS etc., for a description of
soft processes in strong interactions, when the perturbative QCD approach does
not works.

Along with this, It was also recognized that at a large string density in
nucleus-nucleus collisions the strings start overlap in the impact parameter
plain due to finite transverse area of a string, considered as a color flux
tube. So one has to take into account the interaction between strings leading
to the string fusion processes and the formation of string clusters (“color
ropes”) [3, 4]. In the present paper we’ll argue that the string fusion
processes play important role also in pp interactions at LHC energies.

As one of the tools for the investigation of these string fusion effects the
study of the forward-backward (FB) correlations between observables in two
separated rapidity intervals was suggested in [5]. Really, it is known that
the investigations of the long-range rapidity correlations enable to obtain
the information on the very initial stages of a hadronic collision [6], and,
in particular, on the emerging string configuration.

Unfortunately, back in the work [7] it was shown that the traditional FB
correlation coefficient between the charge particle multiplicities, ${n_{F}}$
and ${n_{B}}$, in the observation rapidity windows strongly depends on the
event-by-event variance of the number of cut pomerons (strings) in pp
collisions, i.e. on the so-called “volume” fluctuation \- the trivial
fluctuation in the number of sources. From this point of view it is desirable
to search for another observables, which is not sensitive to the fluctuation
in the number of sources (strings), but is sensitive to the fluctuation of the
properties of sources, e.g. to the formation of string clusters by string
fusion processes.

We can suppress the influence of these trivial “volume” fluctuations going to
the studies of the FB correlation coefficient between the intensive
quantities, such as e.g. event-mean transverse momenta ${p_{F}}$ and ${p_{B}}$
of charge particles in the observation rapidity intervals ${\delta y}_{F}$ and
${\delta y}_{B}$ instead of ${n_{F}}$ and ${n_{B}}$, as in [8, 9], or defining
the more sophisticated correlation observable between ${n_{F}}$ and ${n_{B}}$.
In the present paper we’ll focus on the last approach.

##  2 The strongly intensive observable ${\Sigma({n_{F}},{n_{B}})}$

The general methods of the constructing of variables not affected by the
“volume” fluctuations (the so-called strongly intensive observables) were
developed in the paper [10]. Later [11] it was suggested to study the strongly
intensive quantity ${\Sigma({n_{F}},{n_{B}})}$, defined for the charge
particle multiplicities, ${n_{F}}$ and ${n_{B}}$, in two observation windows:

$\Sigma({n_{F}},{n_{B}})\equiv[\langle{n_{F}}\rangle\,\omega_{n_{B}}+\langle{n_{B}}\rangle\,\omega_{n_{F}}-2\,{\textrm{cov}}({n_{F}}\,{n_{B}})]/[\langle{n_{F}}\rangle+\langle{n_{B}}\rangle]\
,$ (1)

where the $\omega_{n}\equiv D_{n}/\langle n\rangle$ is a scaled variance and
${\textrm{cov}}({n_{F}}\,{n_{B}})\equiv\langle{n_{F}}{n_{B}}\rangle-\langle{n_{F}}\rangle\langle{n_{B}}\rangle$.
For symmetric reaction and symmetric observation windows, ${\delta
y}_{F}$=${\delta y}_{B}$=${\delta y}$, we have
$\langle{n_{F}}\rangle=\langle{n_{B}}\rangle\equiv\langle n\rangle$,
$\omega_{n_{F}}=\omega_{n_{B}}\equiv\omega_{n}$ and the definition (1) can be
simplified to

${\Sigma({n_{F}},{n_{B}})}={\omega_{n}}-{\textrm{cov}}({n_{F}},{n_{B}})/\langle
n\rangle=[\langle n^{2}\rangle-\langle{n_{F}}{n_{B}}\rangle]/{\langle
n\rangle}\ .$ (2)

In paper [13] it was demonstrated that in the model with independent identical
strings the observable ${\Sigma({n_{F}},{n_{B}})}$ can be expressed only
through the parameters of a single string - the one- and two-particle rapidity
distributions arising from decays of a string:

$\lambda(y)\equiv{dN}/{dy}=\langle\mu\rangle/{\delta y}=\mu_{0}\ ,\hskip
14.22636pt\lambda_{2}(y_{1},y_{2})\equiv{d^{2}N}/{dy_{1}\,dy_{2}}=\lambda_{2}(y_{1}\\!-\\!y_{2})\
.$ (3)

Instead of the $\lambda_{2}(y_{1},y_{2})$ one usually exploits the two-
particle correlation function of a string:

$\Lambda(\eta_{1},\eta_{2})\equiv{\lambda_{2}(\eta_{1},\eta_{2})}/{\lambda(\eta_{1})\lambda(\eta_{2})}-1={\lambda_{2}(y_{1}\\!-\\!y_{2})}/{\mu_{0}^{2}}-1=\Lambda(\eta_{1}-\eta_{2})\
,$ (4)

The last transitions in these formulas corresponds to the rapidity invariant
approximation, which is fulfilled at mid-rapidities at LHC energies. It was
shown in [13] that in the model with independent identical strings

${\Sigma({n_{F}},{n_{B}})}={\Sigma({\mu_{F}},{\mu_{B}})}\ .$ (5)

Here the ${\Sigma({\mu_{F}},{\mu_{B}})}$ is the strongly intensive variable
$\Sigma$ defined similarly (1), but for a single source: the ${\mu_{F}}$ and
${\mu_{B}}$ are the number of particles produced in forward and backward
rapidity windows from decays of a single string. It’s connected with the
string parameters as follows

${\Sigma({\mu_{F}},{\mu_{B}})}=1+\mu_{0}\,{\delta
y}\,[J_{FF}-J_{FB}]=\Sigma({\delta y},{\Delta y})\ ,$ (6)

where

${J_{F\\!B}}=\frac{1}{{\delta y_{F}}{\delta y_{B}}}\int_{{\delta
y_{F}}}\\!\\!\\!\\!\\!\\!dy_{1}\int_{{\delta y_{B}}}\\!\\!\\!\\!\\!\\!dy_{2}\
\Lambda(y_{1}\\!-\\!y_{2})\ ,\hskip 28.45274pt{J_{F\\!F}}=\frac{1}{{\delta
y^{2}_{F}}}\int_{{\delta y_{F}}}\\!\\!\\!\\!\\!\\!dy_{1}\int_{{\delta
y_{F}}}\\!\\!\\!\\!\\!\\!dy_{2}\ \Lambda(y_{1}\\!-\\!y_{2})\ .$ (7)

For symmetric case it depends on the acceptance of the observation windows,
${\delta y}$=${\delta y}_{F}$=${\delta y}_{B}$, and the rapidity distance,
${\Delta y}$, between the centers of these windows, so ${\Delta y}\geq{\delta
y}$.

The formula (5) proves that in the framework of the model with independent
identical strings the variable ${\Sigma({n_{F}},{n_{B}})}$ demonstrates the
strongly intensive property. Its value is equal to the value for a single
string ${\Sigma({\mu_{F}},{\mu_{B}})}$ and does not depend on the number of
strings, produced in a given event and their event-by-event fluctuations.

Unfortunately, the string parameters, $\mu_{0}$ and
$\Lambda(\eta_{1}-\eta_{2})$, extracted in papers [12, 13, 14, 15] from the
ALICE data [17] on the FB correlations in pp collisions at initial energies
0.9-7 TeV prove to be dependent on the energy. By (5)–(7) this leads to the
dependence of the ${\Sigma({n_{F}},{n_{B}})}$ on the collision energy, what
can be interpreted as a signature of string fusion processes in pp collisions
at LHC energies [13, 14, 15, 16].

Really, in the paper [16] it was shown, that in the model with different types
of strings the expression (5) for ${\Sigma({n_{F}},{n_{B}})}$ is replaced by

${\Sigma({n_{F}},{n_{B}})}=\sum_{\eta=1}^{\infty}\alpha_{\eta}\,{\Sigma_{\eta}({\mu_{F}},{\mu_{B}})}\
,$ (8)

where ${\Sigma_{\eta}({\mu_{F}},{\mu_{B}})}$ is a $\Sigma$ variable for a
string cluster, arising due to the fusion of the $\eta$ initial strings. For
each such cluster the formulas (6) and (7) are valid, but with the replacement
of the single string parameters to the parameters of the corresponding cluster
formed by the $\eta$ fused initial strings:

$\mu_{0}\to\mu^{\eta}_{0}\ ,\hskip
28.45274pt\Lambda(y_{1}-y_{2})\to\Lambda_{\eta}(y_{1}-y_{2})\ .$ (9)

The weighting factors $\alpha_{\eta}=\langle{n}\rangle_{\eta}/\langle
n\rangle$ in (8) is an average portion of particles produced from decays of
all clusters of a given $\eta$ type. Note that the first term with $\eta=1$
corresponds to the contribution from single strings.

Obviously, with a change in the initial energy and centrality of pp collisions
the portion of the clusters, corresponding to the different number of fused
strings, will changes, leading to the change of the weights $\alpha_{\eta}$
and hence by (8) to the change of the observable ${\Sigma({n_{F}},{n_{B}})}$.
So, given the possibility of the formation of sources of different types the
observable ${\Sigma({n_{F}},{n_{B}})}$ becomes equal to the weighted average
of the ones for different string clusters with weights, $\alpha_{\eta}$, which
depend on the details of the collision - the initial energy and centrality of
the collision.

##  3 Modeling string configurations in pp collision

For the Monte-Carlo generation of the string configurations arising in pp
collision we have used the approach developed earlier in the paper [18]. In
this approach we suppose that at given value of the impact parameter $b$ for
the number of cut pomerons, $N$, we have event-by-event the poissonian
distribution with the parameter
$\overline{N}(b)=N_{0}\,\exp{(-b^{2}/2r_{0}^{2})}$, discarding the cases with
$N=0$. The last condition reflects the fact that we take into account only
inelastic (non-diffractive, ND) pp collisions.

We will consider that each cut pomeron corresponds to the formation of two
initial strings [1, 2], $N_{str}=2N$. To take into account later the string
fusion effects by the MC calculations we need to simulate not only the total
number of strings in each event, but also their distribution in the transverse
plane. According to [18], we suppose that these strings are distributed in the
transverse plane with the probability density:

$w_{str}({\bf s},{\bf b})\sim T({\bf s}-{\bf b}/2)T({\bf s}+{\bf
b}/2)/{\sigma_{\\!pp}(b)}\ ,$ (10)

where the $w_{str}({\bf s},{\bf b})$ is the probability density to find a
string at the point ${\bf s}$ of the transverse plane in the case of the ND pp
collision at the impact parameter $b$, the $\sigma_{\\!pp}(b)$ is the
probability of such interaction. Here the $T({\bf s})$ is the partonic profile
function of a nucleon, for which we use the simplest gaussian distribution
with some parameter $r_{0}$.

As it was shown in [18], in this model along with MC simulations we can
analytically calculate all quantities of physical interest. So, for the cross-
section of the ND $pp$ interaction we have

$\sigma^{\\!N\\!D}=\int\sigma_{\\!pp}(b)\ d^{2}\vec{b}=2\pi
r_{0}^{2}\left[\Phi(N_{0})+\gamma+\ln{N_{0}}\right]\ ,\hskip
28.45274pt\Phi(x)=\int_{x}^{\infty}e^{-t}dt/t\ ,$ (11)

where the $\gamma=0.577...$ is the Euler constant. For the average number of
cut pomerons, $\langle N\rangle$, the scaled variance of the number of cut
pomerons, ${\omega_{N}}\equiv D_{N}/\langle N\rangle$, and the probability,
$P(N)$, to have $N$ cut pomerons in a non-diffractive pp collision at some
given initial energy we have, respectively:

$\langle N\rangle={N_{0}}/[{\Phi(N_{0})+\gamma+\ln{N_{0}}}]\ ,\hskip
28.45274pt{\omega_{N}}=1+{N_{0}}/{2}-\langle N\rangle\ ,$ (12)
$P(N)=\frac{2\pi
r_{0}^{2}}{\sigma^{\\!N\\!D}_{\\!pp}N}\left[1-e^{-N_{0}}\sum_{l=0}^{N-1}N_{0}^{l}/l!\right]\
.$ (13)

In the paper [18] it was also shown that the present approach leads to the
same distribution (13) for the number of cut pomerons in minbias nondifractive
pp collisions as the Gribov-Regge approach (see e.g. [19]). This enables to
connect the model parameters $N_{0}$ and $r_{0}$, with the parameters of the
pomeron trajectory and its couplings to hadrons:

$r_{0}=\sqrt{{2\lambda}/{C}},\
N_{0}=\frac{2\gamma_{pp}C}{\lambda}\exp(\xi\Delta),\
\lambda=R_{pp}^{2}+\alpha^{\prime}\xi,\ \xi=\ln(s/s_{0})\ .$ (14)

Here $\Delta$ and $\alpha^{\prime}$ are the intercept and the slope of the
pomeron trajectory. The parameters $\gamma_{pp}$ and $R_{pp}$ characterize the
coupling of the pomeron trajectory to the initial hadrons. The quasi-eikonal
parameter $C$ is related to the small-mass diffraction dissociation of
incoming hadrons, $s_{0}\simeq 1GeV^{2}$.

The numerical values of the Regge parameters in the formulae (14), obtaining
by the comparison with experimental data, evolve with encrease of the
available initial pp collision energy. The typical modern values are as
follows (see e.g. [19]):

$\Delta=0.139,\ \ \ \alpha^{\prime}=0.21\ GeV^{-2},\ \ \ \gamma_{pp}=1.77\
GeV^{-2},\ \ \ R_{pp}^{2}=3.18\ GeV^{-2},\ \ \ C=1.5\ .$ (15)

However the test calculations with this set of parameters demonstrated that we
can not explaine the rapid increase of the multiplisity at mid rapidities,
experimentally observed at initial energies higher than 900 GeV. This
corresponds with the results obtained in the paper [22], in which for the same
reasons the increased value of the $\Delta=0.1855$ and the reduced value of
the $\alpha^{\prime}=0.175\ GeV^{-2}$ for the soft pomeron were used.

Table 1: The non-diffractive cross section, the multiplicity density at mid-rapidity and the mean number of initial strings in $pp$ collisions at different initial energies. $\sqrt{s}(GeV)$ | $\sigma^{ND}_{th}(mb)$ | $\sigma^{ND}_{MC}(mb)$ | $dN^{ND}/dy$ | $\langle N_{str}\rangle$
---|---|---|---|---
60 | 24.9 | 24.9 | 2.44 | 4.2
900 | 39.9 | 39.9 | 3.76 | 7.8
7000 | 52.5 | 52.4 | 5.44 | 13.4
13000 | 56.5 | 56.6 | 6.03 | 16.0
|  |  |  |

However, within the framework of our model, these changes of the paremeters
are insufficient, since we have additional losses in multiplicity due to the
increase of string fusion effects with increasing energy. So we used in our
calculations the following set of the parameters:

$\Delta=0.2,\ \ \ \alpha^{\prime}=0.05\ GeV^{-2},\ \ \ \gamma_{pp}=1.035\
GeV^{-2},\ \ \ R_{pp}^{2}=3.3\ GeV^{-2},\ \ \ C=1.5\ .$ (16)

which gives an appropriate description the non-diffractive cross section and
the multiplicity density at mid-rapidities in $pp$ collisions (see the Table
1). The problem with the selection of these parameters is that it is difficult
to obtain a rapid increase in multiplicity with increasing energy, while
maintaining a moderate increase in the nondiffractive cross section.

Note that the mean multiplicities in the Table 1 are calculated already with
taking into account the string fusion effects, as it is explained in the next
section with $\mu^{1}_{0}=\mu_{0}=0.7$. To check the MC algorithm we have
calculated the non-diffractive cross section both analytically by formula
(11), $\sigma^{ND}_{th}$, and by MC simulations as
$\sigma^{ND}_{MC}=S_{b}\,{n_{sim}(N=0)}/{n_{sim}(N\geq 0)}$, where the
$n_{sim}(N=0)$ is a number of simulation without cut pomerons, the
$n_{sim}(N\geq 0)$ is a total number of simulation and $S_{b}$ is an area of
the impact parameter simulations.

##  4 Accounting the string fusion effects

In our MC model to account the string fusion effects we introduce a finite
transverse lattice (a grid) [20, 21] with the cells, which area is equal to
the transverse area of a single initial string, $\sigma_{str}=\pi
r^{2}_{str}$. We use $r_{str}=0.2\ fm$ in our calculation [23].

Generating the event we at first generate the impact parameter $b$, then the
number of pomerons, $N$, and hence the number of strings, $N_{str}=2\,N$, as
it was described in the previous section. After that we distribute these
strings in the transverse plane according to the density given by the formula
(10). We consider that all ${\eta_{i}}$ strings, which centers occur within
the same $i$-th cell, are fused into one string cluster.

In the accordance with the string fusion model prescriptions [4] we suppose
the following dependence of the average number of particles,
$\langle{\mu}\rangle_{\eta}$, produced from the decay of this string cluster
on the number of fused strings $\eta$:

$\langle{\mu}\rangle_{\eta}=\mu^{\eta}_{0}\delta y=\mu_{0}\sqrt{\eta}\,\delta
y\ ,$ (17)

where the $\mu^{\eta}_{0}$ and the $\mu^{1}_{0}\equiv\mu_{0}$ is the average
rapidity particle densities produced from the hadronization of the string
cluster and one single string.

Using generated string configurations we can find the weights $\alpha_{\eta}$
in (8) as follows

$\alpha_{\eta}=\frac{\langle{n}\rangle_{\eta}}{\sum_{\eta=1}^{\infty}\langle{n}\rangle_{\eta}}=\frac{\langle{m}\rangle_{\eta}\mu_{0}^{\eta}{\delta
y}}{\sum_{\eta=1}^{\infty}\langle{m}\rangle_{\eta}\mu_{0}^{\eta}{\delta
y}}=\frac{\langle{m}\rangle_{\eta}\sqrt{\eta}}{\sum_{\eta=1}^{\infty}\langle{m}\rangle_{\eta}\sqrt{\eta}}\
,$ (18)

where the $\langle{m}\rangle_{\eta}$ is the mean number of clusters with
$\eta$ fused strings, which we take from our MC simulations of the string
configurations.

To calculate now the ${\Sigma({n_{F}},{n_{B}})}$ by formula (8) we have to
make some assumptions on the two-particle correlation function
$\Lambda_{\eta}(y_{1}-y_{2})$ of a cluster. As in paper [13] we’ll use the
simplest exponential form

${\Lambda_{\eta}}(y_{1}-y_{2})={\Lambda^{\eta}_{0}}\exp(-{|y_{1}-y_{2}|}/{{y^{\eta}_{corr}}})\
,$ (19)

where ${y^{\eta}_{corr}}$ is a characteristic correlation length. In this case
the integrals (7) easily can be calculated analytically and we can find all
${\Sigma_{\eta}({\mu_{F}},{\mu_{B}})}$ by formulas similar to (6).

|
---|---
Figure 1:  The strongly intensive observable ${\Sigma({n_{F}},{n_{B}})}$ (1)
for pp collisions as a function of the rapidity distance ${\Delta y}$ between
the centers of the FB observation windows, for two widths of windows: ${\delta
y}$=0.2 (left panel) and ${\delta y}$=0.4 (right panel), and for two initial
energies: 0.9 TeV (dashed lines) and 7 TeV (solid lines), calculated for
particles with transverse momenta in the interval 0.3-1.5 GeV/c, as in the
experimental analysis in [17].

Since by (17) the multiplicity density from a string cluster decay is
considered proportional to $\sqrt{\eta}$ and the correlations take place only
between the neighbor string cluster segments, then we’ll suppose that the
characteristic correlation length ${y^{\eta}_{corr}}$ decrease with $\eta$ as
$1/\sqrt{\eta}$. Some additional arguments in favor of this assumption were
presented in [13]. So, we suppose the following dependencies on $\eta$ for the
parameters of the cluster correlation function:

${y^{\eta}_{corr}}={y^{1}_{corr}}/\sqrt{\eta}\ ,\hskip
28.45274pt{\Lambda^{\eta}_{0}}=const\ .$ (20)

The values of the parameters ${y^{1}_{corr}}\\!=\\!2.7$ and
${\Lambda^{\eta}_{0}}\\!=\\!0.8$ were chosen so that to obtain a
correspondence with the values of the $\Sigma(n_{F},\ n_{B})$ obtained in [13,
14, 15]. Note that in these papers the $\Sigma(n_{F},\ n_{B})$ was calculated
on the base of the string pair correlation function, $\Lambda(y_{1}-y_{2})$,
extracted in [12] from the ALICE data [17] on the FB correlations in the
approximation of identical (”effective”) strings, which parameters were
depended on initial energy in contrast with present consideration.

The ${\Sigma({n_{F}},{n_{B}})}$ calculated by formula (8) with the listed
parameters is presented in the figure 1 as function of the rapidity distance
${\Delta y}$ between the centers of the observation windows for two values of
the window acceptance, ${\delta y}$=0.2 and 0.4, and for two initial energies,
0.9 and 7 TeV.

It is instructive to compare the parameter dependencies on $\eta$ for string
clusters, given by formulae (17) and (20):

$\mu_{0}^{\eta}=\mu_{0}^{1}\,\sqrt{\eta}\ ,\hskip
28.45274pt\Lambda_{0}^{\eta}=\Lambda_{0}^{1}=const\ ,\hskip
28.45274pty_{corr}^{\eta}=\eta_{corr}^{1}/\sqrt{\eta}\ ,$ (21)

with the ones

$\mu_{0}^{\eta}=\mu_{0}^{1}\,\eta\ ,\hskip
28.45274pt\Lambda_{0}^{\eta}=\Lambda_{0}^{1}/\eta\ ,\hskip
28.45274pty_{corr}^{\eta}=\eta_{corr}^{1}=const\ .$ (22)

for the case, when we have not string fusion in the given transverse cell. In
the last case we find
${\Sigma({n_{F}},{n_{B}})}=\Sigma_{1}({\mu_{F}},{\mu_{B}})$, that in contrast
with (8) does not depends on $\alpha_{\eta}$.

##  5 Conclusion

In present research we performed the MC simulations of string distributions in
the impact parameter plane for pp interactions at few initial energies. We
took into account the string fusion processes leading to the formation of
string clusters, embedding the finite lattice (the grid) in transverse plane.
Making the assumption on the properties of the string cluster correlation
functions we calculate the strongly intensive observable
$\Sigma({n_{F}},{n_{B}})$ for the minbias pp collisions at few initial
energies and analyze its dependence on the distance between the centers of the
observation windows and their width.

We show that the increase of the $\Sigma({n_{F}},{n_{B}})$ with initial energy
of pp collision originates from the growth of the portion of the fused string
clusters in arising string configurations. In the case with different emitting
clusters this observable becomes equal to the combination of the ones for
different clusters with the weights depending on details of the collision -
its energy and centrality. Analyzing these dependencies of the
$\Sigma({n_{F}},{n_{B}})$ we can extract the important information on the
characteristics of string clusters.

Similarly to the considered growth of the observable $\Sigma({n_{F}},{n_{B}})$
with energy, we expect also the enhance of the string fusion effects in the
central pp collisions. The analysis of the arising dependence of this
observable on the pp collision centrality is now in progress (see preliminary
results in [24]).

The research was supported by the Russian Foundation for Basic Research grant
18-02-40075.

## References

## References

  * [1] Kaidalov A B 1982 Phys. Lett. B 116 459
  * [2] Capella A, Sukhatme U P, Tan C I, Tran Thanh Van J 1994 Phys. Rep. 236 225
  * [3] Biro T S, Nielsen H B, Knoll J 1984 Nucl. Phys. B 245 449
  * [4] Braun M A, Pajares C 1993 Nucl. Phys. B 390 542
  * [5] Amelin N S, Armesto N, Braun M A, Ferreiro E G, Pajares C 1994 Phys. Rev. Lett. 73 2813
  * [6] Dumitru A, Gelis F, McLerran L, Venugopalan R 2008 Nucl. Phys. A 810 91
  * [7] Capella A, Krzywicki A 1978 Phys. Rev. D 18 4120
  * [8] Vechernin V 2016 Eur. Phys. J.: Web of Conf. 125 04022
  * [9] Vechernin V, Kovalenko V 2017 J. Phys.: Conf. Ser. 798 012053
  * [10] Gorenstein M I, Gazdzicki M 2011 Phys. Rev. C 84 014904
  * [11] Andronov E V 2015 Theor. Math. Phys. 185 1383
  * [12] Vechernin V 2015 Nucl. Phys. A  939 21
  * [13] Vechernin V 2018 Eur. Phys. J.: Web of Conf. 191 04011
  * [14] Andronov E, Vechernin V 2019 Eur. Phys. J. A 55 14
  * [15] Vechernin V, Andronov E 2019 Universe 5 15
  * [16] Belokurova S N, Vechernin V V 2019 Theor. Math. Phys. 200 1094
  * [17] Adam J, et al. [ALICE Collaboration] 2015 JHEP 05 097
  * [18] Vechernin V, Lakomov I 2013 Proceedings of Science Baldin ISHEPP XXI 072
  * [19] Arakelyan G H, Capella A, Kaidalov A B, Shabelski Yu M 2002 Eur. Phys. J. C 26 81
  * [20] Vechernin V V, Kolevatov R S 2004 Vestn. SPbU: Ser. 4 4 11 (arXiv:0305136 [hep-ph])
  * [21] Braun M A, Kolevatov R S, Pajares C, Vechernin V V 2004 Eur. Phys. J. C  32 535
  * [22] Bleibel J, Bravina L V, Zabrodin E E 2016 Phys. Rev. D 93 114012
  * [23] Kovalenko V N 2017 Phys. Part. Nucl. 48 945
  * [24] S. Belokurova, Study of strongly intense quantities and robust variances in multi-particle production at LHC energies, CERN-STUDENTS-Note-2019-021, CERN, Geneva, 2019 [https://cds.cern.ch/record/2684671]; in Proceedings of the LXX International conference NUCLEUS - 2020, 11-17 October 2020, St.Petersburg, Physics of Particles and Nuclei (in press); arXiv:2011.10434 [hep-ph].