# Heavy flavor as a probe of hot QCD matter produced in proton-proton collisions

Jiaxing Zhao, Joerg Aichelin, Pol-Bernard Gossiaux, Klaus Werner

► To cite this version:

Jiaxing Zhao, Joerg Aichelin, Pol-Bernard Gossiaux, Klaus Werner. Heavy flavor as a probe of hot QCD matter produced in proton-proton collisions. Physical Review D, 2024, 109 (5), pp.054011. 10.1103/PhysRevD.109.054011. hal-04258911

HAL Id: hal-04258911

https://hal.science/hal-04258911v1

Submitted on 16 Jan 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Heavy flavor as a probe of hot QCD matter produced
in proton-proton collisions

Jiaxing Zhao,●, Joerg Aichelin, Pol Bernard Gossiaux, and Klaus Werner
SUBATECH, Nantes University, IMT Atlantique, IN2P3/CNRS,

4 rue Alfred Kastler, 44307 Nantes cedex 3, France

<figure><img src="image_2.png" /></figure>

(Received 30 October 2023; accepted 13 February 2024; published 6 March 2024)

The creation of a quark-gluon plasma (QGP) is expected in heavy-ion collisions. It came as a surprise that proton-proton collisions at ultrarelativistic energies show as well a “QGP-like” behavior and signs of the creation of a fluid, although the corresponding system size is not more than a few cubic femtometers. Even more surprisingly, also heavy-flavor particles seem to be part of the fluid or at least interact with it. In this paper, we will investigate in a quantitative way this “collective behavior” of heavy flavor, by employing the newly developed EPOS4HQ approach, which has proven to be compatible with basic experimental data of light-flavor hadrons. We will investigate all observables, which may manifest collectivity, as particle spectra, elliptic flow, baryon-to-meson ratios, and two-particle correlations, and compare the results with experimental data. We will try to disentangle initial-state effects, those being due to interactions between charm quarks and plasma partons, and final-state effects (hadronization).

DOI: 10.1103/PhysRevD.109.054011

I. INTRODUCTION

There exists ample evidence that a new, deconfined
phase of matter, predicted by lattice QCD calculations [1,2]
and called quark-gluon plasma (QGP), is created in
relativistic heavy-ion collisions. The QGP expands rapidly
leading to a continuous decrease in temperature and
density. Finally, the energy density falls below a critical
value and the QGP hadronizes into hadrons, which are
finally observed. Heavy mesons, which contain either a
heavy quark Q or a heavy antiquark $\bar{Q}$, have turned out to
be an ideal probe to study the QGP due to two reasons:
(a) The heavy quark mass is much larger than the QCD
cutoff, $m_Q \gg \Lambda_{\text{QCD}}$. Therefore their production can be well
described by perturbative QCD (pQCD) and hence their
initial momentum distribution is known; (b) the heavy
quark mass is much larger than the typical temperature of
the hot medium, $m_Q \gg T_{\text{QGP}}$. Their masses are unchanged
in the hot medium and their number is conserved during the
evolution. The transverse momentum spectrum of heavy
mesons in heavy-ion collisions (as compared to that in
proton-proton collisions) can only be understood if a QGP
is created in these collisions. This is confirmed by another
key observable, the elliptic flow $v_2$, which heavy quarks
can only acquire through interactions with the QGP

because initially they are formed in hard processes, which
are azimuthally isotropic.

It came as a surprise that in (high multiplicity) pp
collisions observations have been made, which are con-
sidered in heavy-ion collisions as a proof of the existence of
a QGP. They include the observation of long-range corre-
lations (also called near-side “ridge”) [3–6], of strangeness
enhancement [7], and of a finite elliptic flow of D⁰ mesons
[8]. It was even more astonishing that the multiplicity of
chaired baryons at midrapidity in these collisions [9] is
considerably higher than expected from the analysis of
e⁺e⁻ collisions, questioning the process independence of
fragmentation functions [10–13]. This is also the topic of
several theoretical studies [14–17] which assume that the
hadronization mechanism of heavy quarks in high-energy
pp collisions is quite different from that in e⁺e⁻.

The purpose of this paper is to show that the recently
advanced EPOS4HQ<fnref n="1" /> approach, which allows for the
creation of a QGP in high-energy density regions, inde-
pendent of the system size, describes pp collisions as well
and can reproduce quantitatively the experimental obser-
vations. Therefore pp collisions at ultrarelativistic energies
are just the small system size limit of AA collisions.

II. EPOS4 PRIMARY INTERACTIONS

In the EPOS4 approach, we distinguish “primary inter-
actions” and “secondary interactions.” The former refer to
parallel partonic scatterings, happening at very high

<footnote>Published by the American Physical Society under the terms of
the Creative Commons Attribution 4.0 International license.
Further distribution of this work must maintain attribution to
the author(s) and the published article's title, journal citation,
and DOI. Funded by SCOAP³.</footnote>

