# Interplay between chiral dynamics and repulsive interactions

Michał Marczenko<EMAIL_ADDRESS>Institute of Theoretical
Physics, University of Wrocław, plac Maksa Borna 9, PL-50204 Wrocław, Poland
Krzysztof Redlich  Institute of Theoretical Physics, University of Wrocław,
plac Maksa Borna 9, PL-50204 Wrocław, Poland  Chihiro Sasaki  Institute of
Theoretical Physics, University of Wrocław, plac Maksa Borna 9, PL-50204
Wrocław, Poland

###### Abstract

We investigate fluctuations of the net-baryon number density in hot hadronic
matter. We discuss the interplay between chiral dynamics and repulsive
interactions and their influence on the properties of these fluctuations near
the chiral crossover. The chiral dynamics is modeled by the parity doublet
Lagrangian that incorporates the attractive and repulsive interactions
mediated via the exchange of scalar and vector mesons, respectively. The mean-
field approximation is employed to account for chiral criticality. We focus on
the properties and systematics of the cumulants of the net-baryon number
density up to the sixth order. It is shown that the higher-order cumulants
exhibit a substantial suppression in the vicinity of the chiral phase
transition due to the presence of repulsive interactions. We find, however,
that the characteristic properties of cumulants near the chiral crossover
observed in lattice QCD results are entirely linked to the critical chiral
dynamics and, in general, cannot be reproduced in phenomenological models,
which account only for effective repulsive interactions via excluded-volume
corrections or van-der-Waals type interactions. Consequently, a description of
the higher-order cumulants of the net-baryon density in the chiral crossover
requires a self-consistent treatment of the chiral in-medium effects and
repulsive interactions.

##  I Introduction

Establishing the thermodynamic properties of strongly interacting matter,
described by quantum chromodynamics (QCD), is one of the key directions in
modern high-energy physics. At vanishing density, the first-principles
calculations of lattice QCD (LQCD) provide a reliable description of the
equation of state (EoS) and fluctuations of conserved charges Bazavov  _et
al._ (2014); Borsanyi  _et al._ (2018); Bazavov  _et al._ (2017, 2020, 2020).
There, the EoS exhibits a smooth crossover from confined hadronic matter to a
deconfined quark-gluon plasma, which is linked to the color deconfinement and
the restoration of chiral symmetry Aoki  _et al._ (2006); Bazavov  _et al._
(2019). However, the nature of the EoS at finite density is still not resolved
by LQCD, owing to the sign problem, and remains an open question.

The LQCD results (Aarts  _et al._ , 2015, 2017, 2019) exhibit a clear
manifestation of the parity doubling structure for the low-lying baryons
around the chiral crossover. The masses of the positive-parity ground states
are found to be rather temperature-independent, while the masses of negative-
parity states drop substantially when approaching the chiral crossover
temperature $T_{c}$. The parity doublet states become almost degenerate with a
finite mass in the vicinity of the chiral crossover. Even though these LQCD
results are still not obtained in the physical limit, the observed behavior of
parity partners is likely an imprint of the chiral symmetry restoration in the
baryonic sector of QCD. Such properties of the chiral partners can be
described in the framework of the parity doublet model (De Tar and Kunihiro,
1989; Jido  _et al._ , 2000, 2001). The model has been applied to hot and
dense hadronic matter, neutron stars, as well as the vacuum phenomenology of
QCD (Dexheimer  _et al._ , 2008; Gallas  _et al._ , 2010; Paeng  _et al._ ,
2012; Sasaki  _et al._ , 2011; Gallas  _et al._ , 2011; Zschiesche  _et al._ ,
2007; Benic  _et al._ , 2015; Marczenko and Sasaki, 2018; Marczenko  _et al._
, 2018, 2019; Marczenko, 2020; Marczenko  _et al._ , 2020; Motornenko  _et
al._ , 2020; Mukherjee  _et al._ , 2017a, b; Dexheimer  _et al._ , 2013;
Steinheimer  _et al._ , 2011a; Weyrich  _et al._ , 2015; Sasaki and Mishustin,
2010; Yamazaki and Harada, 2019a, b; Ishikawa  _et al._ , 2019; Steinheimer
_et al._ , 2011b; Giacosa, 2012; Motornenko  _et al._ , 2019; Motohiro  _et
al._ , 2015; Minamikawa  _et al._ , 2020).

It is already confirmed, that at small net-baryon number density, the
thermodynamics of the confined phase of QCD is well-described by the hadron
resonance gas (HRG) model Braun-Munzinger  _et al._ (2004); Andronic  _et al._
(2018). Since the fundamental quarks and gluons are confined, it is to be
expected that at low temperatures the QCD partition function is dominated by
the contribution of hadrons. The HRG model describes well the LQCD data below
the crossover transition to a quark-gluon plasma, as well as the hadron yields
in heavy-ion collisions Andronic  _et al._ (2018). Different extensions of the
HRG model have been proposed to quantify the LQCD EoS and various fluctuation
observable up to near chiral-crossover. They account for consistent
implementation of hadronic interactions within the S-matrix approach
Venugopalan and Prakash (1992); Broniowski  _et al._ (2015); Friman  _et al._
(2015); Huovinen  _et al._ (2017); Lo  _et al._ (2018), a more complete
implementation or a continuously growing exponential mass spectrum and/or
possible repulsive interactions among constituents Majumder and Muller (2010);
Andronic  _et al._ (2012); Albright  _et al._ (2014); Vovchenko  _et al._
(2015); Lo  _et al._ (2015); Man Lo  _et al._ (2016); Andronic  _et al._
(2020). Recently, an interesting suggestion was made that deviations of the
LQCD data on higher-order fluctuations of net-baryon number density from the
HRG baseline in the near vicinity of the chiral transition can be attributed
to repulsive interactions among constituent hadrons Vovchenko  _et al._
(2017).

Fluctuations of conserved charges are known to be auspicious observable for
the search of the chiral-critical behavior at the QCD phase boundary Stephanov
_et al._ (1999); Asakawa  _et al._ (2000); Hatta and Stephanov (2003), and
chemical freeze-out of produced hadrons in heavy-ion collisions Bazavov  _et
al._ (2012); Borsanyi  _et al._ (2014); Karsch and Redlich (2011); Braun-
Munzinger  _et al._ (2015); Vovchenko  _et al._ (2020); Braun-Munzinger  _et
al._ (2020). In particular, fluctuations have been proposed to probe the QCD
critical point in the beam energy scan (BES) programs at the Relativistic
Heavy Ion Collider (RHIC) at BNL and the Super Proton Synchrotron (SPS) at
CERN, as well as the remnants of the $O(4)$ criticality at vanishing and
finite baryon densities Friman  _et al._ (2011); Karsch (2019); Braun-
Munzinger  _et al._ (2020, 2017).

In this work, we analyze the properties and systematics of the fluctuations of
conserved charges in the context of the parity doublet model, which
incorporates the chiral symmetry restoration and repulsive interactions via
the exchange of the scalar and vector mesons, respectively. To account for
critical behavior, the mean-field approximation is employed, which contains
basic features of the $O(4)$ criticality, albeit with different critical
exponents. We study the behavior of the second- and higher-order cumulants of
the net-baryon number density up to the sixth order, as well as the bulk
equation of state. It is systematically examined to what extent the thermal
behaviors are dominated by the chiral criticality and repulsive interactions.

This paper is organized as follows. In Sec. II, we introduce the parity
doublet model. In Sec. III, we discuss the structure of the higher-order
cumulants of the net-baryon number density. In Sec. IV, we present results on
the equation of state and the higher-order cumulants of the net-baryon number
density. Finally, Sec. V is devoted to summary and conclusions.

##  II Parity doublet model

In the conventional Gell-Mann–Levy model of mesons and nucleons Gell-Mann and
Levy (1960), the nucleon mass is entirely generated by the non-vanishing
expectation value of the sigma field. Thus, the nucleon inevitably becomes
massless when the chiral symmetry gets restored. This is led by the particular
chirality-assignment to the nucleon parity doublers, where the nucleons are
assumed to be transformed in the same way as the quarks are under chiral
rotations.

More general allocation of the left- and right-handed chiralities to the
nucleons, the mirror assignment, was proposed in De Tar and Kunihiro (1989).
This allows an explicit mass-term for the nucleons, and consequently, the
nucleons stay massive at the chiral restoration point. For more details, see
Refs. De Tar and Kunihiro (1989); Jido  _et al._ (2000, 2001).

