coalescence, because the fragmentation functions employed are inclusive and
give only the total contribution to $D$ meson including the $D^{*}$ feed-down.
The significantly large $\Lambda_{\mathrm{c}}^{+}/D^{0}$ is due to the
combination of both a lower $D^{0}$ and a larger $\Lambda_{\mathrm{c}}^{+}$
yield with respect to experimental data. However, it has to be noticed that
the value in the plot is the original prediction employing the width
corresponding to the mean square radii of hadrons from the quark model, hence
without any free parameter; a moderate modification of the width parameter
would allow to get quite close to the data.

Figure 15: $p_{\mathrm{T}}$ integrated production cross sections of the
various charm meson and baryon species per unit of rapidity at midrapidity
normalized to that of the $\mathrm{D^{0}}$ meson measured by ALICE in pp
collisions at $\sqrt{s}=5.02~{}\mathrm{TeV}$ [191]. The measurements are
compared with different tunes of PYTHIA 8 calculations. The grey marker shows
the values of $\mathrm{D^{*+}}$/$\mathrm{D^{0}}$ ratio obtained with Monash
2013 PYTHIA8 tune after setting the parameter StringFlav:mesonCvector to 1.5,
which, as shown in the right panel, allows to better describe the
$p_{\mathrm{T}}$-differential $\mathrm{D^{*+}}$/$\mathrm{D^{0}}$ ratio
measured by ALICE in pp collisions at $\sqrt{s}=13$ $\mathrm{TeV}$ [159] (see
text for details).

Figure 16: Left: same as in left panel of Fig. 15 but with additional
comparison to the SHM+RQM [38] and Catania [36] models. Right: charm-quark
fragmentation fractions into charm hadrons measured in pp collisions at
$\sqrt{s}=5.02~{}\mathrm{TeV}$ [191] and at $\sqrt{s}=13~{}\mathrm{TeV}$ [159]
in comparison with experimental measurements performed in
$\mathrm{e^{+}e^{-}}$ collisions at LEP and at B factories, and in ep
collisions at HERA [70].

The charm FF, f($\mathrm{c}\rightarrow H_{c}$), measured by ALICE in pp
collisions at $\sqrt{s}$ = 5.02 TeV [191] and at $\sqrt{s}=13~{}\mathrm{TeV}$
[159], are displayed in Fig. 16 together with values derived from experimental
measurements performed in $\mathrm{e^{+}e^{-}}$ collisions at LEP and B
factories and in ep collisions [70]. An increase by about a factor of 3.3 of
the charm FF to the $\Lambda_{\mathrm{c}}^{+}$ baryon from
$\mathrm{e^{+}e^{-}}$ and ep to pp collisions, with a concomitant decrease by
about a factor 1.4–1.2 for the FF to D mesons is observed. The significance of
the difference is about $5\,\sigma$ for $\Lambda_{\mathrm{c}}^{+}$ baryons and
$6\,\sigma$ for $\mathrm{D^{0}}$ mesons. In previous measurements in
$\mathrm{e^{+}e^{-}}$ and ep collisions, no value for the charm FF to
$\Xi_{\mathrm{c}}^{+,0}$ was obtained. In the analysis performed in Ref. [70]
to estimate charm FF, the $\Xi_{\mathrm{c}}^{0}$ and $\Xi_{\mathrm{c}}^{+}$
yields were assumed to be equal and estimated according to the assumption
$f({\mathrm{c}}\rightarrow\Xi_{\mathrm{c}}^{+}{})/f({\mathrm{c}}\rightarrow\Lambda_{\mathrm{c}}^{+}{})$
=
$f({\mathrm{s}}\rightarrow\Xi^{-})/f({\mathrm{s}}\rightarrow\Lambda^{0})\sim$
0.066 [70]. The fraction $f({\mathrm{c}}\rightarrow\Xi_{\mathrm{c}}^{0}{})$
was measured for the first time in pp collisions by ALICE, which for the
$\Xi_{\mathrm{c}}^{0}/\Lambda_{\mathrm{c}}^{+}$ FF ratio obtained the value
$f({\mathrm{c}}\rightarrow\Xi_{\mathrm{c}}^{0}{})/f({\mathrm{c}}\rightarrow\Lambda_{\mathrm{c}}^{+}{})$
= $0.39\pm 0.07(\mathrm{stat})^{+0.08}_{-0.07}(\mathrm{syst})$ [180]. The
probability of a c quark to hadronize in a $\Xi_{\mathrm{c}}^{0}$ baryon is
measured to be comparable to the probability to hadronize in a
$\mathrm{D_{s}^{+}}$ meson, suggesting that the sum of the
${\mathrm{c}}\rightarrow\Xi_{\mathrm{c}}^{0}$ and
${\mathrm{c}}\rightarrow\Xi_{\mathrm{c}}^{+}$ fragmentation fractions could be
larger than the ${\mathrm{c}}\rightarrow\mathrm{D_{s}^{+}}$ one. Recently, the
ALICE Collaboration measured the charm-quark fragmentation fractions in pp
collisions at $\sqrt{s}$ = 13 TeV [159]. In addition to the charm-quark
fragmentation fractions into $\mathrm{D^{0}}$, $\mathrm{D^{+}}$,
$\mathrm{D_{s}^{+}}$, $\mathrm{D^{*+}}$, $\Lambda_{\mathrm{c}}^{+}$, and
$\Xi_{\mathrm{c}}^{0}$ hadrons, also the results for the fragmentation into
J$/\Psi$ mesons and $\Xi_{\mathrm{c}}^{+}$ and $\Sigma_{\mathrm{c}}^{0,+,++}$
baryons were reported. The results are compared with those in pp collisions at
$\sqrt{s}$ = 5.02 TeV and no significant energy dependence is observed within
the uncertainties. Overall, charm quarks hadronize into baryons approximately
35-40% of the time. These results confirm that the baryon enhancement at the
LHC with respect to $\mathrm{e^{+}e^{-}}$ collisions is caused by different
hadronization mechanisms at play in the parton-rich environment produced in pp
collisions. The fragmentation fractions measured at midrapidity in pp
collisions at the LHC provide strong evidence that the assumption of
universality (colliding-system independence) of parton-to-hadron fragmentation
functions is not valid for charm production.

The LHCb and CMS Collaborations, and earlier the CDF Collaboration, made
similar observations in the beauty sector when computing the
$\Lambda_{\mathrm{b}}^{0}$ to $\mathrm{B}^{0,-}$ yield ratios [170, 169, 192,
193, 173, 194, 195, 196, 13] in pp ($\mathrm{p\overline{p}}$ for CDF)
collisions at various centre-of-mass energies. As shown in the left panel of
Fig. 17, the ratio is observed to strongly depend on $p_{\mathrm{T}}$ and
approaches LEP value at high $p_{\mathrm{T}}$, which can be estimated to be
$\Lambda_{\mathrm{b}}^{0}/\mathrm{B^{0}}=0.22\pm 0.03$ [172]. At low
$p_{\mathrm{T}}$ it is significantly higher than in $\mathrm{e^{+}e^{-}}$
collisions, reaching approximately twice the value measured at LEP for both
the CDF and LHCb experiments. It does not show variations with pseudorapidity,
within the LHCb acceptance of pseudorapidity 2 $<\eta<$ 5\. The ALICE
Collaboration measured the production cross sections of $\mathrm{D^{0}}$ and
$\Lambda_{\mathrm{c}}^{+}$ hadrons originating from beauty-hadron decays at
midrapidity in pp collisions at a centre-of-mass energy $\sqrt{s}$ = 13 TeV
[197]. A trend consistent with the one of
$\Lambda_{\mathrm{b}}^{0}/\mathrm{B^{0,-}}$ measured by LHCb at forward
rapidity is observed. In the left panel of Fig. 17 the
$\mathrm{\Lambda_{\mathrm{b}}^{0}}/B^{0}$ measured by the LHCb Collaboration
[195] is compared with PYTHIA 8 calculations using different tunes. The cyan
line represents the Monash tune computed for the rapidity window 2–4.5, while
the blue solid line shows one of the QCD-based CR-BLC tunes in the same
rapidity range. This comparison confirms that the experimental measurement is
larger than what is expected by models based on fragmentation tuned on
$\mathrm{e^{+}e^{-}}$ collisions. On the other hand, given the inclusion of
junction string topologies and hence junction baryons in the CR-BLC PYTHIA 8
models, though the trend of producing an enhancement of
$\Lambda_{\mathrm{b}}^{0}$ at low $p_{\mathrm{T}}$ and reducing to LEP values
at high $p_{\mathrm{T}}$ is predicted, the magnitude of the enhancement at low
$p_{\mathrm{T}}$ is too large and overpredicts the experimental data. Such a
discrepancy remains an active question in junction string modelling. An
additional point that is important to highlight is the comparison between the
expectations of PYTHIA 8 with CR-BLC in different rapidity ranges shown in the
same figure. Simulations of PYTHIA 8 with CR-BLC for $|y|<0.5$ (dashed blue
line) are compared to the forward rapidity one (2 $<|y|<$ 4.5, solid blue
line). The results indicate that the model does not expect any dependence on
rapidity. Additionally, the expectations of a SHM approach with and without
the set of b-hadron states predicted by the RQM to exist beyond the currently
measured spectrum listed by the particle data group are also shown [51]. The
substantial gap between the LHCb data and the model calculation obtained with
the PDG scenario (dashed dark-red line) is largely overcome by the feeddown of
the “missing” baryons included in the RQM calculation (solid red line),
leading to a better description of the experimental data at intermediate
$p_{\mathrm{T}}$, while overpredicting the measurements at low and high
$p_{\mathrm{T}}$. It would also be interesting to see model calculations
including coalescence, like the Catania model, to verify whether different
baryon enhancements, both in the charm and beauty sector, might be expected in
different rapidity windows. Future measurements of other beauty baryons
($\mathrm{\Xi_{b}}$, $\mathrm{\Omega_{b}}$) are necessary to investigate a
possible enhancement in pp collisions for strange-beauty baryons.

Figure 17: Left panel: ratios
$\mathrm{\Lambda_{\mathrm{b}}^{0}}/\mathrm{B^{0}}$ as a function of
$p_{\mathrm{T}}$ measured at forward rapidity in pp collisions at $\sqrt{s}$ =
13 TeV by the LHCb Collaboration [195]. The data are compared to model
calculations [51, 33, 11, 102] (see text for details). Right panel: comparison
of the $\Lambda_{\mathrm{b}}^{0}$ production asymmetry measured by the LHCb
Collaboration at $\sqrt{s}$ = 7 TeV with PYTHIA 8 predictions [33, 11].

The LHCb Collaboration also published differential measurements of the
asymmetry between $\Lambda_{\mathrm{b}}^{0}$ and $\overline{\Lambda}_{b}^{0}$
baryon production rates in pp collisions at centre-of-mass energies of
$\sqrt{s}$ = 7 TeV and $\sqrt{s}$ = 8 TeV. The production asymmetry is defined
according to the quantity $A_{\mathrm{p}rod}$, which represent the relative
difference between the production rates of $\Lambda_{\mathrm{b}}^{0}$ and
$\overline{\Lambda}_{b}^{0}$. The production asymmetry is measured both in
intervals of rapidity in the range 2.15 $<y<$ 4.10 and transverse momentum in
2 $<p_{\mathrm{T}}<$ 27 GeV$/c$. The results are incompatible with symmetric
production with a significance of 5.8$\sigma$ assuming no CP violation in the
decay. The asymmetry is increasing with rapidity and there is evidence for a
trend with a significance of 4$\sigma$. The comparison between data and
various PYTHIA 8 models is shown in the right panel of Fig. 17. Among the
several tunes of PYTHIA, only that with CR-BLC approximation (red boxes) shows
a good agreement with the measurement. The Monash tune (yellow boxes) and the
gluon move CR (blue boxes) predict asymmetries that are too large, exhibiting
the strongest deviation at low transverse momentum and for the widest rapidity
range investigated [198]. In PYTHIA, a key source of low $p_{\mathrm{T}}$
$\Lambda_{\mathrm{b}}^{0}$ production in pp collisions at large rapidity is
the combination of a b quark with the beam remnants. Obviously, the same does
not apply to the $\overline{\Lambda}_{b}^{0}$ as there are not antiquarks in
proton beam remnants. The CR-BLC adds large amounts of low-$p_{\mathrm{T}}$
junctions, which form $\overline{\Lambda}_{b}^{0}$ and
$\Lambda_{\mathrm{b}}^{0}$ baryons in equal amounts, diluting the asymmetry.

In this context, it is also worth noting that the LHCb Collaboration has
measured the $\mathrm{D_{s}^{+}}$ production asymmetry in pp collisions. The
$\mathrm{D_{s}^{+}}$ meson does not contain any of the proton’s valence
quarks, which means that the production asymmetry is in this case not
influenced by the valence quarks of the colliding protons. PYTHIA simulations
have shown a strong dependence on both $p_{\mathrm{T}}$ and $y$, whereas the
experimental measurements do not [199].

###  III.4 Measurements as a function of multiplicity

The observation of substantially different charm baryon-to-meson ratios in pp
collisions compared to $\mathrm{e^{+}e^{-}}$ collisions evidences differences
in the hadronization process in the two systems. Though at
high-$p_{\mathrm{T}}$ the ratios seem to approach system-universal values,
suggesting the presence of an energy scale delimiting the modification of the
hadronization process, the absence of clear differences between the
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratios measured at LEP and B-factory
energies and the similarity of the trends observed in pp collisions in the
charm and beauty sectors makes it extremely difficult to reconcile pp and
$\mathrm{e^{+}e^{-}}$ data ascribing the difference to a kinematic effect
related to the energy of the initial-parton or of the jet, without implicating
different properties of the systems. One is instinctively tempted to relate
the observed differences to the typically larger number of quarks and gluons
produced in hadronic collisions. The parton density is strictly connected to
the “hadronic activity”, which is encoded in the measurable hadron
multiplicity. Inevitably, the study of the evolution of the heavy-flavour
baryon-to-meson ratios with particle multiplicity was expected to provide
important information to better understand the modification of the
hadronization process and to model it. As it will be discussed later in the
section, it is important to stress that it is not just the number of hadrons
produced that differentiates pp and $\mathrm{e^{+}e^{-}}$ collisions (at LEP
energies the multiplicity can reach values not far from the average ones of pp
collisions at the LHC) but rather the way in which hadrons can be formed. As
discussed in Sect. II.2.1 the concurrence of several initial scatterings, i.e.
of MPI, allow for more and topologically different colour connections, which
can impact the hadronization process. Also the hadronization via coalescence,
which requires multiple quark wavefunctions to overlap in position and
momentum space, is expected to progressively emerge as hadron-formation
process with increasing number of quarks produced in the collision. Moreover,
given the enhancement of multi-strange hadrons observed in high-multiplicity
pp collisions [200], one naturally wonders whether the production of heavy-
flavour hadrons containing a strange quark [200] also increases as particle
multiplicity increases and whether, if present, this increase should related
to heavy-flavour baryon production.

The ALICE and LHCb experiments measured the evolution of the
$\mathrm{D_{s}^{+}}/\mathrm{D^{0}}$, $\mathrm{B^{0}_{s}}/\mathrm{B^{0}}$,
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$, and
$\mathrm{\Lambda_{\mathrm{b}}^{0}}/B^{0}$ yield ratios from very low to high
charged-particle multiplicity in pp collisions at $\sqrt{s}=13$ TeV to provide
additional experimental constraints on the hadronization mechanisms in a high
density partonic environment [201, 202].

Figure 18: The $\mathrm{D_{s}^{+}}/\mathrm{D^{0}}$ (top) and
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ (bottom) ratios measured at
midrapidity by ALICE in pp collisions at $\sqrt{s}=13$ TeV in the lowest
(left) and highest (right) midrapidity multiplicity classes inspected [201].
The measurements are compared to PYTHIA 8 predictions with the Monash [11] and
the CR-BLC [33] tunes, and the CE-SH model [203], calculated in similar
multiplicity classes. Figure from Ref. [201].