<footnote n="1">¹Based on version EPOS4.0.1.s9.</footnote>
energies instantaneously, at $t=0$. Any notion of a sequential ordering makes no sense. The theoretical tool is S-matrix theory, using a particular form of the proton-proton scattering S matrix. The main new development in EPOS4 [18–21] is a way to accommodate simultaneously (i) rigorous parallel scattering, (ii) energy-momentum sharing, and (iii) validity of the Abramovsky-Gribov-Kancheli (AGK) theorem [22], which assures binary scaling (in AA scattering) and factorization [23] (in $pp$) for hard processes, by introducing (in a very particular way) saturation, compatible with recent “low-x-physics” considerations [24–30]. So although energy-momentum sharing makes things complicated, it is not only mandatory for a consistent picture, it also allows one to understand a connection between factorization and saturation.

Validity of AGK means that we can do the same as models based on factorization (defining and using parton distribution functions) to study hard processes (this is needed to prove consistency), but we can do much more. One of the highlights of the past decade in our domain, concerns collective phenomena in small systems. It has been shown that high-multiplicity $pp$ events show very similar collective features as earlier observed in heavy-ion collisions [3]. High multiplicity means automatically “multiple parton scattering,” and the EPOS4 formalism allows perfectly to treat this. Here one does not employ the usual parton distribution functions (representing the partonic structure of a fast nucleon); one treats the different scatterings (happening in parallel) individually, for each one a parton evolution is realized according to some evolution function $E$ (representing the spacelike cascade), as sketched in Fig. 1. One still has Dokshitzer-Gribov-Lipatov-Altarelli-Parisi evolution [31–33], for each of the scatterings, but one introduces saturation scales. But, most importantly, these scales are not constants, they depend on the number of scatterings, and they depend as well on $x^+$ and $x^-$ [18]. In Fig. 1, we show for simplicity only gluons, and we do not show the timelike cascade of further parton emissions from the emitted gluons.

<figure><img src="image_3.png" /><figcaption>FIG. 1. Rigorous parallel scattering scenario, here for three parallel scatterings. The red symbols should remind us that the parts of the diagram representing nonlinear effects (like gluon fusion) are replaced by simply using saturation scales.</figcaption></figure>

### III. HEAVY-FLAVOR PRODUCTION IN THE PRIMARY INTERACTIONS

Concerning heavy flavor, we use the general notation of $Q$ for quarks and $\bar{Q}$ for antiquarks. Heavy flavor may be produced in different ways, as shown in Fig. 2. Starting from a gluon, a $Q-\bar{Q}$ pair may be produced in the spacelike cascade, as shown in Fig. 2(a), provided the virtuality is large enough. The number of allowed flavors is considered to be depending on the virtuality (variable flavor number scheme). It is also possible to create a $Q-\bar{Q}$ in the Born process, via $g+g \to Q+\bar{Q}$ or $q+\bar{q} \to Q+\bar{Q}$ (for light-flavor quarks $q$), as shown in Fig. 2(b), and finally $Q-\bar{Q}$ may be produced in the timelike cascade, via $g \to Q+\bar{Q}$, as shown in Fig. 2(c).

Once the parton emissions are done, one considers the corresponding Feynman diagram and constructs a color flow picture, which defines chains of partons by following the color flow, as shown in Fig. 3. These chains of partons are then mapped (in a unique fashion) to kinky strings, where each parton corresponds to a kink.

For more details see Ref. [19]. The general mapping procedure (chains of partons to kinky strings) as well as the string decay procedures into “string segments” (which finally correspond to hadrons) are described in detail in [34]. Now we use the term “prehadrons” for the string segments, not knowing yet if they eventually become hadrons.

<figure><img src="image_4.png" /><figcaption>FIG. 2. Different possibilities to create heavy flavor, (a) in the spacelike cascade (SLC), (b) in the Born process, and (c) in the timelike cascade (TLC).</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>FIG. 3. The chains 1–2–3–4–5, 6–7–8, and 9–10–11 are mapped to kinky strings (red lines). The black points indicate the kinks, which carry the parton momenta.</figcaption></figure>
#### IV. EPOS4 CORE-CORONA METHOD AND FLUID EVOLUTION