$m_{0}~{}$[MeV] |  $m_{+}~{}$[MeV] |  $m_{-}~{}$[MeV] |  $m_{\pi}~{}$[MeV] |  $f_{\pi}~{}$[MeV] |  $m_{\omega}~{}$[MeV] |  $m_{\rho}~{}$[MeV] | $\lambda_{4}$ | $\lambda_{6}f_{\pi}^{2}$ | $g_{\omega}$ | $g_{1}$ | $g_{2}$
---|---|---|---|---|---|---|---|---|---|---|---
850 | 939 | 1500 | 140 | 93 | 783 | 775 | 13.15 | 4.67 | 4.64 | 12.42 | 6.39
Table 1: Physical vacuum inputs and the parity doublet model parameters used
in this work. See Sec. II for details.

In the mirror assignment, under $SU(2)_{L}\times SU(2)_{R}$ rotation, two
chiral fields $\psi_{1}$ and $\psi_{2}$ are transformed as follows:

$\begin{split}\psi_{1L}\rightarrow L\psi_{1L},\;\;\;\;\psi_{1R}\rightarrow
R\psi_{1R}\textrm{,}\\\ \psi_{2L}\rightarrow
R\psi_{2L},\;\;\;\;\psi_{2R}\rightarrow L\psi_{2R}\textrm{,}\end{split}$ (1)

where $\psi_{i}=\psi_{iL}+\psi_{iR}$, $L\in SU(2)_{L}$ and $R\in SU(2)_{R}$.
The nucleon part of the Lagrangian in the mirror model reads

$\begin{split}\mathcal{L}_{N}&=i\bar{\psi}_{1}\hbox to0.0pt{\hbox
to5.30904pt{\hfil/\hfil}\hss}\partial\psi_{1}+i\bar{\psi}_{2}\hbox
to0.0pt{\hbox
to5.30904pt{\hfil/\hfil}\hss}\partial\psi_{2}+m_{0}\left(\bar{\psi}_{1}\gamma_{5}\psi_{2}-\bar{\psi}_{2}\gamma_{5}\psi_{1}\right)\\\
&+g_{1}\bar{\psi}_{1}\left(\sigma+i\gamma_{5}\boldsymbol{\tau}\cdot\boldsymbol{\pi}\right)\psi_{1}+g_{2}\bar{\psi}_{2}\left(\sigma-i\gamma_{5}\boldsymbol{\tau}\cdot\boldsymbol{\pi}\right)\psi_{2}\\\
&-g_{\omega}\bar{\psi}_{1}\hbox to0.0pt{\hbox
to6.22453pt{\hfil/\hfil}\hss}\omega\psi_{1}-g_{\omega}\bar{\psi}_{2}\hbox
to0.0pt{\hbox
to6.22453pt{\hfil/\hfil}\hss}\omega\psi_{2}\textrm{,}\end{split}$ (2)

where $g_{1}$, $g_{2}$, and $g_{\omega}$ are the baryon-to-meson coupling
constants and $m_{0}$ is a mass parameter.

The mesonic part of the Lagrangian reads

$\begin{split}\mathcal{L}_{M}=\frac{1}{2}\left(\partial_{\mu}\sigma\right)^{2}+\frac{1}{2}\left(\partial_{\mu}\boldsymbol{\pi}\right)^{2}-\frac{1}{4}\left(\omega_{\mu\nu}\right)^{2}-V_{\sigma}-V_{\omega}\textrm{,}\end{split}$
(3)

where $\omega_{\mu\nu}=\partial_{\mu}\omega_{\nu}-\partial_{\nu}\omega_{\mu}$
is the field-strength tensor of the vector field, and the potentials read

$\displaystyle V_{\sigma}$
$\displaystyle=-\frac{\lambda_{2}}{2}\Sigma+\frac{\lambda_{4}}{4}\Sigma^{2}-\frac{\lambda_{6}}{6}\Sigma^{3}-\epsilon\sigma\textrm{,}$
(4a) $\displaystyle V_{\omega}$
$\displaystyle=-\frac{m_{\omega}^{2}}{2}\omega_{\mu}\omega^{\mu}\textrm{.}$
(4b)

where $\Sigma=\sigma^{2}+\boldsymbol{\pi}^{2}$,
$\lambda_{2}=\lambda_{4}f_{\pi}^{2}-\lambda_{6}f_{\pi}^{4}-m_{\pi}^{2}$, and
$\epsilon=m_{\pi}^{2}f_{\pi}$. $m_{\pi}$ and $m_{\omega}$ are the $\pi$ and
$\omega$ meson masses, respectively, and $f_{\pi}$ is the pion decay constant.
Note that the chiral symmetry is explicitly broken by the linear term in
$\sigma$ in Eq. (4a).

The full Lagrangian of the parity doublet model is then

$\mathcal{L}=\mathcal{L}_{N}+\mathcal{L}_{M}\textrm{.}$ (5)

In the diagonal basis, the masses of the chiral partners, $N_{\pm}$, are given
by

$m_{\pm}=\frac{1}{2}\left(\sqrt{\left(g_{1}+g_{2}\right)^{2}\sigma^{2}+4m_{0}^{2}}\mp\left(g_{1}-g_{2}\right)\sigma\right)\textrm{.}$
(6)

From Eq. (6), it is clear that, in contrast to the naive assignment under
chiral symmetry, the chiral symmetry breaking generates only the splitting
between the two masses. When the symmetry is restored, the masses become
degenerate, $m_{\pm}(\sigma=0)=m_{0}$.

To investigate the properties of strongly-interacting matter, we adopt a mean-
field approximation Serot and Walecka (1986). Rotational invariance requires
that the spatial component of the $\omega_{\mu}$ field vanishes, namely
$\langle\boldsymbol{\omega}\rangle=0$#1#1#1Since $\omega_{0}$ is the only non-
zero component in the mean-field approximation, we simply denote it by
$\omega_{0}\equiv\omega$.. Parity conservation on the other hand dictates
$\langle\boldsymbol{\pi}\rangle=0$. The mean-field thermodynamic potential of
the parity doublet model reads

$\Omega=\sum_{x=\pm}\Omega_{x}+V_{\sigma}+V_{\omega}\textrm{,}$ (7)

with

$\Omega_{x}=\gamma_{x}\int\frac{\mathrm{d}^{3}p}{(2\pi)^{3}}\;T\left[\ln\left(1-f_{x}\right)+\ln\left(1-\bar{f}_{x}\right)\right]\textrm{,}$
(8)

where $\gamma_{\pm}=2\times 2$ denotes the spin-isospin degeneracy factor for
both parity partners, and $f_{x}$ $(\bar{f}_{x})$ is the particle
(antiparticle) Fermi-Dirac distribution function,

$\begin{split}f_{x}=\frac{1}{1+e^{\beta\left(E_{x}-\mu^{\ast}\right)}}\textrm{,}\\\
\bar{f}_{x}=\frac{1}{1+e^{\beta\left(E_{x}+\mu^{\ast}\right)}}\textrm{,}\\\
\end{split}$ (9)

with $\beta$ being the inverse temperature, the dispersion relation
$E_{x}=\sqrt{\boldsymbol{p}^{2}+m_{x}^{2}}$ and the effective chemical
potential $\mu^{\ast}=\mu_{B}-g_{\omega}\omega$.

In-medium profiles of the mean fields are obtained by extremizing the
thermodynamic potential in Eq. (7), leading to the following gap equations:

$\begin{split}0=\frac{\partial\Omega}{\partial\sigma}&=\frac{\partial
V_{\sigma}}{\partial\sigma}+s_{+}\frac{\partial
m_{+}}{\partial\sigma}+s_{-}\frac{\partial m_{-}}{\partial\sigma}\textrm{,}\\\
0=\frac{\partial\Omega}{\partial\omega}&=\frac{\partial
V_{\omega}}{\partial\omega}+g_{\omega}\left(n_{+}+n_{-}\right)\textrm{,}\end{split}$
(10)

where the scalar and vector densities are

$s_{\pm}=\gamma_{\pm}\int\frac{\mathrm{d}^{3}p}{\left(2\pi\right)^{3}}\frac{m_{\pm}}{E_{\pm}}\left(f_{\pm}+\bar{f}_{\pm}\right)$
(11)