In the upper panels of Fig. 18, the $\mathrm{D_{s}^{+}}/\mathrm{D^{0}}$ ratio
in the lowest (left) and highest (right) multiplicity classes measured by the
ALICE experiment [201] are shown. The measurements do not indicate a
significant multiplicity dependence of the ratio, which within the
uncertainties is also independent of $p_{\mathrm{T}}$ in the measured
$p_{\mathrm{T}}$ interval. This is consistent with a coalescence picture
because they are both affected in a very similar way by the density of the
medium, while the baryon over meson ratio should show an increase with the
density of the medium. On the other hand, the lower panels of Fig. 18 report
the $p_{\mathrm{T}}$-differential $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$
ratio, which shows evidence of a dependence on multiplicity such that higher
multiplicities result in an enhancement of the
$\Lambda_{\mathrm{c}}^{+}$/$\mathrm{D^{0}}$ ratio. When comparing the highest
multiplicity interval to the lowest one, the effect has a significance of
5.3$\sigma$ for 1 $<p_{\mathrm{T}}<$ 12 GeV$/c$. It is worth noticing that the
measured $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio in the lowest
multiplicity class is still higher (in the measured $p_{\mathrm{T}}$ range)
than the average of corresponding ratios measured in $\mathrm{e^{+}e^{-}}$
collisions at LEP. The ALICE Collaboration performed the measurements
selecting the multiplicity using two different multiplicity estimators located
at different pseudorapidity regions; one overlapping with the pseudorapidity
window in which $\Lambda_{\mathrm{c}}^{+}$ and $\mathrm{D^{0}}$ are
reconstructed, and one with a significant $\eta$ gap. The increasing trend of
the $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio is consistent between the
two measurements. This indicates that the enhancement from low to high
multiplicity is not induced only by the possibly different correlation of
$\Lambda_{\mathrm{c}}^{+}$ and $\mathrm{D^{0}}$ production with the “local”
multiplicity of particles produced close to them in pseudorapidity. Such a
correlation could arise either from a genuine correlation between
hadronization process and the multiplicity, e.g. in case
$\Lambda_{\mathrm{c}}^{+}$ production were favoured compared to
$\mathrm{D^{0}}$ production when the “local” parton density is higher, or from
a correlation dictated by different average local multiplicities accompanying
$\Lambda_{\mathrm{c}}^{+}$ and $\mathrm{D^{0}}$ production independently from
the hadronization process, e.g. if $\Lambda_{\mathrm{c}}^{+}$ and
$\mathrm{D^{0}}$ are produced in jets with different number of constituents.
Correlations induced by a bias on the multiplicity estimate related to signal
reconstruction, like the contribution of a different number of decay particles
to the multiplicity estimate, are instead removed by ALICE in the analysis.

The measured charm-hadron ratios in Fig. 18 for the lowest and highest
multiplicity class are compared to model predictions from MC generators and a
statistical hadronization model. The simulations with the PYTHIA 8 event
generator were performed with the Monash and the CR-BLC tunes. The measured
$\mathrm{D_{s}^{+}}/\mathrm{D^{0}}$ ratios at low and high multiplicity are
compatible with PYTHIA 8 expectations both with Monash and CR-BLC tunes. As
expected by the $\Lambda_{\mathrm{c}}^{+}$/$\mathrm{D^{0}}$ ratio results
shown in Fig. 13, the Monash tune does not reproduce the
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio in either multiplicity
interval, nor does it show a strong multiplicity dependence. By contrast, the
CR-BLC tunes describe the $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ decreasing
trend with $p_{\mathrm{T}}$, and are closer to the overall magnitude in both
multiplicity ranges. All CR-BLC tunes show a clear dependence with
multiplicity, qualitatively reproducing the trend observed in data, but the
evolution predicted by the three “modes” is different, indicating that the
data can be useful for tuning further the model.

The predictions of a canonical-ensemble statistical hadronization (CE-SH)
model [203] describe the $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratios,
reproducing their multiplicity dependence. In the calculation, the underlying
charm-baryon spectrum in the calculations is augmented to include additional
excited baryon states predicted by RQM. The CE-SH model explains the
multiplicity dependence as deriving from the reduced volume size of the
formalism towards smaller multiplicity, where a decrease of the
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio is a consequence of the strict
baryon-number conservation. Such behaviour is also expected for charm-strange
mesons relative to non-strange charm mesons because of strangeness-number
conservation. However, the resulting prediction reproduces the
$\mathrm{D_{s}^{+}}/\mathrm{D^{0}}$ ratio only at low multiplicity, while it
overestimates it in the highest multiplicity interval. For the
$\mathrm{D_{s}^{+}}/\mathrm{D^{0}}$ ratio only the PDG case is shown since the
additional RQM states do not modify significantly the D-meson relative yields
with respect to the PDG set.

Figure 19: The left and middle panel display the ratio of cross sections
$\sigma_{\mathrm{B^{0}_{s}}}/\sigma_{\mathrm{B^{0}}}$ versus the normalized
multiplicity of all VELO tracks and backward VELO tracks, respectively. The
right panel shows the comparison with PYTHIA 8 calculations for the
$0<p_{\mathrm{T}}<6$ $\mathrm{GeV}/c$ interval. The horizontal bands show the
values measured in $\mathrm{e^{+}e^{-}}$ collisions. Figure from Ref. [202].

The LHCb Collaboration measured the production rate of $\mathrm{B^{0}_{s}}$
mesons relative to $\mathrm{B^{0}}$ mesons in pp collisions at a centre-of-
mass energy $\sqrt{s}$ = 13 TeV over the forward rapidity interval 2 $<y<$ 4.5
as a function of the event charged particle multiplicity [202]. The
multiplicity estimators used in the analysis are the total number of charged
tracks reconstructed in the VELO detector ($\mathrm{N_{VELO}}$ tracks), and
the subset of VELO tracks that point in the backward direction, away from the
LHCb spectrometer ($\mathrm{N_{back}}$ tracks). The results reported in Fig.
19 show that the rate of $\mathrm{B^{0}_{s}}$ production relative to
$\mathrm{B^{0}}$ might be enhanced with respect to the one measured in
$\mathrm{e^{+}e^{-}}$. The enhancement is observed when charged-particle
multiplicity is measured in the full VELO detector, while it is not observed
when the charged-particle multiplicity is based on $\mathrm{N_{back}}$ tracks.
In collisions with relatively low charged-particle multiplicity, and for B
mesons with $p_{\mathrm{T}}$ $>$ 6 GeV$/c$, the rate of $\mathrm{B^{0}_{s}}$
production relative to $\mathrm{B^{0}}$ production is consistent with that
measured in $\mathrm{e^{+}e^{-}}$ collisions. The ratio of cross sections
calculated by the PYTHIA event generator with and without colour reconnection
is compared with the experimental measurements for the $p_{\mathrm{T}}$
interval 0–6 GeV$/c$ in the right panel of Fig. 19. Both sets of PYTHIA
calculations show a rise with multiplicity, which is more pronounced when
colour reconnection is included. Note that the source of this rise in
$\mathrm{B^{0}_{s}}/\mathrm{B^{0}}$ with respect to multiplicity in PYTHIA is
not due to a explicit strangeness enhancement model (e.g. Rope hadronization
[47] or close-packing [48]) but instead is already an existing consequence of
the standard Lund fragmentation modelling, however the exact source of this
multiplicity dependence has not yet been studied. Both PYTHIA models reproduce
the data within uncertainties though both scenarios give values that are lower
than the central values of the data at high multiplicity, possibly indicating
a different slope of the rise of the $\mathrm{B^{0}_{s}}/\mathrm{B^{0}}$ ratio
with multiplicity. The enhancement of the rate of $\mathrm{B^{0}_{s}}$
production relative to $\mathrm{B^{0}}$ production might be a consequence of a
bias arising from the coinciding rapidity regions between the multiplicity
estimator and the measurement of interest or it might signal that the effect
seen is related to the “local” multiplicity of particles surrounding the B
mesons. It would be interesting to further investigate these possibilities
with Run 3 data and to perform similar measurements at midrapidity. On the
other hand, the results could show evidence that the enhancement is present in
collisions with high charged-particle multiplicity, while for relatively low
charged-particle multiplicity, and for B mesons with $p_{\mathrm{T}}>$ 6
GeV$/c$, the rate of $\mathrm{B^{0}_{s}}$ production relative to
$\mathrm{B^{0}}$ production is consistent with what is measured in
$\mathrm{e^{+}e^{-}}$ collisions.

Figure 20: Ratio of cross sections
$\sigma_{\Lambda_{\mathrm{b}}^{0}}/\sigma_{\mathrm{B}^{0}}$ measured by LHCb
in $2<|y|<4.5$ in intervals of the multiplicity of all VELO tracks, which is
reported on the $x$ axis after scaling it by the average multiplicity. The
purple point indicates the value measured in $\mathrm{e^{+}e^{-}}$ collisions
at LEP. Figure from [195].

The LHCb Collaboration recently released the measurements of the
$\sigma_{\Lambda_{\mathrm{b}}^{0}}/\sigma_{\mathrm{B}^{0}}$ as a function of
$p_{\mathrm{T}}$ and multiplicity [195]. The $p_{\mathrm{T}}$-integrated
$\sigma_{\Lambda_{\mathrm{b}}^{0}}/\sigma_{\mathrm{B}^{0}}$ ratio is shown as
a function of relative multiplicity in Fig. 20. In the lowest multiplicity
interval, the LHCb data approaches the baryon fraction measured in
$\mathrm{e^{+}e^{-}}$ collisions at the LEP. The distinct rise of the baryon
fraction with multiplicity reaches a plateau at $\sim$ 0.5 for collisions that
produce more than twice the average number of VELO tracks.
These last two measurements show that b-quark hadronization is not universal
between $\mathrm{e^{+}e^{-}}$ and pp collisions. This may provide indication
that quark coalescence plays an important role in forming beauty hadrons at
the LHC, or that one should more carefully examine string models with
multiplicity dependencies. The trend observed for the
$p_{\mathrm{T}}$-integrated
$\sigma_{\Lambda_{\mathrm{b}}^{0}}/\sigma_{\mathrm{B}^{0}}$ ratio resembles a
“turn-on” curve which may mark the transition from two different regimes, with
the onset of a different hadronization mechanism in pp collisions becoming
very quickly effective as soon as a minimal hadronic activity is present. It
would be interesting to try to correlate the observed trend with the number of
MPI. It must be noted that the evolution observed for the
$p_{\mathrm{T}}$-integrated
$\sigma_{\Lambda_{\mathrm{b}}^{0}}/\sigma_{\mathrm{B}^{0}}$ ratio does not
contrast with what previously discussed about the evolution with multiplicity
of the $p_{\mathrm{T}}$-differential $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$
ratio measured by ALICE at midrapidity. This will be further discussed in
Sect. III.6.

###  III.5 Proton–nucleus collisions

Figure 21: Left: $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio as a function
of $p_{\mathrm{T}}$ in pp collisions (at low and high multiplicity and
multiplicity integrated) and in p–Pb collisions [174, 79, 201]. Right:
comparison of $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio in pp and p–Pb
collisions with predictions from QCM model [174]. Figure 22: Nuclear
modification factor of $\Lambda_{\mathrm{c}}^{+}$ baryons as a function of
$p_{\mathrm{T}}$ measured by ALICE [174] and compared with predictions from
POWLANG [204], QCM [93], and POWHEG + PYTHIA 6 + EPPS16 nPDF [205, 206, 10,
207, 208]. See text for details.

The influence of a hadronic environment on the hadronization process is
evidenced by the observed modifications of heavy-flavour baryon-to-meson
ratios from $\mathrm{e^{+}e^{-}}$ to pp collisions in conjunction with the
observed dependence of these ratios on charged particle multiplicity in pp
collisions as discussed in the previous section. From this alone we see that
the study of increasing the collision system size is important, and pA
collision systems provide a tantalising bridge between the small pp and large
AA collision environments. Collisions of protons with nuclei are interesting
for the study of heavy-flavour hadronization in particular for at least two
reasons. The first is that they give access to high multiplicities without the
need of selecting very rare events, unlike pp collisions whereby high-
multiplicity events represent only a small fraction of the pp inelastic cross
section. The second is the possible interplay of hadronization and collective
effects, the presence of which in small collision systems emerged in the last
decade, particularly from the study of particle azimuthal anisotropies and the
observation of relatively large flow coefficients in high-multiplicity p–Pb
collisions [209]. These effects resemble those observed in high-energy
nucleus–nucleus collisions and ascribed to the formation of a QGP, which
expands and cools down following hydrodynamic laws.

The first aspect can be better quantified by the following considerations.
High-multiplicity pp events represent a small fraction of the pp inelastic
cross section. For instance, the highest multiplicity class considered in
ALICE measurement of $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio reported
in Fig. 18 represents a fraction certainly smaller than 0.5% [210, 201] of all
inelastic events333This can be deduced, as a conservative upper limit, from
the numbers reported in Table 1 and Table 4 of Ref. [201] and [210],
respectively.. Thus, it cannot be excluded that the multiplicity trends
observed in these events derive from ’peculiar’ features of the selected
events rather than from a general dependence of the measured particle ratios
from multiplicity. A quite natural next step is to measure these ratios in
collision systems in which similar or higher multiplicities are reached
without selecting rare events. This opportunity is offered by collisions of
nuclei or of protons with nuclei, in which significantly higher multiplicities
are reached with respect to pp collisions. The multiplicity of charged
particles at midrapidity increases with the number of nucleons participating
to the collisions ($N_{\mathrm{part}}$) [211, 212, 213]. In p–Pb collisions at
$\sqrt{s_{\mathrm{NN}}}=5.02$ $\mathrm{TeV}$, the average charged-particle
multiplicity at midrapidity in non-single diffractive events is more than a
factor of 3 higher than the multiplicity measured in non-single diffractive pp
events [214]. The average multiplicity of the highest-multiplicity interval in
which ALICE measured the $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio in pp
collisions is about a factor of two higher and is reached in 0–5% more central
p–Pb collisions. As shown in Fig. 21 (left panel), the
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio measured by ALICE in p–Pb
collisions is higher than that measured in pp collisions at
$\sqrt{s}=13~{}$$\mathrm{TeV}$ without multiplicity selections in the range
$3<p_{\mathrm{T}}<12$ $\mathrm{GeV}/c$, and similar to that measured in high-
multiplicity pp events. At lower $p_{\mathrm{T}}$, there is a hint for a lower
ratio in p–Pb collisions than in pp collisions. Preliminary results from CMS
[215] suggest that the ratio does not evolve significantly with multiplicity
in p–Pb collisions, contrary to what observed for the
$\Lambda/\mathrm{K_{S}^{0}}$ ratio in the light-flavour sector [216]. Overall,
the data suggest that the modification with multiplicity of the
$p_{\mathrm{T}}$-differential $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio
observed in pp collisions is not dictated by peculiar properties of the
selected high-multiplicity events. If CMS preliminary data is confirmed, the
$p_{\mathrm{T}}$-differential $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio
from low-multiplicity pp collisions up to high-multiplicity p–Pb collisions,
would display a stronger evolution at low multiplicities followed by
saturation. Constraining better the low-multiplicity region is a major goal
for the future. New measurements with reduced uncertainty and extended to even
lower multiplicities are necessary to profile the onset of the modification of
the baryon-to-meson ratio from $\mathrm{e^{+}e^{-}}$ collisions to hadronic
collisions. In the right panel of Fig. 21, the
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio measured in pp and p–Pb
collisions are compared with QCM predictions, which are compatible with the
data. In QCM, the ratio has a different $p_{\mathrm{T}}$ trend in p–Pb
collisions and peaks at higher $p_{\mathrm{T}}$ than in pp collisions because
of the larger assumed light-quark velocity, which is needed to reproduce the
hardening of the light-flavour hadron spectra. The evolution of the
$p_{\mathrm{T}}$-differential $\Lambda_{\mathrm{c}}^{+}$ spectrum from pp to
p–Pb collisions is quantified in Fig. 22 in terms of the nuclear modification
factor, $R_{\mathrm{p}Pb}$, which is calculated by dividing the p–Pb
$p_{\mathrm{T}}$-differential production cross section by the pp one scaled by
the Pb mass number. A deviation of $R_{\mathrm{p}Pb}$ from unity is observed
and corresponds to a modification of to the spectra in p–Pb with respect to pp
collisions. The modification of the $\Lambda_{\mathrm{c}}^{+}$-baryon
production spectrum in p–Pb collisions has been additionally confirmed by the
computation of the mean transverse momentum, $\langle p_{\mathrm{T}}\rangle$,
which results larger in p–Pb with respect to pp [174]. For D mesons
$R_{\mathrm{p}Pb}\approx 1$ [217]. Therefore, the $\Lambda_{\mathrm{c}}^{+}$
$R_{\mathrm{p}Pb}$ is almost equivalent to the double ratio
$(\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}})^{\mathrm{p--
Pb}}/(\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}})^{\mathrm{pp}}$. The
POWHEG+PYTHIA6+EPPS16 nPDF model is based on the POWHEG NLO calculation [205,
206] of charm quark production cross section matched with PYTHIA 6 [10] for
parton shower and fragmentation. It uses CT14NLO [207] parton distribution
functions (PDF) and the EPPS16 parameterization [208] of nuclear PDF for the
the Pb nucleus. The latter term causes a departure of the expected
$R_{\mathrm{p}Pb}$ from unity at low $p_{\mathrm{T}}$ due to the suppression
(“nuclear shadowing”) at small Bjorken-x of the gluon PDF. The model, which
can describe D-meson data within uncertainties, cannot describe
$\Lambda_{\mathrm{c}}^{+}$ $R_{\mathrm{p}Pb}$, evidencing that the observed
modification cannot be determined by the modification of PDF of nuclei with
respect to those of protons. For heavy-flavour hadrons at LHC energies,
nuclear shadowing was considered the largest “cold-nuclear matter” effect,
i.e. an effect not related to the formation of a quark-gluon plasma, that
could modify the production cross sections in nuclear collisions to the simple
geometrical scaling of cross sections in pp collisions. In the POWLANG
calculation [204], charm quarks are transported through a “small-size”
expanding QGP, which affects the $p_{\mathrm{T}}$ distributions of charm
hadrons. However, the calculated $R_{\mathrm{p}Pb}$ value is identical for all
charm-hadron species as it does not consider any modifications of the relative
hadron abundances due to quark coalescence. The presence of a QGP phase has a
substantial impact on the predicted $R_{\mathrm{p}Pb}$ but the model cannot
describe either the D-mesons nor the $\Lambda_{\mathrm{c}}^{+}$
$R_{\mathrm{p}Pb}$ in the full $p_{\mathrm{T}}$ range. The QCM model instead
describes the data within uncertainties.