From the above-mentioned primary interactions, we obtain a more or less important number of prehadrons. We employ a core-corona procedure [35–37], where the prehadrons, considered at a given proper time $\tau_0$, are separated into “core” and “corona” prehadrons, depending on the energy loss of each prehadron when traversing the “matter” composed of all the others. Corona prehadrons (per definition) can escape, whereas core prehadrons lose all their energy and constitute what we call “core,” which acts as an initial condition for a hydrodynamic evolution [37,38]. The evolution of the core ends whenever the energy density falls below some critical value $\epsilon_{\text{FO}}$, which marks the point where the fluid “decays” into hadrons. It is not a switch from fluid to particles, it is a sudden decay, called “hadronization.” Let us consider a (randomly chosen, but typical) 7 TeV proton-proton scattering event with six Pomerons (representing roughly 3 times the average). In Fig. 4, we plot the energy density in the transverse plane $(x, y)$. We consider two snapshots, namely at the start time of the hydro evolution $\tau_0 = 0.40 \text{ fm/c}$ (upper plot) and a later time $\tau_1$ close to final freeze-out (lower plot). The initial distribution has an elongated shape (just by accident, due to the random positions of interacting partons). One can clearly see that the final distributions are as well elongated, but perpendicular to the initial ones, as expected in a hydrodynamical expansion. More examples can be found in [21].

In EPOS4, as discussed in detail in [21], we developed a new procedure of energy-momentum flow through the “freeze-out (FO) hypersurface” defined by $\epsilon_{\text{FO}}$, which allows for defining an effective invariant mass. It decays according to microcanonical phase space into hadrons, which are then Lorentz boosted according to the flow velocities computed at the FO hypersurface. We also developed new and very efficient methods for the microcanonical procedure [21]. Also in the full scheme, including primary and secondary interactions, energy-momentum and flavors are conserved.

#### V. HEAVY QUARK ENERGY LOSS IN EPOS4HQ

In EPOS4 heavy quarks do not interact and heavy hadrons are produced by fragmentation functions, which describes $e^+e^-$ collisions. In EPOS4HQ, if in some spatial regions the energy density of the fluid is above the critical value and therefore a QGP is formed, heavy quarks collide with the partons of the QGP. We include in this study both elastic [39] and radiative [40] collisions. To describe these collisions we select in a first step the collision rate, averaging over the thermal distribution of the partons. In a second step we select randomly the momentum of the QGP parton, a light quark or a gluon, from their corresponding thermal distribution assuming that these particles

<figure><img src="image_6.png" /><figcaption>FIG. 4. Energy density in the transverse plane $(x, y)$ for proton-proton scattering involving six Pomerons. The upper plot represents the start time $\tau_0$ (of the hydro evolution), and the lower plot a later time $\tau_1$, close to final freeze-out.</figcaption></figure>

are massless. The thermal distribution is determined by the temperature and the mean velocity at the freeze-out surface. The scattering cross sections of the heavy quark with gluons and light quarks are calculated by pQCD matrix elements with a running coupling constant.

The pQCD elastic scattering cross section diverges for a small momentum transfer in the $t$ and $u$ channels. These infrared divergences are healed by the Debye screening mass $m_D(T)$ of gluons in the hot medium, which is calculated in the hard thermal loop approach. It serves as a regulator of the propagator of the exchanged gluon. Scattering at high momentum transfer is, on the contrary, described by a free gluon propagator. A smooth transition between both regimes can be assured by an effective Debye mass $m_{\text{eff}} = \kappa m_D(T)$ in the gluon propagator [39], with $\kappa = 0.2$. This effective mass is assumed in our calculations whose details may be found in Ref. [39].
The pQCD inelastic scattering cross section has been calculated in [40]. This cross section contains five matrix elements for gluon emission from the heavy quark and the light quark and gluon, respectively. Also for the inelastic cross section the momentum of the plasma particle is chosen by a Monte Carlo approach from the local thermal distribution. As in the elastic cross section, the gluon propagator is regulated by $m_{\text{eff}} = \kappa m_D(T)$. For the gluon emission vertex a constant $\alpha_s = 0.3$ is used. The emitted gluon is considered as massless. The different limits of the pQCD cross section calculations as well as more details of the approach have been discussed in Ref. [40].

Both the elastic as well as the inelastic collisions have been already employed in former EPOS versions to describe heavy meson data in heavy-ion collisions at LHC energies [41,42]. We use this theoretical approach to calculate the energy loss of heavy quarks, without modification, also in this new EPOS4HQ version. In this paper the K factor for elastic as well as for inelastic collisions, which has been varied in the past [43,44], is equal one, so the calculated pQCD cross sections are not modified by an overall factor.

## VI. HEAVY-FLAVOR HADRONIZATION IN EPOS4HQ

When the heavy quark crosses the freeze-out hypersurface, it will be converted into a colorless heavy-flavor hadron. This process is nonperturbative and is usually related to soft gluon radiation. Whereas in the standard EPOS4 approach (called pure EPOS4 in this study) all heavy quarks hadronize by fragmentation, in EPOS4HQ also hadronization by coalescence may contribute if a QGP is formed. In the coalescence process the heavy quark combines with one or two nearby fluid (QGP) partons to form final-state mesons or baryons. Heavy quarks, which do not traverse a QGP, create hadrons by fragmentation.