and

$n_{\pm}=\gamma_{\pm}\int\frac{\mathrm{d}^{3}p}{\left(2\pi\right)^{3}}\left(f_{\pm}-\bar{f}_{\pm}\right)\textrm{,}$
(12)

respectively.

In the grand canonical ensemble, the thermodynamic pressure reads

$P=-\Omega+\Omega_{0}\textrm{,}$ (13)

where $\Omega_{0}$ is the value of the thermodynamic potential in the vacuum,
and the net-baryon number density can be calculated as follows:

$n_{B}=\frac{\partial P(T,\mu_{B})}{\partial\mu_{B}}\textrm{.}$ (14)

The positive-parity state, $N_{+}$, corresponds to the nucleon $N(938)$. Its
negative parity partner is identified with $N(1535)$. Their vacuum masses are
shown in Table 1. The value of the parameter $m_{0}$ has to be chosen so that
a chiral crossover transition is featured at finite temperature and vanishing
chemical potential. The model predicts the chiral transition to be a crossover
for $m_{0}\gtrsim 700~{}$MeV. Following the previous studies of the parity-
doublet-based models (Dexheimer  _et al._ , 2008; Gallas  _et al._ , 2010;
Paeng  _et al._ , 2012; Sasaki  _et al._ , 2011; Gallas  _et al._ , 2011;
Zschiesche  _et al._ , 2007; Benic  _et al._ , 2015; Marczenko and Sasaki,
2018; Marczenko  _et al._ , 2018, 2019; Marczenko, 2020; Marczenko  _et al._ ,
2020; Motornenko  _et al._ , 2020; Mukherjee  _et al._ , 2017a, b; Dexheimer
_et al._ , 2013; Steinheimer  _et al._ , 2011a; Weyrich  _et al._ , 2015;
Sasaki and Mishustin, 2010; Yamazaki and Harada, 2019a, b; Ishikawa  _et al._
, 2019; Steinheimer  _et al._ , 2011b; Giacosa, 2012; Motornenko  _et al._ ,
2019; Motohiro  _et al._ , 2015; Minamikawa  _et al._ , 2020), as well as
recent lattice QCD results (Aarts  _et al._ , 2017, 2019, 2015), we choose a
rather large value, $m_{0}=850$ MeV. We note that, the results presented in
Sec. IV qualitatively do not depend on the choice of $m_{0}$, as long as
chiral crossover transition is featured. The parameters $g_{1}$ and $g_{2}$
are determined by the aforementioned vacuum nucleon masses and the chirally
invariant mass $m_{0}$ via Eq. (6). The parameters $g_{\omega}$, $\lambda_{4}$
and $\lambda_{6}$ are fixed by the properties of the nuclear ground state at
zero temperature, i.e., the saturation density, binding energy and
compressibility parameter at $\mu_{B}=923~{}$MeV. The constrains are as
follows:

$\displaystyle n_{B}$ $\displaystyle=0.16~{}\textrm{fm}^{-3}\textrm{,}$ (15a)
$\displaystyle E/A-m_{+}$ $\displaystyle=-16~{}\textrm{MeV}\textrm{,}$ (15b)
$\displaystyle K=9n^{2}_{B}\frac{\partial^{2}\left(E/A\right)}{\partial
n_{B}^{2}}$ $\displaystyle=240~{}\textrm{MeV}\textrm{.}$ (15c)

We note that the six-point scalar interaction term in Eq. (4a) is essential in
order to reproduce the empirical value of the compressibility in Eq. (15c)
(Motohiro  _et al._ , 2015). The parameters used in this paper are tabulated
in Table 1. For this set of parameters, we obtain the critical temperature of
the chiral crossover transition at vanishing chemical potential,
$T_{c}=208.7~{}$MeV.

In the following, we will also compare our results with the hadron resonance
gas (HRG) Karsch  _et al._ (2003a, b); Karsch (2014); Andronic  _et al._
(2012); Albright  _et al._ (2014, 2015) model formulation of the
thermodynamics of the confined phase of QCD. The thermodynamic potential of
the HRG model is given as a sum of uncorrelated ideal-gas particles:

$\Omega^{\rm HRG}=\sum_{x=\pm}\Omega_{x}\textrm{,}$ (16)

with $\Omega_{x}$ given by Eq. (8), where the masses of $N_{\pm}$ are taken to
be the vacuum masses (see Table 1) and $\mu^{\ast}=\mu_{B}$. We will also
consider a modification of the HRG model, $\sigma$HRG, where the thermodynamic
potential is given as in Eq. (16), but the vacuum masses of $N^{\pm}$ are
substituted by the in-medium masses obtained by solving the parity doublet
model. In both models, the pressure and net-baryon density are obtained
through Eqs. (13) and (14), respectively.

In the next section, we discuss the general structure of the higher-order
cumulants of the net-baryon number.

Figure 1: Thermodynamic pressure (left), and the second-order cumulant of the
net-baryon number density (right), calculated under different schemes, see
text. The temperature is expressed in the units of the chiral-critical
temperature, $T_{c}=T_{c}(\mu_{B}=0)$.

##  III Higher-order cumulants of the net-baryon number

The fluctuations of conserved charges reveal more information about the matter
composition than the equation of state and can be used as probes of a phase
boundary. The critical properties of chiral models, within the functional
renormalization group (FRG) approach Wetterich (1993); Morris (1994);
Ellwanger (1993); Berges  _et al._ (2002), are governed by the same
universality classes as QCD, i.e., the chiral transition belongs to $O(4)$
universality class, which, at large values of the baryon chemical potential,
may develop a $Z(2)$ critical point, followed by the first-order phase
transition Asakawa and Yazaki (1989); Halasz  _et al._ (1998); Berges and
Rajagopal (1999). This criticality is naturally encoded in the hadronic parity
doublet model, as well as in quark-based models Skokov  _et al._ (2010, 2011);
Schaefer and Wambach (2007); Friman  _et al._ (2011); Almasi  _et al._ (2017),
although the mean-field treatment yields different critical exponents.

The main objective of the present studies is to analyze and delineate the
contribution to thermodynamics from chiral dynamics and repulsive baryon-
baryon interactions. To this end, we analyze the fluctuations of the net-
baryon number at finite temperature and vanishing chemical potential.

In the grand canonical ensemble, the cumulants of the net-baryon number,
$\chi_{n}$, are commonly defined as temperature-normalized derivatives w.r.t.
the baryon chemical potential,

$\chi_{n}\left(T,\mu_{B}\right)=T^{n-4}\frac{\partial^{n-1}n_{B}\left(T,\mu_{B}\right)}{\partial\mu_{B}^{n-1}}\textrm{,}$
(17)

where $n_{B}$ is defined in Eq. (14).

In the mean-field approximation, the net-baryon number density, as well as any
other thermodynamic quantity, contains explicit dependence on the mean fields.
Here, we consider only $\sigma$ and $\omega$ mean fields (cf. Eq. (7)), thus
$n_{B}=n_{B}\left(T,\mu_{B},\sigma(T,\mu_{B}),\omega(T,\mu_{B})\right)$.
Consequently, from Eq. (17) one derives the following general form of the
second-order cumulant,

$\chi_{2}=\chi_{2}^{\rm id}\beta_{\rm rep}+\frac{\partial
n_{B}}{\partial\sigma}\frac{\partial\sigma}{\partial\mu_{B}}\textrm{,}$ (18)

where $\chi_{2}^{\rm id}=\chi_{2}^{\rm
id}\left(T,\mu_{B},\sigma(T,\mu_{B}),\omega(T,\mu_{B})\right)$ is the ideal
gas expression for the net-baryon number susceptibility, and

$\beta_{\rm rep}=1-g_{\omega}\frac{\partial\omega}{\partial\mu_{B}}$ (19)

is the suppression factor due to repulsive interactions. From Eq. (18) it is
clear that the non-interacting ideal gas result is retrieved when the mean-
field contribution is neglected.

At vanishing chemical potential, Eq. (18) reduces to

$\chi_{2}=\chi_{2}^{\rm id}\beta_{\rm rep}\textrm{.}$ (20)