Most recently, the ALICE Collaboration released measurements of the
fragmentation fractions of charm quarks to different charm-hadron species,
which are measured for the first time in p–Pb collisions at
$\sqrt{s_{\mathrm{NN}}}$ = 5.02 TeV at midrapidity [218]. The fragmentation
fractions for all charm-hadron species are found to be consistent for pp and
p–Pb collisions at the same collision energy, indicating that there is no
significant modification of the charm-quark hadronisation process due to the
different hadronic collisions’ system sizes.

The LHCb collaboration measured the $\Lambda_{\mathrm{c}}^{+}$ and
$\Xi_{\mathrm{c}}^{+}$ production and their ratios to $\mathrm{D^{0}}$ at
forward and backward rapidity in p–Pb collisions [219, 220]. When compared
with the measurements at midrapidity [221, 174], the baryon-to-meson ratios
are lower for $p_{\mathrm{T}}$ $<$ 8 GeV/$c$, while at higher $p_{\mathrm{T}}$
the measurements are consistent within uncertainties. The origin of this
tension is not yet fully understood, for a more detailed discussion on this
open point see Sect. IV.

###  III.6 Nucleus–nucleus collisions

The study of heavy-quark hadronization in nucleus–nucleus collisions offers a
unique window for addressing hadron formation in an extended system of quarks,
which can be assumed to be in local thermodynamic equilibrium and to move
collectively along with the system expansion. Indeed, since long it was
suggested that in the QGP hadrons may form via recombination of existing
quarks [66], with important implications for the hadron abundances and
momentum distributions, expected to significantly differ from those resulting
from fragmentation and hadronization in a “vacuum-like” system, as that of an
$\mathrm{e^{+}e^{-}}$ collision can be classified. Understanding the
hadronization process in these collisions is also a necessary step to exploit
heavy-quarks for the study of in-medium partonic energy loss and of the medium
transport properties [222]. This was the initial goal motivating heavy-flavour
hadron production measurements in heavy-ion collisions: even at LHC energies,
the medium temperature is significantly lower than the charm-quark mass, so
that charm and beauty quark thermal production is negligible compared to that
in hard-scattering processes occurring in the early stage of the collision
with cross sections known with relatively small uncertainty [191, 158, 157,
223, 224]. Interacting with the medium, constituents heavy quarks exchange
energy but preserve their flavour identity. Thus they can be regarded as
“calibrated” probes of the medium and used as markers to trace the quark
diffusion and the partonic interactions in the QGP, resolving the medium
constituents and allowing to build a connection between the microscopic
“local” partonic interactions and the medium global properties, in particular,
its transport ones. However, the uncertainties in the modelling of the
hadronization process, which significantly influences the momentum
distribution of the final-state hadrons, can reduce the capability of
inspecting the dynamics in the partonic phase. On the other hand, a QGP offers
a (expanding) thermal source of light quarks and gluons whose modelling, e.g.
in coalescence models, might be simpler and less sensitive to event-by-event
fluctuations than that of a pp collision system. Therefore, the measurement of
the production of different heavy-flavour hadron species in heavy-ion
collisions can simultaneously constrain the modelling of hadronization and
parton dynamics in the medium.

Several measurements in Au–Au collisions at RHIC and Pb–Pb collisions at the
LHC demonstrate that charm and beauty quarks interact strongly with the QGP
constituents [225, 226, 227, 228, 229, 162, 230, 231, 232, 233, 234, 235, 236,
237, 238, 239, 240, 241, 242, 243, 244, 245, 246, 247, 248, 160, 249, 250,
251, 252, 253]. The $R_{\mathrm{AA}}$ values observed for both charm and
beauty hadrons are significantly smaller than unity at intermediate
$p_{\mathrm{T}}$ and larger for beauty than charm, indicating a significant
energy loss of both quarks in the medium and evidencing the expected
dependence of energy loss on the quark mass. At moderate/low $p_{\mathrm{T}}$
the $v_{2}$ of prompt D mesons is large, close to that observed for charged
pions [238]. This can indicate that charm quarks are strongly coupled to the
system and are pushed towards kinetic equilibrium via rescattering with the
medium constituents. It can also signal coalescence as a hadronization process
relevant to charm quarks through which charm hadrons partly inherit the flow
of light quarks, a relevant aspect that could be more deeply investigated
through heavy-light correlation of the anisotropic flow $v_{n}$ [254, 85]. At
low $p_{\mathrm{T}}$, where elastic processes are the most important ones,
heavy quark transport in the medium can be described as a diffusion process
[255, 256, 56, 257]. The main transport parameter is the spatial diffusion
coefficient $D_{s}$, which is related to the kinetic equilibration time of a
heavy-quark with mass $m_{\mathrm{Q}}$,
$\tau_{\mathrm{Q}}=D_{s}m_{\mathrm{Q}}/T$ [258, 259]. The comparison of
$R_{\mathrm{AA}}$ and $v_{\mathrm{2}}$ data with model predictions favour
small values of heavy-quark $D_{s}$, in the range $1.5<2\pi
D_{s}T_{\mathrm{c}}<4.5$ at the pseudocritical temperature
$T_{\mathrm{c}}=155$ $\mathrm{MeV}$, in line with other estimates [260, 256,
55, 243, 261]. This interval corresponds to relaxation times for charm quarks
smaller than about 10 fm$/c$, thus of the order of or shorter than the QGP
lifetime, actually suggesting that charm quarks might get close to reaching
kinetic equilibrium in the medium [56, 257]. A recent estimate that combines
lattice QCD to the heavy quark effective theory in the limit $M\gg T$ [257]
implies a very small thermalization time, at least in the $p\rightarrow 0$
limit; namely of about 1-1.5 fm/c for charm quarks and 3-5 fm/c for bottom
which are significantly smaller than those extrapolated from the
phenomenological analysis. This work has been developed in non-quenched QCD
which is a very important advancement w.r.t. previous studies, however, it is
still in the $M\rightarrow\infty$ limit. In Ref. [86, 262] it has been shown
that within a quasi-particle model (QPM-Catania) one expects that $D_{s}$ has
a significant mass dependence from the charm quark mass to the limit
$M\rightarrow\infty$, decreasing by about a factor 1.7. This could to a large
extent compensate the discrepancy between the value of $2\pi TD_{s}$ from the
lQCD driven estimate w.r.t. the one extrapolated from the phenomenology for
charm quarks. However, an even more recent extension of the study in Ref.
[257] has done a first extrapolation to finite mass finding a very small
dependence of $2\pi TD_{s}$ on $M_{Q}$. This, if confirmed by further studies
within lQCD and heavy quark effective field theory, opens the question on how
to reconcile them with the phenomenological estimates. At this stage, we only
mention that a preliminary study have indicated that an initial glasma phase
may lead to infer smaller values of $D_{s}$ [263], and more recently it has
been suggested that also memory effect (non-Markovian dynamics) may affect the
phenomenological extrapolation of $D_{s}$ [264], especially for charm quarks.
Therefore, with more precise data becoming available with the ALICE upgrade
(especially at low $p_{\mathrm{T}}$), it is likely that will be possible to
achieve a significant advancement in the understanding of the heavy quark
dynamics.

However, even within the uncertainties of the mass dependence of the diffusion
or drag coefficient, it remains that because of the larger mass, beauty quarks
are expected to diffuse less than charm quarks and to have a significantly
larger relaxation time that is at least a factor of two [86, 264] if not a
factor of three larger (if $D_{S}$ is nearly mass independent) as suggested by
the ratio $M_{b}/M_{c}$. Energy loss is also expected to be smaller for beauty
than charm quarks. In line with these expectations, experimental measurements
of beauty hadrons or their decay products show higher $R_{\mathrm{AA}}$ and
smaller $v_{2}$ with respect to charm or light-flavour hadrons [265, 225].

Figure 23: Prompt D-meson $R_{\mathrm{AA}}$ in the 0–10% centrality class
(left panel) and $v_{2}$ in the 30–50% centrality class (right panel) compared
with the PHSD [91], POWLANG [266, 267], and Catania [37] predictions obtained
with and without including hadronization via recombination. Figure adapted
from Ref. [233].

The hadronization schemes used in transport models have a significant impact
on the predicted values of $R_{\mathrm{AA}}$ and $v_{2}$ and on their
dependence on $p_{\mathrm{T}}$. The effect of hadronization via recombination
on D-meson $R_{\mathrm{AA}}$ and $v_{2}$ is illustrated in Fig. 23 where ALICE
measurements in Pb–Pb collisions in the centrality intervals 0–10%
($R_{\mathrm{AA}}$, left panel) and 30–50% ($v_{\mathrm{2}}$, right panel) are
compared with the expectations of the PHSD, POWLANG, and Catania models with
and without the inclusion of recombination [233]. In a medium that is radially
expanding with a common velocity, hadronization via quark recombination tends
to equilibrate particle velocities towards the medium expansion velocity. The
combined effect of recombination, diffusion, and energy loss on the
$R_{\mathrm{AA}}$ is the appearance of a maximum at momenta not much larger
than the hadron mass. For what concerns the elliptic flow, the
$v_{\mathrm{2}}$ of heavy-flavour hadrons formed via recombination derives
from a combination of the charm and light quarks flows. From the comparison of
model expectations with and without hadronization via recombination, shown in
the right panel of Fig. 23, it can be appreciated the effect of recombination
on charm-meson $v_{\mathrm{2}}$ in increasing $v_{\mathrm{2}}$ and pushing its
maximum at higher momenta. Recombination is an essential ingredient for
reproducing the data [268, 68, 222, 256, 102].

All models that can reproduce the large observed $v_{2}$ of D mesons [238]
predict that it partly derives from hadronization via coalescence. In order to
single out the contribution of hadronization from the contribution of charm-
quark diffusion in the generation of charm flow, the modification in Pb–Pb
collisions with respect to pp collisions of the $p_{\mathrm{T}}$-differential
yield ratios of different heavy-flavour particle species must be studied. In
particular, models based on statistical hadronization or coalescence were
initially designed to describe hadron formation in AA collisions. They were
expecting an increase of the production of heavy-flavour baryons relative to
heavy-flavour mesons and of hadrons with strange quarks relative to non-
strange ones in heavy-ion with respect to pp collisions, the latter originally
thought to deliver systems with properties closer to $\mathrm{e^{+}e^{-}}$
collisions and in which “in-vacuum” fragmentation should have been the
dominant hadronization mechanism. Such an expectation was motivated by the
standard lore that only the large reservoir of partons and high energy
densities of nucleus–nucleus collision systems could make it likely for quark
coalescence to become a relevant process compared to the mechanism of
fragmentation, which requires the creation of quarks from the vacuum. Hence,
the typical signature of coalescence, that is, an enhancement of baryon-to-
meson ratios, was not expected to occur in pp collisions. However, as
discussed in the previous sections, the data indicates that, particularly in
the heavy-flavor sector, pp collisions exhibit significant differences
compared to smaller collision systems such as $\mathrm{e^{+}e^{-}}$ and
$\mathrm{e^{\pm}p}$. These differences result in modifications to the particle
ratios, providing a sensitive test of changes in the hadronization dynamics.
The coalescence plus fragmentation approach, such as the Catania model,
suggests that even if pp collisions result in the formation of a QGP droplet,
D meson production is likely still dominated by the fragmentation process.
Meanwhile, the coalescence process becomes dominant for baryon production,
significantly increasing it in pp collisions compared to $\mathrm{e^{+}e^{-}}$
collisions. This outcome, though initially unexpected, naturally follows from
the application of coalescence dynamics and effectively explains the large
increase in the measured $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio.
However, the enhancement in pp collisions can also be explained in the context
of PYTHIA, but, also in this case, the main modification affects the baryon
production that has to come from colour reconnections of quarks and gluons in
topologies irrelevant in $\mathrm{e^{+}e^{-}}$ collisions but emerging with an
increasing number of MPI, hence particularly in high-multiplicity events, in
hadronic collisions. This also shows that already in pp (in particular for
high-multiplicity events) more local recombination of quarks (minimizing
string energy) becomes quite relevant for baryon production. The same
coalescence modelling shows that in AA collisions the large-size reservoir of
high-density quarks makes a further difference making the recombination the
dominant mechanism for both heavy quark mesons and baryons.

The measurements of the production of prompt and non-prompt
$\mathrm{D_{s}^{+}}$ mesons [236, 269] and of $\mathrm{B_{s}^{0}}$ mesons
[270] provide first hints supporting the presence of an enhancement of the
production of heavy-strange mesons with respect to non-strange ones. Figure 24
shows in the left panel the comparison of the prompt
$\mathrm{D_{s}^{+}}/\mathrm{D^{0}}$ ratios measured by ALICE in central Pb–Pb
collisions and in pp collisions. A hint of a higher ratio in Pb–Pb collisions
is visible in $3<p_{\mathrm{T}}<10$ $\mathrm{GeV}/c$. The Catania model
reproduces the data within uncertainties, while the TAMU and POWLANG transport
models as well as the GSI/Hd+BW statistical hadronization model, for which the
hadron $p_{\mathrm{T}}$ spectrum is modelled according to a Blast Wave
function, tend to overestimate the data in the low $p_{\mathrm{T}}$ region.
Slightly higher $\mathrm{D_{s}^{+}}/\mathrm{D^{0}}$ values, closer to 0.4 and
without a significance dependence on $p_{\mathrm{T}}$ in $1<p_{\mathrm{T}}<8$
$\mathrm{GeV}/c$, were reported by STAR in semicentral and peripheral Au–Au
collisions at $\sqrt{s_{\mathrm{NN}}}=200$ $\mathrm{GeV}$ [271]. In central
Au–Au collisions, the ratio decreases below 0.3 for $p_{\mathrm{T}}<4$
$\mathrm{GeV}/c$ and, within uncertainties, becomes compatible with PYTHIA
expectation and pp data at LHC energies. New measurements with reduced
uncertainty and lower $p_{\mathrm{T}}$ reach, as well as the measurement of
$\Xi_{\mathrm{c}}^{+,0}$ production (larger than $\mathrm{D_{s}^{+}}$ one in
pp collisions), are necessary for a final assessment and understanding on the
enhancement of char-strange hadron production. In the right panel of Fig. 24,
the $\mathrm{B_{s}^{0}}/\mathrm{B^{+}}$ ratio measured as a function of
$p_{\mathrm{T}}$ by CMS in 0–90% central collisions [270] and normalized to
the $f_{\mathrm{s}}/f_{\mathrm{d}}$ ratio measured at $\sqrt{s}$ =7
$\mathrm{TeV}$ by LHCb in $2<y<4.5$ [272] is shown together with the ratio of
non-prompt $\mathrm{D_{s}^{+}}$ $R_{\mathrm{AA}}$ to $\mathrm{D^{0}}$
$R_{\mathrm{AA}}$ measured by ALICE in the 0–10% centrality interval [269,
225]. In pp collisions about 50% of non-prompt $\mathrm{D_{s}^{+}}$ derive
from the decay of a $\mathrm{B_{s}^{0}}$ meson. Together the data provide a
hint for an enhancement of strange B-meson production relative to non-strange
B-meson production, in line with the expectation of the TAMU model [51].