In the coalescence model, the heavy quarks coalesce with light quarks to a hadron $m$ when the charm passes the hadronization hypersurface determined by the critical energy density. The light (anti)quarks are assumed to be thermalized. The differential yield of the heavy hadron is given by

$$ \frac{dN}{d^3\mathbf{P}} = g_H \sum_{N_c} \int \prod_{i=1}^{k} \frac{d^3\mathbf{p}_i}{(2\pi)^3} f_i(\mathbf{p}_i) W_m(\mathbf{p}_1, \dots, \mathbf{p}_i) \times \delta^{(3)}(\mathbf{P} - \sum_{i=1}^{k} \mathbf{p}_i), \quad (1) $$

where $g_H$ is the degeneracy factor of color and spin. $k=2$ (3) for mesons (baryons). $\mathbf{P}$ and $\mathbf{p}_i$ are the momenta of heavy-flavor hadron and the constituent quarks, respectively. The delta function conserves the momentum. The summation is performed over all heavy quarks in the system.

$f_1(\mathbf{p}_1) = (2\pi)^3 \delta^{(3)}(\mathbf{p}_c - \mathbf{p}_1)$ is the normalized momentum space distribution of the heavy quark and $f_i(\mathbf{p}_i)$ for $i > 1$ is the momentum space distribution of the constituent quarks in the heavy hadron. $W_m(\mathbf{p}_1, \dots, \mathbf{p}_i)$ is the Wigner density of the heavy hadron $m$ in momentum space, which can be constructed from the hadron wave function. The hadron wave function can be approximated by a three-dimensional harmonic oscillator state with the same root mean square radius. For the ground state of charmed meson, the Wigner density in the center-of-mass (c.m.) frame can be expressed as

$$ W(p_r) = (2\sqrt{\pi}\sigma)^3 e^{-\sigma^2 p_r^2}, \quad (2) $$

where $p_r$ is the relative momentum between two constituent quarks in the c.m. frame. It is normalized, $\int W(p_r) d^3\mathbf{p}_r / (2\pi)^3 = 1$. We use $p_r = |E_2\mathbf{p}_1 - E_1\mathbf{p}_2| / (E_1 + E_2)$, where $E_1(\mathbf{p}_1)$ and $E_2(\mathbf{p}_2)$ are the energies and momenta of the quark or antiquark in the c.m. frame, respectively. The width $\sigma$ in the Wigner density is controlled by the root-mean-radius as shown next. Baryons are treated as two two-body systems (baryons are produced by recombining two quarks first and then using their center of mass to recombine with the third one).

The light quarks are assumed to be thermalized. The quark masses are $m_{u/d} = 0.1$ GeV, $m_s = 0.3$ GeV, and $m_c = 1.5$ GeV. The heavy quark coalescence probability in a static hot medium with a temperature $T_{\text{FO}} = 167$ MeV corresponding to a critical energy density of $\epsilon_{\text{FO}} = 0.57$ GeV/fm$^3$ is obtained by integrating Eq. (1). The root-mean-square radius of a heavy hadron is defined as $\langle r^2 \rangle = \sum_{i=1}^k \langle (\mathbf{r}_i - \mathbf{R})^2 \rangle$ with the quark coordinate $\mathbf{r}_i$ and the c.m. coordinate $\mathbf{R}$. For charmed mesons, $\langle r^2 \rangle = \frac{3}{2} \frac{m_c^2+m_q^2}{(m_c+m_q)^2} \sigma^2$. The root-mean-square radius of the ground state charmed meson can be calculated by the two-body Dirac equation [45]. It gives $\sqrt{\langle r^2 \rangle} \approx 0.85$ fm, for in-medium $D^0$. This leads to a width $\sigma = 3.725$ GeV$^{-1}$. We take the same width $\sigma$ for $D_s$ and any two-quark system in charmed baryons.

The coalescence probability of excited states, which can strongly decay into the ground states, is estimated via the statistic model. There the hadron density at the temperature $T_{\text{FO}}$ is given by [46] $n_i = \frac{g_i}{2\pi^2} T_{\text{FO}} m_i^2 K_2(\frac{m_i}{T_{\text{FO}}})$. $g_i$ is the spin isospin degeneracy. $m_i$ is the mass of the hadron. $K_2$ is the second-order Bessel function. In our study, we consider almost all possible excited states, also the missing baryons, which are predicted by the quark model [47] and lattice QCD [48,49]. For each ground state hadron $D$, $D_s$, $\Lambda_c$, $\Xi_c$, and $\Omega_c$ we calculate the density of each excited state $m$ and define the momentum-independent ratio $R^m = n_{\text{excited}}^m / n_{\text{ground}}$. Finally we sum over all excited states $R = \sum R^m$ and multiply the ground state $p_T$-dependent coalescence probability by $1 + R$ to obtain the coalescence
<figure><img src="image_7.png" /><figcaption>FIG. 5. Momentum-dependent coalescence probability of a charm quark into various charmed hadrons in a hot and static medium with temperature $T_{\text{FO}} = 167$ MeV. Here $D$ is the sum of $D^0$ and $D^+$. $\Xi_c$ is the sum of $\Xi^0$ and $\Xi^+$.</figcaption></figure>