We note that, depending on the details of the model, $\chi_{2}^{\rm id}$ in
Eqs. (18) and (20) contains also dependence on the $\sigma$ and $\omega$ mean
fields. However, at vanishing $\mu_{B}$, the expectation value of $\omega$
vanishes as well, i.e., the effective chemical potential is $\mu^{\ast}=0$.
Thus, $\chi_{2}^{\rm id}$ contains only the contribution from the $\sigma$
mean field. Therefore, it encodes the information about attractive
interactions, while the information about repulsive interactions is contained
in the suppression factor $\beta_{\rm rep}$. To some extent, such separation
is qualitatively similar to that of the excluded volume approach. We note that
the expression for the second-order cumulant in Eq. (20) is exact at vanishing
chemical potential.

In similar spirit, one derives the higher-order cumulants as

$\chi_{n}=\chi_{n}^{\rm id}\beta_{\rm rep}^{n-1}+\ldots$ (21)

where $\chi_{n}^{\rm id}$ is the ideal gas expression for the n’th order
cumulant. For $n>2$, Eq. (21) contains extra terms, as explained in Appendix
A. Nevertheless, keeping the first term provides a relatively good
approximation to the full expression #2#2#2In Appendix A, we present the
evaluation of the higher-order cumulants and discuss the comparison of the
full expressions and approximations used in this study.. We note that the
general structure of $\chi_{n}$ derived in Eq. (21) is the same in any kind of
$\sigma-\omega$ model under the mean-field approximation and is independent of
the details of the model.

From Eq. (21), keeping only the first term, we may also estimate the ratio of
the cumulants:

$\frac{\chi_{n}}{\chi_{m}}=\frac{\chi_{n}^{\rm id}}{\chi_{m}^{\rm
id}}\beta_{\rm rep}^{n-m}+\ldots$ (22)

Clearly, the separation of the attractive and repulsive contributions persists
in the approximation of the higher-order cumulants, as well as in their
ratios. This allows to precisely delineate the contribution of chiral symmetry
restoration and repulsive interaction to the critical behavior of the
cumulants in the vicinity of the chiral phase transition.

In the following, we quantify and discuss the properties of the obtained
equation of state, as well as the higher-order cumulants of the net-baryon
number density at vanishing chemical potential in order to identify the
importance of the repulsive interactions near the chiral crossover transition.

##  IV Results

In the left panel of Fig. 1, we show the numerical results on the temperature-
normalized thermodynamic pressure at vanishing baryon chemical potential. The
pressure obtained in the HRG model increases monotonically and does not
resemble any critical behavior. This is expected because $\Omega^{\rm HRG}$ is
just a sum of uncorrelated particles (cf. Eq. (16)) with vacuum hadron masses.
There are clear deviations of the parity doublet model result on thermodynamic
pressure from the corresponding ideal HRG. The increase of the pressure is a
bulk consequence of an interplay between critical chiral dynamics with in-
medium hadron masses and repulsive interactions. The influence of in-medium
hadron masses is identified when considering the pressure $P^{\sigma\rm HRG}$
of the $\sigma$HRG model which increases around $T_{c}$. This is evidently
linked to the in-medium shift of baryon masses due to chiral symmetry
restoration. We note, that $P^{\sigma\rm HRG}$ is systematically higher than
the parity doublet pressure. The reason is that the partition function of
uncorrelated particles does not contain the mean-field potentials which
provide a negative contribution to the pressure, as seen in the parity doublet
model from Eq. (7). All pressures shown in Fig. 1 converge to the Stefan-
Boltzmann limit at high temperatures.

In the right panel of Fig. 1, we show the second-order cumulant, $\chi_{2}$.
We note that, at vanishing chemical potential the expectation value of
$\omega$ is zero, thus $\chi_{n}^{\rm id}$ are equivalent to the $\sigma$HRG
formulation. Similarly to the case of pressure, the result for $\sigma$HRG
lies systematically above the ideal gas result. $\chi_{2}$ in HRG and
$\sigma$HRG models converge to the Stefan-Boltzmann limit at high-
temperatures. In contrast, the parity doublet result saturates above $T_{c}$
and monotonically decreases to zero at high temperature. From Eq. (20), it is
clear that the difference between $\sigma$HRG and parity doublet results is
due to the suppression originating from $\beta_{\rm rep}$.

In Fig. 2, we show the suppression factor, $\beta_{\rm rep}$. It changes
gradually from unity at low temperatures to zero at high temperatures. This
indicates that the repulsive forces become more important with increasing
temperature. At the critical temperature, $\beta_{\rm rep}\simeq 0.8$.
Consequently, $\chi_{2}^{\rm id}$ is reduced by $20\%$ due to repulsive
interactions. From Eq. (21) one may also estimate the suppression of the
higher-order cumulants. At $T\simeq T_{c}$, the suppression of $\chi_{4}$ and
$\chi_{6}$ due to repulsive interactions amounts to $41\%$ and $67\%$,
respectively.

Figure 2: The suppression factor $\beta_{\rm rep}$ from Eq. (21) for different
temperatures, expressed in the units of the chiral-critical temperature,
$T_{c}=T_{c}(\mu_{B}=0)$.

Figure 3: Ratio of different higher-order cumulants of the net-baryon number
density fluctuations, $\chi_{4}/\chi_{2}$ (left) and $\chi_{6}/\chi_{2}$
(right), calculated under different schemes, see text. The temperature is
expressed in the units of the chiral-critical temperature,
$T_{c}=T_{c}(\mu_{B}=0)$.

In the left panel of Fig. 3, we show the net-baryon kurtosis
$\chi_{4}/\chi_{2}$, and the ratio $\chi_{6}/\chi_{2}$. For the ideal HRG
model, these ratios are equal to unity due to the expected Skellam probability
distribution of the net-baryon density Braun-Munzinger  _et al._ (2004). The
chiral dynamics and repulsive interactions implemented in the parity doublet
model imply strong deviations of these fluctuation ratios from the Skellam
baseline. The kurtosis exhibits a peak around the transition temperature,
after which it drastically drops below unity. This is in contrast to the
$\sigma$HRG result, where the peak structure appears as well, however, the
result converges back to the Skellam distribution limit at higher
temperatures. Thus, the appearance of the peak in the kurtosis is attributed
to remnants of the chiral symmetry restoration, whereas strong suppression
around $T\simeq T_{c}$ is due to repulsive interactions between baryons. We
note that at vanishing chemical potential and in the chiral limit $\chi^{n}$
cumulants are non-critical up to $n<6$ order Friman  _et al._ (2011). Thus,
the peak of the kurtosis stays finite at the chiral phase transition and its
strength is strongly model-dependent. It is known that when including quantum
and thermal fluctuations in effective chiral models within the FRG approach
such a peak essentially disappears leaving the kurtosis nearly unaffected by
the chiral symmetry restoration Skokov  _et al._ (2010, 2011); Schaefer and
Wambach (2007); Friman  _et al._ (2011). The influence of repulsive
interactions, however, can imply suppression of kurtosis near the chiral
crossover below the Skellam baseline. In this model calculations, such
suppression at $T\simeq T_{c}$ is of the order of $35\%$. The kurtosis of net
baryon-number density was introduced as an excellent probe of quark
deconfinement, and outside the critical region was shown to be quantified by
the square of the baryon number carried by medium constituents Ejiri  _et al._
(2006); Karsch  _et al._ (2006). Such an interpretation of kurtosis is not
accessible in the present hadronic chiral model since it does not contain
quarks degrees of freedom. Thus, the parity doublet model can be only
applicable up to the near vicinity of the chiral crossover.

The ratio $\chi_{6}/\chi_{2}$ exhibits a strong sensitivity to dynamical
effects related to chiral symmetry restoration, as shown in the right panel of
Fig. 3. The characteristic S-type structure of this ratio obtained in the
parity doublet model with a well-pronounced peak followed by a dip at negative
values in the near vicinity of $T_{c}$ is expected as an imprint of the chiral
criticality Friman  _et al._ (2011). The leading role of the chiral symmetry
restoration on the properties of $\chi_{6}/\chi_{2}$ is also seen by comparing
the full parity doublet and $\sigma$HRG model results in Fig. 3. In both
cases, the S-type structure of this ratio is preserved, albeit with some
quantitative differences which are linked to repulsive interactions. Indeed,
as already discussed in the context of the kurtosis, the presence of repulsive
interactions suppresses $\chi_{6}/\chi_{2}$ when compared to the $\sigma$HRG
results. Nevertheless, the qualitative structure of this ratio remains the
same.