Figure 24: Left (from Ref. [14]): prompt $\mathrm{D_{s}^{+}}$/$\mathrm{D^{0}}$
ratio measured by ALICE in $|y|<0.5$ in pp and 0–10% central Pb–Pb collisions
at $\sqrt{s_{\mathrm{NN}}}$ =5.02 $\mathrm{TeV}$ [236, 163, 158] compared to
model expectations [37, 36, 39, 97, 49]. Right: red (black and grey) markers
show the ratio of strange to non-strange non-prompt D mesons (B mesons) yield
in 0–10% (0–90%) Pb–Pb collisions at $\sqrt{s_{\mathrm{NN}}}$ =5.02
$\mathrm{TeV}$ measured by ALICE [269, 225] (CMS [270]) and normalised to the
analogous ratio measured in pp collisions. In the B-meson case, the
parametrization of the $f_{\mathrm{s}}/f_{\mathrm{d}}$ ratio measured at
$\sqrt{s}$ =7 $\mathrm{TeV}$ by LHCb in $2<y<4.5$ [272] is used as pp
reference, considering the collision-energy and rapidity dependence of the
ratio, as well as the parametrization uncertainty, as negligible compared to
the Pb–Pb data uncertainty. The data are compared to TAMU calculations [273,
39, 274].

The $\Lambda_{\mathrm{c}}^{+}$ is the only heavy-flavour baryon whose
production has been measured in nucleus–nucleus collisions. The left panel of
Fig. 25 shows the evolution of the $p_{\mathrm{T}}$-integrated
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio as a function of the charged-
particle multiplicity at midrapidity measured by ALICE in pp, p–Pb, and Pb–Pb
collisions and by STAR in Au–Au collisions [73]. The ratio is approximately
0.5, thus much higher than the $\approx 0.11$ value measured in
$\mathrm{e^{+}e^{-}}$ and ep collisions, and, within uncertainties, it does
not show a dependence on multiplicity. The same ratio in pp (minimum bias) and
in Pb–Pb collisions is expected by the Catania and TAMU models, which
reproduce also quantitatively the data. Also the SHMc model does not expect a
dependence on centrality in Pb–Pb collisions. However, it predicts a ratio
around 0.2 and underestimates the data. It must be noted that this version of
the SHMc does not include the augmented set of high-mass baryon states
considered by TAMU. In the SHMc charm is not in chemical equilibrium and the
abundance of charm hadron species depends on the charm fugacity factor, which
however cancels in single-charm particle ratios, and by the chemical freezeout
temperature, which results in being independent of centrality. Similarly, in
the TAMU, though the hadronization process is modelled differently in pp and
Pb–Pb collisions, particle ratios mainly depend on statistical factors. For
the Catania model, the similarity between pp and Pb–Pb values derives from the
similar description of the bulk of light-flavour quarks as a QGP medium with
similar temperatures in the two systems and from the applied normalization
condition that coalescence probability equals unity at
$p_{\mathrm{T}}^{\mathrm{hadron}}=0$. There is still a dependence on the
volume of the hadronizing fireball and a subdominat dependence on the radial
flow but the last has a weak effect on the $p_{T}$-integrated yield. In pp
collisions, the PYTHIA 8 Monash tune is even lower than the
$\mathrm{e^{+}e^{-}}$ expectation and largely independent of event
multiplicity. Contrastingly, PYTHIA 8 with CR-BLC expects a significance
increase of the ratio with multiplicity, deriving from the larger probability
of forming junction colour connections in events with a higher number of MPI,
which is not supported by the data. The region of very low multiplicities is
not covered by measurements, thus leaving open the possibility that
$\mathrm{e^{+}e^{-}}$ values are recovered when the hadronic activity in the
event is very small, in line with what the
$\Lambda_{\mathrm{b}}^{0}/\mathrm{B}^{0}$ ratio at forward rapidity, discussed
in Sect. III.4, suggests. It is worth mentioning, that the average
multiplicity in $\mathrm{e^{+}e^{-}}$ collisions at LEP energies, where most
particles are part of high-energy dijets, does not differ much from the
average multiplicity in minimum bias pp collisions at LHC energies (see Fig.
19.6 of Ref. [50]). This suggests that the number of produced particles is not
the only relevant parameter and rather the whole event topology and the impact
of the underlying event and MPI in hadronic collisions, even at the lowest
multiplicity, must be taken into account for explaining the differences
observed with respect to $\mathrm{e^{+}e^{-}}$ collisions.

Figure 25: Left (from [73]): baryon-to-meson
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ cross-section ratio
$p_{\mathrm{T}}$-integrated for $p_{\mathrm{T}}>0$ as a function of charged-
particle multiplicity at midrapidity measured by ALICE in pp, p–Pb, and Pb–Pb
collisions [201, 73] and by STAR [71] in Au–Au collisions. Right (adapted from
[176]): comparison of the $p_{\mathrm{T}}$-differential
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ cross-section ratios measured in pp
collisions by ALICE [78], in Pb–Pb collisions by ALICE [73], CMS [74], and
LHCb [275] in different centrality intervals, and in Au–Au collisions by STAR
[71].

In the right panel of Fig. 25, the evolution of the
$p_{\mathrm{T}}$-differential $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio
with centrality is shown, by comparing the measurements performed in pp
collisions [78] and in different centralities in Pb–Pb collisions at the LHC
[73, 275, 175] and in Au–Au collisions at RHIC [71]. Contrary to the
$p_{\mathrm{T}}$-integrated case, the $p_{\mathrm{T}}$-differential
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio evolves substantially with
centrality, as shown in the right panel of Fig. 25, in which the Pb–Pb ALICE
[73], CMS [74], and LHCb [275] results are reported and compared to ALICE pp
data as well as to the measurement in Au–Au collisions by STAR [71]. At
intermediate $p_{\mathrm{T}}$, the ratio increases with centrality. In
peripheral collisions, LHCb data at forward rapidity overlap with ALICE
minimum-bias pp data at midrapidity. In 30–50% Pb–Pb collisions the ratio is
higher, though the values are compatible with those measured in high-
multiplicity pp and p–Pb data. In central collisions a “radial-flow”-like peak
emerges, possibly shifted at lower $p_{\mathrm{T}}$ at RHIC with respect to
LHC energies. The appearance of this peak can be due to hadronization via
coalescence, possibly reinforced by space-momentum correlations (see Sect.
II.2.4). In Fig. 26 the $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio
measured by ALICE in 0–10% (left panel) and 30–50% (middle panel) central
Pb–Pb collisions is compared to model expectations. Despite the different
values predicted by the models, the large experimental and theoretical
uncertainties do not allow to discriminate between the $\mathrm{S}HM_{c}$
[49], Catania [37], and TAMU [38] models. At $p_{\mathrm{T}}>20$
$\mathrm{GeV}/c$ the ratio is the same at all centralities and consistent with
that measured in $\mathrm{e^{+}e^{-}}$ collisions, possibly signalling the
transition at high-enough momentum to a regime in which “vacuum”-like
fragmentation is an effective description leading to a “universal” value.

Figure 26: Baryon-to-meson $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ cross-
section ratio as a function of $p_{\mathrm{T}}$ in 0-10% (left) and 30-50%
(middle) Pb–Pb collisions, and in pp collisions (right), compared with
expectations from the Catania, SHMc, and TAMU models and from the PYTHIA 8
event generator with the Monash and CR Mode 2 tunes. Figure from Ref. [276,
73].

##  IV Prospects and future developments

The measurements of heavy-flavour baryon production carried out in the last
decade in hadronic collisions have somewhat unexpectedly broken the paradigm
that fragmentation modelled on $\mathrm{e^{+}e^{-}}$ data could be sufficient
in describing hadronization in small collision systems, and that
nucleus–nucleus would be the collision system where one would search for new
phenomena. Rather, a quite sharp distinction between hadronization in
$\mathrm{e^{+}e^{-}}$ and any hadronic collision system emerged, accompanied
by the need of continuity from pp to nucleus–nucleus collisions. Multiplicity,
as an estimator of the event activity (hence MPI?) and system size, provides a
reasonable guide for examining the evolution of the observables from small to
large collision systems. While reasonable explanations and ideas about the
origin of the transition from $\mathrm{e^{+}e^{-}}$ to pp collisions have been
put forward, together with effective modelling capable to reproduce the trends
shown by the measurements, a new theoretical paradigm for hadronization is
still being shaped. Certainly, as extensively discussed in the previous
sections, the measurements performed so far achieved some milestones
delineating features and trends that theoretical modelling of hadronization
will have to consider. They can be recapped in the following list.

  * •

The first is certainly the evidence that heavy-flavour baryon production is
significantly larger in pp than in $\mathrm{e^{+}e^{-}}$ collisions. A direct
consequence is the evidence of the breaking of fragmentation universality
already in pp collisions.

  * •

The second concerns the boundaries out of which fragmentation seems recovered
as a good description of hadronization. The observed evolution of
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ and
$\Lambda_{\mathrm{b}}^{0}/\mathrm{B^{0,-}}$ with $p_{\mathrm{T}}$ and
multiplicity suggests that this occurs at high momentum (above about 15-20
$\mathrm{GeV}/c$) and at extremely low multiplicity.

  * •

The third is that the enhancement from $\mathrm{e^{+}e^{-}}$ to pp of the
baryon-to-meson ratio is significantly larger for charm-strange baryons (e.g.
$\Xi_{\mathrm{c}}^{+,0}$) than for $\Lambda_{\mathrm{c}}^{+}$, despite no
significant difference is observed for the
$\mathrm{D_{s}^{+}}$/$\mathrm{D^{0}}$ ratio. This highlights the problem to
understand strange production.

  * •

The observation for $\Sigma_{\mathrm{c}}^{0,+,++}/\mathrm{D^{0}}$ and the
augmented feed-down to $\Lambda_{\mathrm{c}}^{+}$ from
$\Sigma_{\mathrm{c}}^{0,+,++}$ decay. In general, one may interpret these
observations as a suppression of $\Sigma_{\mathrm{c}}^{0,+,++}$ and charm-
strange baryon production in $\mathrm{e^{+}e^{-}}$ due to “local” dynamical
constraints (e.g. diquark production in string breaking), which are absent or
mitigated by the emergence of other processes in hadronic collisions.

  * •

Finally, remarkable experimental data shows the contrast between the
apparently flat evolution with multiplicity and collision centrality of the
$p_{\mathrm{T}}$-integrated $\Lambda_{\mathrm{c}}^{+}$/$\mathrm{D^{0}}$ ratio
from pp to Pb–Pb collision systems, compared to the strong modification of the
$p_{\mathrm{T}}$-differential ratio.

  * •

On the theoretical side, coalescence and statistical hadronization emerged as
effective models to describe some of the measured yields in pp collisions. For
an extensive model like PYTHIA, which aspire at describing the full event
dynamics, the central role of MPI and colour reconnection in causing and
shaping the observed effects appears as an unavoidable conclusion suggested by
the measurements.

In this section, we discuss prospects and goals for the future, also touching
some topics and open points connected to heavy-flavour hadronization that
could not be deepened in the previous sections. These includes recent
measurements that have not yet conducted to firm conclusions, either because
still preliminary or because affected by large uncertainties, but that by
exploring different observables are pointing to new effects and directions of
search for the future. We also discuss new ideas and measurements, which have
the potential to lead toward an understanding of the comparison between
experimental measurements performed by different collaborations and the
mismatch with some theoretical calculations.

Consolidation of emerging picture

Many of the measurements presented in this review show striking effects, often
with results so far from the original expectation of the “baseline models”
described in Sect. II.1 that even large experimental uncertainties did not
prevent important conclusions. However, the natural first goal for the future
should be to have new measurements of the observables already inspected with a
reduction of the uncertainties and with an extension of the low and high
$p_{\mathrm{T}}$ and multiplicity reach. In particular, profiling the low
$p_{\mathrm{T}}$ and the $p_{\mathrm{T}}$-integrated baryon-to-meson ratios
for other particles than $\Lambda_{\mathrm{b}}^{0}$ (in particular of
$\Lambda_{\mathrm{c}}^{+}$, $\Xi_{\mathrm{c}}^{+,0}$, and
$\Sigma_{\mathrm{c}}^{0,+,++}$) at low multiplicity in pp collisions is of
uttermost importance to investigate the onset of the modification of
hadronization from $\mathrm{e^{+}e^{-}}$ to hadronic collisions. Measurements
of charm and strange baryons in nucleus–nucleus collisions would also be
particularly significant to verify whether the production rate of
$\Xi_{\mathrm{c}}^{+,0}$ and $\Omega_{\mathrm{c}}^{0}$ relative to
$\Lambda_{\mathrm{c}}^{+}$ increases in a system in which strange quarks are
expected to be abundantly produced given that light-flavour strange hadrons
yields match the expectation from chemical equilibrium without canonical
suppression. In the open beauty sector, the experimental study of
hadronization in nucleus–nucleus collisions has substantially only started
with the first measurements of strange and non-strange B mesons and non-prompt
D mesons. Though very important for first assessments, the large
uncertainties, the limited and not granular $p_{\mathrm{T}}$ coverage, and the
reduced sensitivity implied by the decay kinematics for non-prompt signals,
and the lack of measurements of beauty-baryon production, permit only
qualitative comparison with model expectations. With next LHC runs, also
thanks to the LHC detector upgrades, new measurements, including first ones of
$\Lambda_{\mathrm{b}}^{0}$ or non-prompt $\Lambda_{\mathrm{c}}^{+}$, should
allow for firmer conclusions. In pp collisions, LHCb data already provided a
rather detailed picture at forward rapidity, but precise measurements down to
low $p_{\mathrm{T}}$ are missing for both $\Lambda_{\mathrm{b}}^{0}$ and
$\mathrm{B_{s}^{0}}$ at midrapidity. As described in what follows, determining
the dependence on rapidity of heavy-flavour baryon-to-meson ratios must be an
important objective of experimental searches. Finally, measurements of the
production of strange and beauty baryons would also be very important, also
given what was observed for the
$\Xi_{\mathrm{c}}^{+,0}/\Lambda_{\mathrm{c}}^{+}$ ratio, but they are
experimentally extremely challenging and likely difficult to achieve in the
near future.

Rapidity puzzle

The comparison of the measurements performed by ALICE and LHCb at mid and
forward rapidity, respectively, show tantalising signs of a possible rapidity
dependence of the baryon-to-meson ratio in pp, p–Pb, and Pb–Pb collisions
[277, 220, 181, 219, 79, 275, 73, 221] for both the $\Lambda_{\mathrm{c}}^{+}$
and $\Xi_{\mathrm{c}}^{+,0}$ baryons, as displayed in Fig. 27. However, the
current measurements do not allow for conclusions, either because of the large
uncertainties (including those on the BR in the $\Xi_{\mathrm{c}}^{+,0}$ case)
or because of inconsistencies in the accessed $p_{\mathrm{T}}$ intervals or in
the definition and ranges used for multiplicity and centrality. It is
envisioned that the LHC collaborations in the future will, as far as possible,
perform measurements in the same centrality intervals and transverse momentum,
facilitating easier comparisons and reducing potential ambiguities.

In the beauty sector, the LHCb experiment measured a baryon-to-meson cross-
section ratio compatible between p–Pb and pp collisions [278] at both forward
and backward rapidity. However, also in this case, the data uncertainties
prevent to conclude whether $p_{\mathrm{T}}$ modifications of the baryon-to-
meson ratio from pp to p–Pb collisions observed in charmed sector are present
also in the beauty sector. In pp collisions ALICE measured at midrapidity the
ratio of non-prompt $\Lambda_{\mathrm{b}}^{0}$ to non-prompt $\mathrm{D^{0}}$,
which is strictly connected to the $\Lambda_{\mathrm{b}}^{0}/\mathrm{B^{0,-}}$
ratio [197]. However, the comparison with LHCb data at forward rapidity is
quite indirect, requiring a modelling of the $p_{\mathrm{T}}$-differential
cross section of beauty production tuned on LHCb data and of the decay
kinematics, which unavoidably dilutes possible effects over a broader
kinematic range. The CMS collaboration measured $\Lambda_{\mathrm{b}}^{0}$
production at 7 $\mathrm{TeV}$ but only for $p_{\mathrm{T}}>10$
$\mathrm{GeV}/c$ [196]. Measurements of the baryon-to-meson ratios at low
$p_{\mathrm{T}}$ in the beauty sector at midrapidity via fully reconstructed B
and $\Lambda_{\mathrm{b}}^{0}$ decays are envisioned for the next runs at the
LHC: they will be essential to investigate a possible rapidity dependence of
these ratios.