probability for prompt charmed hadrons as shown in Fig. 5. The sum of the effective $p_T$ distributions for all hadrons gives the total coalescence probability $P_{\text{coal}}^{\text{Total}}(p_T)$ as shown with the black line in Fig. 5. Heavy quarks, which do not hadronize via coalescence, will fragment into a heavy-flavor hadron. The fragmentation probability is therefore $1 - P_{\text{coal}}^{\text{Total}}$. The fragmentation function we employ are those from the heavy quark effective theory [50,51] and the fragmentation ratios to various charmed hadrons are taken as the $e^+e^-$ collisions [52]. The evolution of heavy-flavor hadrons in the hadronic phase is controlled by the UrQMD transport model but generally hadronic rescattering is negligible, as we will show.

## VII. RESULTS AND ANALYSIS

When in $pp$ collisions a QGP is created because locally the energy density exceeds the critical energy density, there are two ways in which heavy hadrons are formed: by fragmentation and by coalescence, as described in the last section. It is the coalescence contribution that changes the observables in medium- and high-multiplicity $pp$ collisions as compared to those at low multiplicity and observed in $e^+e^-$ events. We focus in the following on four observables, the transverse momentum spectra, yield ratios, elliptic flow, and correlations.

### A. Transverse momentum spectrum

The EPOS4HQ transverse momentum spectra of $D^0$, $D^+$, $D_s$, $\Lambda_c$, $\Xi_c^0$, and $\Omega_c$ hadrons are shown in Fig. 6 as red lines and compared to the available experimental data from ALICE [10,11] and CMS [12]. For $D^0$ mesons we display as well the distribution of the charm quarks at creation (dotted black line) and before hadronization (dashed magenta line), as well as that of $D^0$ mesons immediately

<figure><img src="image_8.png" /><figcaption>FIG. 6. $p_T$ spectrum of $D^0$, $D^+$, $D_s$, $\Lambda_c$, $\Xi_c^0$, and $\Omega_c$ in 5.02 TeV $pp$ collisions. The experiment data are from the ALICE [10,11] and CMS [12].</figcaption></figure>

after hadronization, before the hadronic rescattering (dashed blue line), which is almost identical to the red line. The black-dashed lines display the results from the pure EPOS4 with fragmentation only. A significant difference can be observed between the pure EPOS4 and EPOS4HQ in the charmed baryon spectra. This reveals the importance of the coalescence in high-energy $pp$ collisions. We see that the energy loss of the charm quarks in the QGP changes the $p_T$ spectrum considerably. The difference in the $p_T$ spectrum of charm quarks at creation and the $D^0$ mesons observed in the detectors is almost exclusively due to the hadronization. Consequently, the transverse momentum spectra are not very sensitive to the presence of a QGP. We observe furthermore that our results reproduce quite well the experimental data.

### B. Yield ratio

The measured yield ratios of charmed baryons and mesons, $\Lambda_c/D^0$ and $\Xi_c^0/D^0$, in $pp$ collisions [3], are shown in Fig. 7. We compare pure EPOS4 and EPOS4HQ results with the experimental data. EPOS4HQ shows a strong enhancement of this ratio at low $p_T$ and describes the experimental data quite reasonably. Pure EPOS4, in which hadrons are exclusively produced via fragmentation, shows a $p_T$ independent ratio. This indicates that the origin of the enhancement by roughly a factor of 2.5 is the rescattering of heavy quarks with QGP partons and the subsequent
<figure><img src="image_9.png" /><figcaption>FIG. 7. Yield ratio of $\Lambda_c/D^0$ and $\Xi_c^0/D^0$ as a function of the transverse momentum $p_T$ in 5.02 TeV pp collisions. The experiment data are from the ALICE [10,11] and CMS [12].</figcaption></figure>

hadronization via coalescence. The difference of this ratio in $e^+e^-$ and $pp$ collisions has been interpreted as a sign of process-dependent fragmentation functions [53]. Our calculation does not point into this direction. The origin is rather a more complex reaction scenario in $pp$ as compared to $e^+e^-$, in particular the creation of regions with a high energy density what can be interpreted as the creation of a QGP. In Fig. 8 we study the momentum dependence of $\Lambda_c/D^0$ in four multiplicity bins, left for EPOS4HQ and right for pure EPOS4. In pure EPOS4 this ratio is, as expected, constant and independent of the multiplicity. In EPOS4HQ this ratio is considerably higher and approaches only for decreasing multiplicity the value of EPOS4. With decreasing multiplicity, the energy density of the produced partons gets lower. As a result, both the size of the QGP and the charm quark fraction in the QGP are decreasing. In Fig. 9 we display the multiplicity-dependent charm fraction. The charm quarks, which do not pass a QGP, do not interact and hadronize only via fragmentation, while those which pass the QGP will interact with the thermal partons and hadronize finally via coalescence plus fragmentation. The increase of the fraction of charm quarks, which pass a QGP, leads to an enhancement of the coalescence contribution and as a result of the $\Lambda_c/D^0$ ratio.