The magnitude of the peak-dip structure observed in $\chi_{6}/\chi_{2}$ is a
direct consequence of the mean-field approximation employed in our
calculations. We note that the inclusion of mesonic fluctuations weakens the
critical behavior. This was presented in other models exhibiting $O(4)$ chiral
criticality, e.g., quark-meson (QM) Morita  _et al._ (2015) and Polyakov loop-
extended quark-meson (PQM) Skokov  _et al._ (2010, 2011); Friman  _et al._
(2011); Almasi  _et al._ (2017) models within FRG approach. In these models,
the $\chi_{4}/\chi_{2}$ ratio decreases monotonously with temperature and
practically no peak structure is exhibited neat $T_{c}$. In contrast, the
general peak-dip structure of $\chi_{6}/\chi_{2}$ obtained in the mean-field
approximation prevails when mesonic fluctuations are included. This highly
non-monotonic behavior is also seen in lattice QCD simulations. In particular,
the peak in $\chi_{6}/\chi_{2}$ is seen, despite huge systematic error at low
temperatures Bazavov  _et al._ (2017); Borsanyi  _et al._ (2018).

Lastly, we compare the properties of $\chi_{4}/\chi_{2}$ and
$\chi_{6}/\chi_{2}$ fluctuation ratios with the excluded volume formulation of
the repulsive interactions. The effect of excluded volume on thermodynamic
properties of a hadronic medium was extensively studied Andronic  _et al._
(2012); Rischke  _et al._ (1991); Albright  _et al._ (2014, 2015); Redlich and
Zalewski (2016); Lo  _et al._ (2017); Venugopalan and Prakash (1992). We
consider the common formulation of the excluded volume effect, in which it is
considered for the bulk pressure of the system. In Ref. Lo  _et al._ (2017),
it was pointed out that this may not be a robust gauge of the repulsion in the
individual interaction channels.

The thermodynamic pressure in the excluded volume approach is given by a
thermodynamically self-consistent equation:

$P^{\rm ex}(T,\mu^{\star})=\sum_{k}P^{\rm id}\left(T,\mu_{k}\right)\textrm{,}$
(23)

with $\mu_{k}=\pm\mu_{B}-v_{0}P^{\rm ex}(T,\mu^{\star})$, and $P^{\rm id}$ is
the ideal-gas pressure. The ($\pm$) sign applies to particles and
antiparticles, respectively. The summation over $k$ goes over particles and
anti-particles. The constant eigenvolume,

$v_{0}=\frac{16}{3}\pi r_{0}^{3}\textrm{,}$ (24)

where $r_{0}=0.3~{}$fm. Once the excluded volume, $v_{0}$, is specified, the
pressure can be obtained by solving Eq. (23) self-consistently. From this,
higher-order cumulants are obtained numerically.

In Fig. 3, we show the ratios of cumulants, $\chi_{4}/\chi_{2}$ and
$\chi_{6}/\chi_{2}$, obtained in the excluded volume approach labeled as
exHRG. For consistency, the temperature is normalized to the critical
temperature, $T_{c}$, obtained in the parity doublet model. The excluded
volume model shows a swift decrease from the ideal HRG behavior at low
temperature and turns negative at around $1.35~{}T_{c}$. Recently, it was
suggested that this behavior may call into question the connection to
deconfinement transition in QCD Vovchenko  _et al._ (2017). For
$\chi_{6}/\chi_{2}$ the excluded volume approach deviates from the ideal HRG
result, turns negative, and predicts a dip above $T_{c}$. Very similar
behavior is also observed in models incorporating Van-der-Waals type
formulation of repulsive and attractive interactions between hadrons Vovchenko
_et al._ (2017).

Clearly, the results of the excluded volume and the parity doublet model are
qualitatively different. However, our consistent chiral model calculations
confirmed that indeed repulsive interactions between hadrons imply suppression
of $\chi_{4}/\chi_{2}$ and $\chi_{6}/\chi_{2}$ fluctuation ratios near the
chiral crossover as pointed out in Ref. Vovchenko  _et al._ (2017). Thus, the
phenomenological hadronic models that account for repulsive and attractive
interactions between constituents can be successful in describing some
deviations of net-charge fluctuations from the HRG baseline observed in LQCD.
However, this is not the case if such fluctuation observables are affected by
the chiral criticality. This is very transparent when considering
$\chi_{6}/\chi_{2}$ ratio shown in Fig. 3. The lack of in-medium effects due
to the chiral symmetry restoration in the excluded volume approach directly
implies that such a model is not capable of reproducing a characteristic
structure of this ratio near the chiral crossover. This indicates that in
order to fully describe the properties of cumulants of net-baryon number
fluctuations near the chiral crossover, it is not sufficient to account only
for repulsive interactions, but it is essential to formulate a consistent
framework that implements the chiral in-medium effects and repulsive
interactions simultaneously.

##  V Summary

We have discussed the role of attractive and repulsive nucleon-nucleon
interactions on the thermodynamic and chiral-critical properties of a strongly
interacting medium at finite temperature. To this end, we have used the parity
doublet model in the mean-field approximation. We have analyzed the
thermodynamic pressure and the cumulants of the net-baryon number density,
$\chi_{n}$, up to the sixth order, as well as their ratios.

We have shown that, at vanishing chemical potential, the second-order cumulant
factorizes as a product of a term that is directly linked to attractive scalar
interactions and a suppression factor due to repulsive interactions.
Furthermore, we have found that to a good approximation, a similar separation
also holds for higher-order cumulants. This allowed to consistently delineate
different in-medium effects and to identify the role of repulsive interactions
near the chiral crossover transition.

We have pointed out that even a moderate influence of the repulsive
interactions between hadrons on the equation of state becomes more readily
exposed in the quantitative structure of the cumulants with increasing their
order. Consequently, in the phenomenological description of cumulants of net
baryon density fluctuations calculated in LQCD or measured in heavy-ion
collisions attention is to be given to account for possible repulsive
interactions between baryons.

A frequently used approach to account for repulsive hadronic interactions is
the hard-core repulsion or van-der-Waals type interaction model. We have
compared our results for the nth-order cumulants of net-baryon number
fluctuations with an excluded volume formulation of the repulsive
interactions. Such phenomenological model provides suppression of cumulants
with increasing temperature due to hadronic repulsion. In particular, the
kurtosis $\chi_{4}/\chi_{2}$ in this model is reduced from the Skellam limit
towards the chiral crossover, as observed in LQCD results. However, when
considering the $\chi_{6}/\chi_{2}$ fluctuation ratio, which exhibits a
dominant contribution from the chiral criticality, such phenomenological model
fails to capture the characteristic properties of this ratio. Consequently, a
description of the higher-order cumulants of the net-baryon density in the
chiral crossover requires a consistent framework that accounts for a self-
consistent treatment of the chiral in-medium effects and repulsive
interactions.

At low temperature, the parity doublet model predicts sequential liquid-gas
and chiral phase transitions in the baryon-rich matter, with critical
endpoints of both transitions at moderate temperatures of tens of MeV
Zschiesche  _et al._ (2007); Sasaki and Mishustin (2010). It is challenging to
identify the role of in-medium effects and hadronic interactions on the
properties of higher-order cumulants near these distinct phase transitions.
Work in this direction is in progress.

###### Acknowledgements.

This work was partly supported by the Polish National Science Center (NCN),
under OPUS Grant No. 2018/31/B/ST2/01663 (K.R. and C.S.) and Preludium Grant
No. UMO-2017/27/N/ST2/01973 (M.M.). K.R. also acknowledges the support of the
Polish Ministry of Science and Higher Education. M.M. acknowledges helpful
comments from Pok Man Lo.

##  Appendix A Evaluation of the higher-order cumulants

Figure 4: Fourth- and sixth-order cumulants of the net-baryon number obtained
in the parity doublet model and their approximations at $\mu_{B}=0$, see text.
The temperature is expressed in the units of the chiral-critical temperature,
$T_{c}=T_{c}(\mu_{B}=0)$.

From the definition of the nth order cumulant of the net-baryon number
density,

$\chi_{n}=T^{n-4}\frac{\partial^{n-1}n_{B}}{\partial\mu_{B}^{n-1}}\textrm{,}$
(25)

one can derive the following relation:

$\begin{split}\chi_{n+1}&=\frac{\partial\chi_{n}\left(T,\mu_{B},\sigma(T,\mu_{B}),\omega(T,\mu_{B})\right)}{\partial\mu_{B}}=\frac{\partial\chi_{n}}{\partial\mu_{B}}+\frac{\partial\chi_{n}}{\partial\omega}\frac{\partial\omega}{\partial\mu_{B}}+\frac{\partial\chi_{n}}{\partial\sigma}\frac{\partial\sigma}{\partial\mu_{B}}=\\\
&\frac{\partial\chi_{n}}{\partial\mu^{\ast}}\frac{\partial\mu^{\ast}}{\partial\mu_{B}}+\frac{\partial\chi_{n}}{\partial\mu^{\ast}}\frac{\partial\mu^{\ast}}{\partial\omega}\frac{\partial\omega}{\partial\mu_{B}}+\frac{\partial\chi_{n}}{\partial\sigma}\frac{\partial\sigma}{\partial\mu_{B}}=\chi_{n+1}^{\rm
id}-g_{\omega}\chi_{n+1}^{\rm
id}\frac{\partial\omega}{\partial\mu_{B}}+\frac{\partial\chi_{n}}{\partial\sigma}\frac{\partial\sigma}{\partial\mu_{B}}=\chi_{n+1}^{\rm
id}\beta_{\rm
rep}+\frac{\partial\chi_{n}}{\partial\sigma}\frac{\partial\sigma}{\partial\mu_{B}}\textrm{,}\end{split}$
(26)

where $\chi_{n+1}^{\rm id}$ is the ideal gas expression for the (n+1)th order
cumulant, $\beta_{\rm rep}=1-g_{\omega}\partial\omega/\partial\mu_{B}$, and
$\mu^{\ast}=\mu_{B}-g_{\omega}\omega$ is the effective chemical potential.
Applying Eq. (26) to the net-baryon density, one obtains the second-order
cumulant as

$\chi_{2}=\chi_{2}^{\rm id}\beta_{\rm rep}+\frac{\partial
n_{B}}{\partial\sigma}\frac{\partial\sigma}{\partial\mu_{B}}\textrm{.}$ (27)

At $\mu_{B}=0$, $n_{B}$ vanishes, owing to particle-antiparticle symmetry.
Hence, the last term in Eq. (27) vanishes and one obtains Eq. (20).

Higher-order cumulants are obtained by applying Eq. (26) iteratively. In
particular, the third- and fourth-order cumulants, read

$\displaystyle\chi_{3}=$ $\displaystyle\chi_{3}^{\rm id}\beta^{2}_{\rm
rep}+2\frac{\partial\chi_{2}^{\rm
id}}{\partial\sigma}\frac{\partial\sigma}{\partial\mu_{B}}\beta_{\rm
rep}-g_{\omega}\chi_{2}^{\rm
id}\frac{\partial^{2}\omega}{\partial\mu_{B}^{2}}+\frac{\partial^{2}n_{B}}{\partial\sigma^{2}}\left(\frac{\partial\sigma}{\partial\mu_{B}}\right)^{2}+\frac{\partial
n_{B}}{\partial\sigma}\frac{\partial^{2}\sigma}{\partial\mu_{B}^{2}}\textrm{,}$
(28a) $\displaystyle\chi_{4}=$ $\displaystyle\chi_{4}^{\rm id}\beta_{\rm
rep}^{3}+3\frac{\partial\chi_{3}^{\rm
id}}{\partial\sigma}\frac{\partial\sigma}{\partial\mu_{B}}\beta_{\rm
rep}^{2}-3g_{\omega}\chi_{3}^{\rm
id}\frac{\partial^{2}\omega}{\partial\mu_{B}^{2}}\beta_{\rm
rep}+3\frac{\partial^{2}\chi_{2}^{\rm
id}}{\partial\sigma^{2}}\left(\frac{\partial\sigma}{\partial\mu_{B}}\right)^{2}\beta_{\rm
rep}$ $\displaystyle+3\frac{\partial\chi_{2}^{\rm
id}}{\partial\sigma}\left(\frac{\partial^{2}\sigma}{\partial\mu_{B}^{2}}\beta_{\rm
rep}-g_{\omega}\frac{\partial\sigma}{\partial\mu_{B}}\frac{\partial^{2}\omega}{\partial\mu_{B}^{2}}\right)-g_{\omega}\chi_{2}^{\rm
id}\frac{\partial^{3}\omega}{\partial\mu_{B}^{3}}+\frac{\partial^{3}n_{B}}{\partial\sigma^{3}}\left(\frac{\partial\sigma}{\partial\mu_{B}}\right)^{3}+3\frac{\partial^{2}n_{B}}{\partial\sigma^{2}}\frac{\partial\sigma}{\partial\mu_{B}}\frac{\partial^{2}\sigma}{\partial\mu_{B}^{2}}+\frac{\partial
n_{B}}{\partial\sigma}\frac{\partial^{3}\sigma}{\partial\mu_{B}^{3}}.$ (28b)

We note, that $\sigma(T,\mu_{B})=\sigma(T,-\mu_{B})$ and
$\omega(T,\mu_{B})=-\omega(T,-\mu_{B})$. Thus, the odd derivatives of $\sigma$
and even derivatives of $\omega$ w.r.t $\mu_{B}$ vanish, i.e.,

$\frac{\partial^{2k+1}\sigma}{\partial\mu_{B}^{2k+1}}\Bigg{|}_{\mu_{B}=0}=0\textrm{,}~{}~{}~{}~{}~{}\frac{\partial^{2k}\omega}{\partial\mu_{B}^{2k}}\Bigg{|}_{\mu_{B}=0}=0\textrm{.}$
(29)

Consequently, at $\mu_{B}=0$, each term in $\chi_{3}$ vanishes separately,
thus $\chi_{3}(T,\mu_{B}=0)=0$. On the other hand, the even-order cumulants,
$\chi_{2}$ and $\chi_{4}$, are reduced to

$\displaystyle\chi_{2}$ $\displaystyle=\chi_{2}^{\rm id}\beta_{\rm
rep}\textrm{,}$ (30a) $\displaystyle\chi_{4}$ $\displaystyle=\chi_{4}^{\rm
id}\beta_{\rm rep}^{3}+3\frac{\partial\chi_{2}^{\rm
id}}{\partial\sigma}\frac{\partial^{2}\sigma}{\partial\mu_{B}^{2}}\beta_{\rm
rep}-g_{\omega}\chi_{2}^{\rm
id}\frac{\partial^{3}\omega}{\partial\mu_{B}^{3}}\textrm{.}$ (30b)

The higher-order cumulants are calculated following a similar method as
explained above. The approximations used in Eqs. (20) and (21) are obtained by
dropping all terms except the first one. In Fig. 4, we present a direct
comparison of the exact results for the fourth- and sixth-order cumulants and
their approximations given in Eq. (21). The agreement is quite satisfactory
for temperatures $T\leq 1.1~{}T_{c}$, keeping the basic properties of the
cumulants.