From a theoretical point of view, it was already briefly mentioned that PYTHIA
8 does not predict any difference with rapidity in the charm or beauty baryon-
to-meson ratio in pp collisions, and expects larger baryon-to-meson ratios at
forward rapidity than those measured by the LHCb Collaboration. It would be
extremely interesting to systematically investigate models that include
coalescence to see if any rapidity dependence in the baryon-to-meson emerges
and to understand if a possible effect would be driven by the medium density
and/or by possible differences in the heavy-quark density, which would be
important for multi-charm hadron states as well as for quarkonia and the
$\mathrm{B^{+}_{c}}$ meson, between mid- and forward/backward rapidity.
However, within a coalescence approach, one would need an appropriate
description of the evolution of both the bulk medium and the charm quark
momentum distribution from central rapidities to forward ones. In fact, the
fireball degree of thermalization and, in general, the time evolution at
forward rapidity already for the bulk matter have been less investigated in
hydrodynamics or transport theory, and even much less the dynamics of charm
quarks. Before making predictions of the
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio one should check that
$R_{\mathrm{A}A}(p_{\mathrm{T}})$ and $v_{2}(p_{\mathrm{T}})$ of D mesons at
forward rapidity can be correctly predicted, as done at midrapidity. This
would give good confidence on the level of charm-quark equilibration at
forward rapidity and would provide a reasonable knowledge of the charm-quark
$p_{\mathrm{T}}$ distribution.

Furthermore, the study of the directed flow of D mesons in nucleus–nucleus
collisions at RHIC and LHC has revealed that the bulk has a tilt with respect
to the longitudinal direction and has, therefore, a displacement with respect
to the charm quark at forward rapidity [279, 280, 281]. This could also imply
a different impact of the corona effect, or more generally of the amount of
non-thermalized heavy quarks that escape from the fireball. Summarizing, a
solid prediction for $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ within the
coalescence approach should be based on a solid description of the bulk and
charm distribution function in space and momentum as already done at
midrapidity.

Figure 27: Top left panel: the $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio
measured in p–Pb collisions at $\sqrt{s_{\mathrm{NN}}}$ = 5.02 $\mathrm{TeV}$
at mid, forward, and backward rapidity [277, 79, 219]. Top right panel (from
Refs. [221, 220]): comparison of the $\Xi_{\mathrm{c}}^{0}/\mathrm{D^{0}}$
ratio measured at midrapidity in p–Pb collisions at $\sqrt{s}$ = 5.02 TeV with
the measurements at forward and backward rapidity at
$\sqrt{s_{\mathrm{NN}}}=8.16$ $\mathrm{TeV}$. Bottom panel: the
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio measured in wide
$p_{\mathrm{T}}$ intervals shown as a function of rapidity in central,
semicentral collisions [73], and peripheral [275] Pb–Pb collisions. The
vertical lines represent the total uncertainties. ALICE data derive from the
extrapolation to $p_{\mathrm{T}}>0$ of the ratios measured in the interval
$1<p_{\mathrm{T}}<24$ $\mathrm{GeV}/c$. In the LHCb case, the centrality
interval indicated in Ref. [275] is $\langle N_{\mathrm{part}\rangle}=15.8\pm
10$. The percentile range corresponding to this interval can be estimated to
be approximately 65-90% based on Ref. [282, 283].

The role of diquarks

Coalescence models often calculate the formation of a baryon not as a
$3\rightarrow 1$ process but rather as a two-step process, in which first two
quarks form a diquark and then the diquark combines with a third quark.
However, this does not necessarily correspond to requiring the formation of
diquarks as real bound states nor as effective degrees of freedom of the
system. Whether this occurs or not is an open question. The suppression of the
production of $\Sigma_{\mathrm{c}}^{0,+,++}$ states relative to
$\Lambda_{\mathrm{c}}^{+}$ ones in $\mathrm{e^{+}e^{-}}$ collisions and its
removal in hadronic collisions is striking and appears to support the picture
of diquark formation in string breaking.

While it is very hard to design what could be a “smoking gun” observation
supporting the formation of diquarks is worth noticing that in [284] it was
suggested that a measurement of $\Lambda_{c}$ elliptic flow in Pb–Pb and its
difference w.r.t. the D-meson one would allow to infer the relevance of
diquark degrees of freedom in the QGP. This possibility is motivated by the
comparison of the $p_{\mathrm{T}}$ dependence of the difference between
$v_{2}(\Lambda_{c})$ and $v_{2}(D^{0})$ predicted in a coalescence approach
and in the POWLANG mechanism. In the coalescence model the $v_{2}$ of
$\Lambda_{c}$ receives a contribution from two light quarks while the D meson
only from one. In the POWLANG approach instead $\Lambda_{\mathrm{c}}^{+}$
combines to a diquark evolving hydrodynamically with the bulk medium implying
a milder difference w.r.t. D, in particular at intermediate $p_{\mathrm{T}}$.

String-related investigations

There are several areas of interest in extending the current modelling of
string hadronization. Some have already been mentioned above, including
modelling of exotica and the heavy-ion limit of strings. Further areas of
interest include considering colour reconnections in the initial state, the
collective behaviour of strings, interactions between final state hadrons,
studies into baryon production mechanisms, and even looking beyond Lund
strings. Here we elaborate on just a few of these areas of interest.

One noticeable flaw as shown in Fig. 17 is the overprediction of the
$\mathrm{\Lambda_{\mathrm{b}}^{0}}/B^{0}$ ratio by PYTHIA despite well
predicting the $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio due to the
inclusion of junction topologies. Exploration into other ratios involving
these heavy flavour hadrons such as
$\Lambda_{\mathrm{b}}^{0}/\Lambda_{\mathrm{c}}^{+}$ and $\mathrm{B}^{0}/D^{0}$
would be particularly insightful into attempting to find the root of the
discrepancy between model and experiment. Another set of observables that
would be particularly interesting for studying junction modelling are baryon
correlations. In the string model, baryons produced via standard diquark-
antidiquark pair creation from string breaks would be highly correlated. This
correlation would be reduced given the popcorn model for diquark creation,
which allows meson production between diquark-antidiquark pairs. Comparatively
there is generally no expectation of a correlation between junction baryons
and antijunction baryons.

As presented throughout this paper, there have been several already
experimentally observed behaviours correlated to the size of the system,
whether it be observables as a function of charged particle multiplicities or
contrasting $\mathrm{e^{+}e^{-}}$, pp and AA results. Many measurements have
already shown clear indications of collective-like effects in pp collision
environments, which prompts two key questions; are these effects present in
$\mathrm{e^{+}e^{-}}$ collisions (and/or in
$\mathrm{e^{\pm}p}$/$\mathrm{e^{\pm}}$A ones), and how can we model such
collective effects in the string picture? With regards to the former, many of
these effects were not studied at LEP and hence there is little indication
whether they are present in $\mathrm{e^{+}e^{-}}$ or not, thus it would be
intriguing to study them at future high energy $\mathrm{e^{+}e^{-}}$,
$\mathrm{e^{\pm}p}$, and $\mathrm{e^{\pm}}$A experiments [285].

Many collective effects in the string picture have already been studied in
PYTHIA, including string tension modifications (e.g. Rope hadronization [47]
and close-packing [48]), string interactions (e.g. string shoving [286]) and
final state hadron interactions (e.g. hadron rescattering [287]). Such models
have already shown interesting effects and compelling agreement to data. For
example, the string shoving model [286] has already been shown to produce the
near-sided ridge experimentally observed in pp collision events which is
unexpected with the baseline string model. Studies of such models in AA
environments will be particularly insightful to help determine whether a
string model can reproduce other behaviours also typically considered QGP
signatures such as flow. These models would also have particularly interesting
effects on $p_{\mathrm{T}}$ spectra of hadrons, and further phenomenology
studies would be useful to determine any signature behaviours of such models
that experiments could attempt to identify. One can also consider interactions
between final state hadrons via hadron rescattering [287], which would have
effects particularly for high multiplicity events. For heavy flavour hadrons
in particular, one would expect effects on the shape of the $p_{\mathrm{T}}$
spectra. It is interesting to note that such rescattering effects would also
be interesting for studying the $p/\pi$ ratio with respect to charged
multiplicity in pp collision events, which PYTHIA currently overpredicts even
without the inclusion of junctions.

Perhaps one of the most intriguing questions for the string model is whether a
string in a densely packed string system hadronizes differently from a string
in vacuum. There are already indications that this is the case given
enhancement of the light strange hadron-to-pion ratio as a function of charged
multiplicity as seen by ALICE in Ref. [44]. Higher multiplicity events either
correspond the hadronization of longer strings or systems with a higher
quantity of strings, and hence generally more densely packed string
environments. There are already existing models such as Rope hadronization
[47] and close-packing [48] which attempt to implement a strangeness
enhancement mechanism by enhancing the string tension in more dense string
systems. An increased string tension results in the reduction of the
strangeness suppression via the Schwinger mechanism due to the strange quark
mass. Both Rope hadronization and close-packing models have already been able
to describe the general behaviour seen in the ALICE [44] data for the strange
hadron-to-pion ratios as a function of multiplicity. Additionally, in pp
collisions the majority of strings are aligned along the beam axis, which in
turn would result in strangeness enhancement primarily for strings at low
$p_{\mathrm{T}}$, which would lead to strangeness enhancement for low
$p_{\mathrm{T}}$ hadrons. One may be hopeful that this low $p_{\mathrm{T}}$
strangeness enhancement could describe the ratios involving heavy strange
hadrons such as $\Xi_{\mathrm{c}}^{+,0}$, however as mention in Sec. III.3
these strangeness enhancement mechanisms alone are unable to describe the
level of strangeness enhancement needed to fit the
$\Xi_{\mathrm{c}}^{+,0}$/$\Lambda_{\mathrm{c}}^{+}$ data whilst maintaining
reasonable model parameters. Additional measurements involving strange
b-hadrons would be vital in examining such mechanisms, so that one can make
comparisons between distributions such as the
$\Xi_{\mathrm{c}}/\Lambda_{\mathrm{c}}$ and
$\Xi_{\mathrm{b}}/\Lambda_{\mathrm{b}}^{0}$ ratios alongside the
$\Xi_{\mathrm{c}}/\mathrm{D^{0}}$ and $\Xi_{\mathrm{b}}/\mathrm{B}^{0}$
ratios. These ratios would provide particularly insightful indications to the
source of strangeness enhancement and whether the behaviour is different in
the beauty and charm baryon sectors. Further data for ratios involving double
strange heavy baryons such as $\Omega_{\mathrm{c}}$ and $\Omega_{\mathrm{b}}$
would also provide very useful insight. One could also explore the enhancement
of strangeness around the junction itself which is not an unreasonable
assumption as we know little about what the field around a junction itself
should look like, and such a model would have significant impact as the
majority of heavy flavour baryons are sourced from junctions. Likewise due to
the sensitivity of heavy baryons on junction formation, further data for
ratios involving $\Xi_{\mathrm{b}}$ and $\Omega_{\mathrm{b}}$ may even help
identify the source of the overprediction of the
$\Lambda_{\mathrm{b}}^{0}$/$\mathrm{B}^{0}$ ratio using PYTHIA CR-BLC.

To study in further detail the hadronization of a string in vacuum compared to
a densely packed string system, phenomenological studies would be insightful
to examine hadron distributions for a single string compared to higher
multiplets, given all the aforementioned models, in order to check if there
are any identifiable signatures experiments could attempt to find. To directly
compare the fragmentation of a triplet vs. octet string, it would be
insightful to examine hadron spectra in the tips of gluon jets compared to
quark jets, LEP hairpin configurations, and perhaps even diffractive LHC
events may allow for relatively clean environments to search for such
behaviours.

Strings in heavy-ion collisions?

Through this review, models in PYTHIA useful for $\mathrm{e^{+}e^{-}}$ up to
pp collisions were discussed. PYTHIA itself does not have the ability to
perform pA and AA collisions. Instead, to simulate pA and AA collisions one
can use the Angantyr model [288], which uses PYTHIA 8 as its base, taking
inspiration from the Fritiof model [289, 290] and using the idea of wounded
nucleons. A key feature of Angantyr is testing the assumption of continuing
the usage of PYTHIA’s string fragmentation mechanisms to model the
hadronization process also in pA and AA collisions, rather than assuming the
formation of a QGP. This includes the implementation of the QCD-CR mechanisms
described in Sect. II.2.1, which allows CR between different nucleon-nucleon
collisions [291, 292]. In principle one can also use Rope hadronization [47]
in conjunction with Angantyr, however the combination of Rope hadronization
and Angantyr is yet to be validated to work properly, nor have the model
parameters been tuned. Though Angantyr in its current state does not fully
describe the collective behaviour of high string densities, nonetheless the
existing implementation already provides insight into the non-collective
background of such observables and provides an important comparison between
the QGP and the string model for hadronization.

This also prompts the interesting question of whether there is a some
quenching limit for strings. That is, within this type of modelling paradigm,
at high enough string densities, do the string states cease to exist and
instead do we get a transition into some deconfined QGP state? This then begs
the question when would this transition occur, what would this transitional
state look like, and how would one model it. These are largely unexplored
questions within the string picture, and the study of AA collisions with
string hadronization using Angantyr provides the first step to exploring such
questions.

Heavy-flavour jets and correlations

The study of beauty- and charm-tagged jets provides additional information for
interpreting the modification of the baryon-to-meson ratios from
$\mathrm{e^{+}e^{-}}$ to hadronic collisions, giving a handle to connect it to
possible modification of heavy-quark fragmentation and jet properties.
Measurements of jets tagged by the presence of a fully reconstructed heavy-
flavour hadron can be particularly informative. This tagging technique, with
respect to others based on the recognition of a displaced secondary vertex or
on lepton tagging, has the advantage of providing the momentum vector of the
heavy-flavour hadron and its correlation to the jet momentum. Though used
since UA1 times, in recent years it was refined and exploited to study the
production and sub-structure properties of D-meson tagged jets [293, 294, 295,
165, 296, 297, 298], as well as for the first direct observation of the “dead-
cone” effect [299]. These measurements allow testing the modelling of jet
fragmentation, as well as of charm-production processes, in event generators
like PYTHIA and HERWIG [300, 301, 17], possibly interfaced with POWHEG [205,
302, 303]. Relevant variables include the fraction of jet momentum carried by
the charm hadron, calculated along the jet axis direction
($z_{\mathrm{||}}^{\mathrm{ch}}$), and the transverse momenta and
(pseudo)rapidities of hadrons measured with respect to the jet axis.

As extensively discussed in Ref. [296], the $z_{\mathrm{||}}^{\mathrm{ch}}$
distribution of D-meson tagged jets measured by the ALICE Collaboration is
quite well described by Monte Carlo event generators. Recently, ALICE reported
a first measurement of the same distribution for
$\Lambda_{\mathrm{c}}^{+}$-tagged jets. The comparison with D-meson tagged jet
data done with the same jet and hadron kinematic selection is reported in Fig.
28. Though the uncertainties prevent firm conclusions, the data suggest that
$\Lambda_{\mathrm{c}}^{+}$ carry a smaller fraction of jet momentum along the
jet axis with respect to D mesons. This trend contrasts with that expected by
the default version of PYTHIA 8 (Monash tune), while it agrees with that from
PYTHIA 8 CR-BLC. The trend expected by Monash reflects the tendency in jet
fragmentation that heavier particles carry a larger fraction of jet momentum,
which in a collinear fragmentation picture can be related to energy
conservation. The data therefore suggest that either the
$\Lambda_{\mathrm{c}}^{+}$ production is favoured in jets with higher momentum
and softer fragmentation, either intrinsically or because modified by the
interplay with the underlying event, or that the $\Lambda_{\mathrm{c}}^{+}$
momentum vector is less aligned with the jet one, either as an intrinsic
feature of charm-jet fragmentation to $\Lambda_{\mathrm{c}}^{+}$ baryons (not
expected however in Monash) or because of a angular decorrelation at
hadronization time. A partly related study is that of the azimuthal
correlation of $\Lambda_{\mathrm{c}}^{+}$ baryons with charged particles
produced in the event.