<figure><img src="image_10.png" /><figcaption>FIG. 8. $\Lambda_c/D^0$ ratio as a function of $p_T$ in 13 TeV pp collisions. The results are shown with the EPOS4HQ (left panel) and pure EPOS4 (right panel) in different multiplicity range. The experiment data are from the ALICE [54]. Here the multiplicity bin is selected in midrapidity $|\eta| < 1$.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>FIG. 9. Fraction of charm quarks which interact with the QGP in $pp$ collisions at 13 TeV as a function of the multiplicity of charged particles with $|\eta| < 1$.</figcaption></figure>

Figure 10 shows the $p_T$ integrated ratio $\Lambda_c/D^0$ as a function of the multiplicity in comparison with experimental data [54]. This enhancement increases smoothly from low (where the results agree with those measured in $e^+e^-$ collisions) to high multiplicity and finally reaches saturation as the fraction of the charm quarks, which pass a QGP (Fig. 9).

### C. Elliptic flow $v_2$ of $D^0$ meson

Another observable, which is sensitive to collectivity, is the elliptic flow $v_2$. It is displayed, as a function of $p_T$, in Fig. 11. Initially heavy quarks show an isotropic distribution (black dotted line). They gain $v_2$ by interactions with the QGP partons before hadronization (dashed magenta line). $v_2$ of $D$ mesons immediately after hadronization is shown as a short dashed blue line, whereas the result after

<figure><img src="image_12.png" /><figcaption>FIG. 10. Multiplicity dependent of $\Lambda_c/D^0$. The results are shown for the EPOS4HQ (solid lines) and pure EPOS4 (dashed lines) in 5.02 TeV, 13 TeV $pp$ and 30%-50% central 5.02 TeV PbPb collisions. The experimental data are from Ref. [54].</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>FIG. 11. Elliptic flow of charm quarks at creation and after passing through the QGP as well as that of D mesons immediately after production and after hadronic rescattering using the UrQMD model. $v_2$ of the K meson in comparison with the experimental data ($v_2\{2\}$) from CMS [6] is displayed as well.</figcaption></figure>

hadronic rescattering, employing the UrQMD model, is shown as a red line. It is evident that the $v_2$ of heavy mesons is created while the c quark passes the QGP. Hadronization and hadronic rescatter modify $v_2$ only marginally. We plot as well the $v_2\{4\}$ of kaons and compare the result to $v_2\{2\}$ CMS data to show that $v_2$ in the light sector is described by EPOS4. To be sure that this is not an artifact due to a too high $v_2$ value of the QGP partons, we compare in Fig. 11 as well $v_2(p_T)$ of kaons (orange line) with experiment and find that our calculation agrees with the experimental finding at least at low $p_T$. Hence the $v_2$ of the QGP partons is correctly reproduced. We note that $v_2(p_T)$ of D mesons is well below the $v_2(p_T)$ of the kaons. In pure EPOS4 calculations, where heavy hadrons are produced by fragmentation, $v_2$ of heavy hadrons is zero, as expected.

Selecting for the EPOS4HQ results the same multiplicity bins as in the experiments, we can compare our calculations with the experimental data of the ATLAS and CMS Collaborations. This comparison is shown in Fig. 12. The calculated (as compared to the kaons lower) $v_2$ values (see Fig. 11) agree quantitatively with the experimental data for different multiplicities as measured by the two collaborations.

The multiplicity dependence of $v_2$ is shown in Fig. 13. Also here one can see that $v_2$ increases with multiplicity. As expected from Fig. 9 we reach an asymptotic value at $N_{ch} \approx 50$.

## D. Correlations

Pure EPOS4 reproduces well DD and D$\bar{D}$ correlations [56]. The correlations of $c\bar{c}$ and D$\bar{D}$ in EPOS4HQ are displayed in Fig. 14 for $p_T > 3$ GeV (for both heavy quarks or mesons) and are more pronounced with this $p_T$ cut than without. Neither the interaction of heavy quarks