## References

  * Bazavov  _et al._ (2014) A. Bazavov _et al._ (HotQCD), Phys. Rev. D90, 094503 (2014), arXiv:1407.6387 [hep-lat] .
  * Borsanyi  _et al._ (2018) S. Borsanyi, Z. Fodor, J. N. Guenther, S. K. Katz, K. K. Szabo, A. Pasztor, I. Portillo, and C. Ratti, JHEP 10, 205 (2018), arXiv:1805.04445 [hep-lat] .
  * Bazavov  _et al._ (2017) A. Bazavov _et al._ , Phys. Rev. D 95, 054504 (2017), arXiv:1701.04325 [hep-lat] .
  * Bazavov  _et al._ (2020) A. Bazavov _et al._ , Phys. Rev. D 101, 074502 (2020), arXiv:2001.08530 [hep-lat] .
  * Aoki  _et al._ (2006) Y. Aoki, G. Endrodi, Z. Fodor, S. Katz, and K. Szabo, Nature 443, 675 (2006), arXiv:hep-lat/0611014 .
  * Bazavov  _et al._ (2019) A. Bazavov _et al._ (HotQCD), Phys. Lett. B 795, 15 (2019), arXiv:1812.08235 [hep-lat] .
  * Aarts  _et al._ (2015) G. Aarts, C. Allton, S. Hands, B. Jäger, C. Praki, and J.-I. Skullerud, Phys. Rev. D92, 014503 (2015), arXiv:1502.03603 [hep-lat] .
  * Aarts  _et al._ (2017) G. Aarts, C. Allton, D. De Boni, S. Hands, B. Jäger, C. Praki, and J.-I. Skullerud, JHEP 06, 034 (2017), arXiv:1703.09246 [hep-lat] .
  * Aarts  _et al._ (2019) G. Aarts, C. Allton, D. De Boni, and B. Jäger, Phys. Rev. D99, 074503 (2019), arXiv:1812.07393 [hep-lat] .
  * De Tar and Kunihiro (1989) C. E. De Tar and T. Kunihiro, Phys. Rev. D39, 2805 (1989).
  * Jido  _et al._ (2000) D. Jido, T. Hatsuda, and T. Kunihiro, Phys. Rev. Lett. 84, 3252 (2000), arXiv:hep-ph/9910375 [hep-ph] .
  * Jido  _et al._ (2001) D. Jido, M. Oka, and A. Hosaka, Prog. Theor. Phys. 106, 873 (2001), arXiv:hep-ph/0110005 [hep-ph] .
  * Dexheimer  _et al._ (2008) V. Dexheimer, S. Schramm, and D. Zschiesche, Phys. Rev. C77, 025803 (2008), arXiv:0710.4192 [nucl-th] .
  * Gallas  _et al._ (2010) S. Gallas, F. Giacosa, and D. H. Rischke, Phys. Rev. D82, 014004 (2010), arXiv:0907.5084 [hep-ph] .
  * Paeng  _et al._ (2012) W.-G. Paeng, H. K. Lee, M. Rho, and C. Sasaki, Phys. Rev. D85, 054022 (2012), arXiv:1109.5431 [hep-ph] .
  * Sasaki  _et al._ (2011) C. Sasaki, H. K. Lee, W.-G. Paeng, and M. Rho, Phys. Rev. D84, 034011 (2011), arXiv:1103.0184 [hep-ph] .
  * Gallas  _et al._ (2011) S. Gallas, F. Giacosa, and G. Pagliara, Nucl. Phys. A 872, 13 (2011), arXiv:1105.5003 [hep-ph] .
  * Zschiesche  _et al._ (2007) D. Zschiesche, L. Tolos, J. Schaffner-Bielich, and R. D. Pisarski, Phys. Rev. C75, 055202 (2007), arXiv:nucl-th/0608044 [nucl-th] .
  * Benic  _et al._ (2015) S. Benic, I. Mishustin, and C. Sasaki, Phys. Rev. D91, 125034 (2015), arXiv:1502.05969 [hep-ph] .
  * Marczenko and Sasaki (2018) M. Marczenko and C. Sasaki, Phys. Rev. D97, 036011 (2018), arXiv:1711.05521 [hep-ph] .
  * Marczenko  _et al._ (2018) M. Marczenko, D. Blaschke, K. Redlich, and C. Sasaki, Phys. Rev. D98, 103021 (2018), arXiv:1805.06886 .
  * Marczenko  _et al._ (2019) M. Marczenko, D. Blaschke, K. Redlich, and C. Sasaki, Universe 5, 180 (2019), arXiv:1905.04974 [nucl-th] .
  * Marczenko (2020) M. Marczenko, Eur. Phys. J. ST 229, 3651 (2020), arXiv:2005.14535 [nucl-th] .
  * Marczenko  _et al._ (2020) M. Marczenko, D. Blaschke, K. Redlich, and C. Sasaki, Astron. Astrophys. 643, A82 (2020), arXiv:2004.09566 [astro-ph.HE] .
  * Motornenko  _et al._ (2020) A. Motornenko, J. Steinheimer, V. Vovchenko, S. Schramm, and H. Stoecker, Phys. Rev. C101, 034904 (2020), arXiv:1905.00866 [hep-ph] .
  * Mukherjee  _et al._ (2017a) A. Mukherjee, S. Schramm, J. Steinheimer, and V. Dexheimer, Astron. Astrophys. 608, A110 (2017a), arXiv:1706.09191 [nucl-th] .
  * Mukherjee  _et al._ (2017b) A. Mukherjee, J. Steinheimer, and S. Schramm, Phys. Rev. C96, 025205 (2017b), arXiv:1611.10144 [nucl-th] .
  * Dexheimer  _et al._ (2013) V. Dexheimer, J. Steinheimer, R. Negreiros, and S. Schramm, Phys. Rev. C87, 015804 (2013), arXiv:1206.3086 [astro-ph.HE] .
  * Steinheimer  _et al._ (2011a) J. Steinheimer, S. Schramm, and H. Stocker, Phys. Rev. C84, 045208 (2011a), arXiv:1108.2596 [hep-ph] .
  * Weyrich  _et al._ (2015) J. Weyrich, N. Strodthoff, and L. von Smekal, Phys. Rev. C92, 015214 (2015), arXiv:1504.02697 [nucl-th] .
  * Sasaki and Mishustin (2010) C. Sasaki and I. Mishustin, Phys. Rev. C82, 035204 (2010), arXiv:1005.4811 [hep-ph] .
  * Yamazaki and Harada (2019a) T. Yamazaki and M. Harada, Phys. Rev. D 99, 034012 (2019a), arXiv:1809.02359 [hep-ph] .
  * Yamazaki and Harada (2019b) T. Yamazaki and M. Harada, Phys. Rev. C100, 025205 (2019b), arXiv:1901.02167 [nucl-th] .
  * Ishikawa  _et al._ (2019) T. Ishikawa, K. Nakayama, and K. Suzuki, Phys. Rev. D99, 054010 (2019), arXiv:1812.10964 [hep-ph] .
  * Steinheimer  _et al._ (2011b) J. Steinheimer, S. Schramm, and H. Stocker, J. Phys. G38, 035001 (2011b), arXiv:1009.5239 [hep-ph] .
  * Giacosa (2012) F. Giacosa, Prog. Part. Nucl. Phys. 67, 332 (2012), arXiv:1111.4944 [hep-ph] .
  * Motornenko  _et al._ (2019) A. Motornenko, V. Vovchenko, J. Steinheimer, S. Schramm, and H. Stoecker,  _Proceedings, 27th International Conference on Ultrarelativistic Nucleus-Nucleus Collisions (Quark Matter 2018): Venice, Italy, May 14-19, 2018_ , Nucl. Phys. A982, 891 (2019), arXiv:1809.02000 [hep-ph] .
  * Motohiro  _et al._ (2015) Y. Motohiro, Y. Kim, and M. Harada, Phys. Rev. C92, 025201 (2015), [Erratum: Phys. Rev.C95,no.5,059903(2017)], arXiv:1505.00988 [nucl-th] .
  * Minamikawa  _et al._ (2020) T. Minamikawa, T. Kojo, and M. Harada, (2020), arXiv:2011.13684 [nucl-th] .
  * Braun-Munzinger  _et al._ (2004) P. Braun-Munzinger, K. Redlich, and J. Stachel,  _Quark-Gluon Plasma 3. Edited by Hwa Rudolph & Wang Xin-Nian. Published by World Scientific Publishing Co. Pte. Ltd_ (2004) pp. 491–599, arXiv:nucl-th/0304013 .
  * Andronic  _et al._ (2018) A. Andronic, P. Braun-Munzinger, K. Redlich, and J. Stachel, Nature 561, 321 (2018), arXiv:1710.09425 [nucl-th] .
  * Venugopalan and Prakash (1992) R. Venugopalan and M. Prakash, Nucl. Phys. A 546, 718 (1992).
  * Broniowski  _et al._ (2015) W. Broniowski, F. Giacosa, and V. Begun, Phys. Rev. C 92, 034905 (2015), arXiv:1506.01260 [nucl-th] .
  * Friman  _et al._ (2015) B. Friman, P. M. Lo, M. Marczenko, K. Redlich, and C. Sasaki, Phys. Rev. D 92, 074003 (2015), arXiv:1507.04183 [hep-ph] .
  * Huovinen  _et al._ (2017) P. Huovinen, P. M. Lo, M. Marczenko, K. Morita, K. Redlich, and C. Sasaki, Phys. Lett. B 769, 509 (2017), arXiv:1608.06817 [hep-ph] .
  * Lo  _et al._ (2018) P. M. Lo, B. Friman, K. Redlich, and C. Sasaki, Phys. Lett. B 778, 454 (2018), arXiv:1710.02711 [hep-ph] .
  * Majumder and Muller (2010) A. Majumder and B. Muller, Phys. Rev. Lett. 105, 252002 (2010), arXiv:1008.1747 [hep-ph] .
  * Andronic  _et al._ (2012) A. Andronic, P. Braun-Munzinger, J. Stachel, and M. Winn, Phys. Lett. B 718, 80 (2012), arXiv:1201.0693 [nucl-th] .
  * Albright  _et al._ (2014) M. Albright, J. Kapusta, and C. Young, Phys. Rev. C 90, 024915 (2014), arXiv:1404.7540 [nucl-th] .
  * Vovchenko  _et al._ (2015) V. Vovchenko, D. Anchishkin, and M. Gorenstein, Phys. Rev. C 91, 024905 (2015), arXiv:1412.5478 [nucl-th] .
  * Lo  _et al._ (2015) P. M. Lo, M. Marczenko, K. Redlich, and C. Sasaki, Phys. Rev. C 92, 055206 (2015), arXiv:1507.06398 [nucl-th] .
  * Man Lo  _et al._ (2016) P. Man Lo, M. Marczenko, K. Redlich, and C. Sasaki, Eur. Phys. J. A 52, 235 (2016).
  * Andronic  _et al._ (2020) A. Andronic, P. Braun-Munzinger, D. Gündüz, Y. Kirchhoff, M. Köhler, J. Stachel, and M. Winn, (2020), arXiv:2011.03826 [nucl-th] .
  * Vovchenko  _et al._ (2017) V. Vovchenko, M. I. Gorenstein, and H. Stoecker, Phys. Rev. Lett. 118, 182301 (2017), arXiv:1609.03975 [hep-ph] .
  * Stephanov  _et al._ (1999) M. A. Stephanov, K. Rajagopal, and E. V. Shuryak, Phys. Rev. D 60, 114028 (1999), arXiv:hep-ph/9903292 .
  * Asakawa  _et al._ (2000) M. Asakawa, U. W. Heinz, and B. Muller, Phys. Rev. Lett. 85, 2072 (2000), arXiv:hep-ph/0003169 .
  * Hatta and Stephanov (2003) Y. Hatta and M. Stephanov, Phys. Rev. Lett. 91, 102003 (2003), [Erratum: Phys.Rev.Lett. 91, 129901 (2003)], arXiv:hep-ph/0302002 .
  * Bazavov  _et al._ (2012) A. Bazavov _et al._ , Phys. Rev. Lett. 109, 192302 (2012), arXiv:1208.1220 [hep-lat] .
  * Borsanyi  _et al._ (2014) S. Borsanyi, Z. Fodor, S. Katz, S. Krieg, C. Ratti, and K. Szabo, Phys. Rev. Lett. 113, 052301 (2014), arXiv:1403.4576 [hep-lat] .
  * Karsch and Redlich (2011) F. Karsch and K. Redlich, Phys. Lett. B 695, 136 (2011), arXiv:1007.2581 [hep-ph] .
  * Braun-Munzinger  _et al._ (2015) P. Braun-Munzinger, A. Kalweit, K. Redlich, and J. Stachel, Phys. Lett. B 747, 292 (2015), arXiv:1412.8614 [hep-ph] .
  * Vovchenko  _et al._ (2020) V. Vovchenko, O. Savchuk, R. V. Poberezhnyuk, M. I. Gorenstein, and V. Koch, Phys. Lett. B 811, 135868 (2020), arXiv:2003.13905 [hep-ph] .
  * Braun-Munzinger  _et al._ (2020) P. Braun-Munzinger, B. Friman, K. Redlich, A. Rustamov, and J. Stachel, (2020), arXiv:2007.02463 [nucl-th] .
  * Friman  _et al._ (2011) B. Friman, F. Karsch, K. Redlich, and V. Skokov, Eur. Phys. J. C 71, 1694 (2011), arXiv:1103.3511 [hep-ph] .
  * Karsch (2019) F. Karsch, PoS CORFU2018, 163 (2019), arXiv:1905.03936 [hep-lat] .
  * Braun-Munzinger  _et al._ (2017) P. Braun-Munzinger, A. Rustamov, and J. Stachel, Nucl. Phys. A 960, 114 (2017), arXiv:1612.00702 [nucl-th] .
  * Gell-Mann and Levy (1960) M. Gell-Mann and M. Levy, Nuovo Cim. 16, 705 (1960).
  * Serot and Walecka (1986) B. D. Serot and J. D. Walecka, Adv. Nucl. Phys. 16, 1 (1986).
  * Karsch  _et al._ (2003a) F. Karsch, K. Redlich, and A. Tawfik, Eur. Phys. J. C 29, 549 (2003a), arXiv:hep-ph/0303108 .
  * Karsch  _et al._ (2003b) F. Karsch, K. Redlich, and A. Tawfik, Phys. Lett. B 571, 67 (2003b), arXiv:hep-ph/0306208 .
  * Karsch (2014) F. Karsch, Acta Phys. Polon. Supp. 7, 117 (2014), arXiv:1312.7659 [hep-lat] .
  * Albright  _et al._ (2015) M. Albright, J. Kapusta, and C. Young, Phys. Rev. C 92, 044904 (2015), arXiv:1506.03408 [nucl-th] .
  * Wetterich (1993) C. Wetterich, Phys. Lett. B 301, 90 (1993), arXiv:1710.05815 [hep-th] .
  * Morris (1994) T. R. Morris, Int. J. Mod. Phys. A 9, 2411 (1994), arXiv:hep-ph/9308265 .
  * Ellwanger (1993) U. Ellwanger,  , 206 (1993), arXiv:hep-ph/9308260 .
  * Berges  _et al._ (2002) J. Berges, N. Tetradis, and C. Wetterich, Phys. Rept. 363, 223 (2002), arXiv:hep-ph/0005122 .
  * Asakawa and Yazaki (1989) M. Asakawa and K. Yazaki, Nucl. Phys. A 504, 668 (1989).
  * Halasz  _et al._ (1998) A. M. Halasz, A. Jackson, R. Shrock, M. A. Stephanov, and J. Verbaarschot, Phys. Rev. D 58, 096007 (1998), arXiv:hep-ph/9804290 .
  * Berges and Rajagopal (1999) J. Berges and K. Rajagopal, Nucl. Phys. B 538, 215 (1999), arXiv:hep-ph/9804233 .
  * Skokov  _et al._ (2010) V. Skokov, B. Stokic, B. Friman, and K. Redlich, Phys. Rev. C82, 015206 (2010), arXiv:1004.2665 [hep-ph] .
  * Skokov  _et al._ (2011) V. Skokov, B. Friman, and K. Redlich, Phys. Rev. C 83, 054904 (2011), arXiv:1008.4570 [hep-ph] .
  * Schaefer and Wambach (2007) B.-J. Schaefer and J. Wambach, Phys. Rev. D 75, 085015 (2007), arXiv:hep-ph/0603256 .
  * Almasi  _et al._ (2017) G. A. Almasi, B. Friman, and K. Redlich, Phys. Rev. D 96, 014027 (2017), arXiv:1703.05947 [hep-ph] .
  * Ejiri  _et al._ (2006) S. Ejiri, F. Karsch, and K. Redlich, Phys. Lett. B 633, 275 (2006), arXiv:hep-ph/0509051 .
  * Karsch  _et al._ (2006) F. Karsch, S. Ejiri, and K. Redlich, Nucl. Phys. A 774, 619 (2006), arXiv:hep-ph/0510126 .
  * Morita  _et al._ (2015) K. Morita, B. Friman, and K. Redlich, Phys. Lett. B 741, 178 (2015), arXiv:1402.5982 [hep-ph] .
  * Rischke  _et al._ (1991) D. H. Rischke, M. I. Gorenstein, H. Stoecker, and W. Greiner, Z. Phys. C 51, 485 (1991).
  * Redlich and Zalewski (2016) K. Redlich and K. Zalewski, Phys. Rev. C 93, 014910 (2016), arXiv:1507.05433 [hep-ph] .
  * Lo  _et al._ (2017) P. M. Lo, B. Friman, M. Marczenko, K. Redlich, and C. Sasaki, Phys. Rev. C 96, 015207 (2017), arXiv:1703.00306 [nucl-th] .