A first preliminary measurement, performed by ALICE [304], indicates that the
yields of low-$p_{\mathrm{T}}$ associated particles are large, both in the
near-side and in the opposite direction to the $\Lambda_{\mathrm{c}}^{+}$
(away-side jet region), and resemble those observed only at much higher
D-meson $p_{\mathrm{T}}$ ($p_{\mathrm{T}}\gtrsim 8$ $\mathrm{GeV}/c$) in
correlations of D mesons with charged particles. Naively one could expect
that, except for the initial-parton energy correlation, the fragmentation
(parton-shower) and hadronization of the near- and away-side jets should
proceed almost independently. It is thus surprising that the larger near-side
yield, which is consistent with the softer jet fragmentation suggested by the
$z_{\mathrm{||}}^{\mathrm{ch}}$ analysis, is accompanied by a difference in
the away-side jet.

Taken together, the $z_{\mathrm{||}}^{\mathrm{ch}}$ and azimuthal correlation
data may suggest that, with respect to D mesons, low-$p_{\mathrm{T}}$
$\Lambda_{\mathrm{c}}^{+}$ are produced on average in jets with higher
$p_{\mathrm{T}}$ and softer fragmentation, with a sharing of the jet energy
among several soft particles. Whether this could be due to an interplay with
the underlying event is not clear, though possible, considering also the
multiplicity dependence of $\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ ratio.

More precise measurements, possibly more differential in jet and hadron
kinematics, and a detailed analysis of the jet radial profiles and
substructure and of the correlation-function widths could in the future shed
light on these intriguing observations. Similar measurements in the beauty
sector would also be very useful.

On top of new measurements, it would be desirable in the future to understand
what hadron-formation via coalescence implies for these observables. In the
coalescence process, the $\Lambda_{\mathrm{c}}^{+}$ at a given
$p_{\mathrm{T}}$ gets on average more of its momentum from the light quarks
with respect to the D meson (because two light quarks are involved for forming
$\Lambda_{\mathrm{c}}^{+}$). This can decorrelate the resulting hadron
momentum vector from the jet axis. Therefore, though explicit calculations are
missing, it would be in principle compatible with the observation to assume
that in the jet a $\Lambda_{\mathrm{c}}^{+}$ formed by coalescence of the
charm quark with light quarks from the underlying event carries a smaller
longitudinal momentum fraction (along the jet axis) relative to that carried
by D mesons.

Figure 28: Comparison (top panel) and ratio (bottom panel) of the
distributions of the $\Lambda_{\mathrm{c}}^{+}$ and $\mathrm{D^{0}}$
longitudinal track-based jet-momentum fraction
($z_{\mathrm{||}}^{\mathrm{ch}}$). Track-based jets, defined by their charged-
particle constituents only, with $7<p_{\mathrm{T}}<15$ $\mathrm{GeV}/c$ are
tagged with fully reconstructed $\Lambda_{\mathrm{c}}^{+}$ and
$\mathrm{D^{0}}$ with $3<p_{\mathrm{T}}<15$ $\mathrm{GeV}/c$. In the bottom
panel the ratio is compared to expectations from PYTHIA 8 default tune
(Monash) and PYTHIA 8 tune with colour-reconnection implemented beyond leading
colour. Figure from Ref. [298]

High-mass excited baryon states

Experimental search of not yet observed charm and beauty hadron resonances is
of primary importance in the context of hadronization. It is needed to
determine whether the additional baryon resonance states predicted by the RQM,
which are largely supported by lattice-QCD computations of the vacuum
spectroscopy, exist and to determine their spectral properties. As discussed
extensively in this review, the existence of these states is a necessary
assumption in many models to reproduce the
$\Lambda_{\mathrm{c}}^{+}/\mathrm{D^{0}}$ yield ratio in pp and heavy-ion
collisions, in particular for models based on statistical hadronization.
Though it remains an open question why these high-mass resonance states have
not been observed so far in $\mathrm{e^{+}e^{-}}$ collisions, proving their
existence and measuring the branching ratio to relevant decay channels could
be pursued at $\mathrm{e^{+}e^{-}}$ colliders, e.g. by the Belle II
Collaboration, as well as in pp collisions by the ATLAS, CMS, and LHCb
Collaborations, leveraging the large LHC luminosities.

Actually, on top of verifying the number and spectral properties of yet
unobserved states, a fundamental test for SHM would be given by a measurement
of the production yields of excited mesons and baryons, both known and unknown
ones. ALICE recently performed the first measurement of
$\mathrm{D}^{+}_{s1}(2536)$ and $\mathrm{D}^{*+}_{s2}(2573)$ resonance yields
in pp collisions, which provide further pieces of information to investigate
the charm-quark hadronization [305]. Cross-section measurements of excited
$\Lambda_{\mathrm{c}}^{+}$, $\Sigma_{\mathrm{c}}^{0,+,++}$,
$\Xi_{\mathrm{c}}^{+,0}$ states should be feasible with the Run 3 data at the
LHC.

$\mathrm{B^{+}_{c}}$ meson

The $\mathrm{B^{+}_{c}}$ meson contains a b and a c quark, making it (and its
excited states) the only hadron observed so far that contains heavy quarks of
two different flavours. The $\mathrm{B^{+}_{c}}$ meson has a net flavour
number, which in principle secludes it from the quarkonium family. However,
one could support that the proximity of its quark content to the charmonium
and bottomonium states (“hidden” heavy-flavour mesons) makes it an exotic
heavy quarkonium.

The $\mathrm{B^{+}_{c}}$ meson has a particularly low production cross section
in pp collisions, because two heavy quark pairs (b$\mathrm{\overline{}}{b}$
and c$\mathrm{\overline{}}{c}$) need to be produced. The obvious mechanism is
for both of these heavy-quark pairs to be produced in a single hard-scattering
process, but since this is quite suppressed, additional production mechanisms
may play a role such as multi-parton interactions in high-multiplicity pp
collisions [306] and/or quark recombination in heavy-ion ones. That is, a b
quark produced in a given hard scattering could recombine with a c quark from
another parton-parton or nucleon-nucleon collision (presumably with some
space-time suppression since the pair needs to be sufficiently close to form a
$\mathrm{B^{+}_{c}}$).

If this is a significant mechanism in heavy-ion collisions, the
$\mathrm{B^{+}_{c}}$ yields might be dramatically augmented compared to the
expectations from $N_{\mathrm{{coll}}}$-scaling, considering the small primary
production expected from the pp cross section. The CMS Collaboration has
measured for the first time the $\mathrm{B^{+}_{c}}$ production in heavy-ion
collisions and reported that the $\mathrm{B^{+}_{c}}$ meson is shown to be
less suppressed than quarkonia and most of the open heavy-flavour mesons,
suggesting that effects of the hot and dense nuclear matter created in heavy-
ion collisions contribute to its production [307]. However, the acceptance of
the CMS detector only allows to reconstruct $\mathrm{B^{+}_{c}}$ mesons of
relatively high transverse momentum, typically $p_{\mathrm{T}}$ $>$ 6 GeV/$c$,
while the majority of the $\mathrm{B^{+}_{c}}$ mesons produced by
recombination should have a $p_{\mathrm{T}}$ lower than the acceptance
threshold of CMS. This is because most of the charm quarks present in the QGP
have low $p_{\mathrm{T}}$ and because the recombining quarks need to have a
small relative momentum (typically not larger than the binding energy) to be
able to recombine.

Hence, it would be extremely beneficial to measure low transverse momentum
$\mathrm{B^{+}_{c}}$ in heavy-ion collisions and in high-multiplicity pp ones.
Recently, Egede et al. [306] predicted a large fraction of the produced yield
in high-multiplicity pp collisions to come from MPI (albeit without accounting
for possible space-time suppression effects so this should probably be taken
as a maximum), and Rapp et al. [308] predict a similar trend in heavy-ion
collisions from regenerated $\mathrm{B^{+}_{c}}$ mesons. In this latter study
they observed that the large regeneration contributions cause a markedly
rising $R_{\mathrm{AA}}$ with centrality, reaching values of up to 4-6 in
central collisions. In this model, the inclusive $p_{\mathrm{T}}$-dependent
$R_{\mathrm{AA}}$ for the $\mathrm{B^{+}_{c}}$ is dominated by regeneration
contributions for $p_{\mathrm{T}}$ $<$ 10-15 GeV/$c$ in semi-/central
collisions, reaching values of around 10 at low $p_{\mathrm{T}}$. It is
pointed out in the publication that the results for the $R_{\mathrm{AA}}$ are
rather sensitive to the $\mathrm{B^{+}_{c}}$ production cross section in pp
collisions. A more precise measurement of this quantity, and a significant
reduction of the statistical and systematic uncertainties of
$\mathrm{B^{+}_{c}}$ in heavy-ion collisions will provide unprecedented
insights into heavy-flavour in-medium and is a promising new probe to study
the interplay between the suppression and enhancement mechanisms in the
production of heavy-flavour mesons in the QGP.

Multi-charm hadron states and charm nuclei

Analogously to the $\mathrm{B^{+}_{c}}$ case, another window to investigate
colour-reconnection and coalescence are the measurements of multi-charm
baryons and of charm nuclei. Measurements of multi-charm hadrons, such as the
$\Xi^{+}_{\mathrm{c}c}$, $\Xi^{++}_{\mathrm{c}c}$, $\Omega^{+}_{\mathrm{c}c}$,
$\Omega^{++}_{\mathrm{c}cc}$ would provide a direct window on hadron formation
from a deconfined quark–gluon plasma. In fact, the yields of multicharm
baryons relative to the number of produced charm quarks are predicted to be
significantly enhanced in AA relative to pp collisions [49, 75, 309, 310, 311,
312, 313]. Enhancements are expected by as much as a factor $10^{2}$ for
$\mathrm{\Xi}_{cc}$ baryon and even by as much as a factor $10^{3}$ for the as
yet undiscovered $\mathrm{\Omega}_{ccc}$ baryon. Most recently, the SHMc and
the Catania models [49, 75] provided predictions for the production of multi-
charmed hadrons in Pb–Pb collisions and they explored the system size
dependence through Kr–Kr, to Ar–Ar and O–O collisions. The yields are
predicted to be quite similar in the case of full thermalization of the charm
quark, a condition that might not be achieved in light systems like O-O. The
coalescence picture coupled to a realistic dynamical simulation for charm
quarks reveals that $\mathrm{\Omega}_{ccc}^{++}$ yield, and even more its
$p_{\mathrm{T}}$ spectrum, would provide a sensitive probe of charm
thermalization and its evolution with the system size of the colliding system
[75]. Such studies are currently far beyond the experimental reach and a
comprehensive program of measurements in the multi-charm sector is discussed
to be carried out with the proposed future ALICE 3 detector [76].

The unique combination of high rates, wide acceptance, strong particle
identification capabilities, and excellent secondary vertex determination of
the proposed ALICE 3 detector, also provide a toolbox for the search for light
nuclei with charm. At LHC energies, the most promising candidates are the
c-deuteron, c-triton and c-3He. The existence of such weakly decaying bound
states with lifetimes similar to those of other unbound charmed baryons is
debated in the literature. In case these states are bound and under the
assumption that their abundance in collisions of nuclei is described by the
SHM, the yields are large enough to bring their experimental discovery within
reach in the next decades and to additionally test different hadronization
models.

Exotica

The $\chi_{c1}(3872)$ is an exotic particle that was first observed by the
Belle Collaboration [314], and then confirmed and studied by other experiments
at electron-positron [315] and hadron colliders [316, 317, 318, 319, 320]. The
nature of this particle is still not fully understood and interpretations in
terms of conventional charmonium (a bound state of charm-anticharm quarks),
D-meson molecule, tetraquark, or an admixture these states have been proposed
[321]. Heavy-ion collisions, as well as high multiplicity pp collisions, offer
a unique window into the properties of this poorly understood hadron. In these
systems, promptly produced $\chi_{c1}(3872)$ hadrons can interact with other
particles in the nucleus and/or those produced in the collision. The influence
of these interactions on the observed $\chi_{c1}(3872)$ yields provides
information that can help discriminating between the various models of its
structure, as well as give insight into the dynamics of the bulk of the
particles produced in these collisions. The LHCb Collaboration has measured a
significant decrease in the ratio of prompt $\chi_{c1}(3872)$ to $\psi$(2S)
cross sections in pp collisions with increasing charged-particle multiplicity
[322]. These data were interpreted in terms of breakup of the
$\chi_{c1}(3872)$ hadron due to interactions with comoving particles produced
in the event, for both compact and molecular models of $\chi_{c1}(3872)$
structure [323, 324]. A relatively-compact tetraquark is preferred in this
case, but it is currently discussed how conclusive this comparison is [325,
326, 327, 328]. The LHCb Collaboration has also recently performed the same
measurement in p–Pb collisions [329]. In this collision system the ratio is
measured to be higher than the one measured in pp collisions. It must be noted
that the evolution of $\psi$(2S) production with particle multiplicity is not
fully understood yet. This complicates the interpretation of $\chi_{c1}(3872)$
data, when $\psi$(2S) is used as a reference. The CMS collaboration has
measured the $\mathrm{\sigma}_{\chi_{c1}(3872)}/\sigma_{\psi(2S)}$ ratio in
Pb–Pb collisions at high $p_{\mathrm{T}}$, and found that the ratio is further
enhanced relative to pp and p–Pb collisions, although the uncertainties and
the significantly different $p_{\mathrm{T}}$ intervals of the measurements
preclude drawing firm conclusions [330]. Calculations based on quark
coalescence show that production rates of $\chi_{c1}(3872)$ hadrons in AA
collisions are sensitive to its structure. In these models, the production of
compact tetraquarks is expected to be enhanced over hadronic molecules [331,
332]. However, a recent transport calculation reaches the opposite conclusion
[325]. Late-stage interactions in the hadron gas phase of a heavy-ion
collision can also modify the yields [333]. It is currently unknown where the
crossover between the suppressing effect of breakup and the enhancing effect
of coalescence may occur. The increase of the ratio of cross sections
$\mathrm{\sigma}_{\chi_{c1}(3872)}/\sigma_{\psi(2S)}$ from pp to p–Pb to Pb–Pb
collisions may indicate that the exotic $\chi_{c1}(3872)$ hadron experiences
different dynamics in the medium than those relevant for a “conventional”
charmonium state like $\psi$(2S). If the experiments could in the future cover
similar $p_{\mathrm{T}}$ intervals in the different collision systems and
multiplicities/centralities, possibly better assessing the lower
$p_{\mathrm{T}}$ region, these measurements will provide new constraints to
determine the internal $\chi_{c1}(3872)$ structure as well as to model parton
transport and hadronization in nuclear collisions. In this context, it will
also be important to precisely measure the production cross section of other
exotic states, like the $\mathrm{T}^{+}_{cc}$.

A related possible future implementation in PYTHIA is the inclusion of
exotica. Not only does QCD-CR allow beyond LC topologies such as junctions,
but it also allows for more complicated structures called junction networks.
These junction networks are string configurations that consist of multiple
junctions within the structure and provide a natural way within the string
picture to represent tetraquarks, pentaquarks, etc. For example, the simplest
of these structures are junction-antijunction configurations which consist of
two colour-charged endpoints and two anticolour-charged endpoints as shown in
Fig. 29. When these junction-antijunction configurations are encountered as
high mass systems, PYTHIA currently handles these structures by making breaks
on the string piece connecting the junction to the antijunction, and then
fragmenting each junction system separately. However, should a junction-
antijunction system have a low invariant mass, one could propose that the
system does not fragment and instead can be treated as a tetraquark. The same
argument holds for larger junction networks. The inclusion of exotica is not
yet implemented in PYTHIA and it would not be straightforward to implement the
required significant modifications to the existing code, however junction
networks are already present and could provide a sensible starting point for
modelling such structures.

Figure 29: Junction-antijunction formation in the string model. The left
images show the initial string configuration in the LC limit, and the right
images depict the possible alternative configuration given the respective
colour reconnection models. Figure 30:  Elliptic flow of prompt
$\mathrm{J}/\psi$ at forward+backward rapidity ($-2.86<y_{\mathrm{cm}}<-1.86$
or $0.94<y_{\mathrm{cm}}<1.94$), and of $\mathrm{D^{0}}$, $\Lambda$ and
$\mathrm{K^{0}_{s}}$ at midrapidity shown as a function of $p_{\mathrm{T}}$ in
high-multiplicity p–Pb collisions at $\sqrt{s}$ =8.16$\mathrm{TeV}$ [334,
335].