<figure><img src="image_14.png" /><figcaption>FIG. 12. $D^0$ meson $v_2$ as a function of transverse momentum $p_T$ in 13 TeV pp collisions. The experiment data are from ATLAS [55] and CMS [8].</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>FIG. 13. Dependence of $v_2$ of $D^0$ mesons on the charged particle multiplicity with $|\eta| < 1$. The EPOS4HQ results are shown for $\sqrt{s_{NN}} = 5.02$ and 13 TeV.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>FIG. 14. $c\bar{c}$ and D$\bar{D}$ correlations in 5.02 TeV pp collisions. A selection of $p_T > 3$ GeV/c for both trigger and associated particles is used in all cases.</figcaption></figure>
with QGP partons (blue dashed line) nor the subsequent hadronic rescattering (red line) has an influence on these correlations and therefore EPOS4 and the EPOS4HQ give almost same result, which agrees nicely with experiment [56]. So we have to conclude that the creation of a QGP does not affect the $D\bar{D}$ correlations for $p_T > 3$ GeV. This gives the chance to probe the initial charm production via the final correlation of charmed hadron.

## VIII. DISCUSSION AND CONCLUSION

We have investigated particle spectra, baryon-to-meson ratios, the elliptic flow $v_2$, correlations of heavy mesons, employing the new EPOS4HQ approach, which is compatible with the experimental light hadron observables, in $pp$ as well as in heavy-ion collisions at LHC and RHIC energies. Employing a system-size-independent critical energy density for the formation of a QGP we observe that also in a $pp$ collisions a QGP can be formed.

If a QGP is formed, the heavy quarks interact with the QGP constituents. These interactions are the origin of the finite $v_2$ values. Our result agrees with the experimental findings. The hadronization itself and the later hadronic rescattering change the elliptic flow of D mesons only little. Due to the dominance of the t-channel elastic scattering, we can obtain a finite $v_2$, even if the energy loss of the heavy quarks in the QGP is small and the change of the transverse

momentum spectrum only marginal. The difference between the $p_T$ spectrum of heavy quarks at production and the $p_T$ spectrum of the observed D mesons is almost exclusively due to the hadronization process.

In EPOS4 hadrons are formed by string fragmentation. The hadronization of heavy quarks, which pass a QGP, is described in EPOS4HQ in addition by coalescence. We observe that coalescence enhances strongly the heavy baryon yield at low $p_T$ reproducing the observed experimental enhancement as compared to $e^+e^-$ collisions. That the heavy baryon to heavy meson ratio in EPOS4HQ increases with the charged particle multiplicity for low multiplicities and saturates at higher multiplicities, in accordance with the experimental findings, presents further evidence that the formation of a QGP is at the origin of this enhancement.

We can therefore conclude that those observables, which are sensitive to the formation of a QGP, suggest that $pp$ collisions are not elementary reactions but the low system size limit of heavy-ion reactions in which the formation of a QGP is observed if the energy density exceeds a system-size-independent critical value.

## ACKNOWLEDGMENTS

This work is funded by the European Union's Horizon 2020 research and innovation program under Grant Agreement No. 824093 (STRONG-2020).

[1] C. Bernard, T. Burch, E. B. Gregory, D. Toussaint, C. E. DeTar, J. Osborn, S. Gottlieb, U. M. Heller, and R. Sugar (MILC Collaboration), Phys. Rev. D **71**, 034504 (2005).

[2] A. Bazavov et al., Phys. Rev. D **85**, 054503 (2012).

[3] V. Khachatryan et al. (CMS Collaboration), J. High Energy Phys. 09 (2010) 091.

[4] G. Aad et al. (ATLAS Collaboration), Phys. Rev. Lett. **116**, 172301 (2016).

[5] V. Khachatryan et al. (CMS Collaboration), Phys. Rev. Lett. **116**, 172302 (2016).

[6] V. Khachatryan et al. (CMS Collaboration), Phys. Lett. B **765**, 193 (2017).

[7] J. Adam et al. (ALICE Collaboration), Nat. Phys. **13**, 535 (2017).

[8] A. M. Sirunyan et al. (CMS Collaboration), Phys. Lett. B **813**, 136036 (2021).

[9] S. Acharya et al. (ALICE Collaboration), Phys. Rev. C **99**, 024906 (2019).

[10] S. Acharya et al. (ALICE Collaboration), Phys. Rev. Lett. **127**, 202301 (2021).

[11] S. Acharya et al. (ALICE Collaboration), J. High Energy Phys. 10 (2021) 159.

[12] A. Tumasyan et al. (CMS Collaboration), J. High Energy Phys. 01 (2024) 128.

[13] S. Acharya et al. (ALICE Collaboration), Phys. Rev. D **108**, 112003 (2023).

[14] M. He and R. Rapp, Phys. Lett. B **795**, 117 (2019).

[15] V. Minissale, S. Plumari, and V. Greco, Phys. Lett. B **821**, 136622 (2021).

[16] H.-h. Li, F.-l. Shao, and J. Song, Chin. Phys. C **45**, 113105 (2021).

[17] A. Beraudo, A. De Pace, D. Pablos, F. Prino, M. Monteno, and M. Nardi, Phys. Rev. D **109**, L011501 (2024).

[18] K. Werner, Phys. Rev. C **108**, 064903 (2023).

[19] K. Werner and B. Guiot, Phys. Rev. C **108**, 034904 (2023).

[20] K. Werner, arXiv:2310.09380.

[21] K. Werner, Phys. Rev. C **109**, 014910 (2024).

[22] V. A. Abramovsky, V. N. Gribov, and O. V. Kancheli, Yad. Fiz. **18**, 595 (1973).

[23] J. Collins, D. Soper, and G. Sterman, in *Perturbative Quantum Chromodynamics*, edited by A. H. Mueller (World Scientific, Singapore, 1989).

[24] L. V. Gribov, E. M. Levin, and M. G. Ryskin, Phys. Rep. **100**, 1 (1983).

[25] L. D. McLerran and R. Venugopalan, Phys. Rev. D **49**, 3352 (1994).

[26] A. Kovner, L. McLerran, and H. Weigert, Phys. Rev. D **52**, 3809 (1995).
[27] Y. V. Kovchegov, Phys. Rev. D **54**, 5463 (1996).

[28] J. Jalilian-Marian, A. Kovner, L. McLerran, and H. Weigert, Phys. Rev. D **55**, 5414 (1997).

[29] J. Jalilian-Marian, A. Kovner, A. Leonidov, and H. Weigert, Nucl. Phys. **B504**, 415 (1997).

[30] J. Jalilian-Marian, A. Kovner, A. Leonidov, and H. Weigert, Phys. Rev. D **59**, 034007 (1999).

[31] V. N. Gribov and L. N. Lipatov, Sov. J. Nucl. Phys. **15**, 438 (1972).

[32] G. Altarelli and G. Parisi, Nucl. Phys. **B126**, 298 (1977).

[33] Y. L. Dokshitzer, Sov. Phys. JETP **46**, 641 (1977).

[34] H. J. Drescher, M. Hladik, S. Ostapchenko, T. Pierog, and K. Werner, Phys. Rep. **350**, 93 (2001).

[35] K. Werner, Phys. Rev. Lett. **98**, 152301 (2007).

[36] K. Werner, I. Karpenko, T. Pierog, M. Bleicher, and K. Mikhailov, Phys. Rev. C **82**, 044904 (2010).

[37] K. Werner, B. Guiot, I. Karpenko, and T. Pierog, Phys. Rev. C **89**, 064903 (2014).

[38] I. Karpenko, P. Huovinen, and M. Bleicher, Comput. Phys. Commun. **185**, 3016 (2014).

[39] P. B. Gossiaux and J. Aichelin, Phys. Rev. C **78**, 014904 (2008).

[40] J. Aichelin, P. B. Gossiaux, and T. Gousset, Phys. Rev. D **89**, 074018 (2014).

[41] P. B. Gossiaux, R. Bierkandt, and J. Aichelin, Phys. Rev. C **79**, 044906 (2009).

[42] P. B. Gossiaux, J. Aichelin, T. Gousset, and V. Guiho, J. Phys. G **37**, 094019 (2010).

[43] M. Nahrgang, J. Aichelin, P. B. Gossiaux, and K. Werner, Phys. Rev. C **90**, 024907 (2014).

[44] M. Nahrgang, J. Aichelin, S. Bass, P. B. Gossiaux, and K. Werner, Phys. Rev. C **91**, 014904 (2015).

[45] J. Zhao, S. Shi, N. Xu, and P. Zhuang, arXiv:1805.10858.

[46] A. Andronic, P. Braun-Munzinger, K. Redlich, and J. Stachel, Phys. Lett. B **659**, 149 (2008).

[47] D. Ebert, R. N. Faustov, and V. O. Galkin, Phys. Rev. D **84**, 014025 (2011).

[48] A. Bazavov et al., Phys. Lett. B **737**, 210 (2014).

[49] M. Padmanath, R. G. Edwards, N. Mathur, and M. J. Peardon, Proc. Sci. LATTICE2014 (2015) 084 [arXiv: 1410.8791].

[50] E. Braaten, K.-m. Cheung, S. Fleming, and T. C. Yuan, Phys. Rev. D **51**, 4819 (1995).

[51] M. Cacciari, P. Nason, and R. Vogt, Phys. Rev. Lett. **95**, 122001 (2005).

[52] M. Lisovyi, A. Verbytskyi, and O. Zenaiev, Eur. Phys. J. C **76**, 397 (2016).

[53] J. R. Christiansen and P. Z. Skands, J. High Energy Phys. 08 (2015) 003.

[54] S. Acharya et al. (ALICE Collaboration), Phys. Lett. B **829**, 137065 (2022).

[55] G. Aad et al. (ATLAS Collaboration), Phys. Rev. Lett. **124**, 082301 (2020).

[56] K. Werner and B. Guiot, Phys. Rev. C **108**, 034904 (2023).