Collectivity in small systems

Positive $v_{\mathrm{2}}$ values have been observed for both open charm [335,
336, 337] and $\mathrm{J}/\psi$ [338, 334] also in high-multiplicity p–Pb
collisions, as well as in high-multiplicity pp collisions [339]. As shown in
Fig. 30, in high-multiplicity p–Pb collisions $\mathrm{D^{0}}$ and
$\mathrm{J}/\psi$ $v_{\mathrm{2}}$ reach similar maximum values, which are
about a factor of two smaller than the maximum $v_{\mathrm{2}}$ of light
hadrons [340] and significantly smaller than the maximum observed in semi-
central Pb–Pb collisions for D mesons. Though the maximum $\mathrm{J}/\psi$
$v_{\mathrm{2}}$ values are similar in high-multiplicity p–Pb collisions and
semi-central Pb–Pb collisions, in the latter case positive values, close to
open-charm and light-flavour ones, are observed up to higher $p_{\mathrm{T}}$.
The origin of flow in small collisions systems are still debated, even for
light-flavour particles [341]. Hadronization via recombination could induce a
positive $v_{\mathrm{2}}$ for open-charm hadrons even in the case of a
small/negligible $v_{\mathrm{2}}$ of charm quarks [66]. This is of course not
the case for charmonia. The absence of clear experimental evidences of charm
energy loss, makes the observation of open charm and charmonia flow
particularly interesting to understand the interplay of initial-state and
final-state effects for the generation of azimuthal anisotropy in momentum
space. A possible dynamic not yet firmly identified could lead to an
$R_{\mathrm{A}A}(p_{\mathrm{T}})$ close to unity together with a sizeable
$v_{2}(p_{\mathrm{T}})$ is the super-diffusion in the initial stage glasma
phase that would lead to an increase of the $R_{\mathrm{A}A}(p_{\mathrm{T}})$
[263, 342, 343] followed by a moderate energy loss that generate a sizeable
elliptic flow.

## Acknowledgements

J.A. and P.S. acknowledge support from The Royal Society and the Wolfson
Foundation (UK), from the Australian Research Council Discovery Project
DP230103014, and from the Monash-Warwick Alliance for Particle Physics. V.G.
acknowledge support PRIN2022 under fundings EU Next Generation (project code
2022SM5YAS).

## References

  * Schwarz [2003] D. J. Schwarz, The first second of the universe, Annalen Phys. 12, 220 (2003), arXiv:astro-ph/0303574 .
  * Castorina  _et al._ [2015] P. Castorina, V. Greco, and S. Plumari, QCD equation of state and cosmological parameters in the early universe, Phys. Rev. D 92, 063530 (2015), arXiv:1505.07655 [hep-ph] .
  * Cacciari  _et al._ [1998] M. Cacciari, M. Greco, and P. Nason, The P(T) spectrum in heavy flavor hadroproduction, JHEP 05, 007, arXiv:hep-ph/9803400 .
  * Cacciari  _et al._ [2012] M. Cacciari, S. Frixione, N. Houdeau, M. L. Mangano, P. Nason, and G. Ridolfi, Theoretical predictions for charm and bottom production at the LHC, JHEP 10, 137, arXiv:1205.6344 [hep-ph] .
  * Kniehl  _et al._ [2005a] B. A. Kniehl, G. Kramer, I. Schienbein, and H. Spiesberger, Inclusive D*+- production in p anti-p collisions with massive charm quarks, Phys. Rev. D 71, 014018 (2005a), arXiv:hep-ph/0410289 .
  * Kniehl  _et al._ [2005b] B. A. Kniehl, G. Kramer, I. Schienbein, and H. Spiesberger, Collinear subtractions in hadroproduction of heavy quarks, Eur. Phys. J. C41, 199 (2005b), arXiv:hep-ph/0502194 [hep-ph] .
  * Kniehl  _et al._ [2012] B. Kniehl, G. Kramer, I. Schienbein, and H. Spiesberger, Inclusive Charmed-Meson Production at the CERN LHC, Eur. Phys. J. C 72, 2082 (2012), arXiv:1202.0439 [hep-ph] .
  * Helenius and Paukkunen [2018] I. Helenius and H. Paukkunen, Revisiting the D-meson hadroproduction in general-mass variable flavour number scheme, JHEP 05, 196, arXiv:1804.03557 [hep-ph] .
  * Braaten  _et al._ [1995] E. Braaten, K.-M. Cheung, S. Fleming, and T. C. Yuan, Perturbative QCD fragmentation functions as a model for heavy quark fragmentation, Phys. Rev. D 51, 4819 (1995), arXiv:hep-ph/9409316 .
  * Sjöstrand  _et al._ [2006] T. Sjöstrand, S. Mrenna, and P. Z. Skands, PYTHIA 6.4 Physics and Manual, JHEP 05, 026, arXiv:hep-ph/0603175 [hep-ph] .
  * Skands  _et al._ [2014] P. Skands, S. Carrazza, and J. Rojo, Tuning PYTHIA 8.1: the Monash 2013 Tune, Eur. Phys. J. C 74, 3024 (2014), arXiv:1404.5630 [hep-ph] .
  * Bahr  _et al._ [2008] M. Bahr _et al._ , Herwig++ Physics and Manual, Eur. Phys. J. C 58, 639 (2008), arXiv:0803.0883 [hep-ph] .
  * Aaltonen  _et al._ [2008] T. Aaltonen _et al._ (CDF), Measurement of Ratios of Fragmentation Fractions for Bottom Hadrons in $p\bar{p}$ Collisions at $\sqrt{s}$ = 1.96-TeV, Phys. Rev. D 77, 072003 (2008), arXiv:0801.4375 [hep-ex] .
  * ALI [2022a] The ALICE experiment – A journey through QCD, (2022a), arXiv:2211.04384 [nucl-ex] .
  * d’Enterria and Loizides [2021] D. d’Enterria and C. Loizides, Progress in the Glauber Model at Collider Energies, Ann. Rev. Nucl. Part. Sci. 71, 315 (2021), arXiv:2011.14909 [hep-ph] .
  * Bierlich  _et al._ [2022a] C. Bierlich _et al._ , A comprehensive guide to the physics and usage of PYTHIA 8.3 10.21468/SciPostPhysCodeb.8 (2022a), arXiv:2203.11601 [hep-ph] .
  * Bellm  _et al._ [2016] J. Bellm _et al._ , Herwig 7.0/Herwig++ 3.0 release note, Eur. Phys. J. C 76, 196 (2016), arXiv:1512.01178 [hep-ph] .
  * Andersson  _et al._ [1983] B. Andersson, G. Gustafson, G. Ingelman, and T. Sjöstrand, Parton Fragmentation and String Dynamics, Phys. Rept. 97, 31 (1983).
  * Eichten  _et al._ [1978] E. Eichten, K. Gottfried, T. Kinoshita, K. D. Lane, and T.-M. Yan, Charmonium: The Model, Phys. Rev. D 17, 3090 (1978), [Erratum: Phys.Rev.D 21, 313 (1980)].
  * Bali and Schilling [1992] G. S. Bali and K. Schilling, Static quark - anti-quark potential: Scaling behavior and finite size effects in SU(3) lattice gauge theory, Phys. Rev. D 46, 2636 (1992).
  * Andersson  _et al._ [2001] B. Andersson, S. Mohanty, and F. Söderberg, The Lund fragmentation process for a multigluon string according to the area law, Eur. Phys. J. C 21, 631 (2001), arXiv:hep-ph/0106185 .
  * Schwinger [1951] J. S. Schwinger, On gauge invariance and vacuum polarization, Phys. Rev. 82, 664 (1951).
  * Bierlich  _et al._ [2022b] C. Bierlich, S. Chakraborty, G. Gustafson, and L. Lönnblad, Hyperfine splitting effects in string hadronization, Eur. Phys. J. C 82, 228 (2022b), arXiv:2201.06316 [hep-ph] .
  * Hamacher and Weierstall [1995] K. Hamacher and M. Weierstall (DELPHI), The Next round of hadronic generator tuning heavily based on identified particle data, (1995), arXiv:hep-ex/9511011 .
  * Barate  _et al._ [1998] R. Barate _et al._ (ALEPH), Studies of quantum chromodynamics with the ALEPH detector, Phys. Rept. 294, 1 (1998).
  * Amoroso  _et al._ [2019] S. Amoroso, S. Caron, A. Jueid, R. Ruiz de Austri, and P. Skands, Estimating QCD uncertainties in Monte Carlo event generators for gamma-ray dark matter searches, JCAP 05, 007, arXiv:1812.07424 [hep-ph] .
  * Bellm and Gellersen [2020] J. Bellm and L. Gellersen, High dimensional parameter tuning for event generators, Eur. Phys. J. C 80, 54 (2020), arXiv:1908.10811 [hep-ph] .
  * Jueid  _et al._ [2023] A. Jueid, J. Kip, R. R. de Austri, and P. Skands, The Strong Force meets the Dark Sector: a robust estimate of QCD uncertainties for anti-matter dark matter searches, (2023), arXiv:2303.11363 [hep-ph] .
  * Bowler [1981] M. G. Bowler, e+ e- Production of Heavy Quarks in the String Model, Z. Phys. C 11, 169 (1981).
  * Webber [2000] B. R. Webber, Fragmentation and hadronization, Int. J. Mod. Phys. A 15S1, 577 (2000), arXiv:hep-ph/9912292 .
  * Gieseke  _et al._ [2018] S. Gieseke, P. Kirchgaeßer, and S. Plätzer, Baryon production from cluster hadronisation, Eur. Phys. J. C 78, 99 (2018), arXiv:1710.10906 [hep-ph] .
  * Bellm  _et al._ [2019] J. Bellm, C. B. Duncan, S. Gieseke, M. Myska, and A. Siódmok, Spacetime colour reconnection in Herwig 7, Eur. Phys. J. C 79, 1003 (2019), arXiv:1909.08850 [hep-ph] .
  * Christiansen and Skands [2015] J. R. Christiansen and P. Z. Skands, String Formation Beyond Leading Colour, JHEP 08, 003, arXiv:1505.01681 [hep-ph] .
  * Oh  _et al._ [2009] Y. Oh, C. M. Ko, S. H. Lee, and S. Yasui, Heavy baryon/meson ratios in relativistic heavy ion collisions, Phys. Rev. C 79, 044905 (2009), arXiv:0901.1382 [nucl-th] .
  * Song  _et al._ [2018] J. Song, H.-h. Li, and F.-l. Shao, New feature of low $p_{T}$ charm quark hadronization in $pp$ collisions at $\sqrt{s}=7$ TeV, Eur. Phys. J. C 78, 344 (2018), arXiv:1801.09402 [hep-ph] .
  * Minissale  _et al._ [2021] V. Minissale, S. Plumari, and V. Greco, Charm hadrons in pp collisions at LHC energy within a coalescence plus fragmentation approach, Phys. Lett. B 821, 136622 (2021), arXiv:2012.12001 [hep-ph] .
  * Plumari  _et al._ [2018] S. Plumari, V. Minissale, S. K. Das, G. Coci, and V. Greco, Charmed Hadrons from Coalescence plus Fragmentation in relativistic nucleus-nucleus collisions at RHIC and LHC, Eur. Phys. J. C 78, 348 (2018), arXiv:1712.00730 [hep-ph] .
  * He and Rapp [2019] M. He and R. Rapp, Charm-Baryon Production in Proton-Proton Collisions, Phys. Lett. B 795, 117 (2019), arXiv:1902.08889 [nucl-th] .
  * He and Rapp [2020] M. He and R. Rapp, Hadronization and Charm-Hadron Ratios in Heavy-Ion Collisions, Phys. Rev. Lett. 124, 042301 (2020), arXiv:1905.09216 [nucl-th] .
  * Andronic  _et al._ [2021a] A. Andronic, P. Braun-Munzinger, M. K. Köhler, A. Mazeliauskas, K. Redlich, J. Stachel, and V. Vislavicius, The multiple-charm hierarchy in the statistical hadronization model, JHEP 07, 035, arXiv:2104.12754 [hep-ph] .
  * Andersson [2005] B. Andersson,  _The Lund model_ , Vol. 7 (Cambridge University Press, 2005).
  * Andersson  _et al._ [1985] B. Andersson, G. Gustafson, and T. Sjöstrand, Baryon Production in Jet Fragmentation and $\Upsilon$ Decay, Phys. Scripta 32, 574 (1985).
  * Eden and Gustafson [1997] P. Eden and G. Gustafson, Baryon production in the string fragmentation picture, Z. Phys. C 75, 41 (1997), arXiv:hep-ph/9606454 .
  * Adam  _et al._ [2017a] J. Adam _et al._ (ALICE), Enhanced production of multi-strange hadrons in high-multiplicity proton-proton collisions, Nature Phys. 13, 535 (2017a), arXiv:1606.07424 [nucl-ex] .
  * Acharya  _et al._ [2020a] S. Acharya _et al._ (ALICE), Multiplicity dependence of (multi-)strange hadron production in proton-proton collisions at $\sqrt{s}$ = 13 TeV, Eur. Phys. J. C 80, 167 (2020a), arXiv:1908.01861 [nucl-ex] .
  * Bierlich [2016] C. Bierlich, Hadronisation Models and Colour Reconnection, PoS DIS2016, 051 (2016), arXiv:1606.09456 [hep-ph] .
  * Bierlich [2018] C. Bierlich, Rope Hadronization and Strange Particle Production, EPJ Web Conf. 171, 14003 (2018), arXiv:1710.04464 [nucl-th] .
  * Fischer and Sjöstrand [2017] N. Fischer and T. Sjöstrand, Thermodynamical String Fragmentation, JHEP 01, 140, arXiv:1610.09818 [hep-ph] .
  * Andronic  _et al._ [2021b] A. Andronic, P. Braun-Munzinger, M. K. Köhler, A. Mazeliauskas, K. Redlich, J. Stachel, and V. Vislavicius, The multiple-charm hierarchy in the statistical hadronization model, JHEP 07, 035, arXiv:2104.12754 [hep-ph] .
  * Workman and Others [2022] R. L. Workman and Others (Particle Data Group), Review of Particle Physics, PTEP 2022, 083C01 (2022).
  * He and Rapp [2023] M. He and R. Rapp, Bottom Hadrochemistry in High-Energy Hadronic Collisions, Phys. Rev. Lett. 131, 012301 (2023), arXiv:2209.13419 [hep-ph] .
  * Ebert  _et al._ [2011] D. Ebert, R. Faustov, and V. Galkin, Spectroscopy and Regge trajectories of heavy baryons in the relativistic quark-diquark picture, Phys. Rev. D 84, 014025 (2011), arXiv:1105.0583 [hep-ph] .
  * Padmanath  _et al._ [2015] M. Padmanath, R. G. Edwards, N. Mathur, and M. J. Peardon, Spectroscopy of charmed baryons from lattice QCD, PoS LATTICE2014, 084 (2015), arXiv:1410.8791 [hep-lat] .
  * Andronic  _et al._ [2023] A. Andronic, P. Braun-Munzinger, H. Brunßen, J. Crkovská, J. Stachel, V. Vislavicius, and M. Völkl, Transverse dynamics of charmed hadrons in ultra-relativistic nuclear collisions, (2023), arXiv:2308.14821 [hep-ph] .
  * Capellino  _et al._ [2023] F. Capellino, A. Dubla, S. Floerchinger, E. Grossi, A. Kirchner, and S. Masciocchi, Fluid-dynamics of charm quarks in the quark–gluon plasma, (2023), arXiv:2307.14449 [hep-ph] .
  * Capellino  _et al._ [2022] F. Capellino, A. Beraudo, A. Dubla, S. Floerchinger, S. Masciocchi, J. Pawlowski, and I. Selyuzhenkov, Fluid-dynamic approach to heavy-quark diffusion in the quark-gluon plasma, Phys. Rev. D 106, 034021 (2022), arXiv:2205.07692 [nucl-th] .
  * Greco  _et al._ [2003a] V. Greco, C. M. Ko, and P. Levai, Parton coalescence at RHIC, Phys. Rev. C 68, 034904 (2003a), arXiv:nucl-th/0305024 .
  * Greco  _et al._ [2003b] V. Greco, C. M. Ko, and P. Levai, Parton coalescence and anti-proton / pion anomaly at RHIC, Phys. Rev. Lett. 90, 202302 (2003b), arXiv:nucl-th/0301093 .
  * Fries  _et al._ [2003a] R. J. Fries, B. Muller, C. Nonaka, and S. A. Bass, Hadron production in heavy ion collisions: Fragmentation and recombination from a dense parton phase, Phys. Rev. C 68, 044902 (2003a), arXiv:nucl-th/0306027 .
  * Fries  _et al._ [2003b] R. J. Fries, B. Muller, C. Nonaka, and S. A. Bass, Hadronization in heavy ion collisions: Recombination and fragmentation of partons, Phys. Rev. Lett. 90, 202303 (2003b), arXiv:nucl-th/0301087 .
  * Molnar and Voloshin [2003] D. Molnar and S. A. Voloshin, Elliptic flow at large transverse momenta from quark coalescence, Phys. Rev. Lett. 91, 092301 (2003), arXiv:nucl-th/0302014 .
  * Fries  _et al._ [2008] R. J. Fries, V. Greco, and P. Sorensen, Coalescence Models For Hadron Formation From Quark Gluon Plasma, Ann. Rev. Nucl. Part. Sci. 58, 177 (2008), arXiv:0807.4939 [nucl-th] .
  * Greco and Ko [2004] V. Greco and C. M. Ko, Effect of resonance decays on hadron elliptic flows, Phys. Rev. C 70, 024901 (2004), arXiv:nucl-th/0402020 .
  * Ravagli and Rapp [2007] L. Ravagli and R. Rapp, Quark Coalescence based on a Transport Equation, Phys. Lett. B 655, 126 (2007), arXiv:0705.0021 [hep-ph] .
  * Minissale  _et al._ [2015] V. Minissale, F. Scardina, and V. Greco, Hadrons from coalescence plus fragmentation in AA collisions at energies available at the BNL Relativistic Heavy Ion Collider to the CERN Large Hadron Collider, Phys. Rev. C 92, 054904 (2015), arXiv:1502.06213 [nucl-th] .
  * Greco  _et al._ [2004] V. Greco, C. M. Ko, and R. Rapp, Quark coalescence for charmed mesons in ultrarelativistic heavy ion collisions, Phys. Lett. B 595, 202 (2004), arXiv:nucl-th/0312100 [nucl-th] .
  * Das  _et al._ [2015] S. K. Das, F. Scardina, S. Plumari, and V. Greco, Toward a solution to the $R_{AA}$ and $v_{2}$ puzzle for heavy quarks, Phys. Lett. B 747, 260 (2015), arXiv:1502.03757 [nucl-th] .
  * Scardina  _et al._ [2017] F. Scardina, S. K. Das, V. Minissale, S. Plumari, and V. Greco, Estimating the charm quark diffusion coefficient and thermalization time from D meson spectra at energies available at the BNL Relativistic Heavy Ion Collider and the CERN Large Hadron Collider, Phys. Rev. C 96, 044905 (2017), arXiv:1707.05452 [nucl-th] .
  * Lin and Molnar [2003] Z.-w. Lin and D. Molnar, Quark coalescence and elliptic flow of charm hadrons, Phys. Rev. C 68, 044901 (2003), arXiv:nucl-th/0304045 .
  * Lisovyi  _et al._ [2016] M. Lisovyi, A. Verbytskyi, and O. Zenaiev, Combined analysis of charm-quark fragmentation-fraction measurements, Eur. Phys. J. C 76, 397 (2016), arXiv:1509.01061 [hep-ex] .
  * Adam  _et al._ [2020] J. Adam _et al._ (STAR), First measurement of $\Lambda_{c}$ baryon production in Au+Au collisions at $\sqrt{s_{\mathrm{N}N}}$ = 200 GeV, Phys. Rev. Lett. 124, 172301 (2020), arXiv:1910.14628 [nucl-ex] .
  * Acharya  _et al._ [2019a] S. Acharya _et al._ (ALICE), $\Lambda_{\mathrm{c}}^{+}$ production in Pb-Pb collisions at $\sqrt{s_{\mathrm{N}N}}=5.02$ TeV, Phys. Lett. B 793, 212 (2019a), arXiv:1809.10922 [nucl-ex] .
  * Acharya  _et al._ [2023a] S. Acharya _et al._ (ALICE), Constraining hadronization mechanisms with $\Lambda_{\mathrm{c}}^{+}/\mathrm{D}^{0}$ production ratios in Pb–Pb collisions at $\sqrt{s_{\mathrm{N}N}}$ = 5.02 TeV, Phys. Lett. B 839, 137796 (2023a), arXiv:2112.08156 [nucl-ex] .
  * Tumasyan  _et al._ [2024] A. Tumasyan _et al._ (CMS), Study of charm hadronization with prompt ${\Lambda}_{\textrm{c}}^{+}$ baryons in proton-proton and lead-lead collisions at $\sqrt{s_{\textrm{NN}}}$ = 5.02 TeV, JHEP 01, 128, arXiv:2307.11186 [nucl-ex] .
  * Minissale  _et al._ [2024] V. Minissale, S. Plumari, Y. Sun, and V. Greco, Multi-charmed and singled charmed hadrons from coalescence: yields and ratios in different collision systems at LHC, Eur. Phys. J. C 84, 228 (2024), arXiv:2305.03687 [hep-ph] .
  * ALI [2022b] Letter of intent for ALICE 3: A next-generation heavy-ion experiment at the LHC, (2022b), arXiv:2211.02491 [physics.ins-det] .
  * Khachatryan  _et al._ [2016a] V. Khachatryan _et al._ (CMS), Measurement of long-range near-side two-particle angular correlations in pp collisions at $\sqrt{s}=$13 TeV, Phys. Rev. Lett. 116, 172302 (2016a), arXiv:1510.03068 [nucl-ex] .
  * Acharya  _et al._ [2021a] S. Acharya _et al._ (ALICE), $\Lambda^{+}_{c}$ Production and Baryon-to-Meson Ratios in pp and p-Pb Collisions at $\sqrt{s_{NN}}$=5.02 TeV at the LHC, Phys. Rev. Lett. 127, 202301 (2021a), arXiv:2011.06078 [nucl-ex] .
  * Acharya  _et al._ [2021b] S. Acharya _et al._ (ALICE), $\Lambda^{+}_{c}$ production in $pp$ and in $p$-Pb collisions at $\sqrt{s_{NN}}$=5.02 TeV, Phys. Rev. C 104, 054905 (2021b), arXiv:2011.06079 [nucl-ex] .
  * Weller and Romatschke [2017] R. D. Weller and P. Romatschke, One fluid to rule them all: viscous hydrodynamic description of event-by-event central p+p, p+Pb and Pb+Pb collisions at $\sqrt{s}=5.02$ TeV, Phys. Lett. B 774, 351 (2017), arXiv:1701.07145 [nucl-th] .
  * Greif  _et al._ [2017] M. Greif, C. Greiner, B. Schenke, S. Schlichting, and Z. Xu, Importance of initial and final state effects for azimuthal correlations in p+Pb collisions, Phys. Rev. D 96, 091504 (2017), arXiv:1708.02076 [hep-ph] .
  * Sun  _et al._ [2020] Y. Sun, S. Plumari, and V. Greco, Study of collective anisotropies $v_{2}$ and $v_{3}$ and their fluctuations in $pA$ collisions at LHC within a relativistic transport approach, Eur. Phys. J. C 80, 16 (2020), arXiv:1907.11287 [nucl-th] .
  * Zhao  _et al._ [2020] W. Zhao, C. M. Ko, Y.-X. Liu, G.-Y. Qin, and H. Song, Probing the Partonic Degrees of Freedom in High-Multiplicity $p-Pb$ collisions at $\sqrt{s_{NN}}$ = 5.02 TeV, Phys. Rev. Lett. 125, 072301 (2020), arXiv:1911.00826 [nucl-th] .
  * Retiere and Lisa [2004] F. Retiere and M. A. Lisa, Observable implications of geometrical and dynamical aspects of freeze out in heavy ion collisions, Phys. Rev. C 70, 044907 (2004), arXiv:nucl-th/0312024 .
  * Sambataro  _et al._ [2022] M. L. Sambataro, Y. Sun, V. Minissale, S. Plumari, and V. Greco, Event-shape engineering analysis of D meson in ultrarelativistic heavy-ion collisions, Eur. Phys. J. C 82, 833 (2022), arXiv:2206.03160 [hep-ph] .
  * Sambataro  _et al._ [2024a] M. L. Sambataro, V. Minissale, S. Plumari, and V. Greco, B meson production in Pb+Pb at 5.02 ATeV at LHC: Estimating the diffusion coefficient in the infinite mass limit, Phys. Lett. B 849, 138480 (2024a), arXiv:2304.02953 [hep-ph] .
  * Hwang [2002] C.-W. Hwang, Charge radii of light and heavy mesons, Eur. Phys. J. C 23, 585 (2002), arXiv:hep-ph/0112237 .
  * Albertus  _et al._ [2004] C. Albertus, J. E. Amaro, E. Hernandez, and J. Nieves, Charmed and bottom baryons: A Variational approach based on heavy quark symmetry, Nucl. Phys. A 740, 333 (2004), arXiv:nucl-th/0311100 .
  * Cao  _et al._ [2020] S. Cao, K.-J. Sun, S.-Q. Li, S. Y. F. Liu, W.-J. Xing, G.-Y. Qin, and C. M. Ko, Charmed hadron chemistry in relativistic heavy-ion collisions, Phys. Lett. B 807, 135561 (2020), arXiv:1911.00456 [nucl-th] .
  * Cassing and Bratkovskaya [2008] W. Cassing and E. L. Bratkovskaya, Parton transport and hadronization from the dynamical quasiparticle point of view, Phys. Rev. C 78, 034919 (2008), arXiv:0808.0022 [hep-ph] .
  * Song  _et al._ [2015] T. Song, H. Berrehrah, D. Cabrera, J. M. Torres-Rincon, L. Tolos, W. Cassing, and E. Bratkovskaya, Tomography of the Quark-Gluon-Plasma by Charm Quarks, Phys. Rev. C 92, 014910 (2015), arXiv:1503.03039 [nucl-th] .
  * Song  _et al._ [2023] T. Song, J. Aichelin, J. Zhao, P. B. Gossiaux, and E. Bratkovskaya, Bottomonium production in pp and heavy-ion collisions, Phys. Rev. C 108, 054908 (2023), arXiv:2305.10750 [nucl-th] .
  * Li  _et al._ [2018] H.-H. Li, F.-L. Shao, J. Song, and R.-Q. Wang, Production of single-charm hadrons by quark combination mechanism in $p$-Pb collisions at $\sqrt{s_{NN}}=5.02$ TeV, Phys. Rev. C 97, 064915 (2018), arXiv:1712.08921 [hep-ph] .
  * Rapp and Shuryak [2003] R. Rapp and E. V. Shuryak, D meson production from recombination in hadronic collisions, Phys. Rev. D 67, 074036 (2003), arXiv:hep-ph/0301245 .
  * Zhao  _et al._ [2023] J. Zhao _et al._ , Hadronization of Heavy Quarks, (2023), arXiv:2311.10621 [hep-ph] .
  * Alberico  _et al._ [2013] W. Alberico, A. Beraudo, A. De Pace, A. Molinari, M. Monteno, _et al._ , Heavy flavors in $AA$ collisions: production, transport and final spectra, Eur.Phys.J. C73, 2481 (2013), arXiv:1305.7421 [hep-ph] .
  * Beraudo  _et al._ [2022] A. Beraudo, A. De Pace, M. Monteno, M. Nardi, and F. Prino, In-medium hadronization of heavy quarks and its effect on charmed meson and baryon distributions in heavy-ion collisions, Eur. Phys. J. C 82, 607 (2022), arXiv:2202.08732 [hep-ph] .
  * Beraudo  _et al._ [2024] A. Beraudo, A. De Pace, D. Pablos, F. Prino, M. Monteno, and M. Nardi, Heavy-flavor transport and hadronization in pp collisions, Phys. Rev. D 109, L011501 (2024), arXiv:2306.02152 [hep-ph] .
  * Pierog  _et al._ [2015] T. Pierog, I. Karpenko, J. M. Katzy, E. Yatsenko, and K. Werner, EPOS LHC: Test of collective hadronization with data measured at the CERN Large Hadron Collider, Phys. Rev. C 92, 034906 (2015), arXiv:1306.0121 [hep-ph] .
  * Werner [2024] K. Werner, Core-corona procedure and microcanonical hadronization to understand strangeness enhancement in proton-proton and heavy ion collisions in the EPOS4 framework, Phys. Rev. C 109, 014910 (2024), arXiv:2306.10277 [hep-ph] .
  * Werner [2007] K. Werner, Core-corona separation in ultra-relativistic heavy ion collisions, Phys. Rev. Lett. 98, 152301 (2007), arXiv:0704.1270 [nucl-th] .
  * Zhao  _et al._ [2024a] J. Zhao, J. Aichelin, P. B. Gossiaux, and K. Werner, Heavy flavor as a probe of hot QCD matter produced in proton-proton collisions, Phys. Rev. D 109, 054011 (2024a), arXiv:2310.08684 [hep-ph] .
  * Zhao  _et al._ [2024b] J. Zhao, J. Aichelin, P. B. Gossiaux, V. Ozvenchuk, and K. Werner, Heavy flavour hadron production in relativistic heavy ion collisions at RHIC and LHC in EPOS4HQ, (2024b), arXiv:2401.17096 [hep-ph] .
  * Niiyama  _et al._ [2018] M. Niiyama _et al._ (Belle), Production cross sections of hyperons and charmed baryons from $e^{+}e^{-}$ annihilation near $\sqrt{s}=10.52$ GeV, Phys. Rev. D 97, 072005 (2018), arXiv:1706.06791 [hep-ex] .
  * Dokshitzer  _et al._ [1991] Y. L. Dokshitzer, V. A. Khoze, A. H. Mueller, and S. I. Troian,  _Basics of perturbative QCD_ (1991).
  * Barate  _et al._ [2000] R. Barate _et al._ (ALEPH), Study of charm production in Z decays, Eur. Phys. J. C 16, 597 (2000), arXiv:hep-ex/9909032 .
  * Alexander  _et al._ [1996] G. Alexander _et al._ (OPAL), A Study of charm hadron production in Z0 —$>$ c anti-c and Z0 —$>$ b anti-b decays at LEP, Z. Phys. C 72, 1 (1996).
  * Abreu  _et al._ [2000a] P. Abreu _et al._ (DELPHI), Measurements of the Z partial decay width into c anti-c and multiplicity of charm quarks per b decay, Eur. Phys. J. C 12, 225 (2000a).
  * Gladilin [2015] L. Gladilin, Fragmentation fractions of $c$ and $b$ quarks into charmed hadrons at LEP, Eur. Phys. J. C 75, 19 (2015), arXiv:1404.3888 [hep-ex] .
  * Sjöstrand [1994] T. Sjöstrand, High-energy physics event generation with PYTHIA 5.7 and JETSET 7.4, Comput. Phys. Commun. 82, 74 (1994).
  * Bortoletto  _et al._ [1988] D. Bortoletto _et al._ (CLEO), Charm Production in Nonresonant $\mathrm{e}^{+}e^{-}$ Annihilations at $\sqrt{s}$ = 10.55 GeV, Phys. Rev. D 37, 1719 (1988), [Erratum: Phys.Rev.D 39, 1471 (1989)].
  * Avery  _et al._ [1991] P. Avery _et al._ (CLEO), Inclusive production of the charmed baryon Lambda(c) from $\mathrm{e}^{+}e^{-}$ annihilations at $\sqrt{s}$ = 10.55-GeV, Phys. Rev. D 43, 3599 (1991).
  * Albrecht  _et al._ [1991] H. Albrecht _et al._ (ARGUS), Inclusive production of D0, D+ and D*+ (2010) mesons in B decays and nonresonant e+ e- annihilation at 10.6-GeV, Z. Phys. C 52, 353 (1991).
  * Albrecht  _et al._ [1992] H. Albrecht _et al._ (ARGUS), Production of D(s)+ mesons in B decays and determination of f(D(s)), Z. Phys. C 54, 1 (1992).
  * Albrecht  _et al._ [1988] H. Albrecht _et al._ (ARGUS), Observation of the Charmed Baryon $\Lambda(c$) in $e^{+}e^{-}$ Annihilation at 10-GeV, Phys. Lett. B 207, 109 (1988).